Commit 399108c7 authored by Mario Chirinos Colunga's avatar Mario Chirinos Colunga

Merge branch 'dev' of gitlab.geoint.mx:adan.salazar/GeoSentinel

parents 8ce9c577 b6d7b956
{
"PATHS": {
"PATH_GEOSENTINEL": "/home/david/centroGEO/repsat/GeoSentinel",
"PATH_NAS": "/home/david/NAS/"
},
"API_SENTINEL": {
"SENTINEL_USER" : "emmhp",
"SENTINEL_PASS" : "geoemm29"
}
}
...@@ -94,6 +94,7 @@ class APISentinel(object): ...@@ -94,6 +94,7 @@ class APISentinel(object):
os.chdir(dir) os.chdir(dir)
# self.api.download_all(products) # self.api.download_all(products)
for p in products: for p in products:
print(products[p]['filename'])
self.api.download(p) self.api.download(p)
def filterProducts(self, productList): def filterProducts(self, productList):
......
...@@ -5,32 +5,45 @@ import osgeo.osr as osr ...@@ -5,32 +5,45 @@ import osgeo.osr as osr
import numpy as np import numpy as np
import gdal import gdal
import sys import sys
import os
import json import json
from collections import Counter from collections import Counter
sys.path.append('/home/mario/git/GeoSentinel')
########### lee archivo de configuración ################
dirname = os.path.dirname(__file__)
configfile = os.path.join(dirname, '../config/config.json')
with open(configfile, 'r') as f:
config = json.load(f)
#print(config['PATHS']['PATH_GEOSENTINEL'])
SENTINEL_PATH = config['PATHS']['PATH_GEOSENTINEL']
###########################################################
sys.path.append(SENTINEL_PATH)
from geosentinel import rasterWkt from geosentinel import rasterWkt
CLASS_MAP = {0:0, 1:1, 2:2, 3:3, 4:4, 5:5, 6:6, 7:7, 8:8, 9:8, 10:8, 11:9}
#------------------------------------------------------------------------------- #-------------------------------------------------------------------------------
def L2ASCLtoDict(filename, wkt): def L2ASCLtoDict(filename, wkt):
''' '''
L2ASCLtoDict L2ASCLtoDict
''' '''
classMap = lambda x: [ CLASS_MAP[xx] for xx in x ]
inputImage = gdal.Open(filename) inputImage = gdal.Open(filename)
rows, cols, geotransform = inputImage.RasterYSize, inputImage.RasterXSize, inputImage.GetGeoTransform() rows, cols, geotransform = inputImage.RasterYSize, inputImage.RasterXSize, inputImage.GetGeoTransform()
data = inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows) data = inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows)
# count = Counter([])
# for r in data:
# count = count + Counter(r)
# sclDict = {str(k):v for k,v in dict(count).items()}
# print ("data")
# print (sclDict)
data = np.bitwise_and(data, rasterWkt.getPolygonArray(inputImage, wkt)) data = np.bitwise_and(data, rasterWkt.getPolygonArray(inputImage, wkt))
count = Counter([]) count = Counter([])
for r in data: for r in data:
count = count + Counter(r) count = count + Counter(classMap(r))
sclDict = {str(k):v for k,v in dict(count).items()} sclDict = {str(k):v for k,v in dict(count).items()}
# print ("data & POLY")
# print (sclDict)
return sclDict return sclDict
#------------------------------------------------------------------------------- #-------------------------------------------------------------------------------
def main(argv): def main(argv):
......
...@@ -4,20 +4,32 @@ ...@@ -4,20 +4,32 @@
import sys, os import sys, os
import json import json
sys.path.append('/home/mario/git/GeoSentinel') ########### lee archivo de configuración ################
dirname = os.path.dirname(__file__)
configfile = os.path.join(dirname, '../config/config.json')
with open(configfile, 'r') as f:
config = json.load(f)
#print(config['PATHS']['PATH_GEOSENTINEL'])
SENTINEL_PATH = config['PATHS']['PATH_GEOSENTINEL']
###########################################################
sys.path.append(SENTINEL_PATH)
from geosentinel import APISentinel from geosentinel import APISentinel
from geosentinel import polygonToBox from geosentinel import polygonToBox
from datetime import date from datetime import date
#def findSentinelProducts(wkt, startDate, endDate, platform, cloud): #def findSentinelProducts(wkt, startDate, endDate, platform, cloud):
# sentinel = APISentinel.APISentinel('asalazarg', 'geo135asg') # sentinel = APISentinel.APISentinel('asalazarg', 'geo135asg')
# products = sentinel.getProducts(wkt, (startDate, endDate), {"platformname":platform, "cloudcoverpercentage":"[0 TO "+str(cloud)+"]"}) # products = sentinel.getProducts(wkt, (startDate, endDate), {"platformname":platform, "cloudcoverpercentage":"[0 TO "+str(cloud)+"]"})
# return products # return products
#productsCodes = {"L1C":"L1C","L2A":"L1C"}
productType = {"L1C":"S2MSI1C", "L2A":"S2MSI2A"}
def main(argv): def main(argv):
if len(sys.argv) != 2: if len(sys.argv) != 2:
...@@ -29,28 +41,72 @@ def main(argv): ...@@ -29,28 +41,72 @@ def main(argv):
sentinel = APISentinel.APISentinel(cfg['username'], cfg['password']) sentinel = APISentinel.APISentinel(cfg['username'], cfg['password'])
wkt = polygonToBox.getWKTPolygonBoundingBox(cfg['wkt'], True) wkt = polygonToBox.getWKTPolygonBoundingBox(cfg['wkt'], True)
productList = sentinel.getProducts(wkt, (cfg['startDate'], cfg['endDate']), {"platformname":cfg['platform'], "cloudcoverpercentage":"[0 TO "+str(cfg['clouds'])+"]"}, "Contains") productList = sentinel.getProducts(wkt, (cfg['startDate'], cfg['endDate']), {"producttype":productType[cfg["productLevel"]], "platformname":cfg['platform'], "cloudcoverpercentage":"[0 TO "+str(cfg['clouds'])+"]"}, "Contains")
if len(productList)<=0: if len(productList)<=0:
print ("No products found with 'Contains' trying with 'Intersects'") print ("No products found with 'Contains' trying with 'Intersects'")
productList = sentinel.getProducts(wkt, (cfg['startDate'], cfg['endDate']), {"platformname":cfg['platform'], "cloudcoverpercentage":"[0 TO "+str(cfg['clouds'])+"]"}) productList = sentinel.getProducts(wkt, (cfg['startDate'], cfg['endDate']), {"producttype":productType[cfg["productLevel"]], "platformname":cfg['platform'], "cloudcoverpercentage":"[0 TO "+str(cfg['clouds'])+"]"})
productsCodes = {"L1C":"L1C","L2A":"L1C"}
fileNames = [productList[k]['filename'].replace("SAFE", "zip").replace(productsCodes[cfg['productLevel']], cfg['productLevel']) for k in productList.keys() ] # fileNames = [productList[k]['filename'].replace("SAFE", "zip").replace(productsCodes[cfg['productLevel']], cfg['productLevel']) for k in productList.keys() ]
inDir = cfg["productsDir"] fileNames = [productList[k]['filename'].replace("SAFE", "zip") for k in productList.keys() ]
outdir = cfg["linksDir"]
dirList = os.listdir(inDir) downloadDir = cfg["productsDir"]+cfg["productLevel"]+"/"
dirList = os.listdir(downloadDir)
matchingProducts = set(fileNames).intersection(set(dirList)) matchingProducts = set(fileNames).intersection(set(dirList))
print ( str(len(matchingProducts))+" of " + str(len(fileNames)) +" products found.") print ( str(len(matchingProducts))+" of " + str(len(fileNames)) +" "+ cfg["productLevel"] + " products found.")
# rawDir = inDir.replace(cfg['productLevel'],productsCodes[cfg['productLevel']])
# print(rawDir)
text =""
while text != "yes" and text != "no":
text = raw_input("Do you want to download this products to "+downloadDir+ " ? (yes, no)")
if text=="yes":
sentinel.downloadProducts(productList,downloadDir)
L2ADir = cfg["productsDir"]+"L2A/"
if cfg["productLevel"]=="L1C":
fileNamesL2A = [productList[k]['filename'].replace("SAFE", "zip").replace("L1C", "L2A") for k in productList.keys() ]
L1CDir = cfg["productsDir"]+"L1C/"
linksDir=cfg["projectDir"]+"L1C/"
if not os.path.exists(linksDir):
os.mkdir(linksDir)
L1CDirList = os.listdir(L1CDir)
#LinkProducts
fileNames = [productList[k]['filename'].replace("SAFE", "zip") for k in productList.keys() ]
matchingProducts = set(fileNames).intersection(set(L1CDirList))
for f in matchingProducts:
print (f)
print("ln -s " + L1CDir+f + " " + linksDir+f)
os.system("ln -s " + L1CDir+f + " " + linksDir+f)
print (str(len(matchingProducts)) + " Linked to " + linksDir)
print ( str(len(matchingProducts))+" of " + str(len(fileNames)) +" L1C products found in L2A.")
text ="" text =""
while text != "yes" and text != "no": while text != "yes" and text != "no":
text = raw_input("Do you want to link this products to "+outdir+ " ? (yes, no)") text = raw_input("Do you want to convert this products and save them to "+downloadDir+ " ? (yes, no)")
if text=="yes":
#print("L1CProductListToL2A.sh "+rawDir+" "+inDir+" 1")
os.system("L1CProductListToL2A.sh "+linksDir+" "+L2ADir+" 1")
text =""
linksDir=cfg["projectDir"]+"L2A/"
while text != "yes" and text != "no":
text = raw_input("Do you want to link this products to "+linksDir+ " ? (yes, no)")
if text=="yes": if text=="yes":
fileNames = [productList[k]['filename'].replace("SAFE", "zip").replace("L1C", "L2A") for k in productList.keys() ]
dirList = os.listdir(L2ADir)
matchingProducts = set(fileNames).intersection(set(dirList))
for f in matchingProducts: for f in matchingProducts:
print (f) print (f)
os.system("ln -s " + inDir+f + " " + outdir+f) os.system("ln -s " + L2ADir+f + " " + linksDir+f)
print (str(len(matchingProducts)) + " Linked to " + outdir) print (str(len(matchingProducts)) + " Linked to " + linksDir)
if __name__ == "__main__": if __name__ == "__main__":
main(sys.argv) main(sys.argv)
......
...@@ -5,25 +5,49 @@ import osgeo.osr as osr ...@@ -5,25 +5,49 @@ import osgeo.osr as osr
import numpy as np import numpy as np
import gdal import gdal
import sys, os import sys, os
import json
########### lee archivo de configuración ################
dirname = os.path.dirname(__file__)
configfile = os.path.join(dirname, '../config/config.json')
with open(configfile, 'r') as f:
config = json.load(f)
#print(config['PATHS']['PATH_GEOSENTINEL'])
SENTINEL_PATH = config['PATHS']['PATH_GEOSENTINEL']
###########################################################
sys.path.append(SENTINEL_PATH)
from geosentinel import rasterWkt
#------------------------------------------------------------------------------- #-------------------------------------------------------------------------------
def imageMissingData(filename): def imageMissingData(filename, wkt):
inputImage = gdal.Open(filename) inputImage = gdal.Open(filename)
rows, cols, geotransform = inputImage.RasterYSize, inputImage.RasterXSize, inputImage.GetGeoTransform() rows, cols, geotransform = inputImage.RasterYSize, inputImage.RasterXSize, inputImage.GetGeoTransform()
data = inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows) data = inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows)
inputImage = None rasterwkt = rasterWkt.getPolygonArray(inputImage, wkt)
data = np.bitwise_and(data, rasterwkt)
return (data>0).sum()/float((rows*cols)) inputImage = None
return (data>0).sum()/float((rasterwkt>0).sum())
#------------------------------------------------------------------------------- #-------------------------------------------------------------------------------
def main(argv): def main(argv):
if len(sys.argv) != 3: if len(sys.argv) != 4:
print("Usage: " + argv[0] + " <File> <%>") print("Usage: " + argv[0] + " <File> <JSON> <%>")
else: else:
dataArea = imageMissingData(argv[1]) jsonFile=open(argv[2]).read()
if dataArea<float(argv[2])/100.0 : cfg = json.loads(jsonFile)
print ("Deleting " + argv[1] + " " + str(dataArea))
dataArea = imageMissingData(argv[1], cfg["wkt"])
th = float(argv[3])/100.0
if dataArea<th :
print ("Deleting " + argv[1] + " " + str(dataArea)+" < " + str(th))
os.system("rm " + argv[1]) os.system("rm " + argv[1])
else: else:
print ("Keeping " + argv[1] + " " + str(dataArea)) print ("Keeping " + argv[1] + " " + str(dataArea))
......
#!/usr/bin/python
# -*- coding: utf-8 -*-
import osgeo.ogr as ogr
import osgeo.osr as osr
import numpy as np
import gdal
import sys
SCL_COLOR = [
(0,0,0,255), #NO_DATA
(255,0,0,255), #SATURATED_OR_DEFECTIVE
(64,64,64,255), #DARK_AREA_PIXELS
(102,51,0,255), #CLOUD_SHADOWS
(0,255,0,255), #VEGETATION
(255,255,0,255), #NOT_VEGETATED
(0,0,255,255), #WATER
(128,128,128,255), #UNCLASSIFIED
(192,192,192,255), #CLOUD_MEDIUM_PROBABILITY
(255,255,255,255), #CLOUD_HIGH_PROBABILITY
(102,204,255,255), #THIN_CIRRUS
(255,153,255,255), #SNOW
]
CLASS_MAP = {0:0, 1:1, 2:2, 3:3, 4:4, 5:5, 6:6, 7:7, 8:8, 9:8, 10:8, 11:9}
CLASS_GROUP = [[0], [1], [2], [3], [4], [5], [6], [7], [8,9,10], [11]]
#-------------------------------------------------------------------------------
def splitSCL(filename, outdir):
inputImage = gdal.Open(filename)
rows = inputImage.RasterYSize
cols = inputImage.RasterXSize
geotransform = inputImage.GetGeoTransform()
projection = inputImage.GetProjection()
sr = osr.SpatialReference()
sr.ImportFromWkt(projection)
driverTiff = gdal.GetDriverByName('GTiff')
# for c in range(len(SCL_COLOR)):
for c in range(len(CLASS_GROUP)):
outputfile=outdir+filename[:-4]+"_"+str(c)+".tif"
output = driverTiff.Create(outputfile, cols, rows, 4, gdal.GDT_Byte)
# classMask = inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows)==c
classMask = np.isin(inputImage.GetRasterBand(1).ReadAsArray(0,0,cols,rows), CLASS_GROUP[c])
for b in range(1, output.RasterCount+1):
output.GetRasterBand(b).SetNoDataValue(0)
output.GetRasterBand(b).WriteArray(classMask*SCL_COLOR[c][b-1])
output.GetRasterBand(b).FlushCache()
output.SetGeoTransform(geotransform)
output.SetProjection(projection)
output = None
inputImage = None
#-------------------------------------------------------------------------------
def main(argv):
splitSCL(argv[1], argv[2])
if __name__ == "__main__":
main(sys.argv)
...@@ -29,6 +29,8 @@ fi ...@@ -29,6 +29,8 @@ fi
if [ ! -e $filePrefix"TCI_10m.jp2" -o ! -e $filePrefix"TCI_60m.jp2" -o ! -e $filePrefix"SCL_60m.jp2" ]; then if [ ! -e $filePrefix"TCI_10m.jp2" -o ! -e $filePrefix"TCI_60m.jp2" -o ! -e $filePrefix"SCL_60m.jp2" ]; then
unzip -n -j $FILE *TCI_10m.jp2 *TCI_60m.jp2 *SCL_60m.jp2 -d $outDir unzip -n -j $FILE *TCI_10m.jp2 *TCI_60m.jp2 *SCL_60m.jp2 -d $outDir
# unzip -n -j $FILE *TCI_10m.jp2 *SCL_60m.jp2 -d $outDir
TCI10m=$outDir"/L2A_"$filepattern"TCI_10m.jp2" TCI10m=$outDir"/L2A_"$filepattern"TCI_10m.jp2"
TCI60m=$outDir"/L2A_"$filepattern"TCI_60m.jp2" TCI60m=$outDir"/L2A_"$filepattern"TCI_60m.jp2"
SCL60m=$outDir"/L2A_"$filepattern"SCL_60m.jp2" SCL60m=$outDir"/L2A_"$filepattern"SCL_60m.jp2"
......
...@@ -5,10 +5,10 @@ filename=$1 #L1C Filename /home/geoint/NAS/meridaL1C_2c/S2B_MSIL1C_20181005T1631 ...@@ -5,10 +5,10 @@ filename=$1 #L1C Filename /home/geoint/NAS/meridaL1C_2c/S2B_MSIL1C_20181005T1631
basename=${filename##*/} #S2B_MSIL1C_20181005T163129_N0206_R083_T16QBH_20181005T214055.zip basename=${filename##*/} #S2B_MSIL1C_20181005T163129_N0206_R083_T16QBH_20181005T214055.zip
#echo $filename #echo $filename
#echo $basename #echo $basename
Z1=$(echo $(echo $basename | cut -d "_" -f1) | cut -d "/" -f3) Z1=$(echo $(echo $basename | cut -d "_" -f1) | cut -d "/" -f3) #S2B
Z2=$(echo $basename | cut -d "_" -f3,4,5,6,7 | cut -d "." -f1) Z2=$(echo $basename | cut -d "_" -f3,4,5,6,7 | cut -d "." -f1) #20181005T163129_N0206_R083_T16QBH_20181005T214055
L2A=$Z1"_MSIL2A_"$Z2 L2A=$Z1"_MSIL2A_"$Z2 #S2B_MSIL2A_20181005T163129_N0206_R083_T16QBH_20181005T214055
#echo $L2A #echo $L2A
if [ ! -e $L2A.zip ]; then if [ ! -e $L2A.zip ]; then
echo "Processing: " $filename " into " $L2A.zip echo "Processing: " $filename " into " $L2A.zip
......
...@@ -2,10 +2,10 @@ ...@@ -2,10 +2,10 @@
PRODCUTSDIR=$1 #Products Intput directory PRODCUTSDIR=$1 #Products Intput directory
JP2DIR=$2 #JP2 Output Directory JP2DIR=$2 #JP2 Output Directory
JOBS=${3:-1} JOBS=${3:-1}
BYTILE=${4:-0}
cd $PRODCUTSDIR cd $PRODCUTSDIR
ls *.zip | parallel --jobs $JOBS ExtractData.sh {} $JP2DIR 1 ls *.zip | parallel --jobs $JOBS ExtractData.sh {} $JP2DIR $BYTILE
......
#!/bin/sh
#Convert L2A Products to L3A using L2AtoL3C Process
DATE=$1
TILE=$2
INDIR=$3
OUTDIR=$4
cd $INDIR
for f in $(ls | egrep $DATE.{21}$TILE)
echo $f
done
...@@ -8,6 +8,8 @@ BYTILE=${2:-0} ...@@ -8,6 +8,8 @@ BYTILE=${2:-0}
# exit # exit
#fi #fi
echo "Starting Process..."
echo $USERDIR"findProducts.json"
wkt=$(getJSONparameter.py $USERDIR"findProducts.json" wkt) wkt=$(getJSONparameter.py $USERDIR"findProducts.json" wkt)
RED='\033[0;31m' RED='\033[0;31m'
NC='\033[0m' # No Color NC='\033[0m' # No Color
...@@ -17,27 +19,39 @@ echo $wkt ...@@ -17,27 +19,39 @@ echo $wkt
echo $USERDIR echo $USERDIR
cd $USERDIR cd $USERDIR
#0.- Create Shape File From WKT
echo "${RED}Creating Shape file...${NC}"
rm -r myshape
wktToShape.py "$wkt" "myshape"
#1.- Link L2A products #1.- Link L2A products
echo "Linking Products..." echo ${RED}"Linking Products..."${NC}
if [ ! -d "L2A" ]; then if [ ! -d "L2A" ]; then
mkdir L2A mkdir L2A
fi fi
cd L2A cd L2A
rm * rm *
cd .. cd ..
echo yes | findProducts.py $USERDIR"findProducts.json" yes yes | findProducts.py $USERDIR"findProducts.json"
#2.- Extract Images #2.- Extract Images
echo "Extracting JP2 Images..." echo ${RED}"Extracting JP2 Images..."${NC}
if [ ! -d "jp2" ]; then if [ ! -d "jp2" ]; then
mkdir jp2 mkdir jp2
fi fi
L2AProductListExtractData.sh $USERDIR"L2A/" $USERDIR"jp2/" 4 if [ $BYTILE -ne "0" ]; then
echo "BY TILE"
L2AProductListExtractData.sh $USERDIR"L2A/" $USERDIR"jp2/" 4 1
else
echo "BY DATE"
L2AProductListExtractData.sh $USERDIR"L2A/" $USERDIR"jp2/" 4 0
fi
#3.- Merge Images #3.- Merge Images
echo "Merging Images..." echo ${RED}"\nMerging Images..."${NC}
if [ $BYTILE -ne "0" ]; then if [ $BYTILE -ne "0" ]; then
echo "BY TILE" echo "BY TILE"
cd $USERDIR"jp2/" cd $USERDIR"jp2/"
...@@ -49,51 +63,92 @@ if [ $BYTILE -ne "0" ]; then ...@@ -49,51 +63,92 @@ if [ $BYTILE -ne "0" ]; then
if [ ! -d "../out/"$i ]; then if [ ! -d "../out/"$i ]; then
mkdir "../out/"$i mkdir "../out/"$i
fi fi
"mergeImagesByDirectory.sh" $USERDIR"jp2/"$i ../../../out/$i "$wkt" mergeImagesByDirectory.sh $USERDIR"jp2/"$i ../../../out/$i "$wkt"
done done
else else
echo "BY DATE" echo "BY DATE"
# mergeImagesByDirectory.sh $USERDIR"jp2/" "$wtk" 4 cd $USERDIR"jp2/"
if [ ! -d "../out" ]; then
mkdir "../out"
fi
mergeImagesByDirectory.sh $USERDIR"jp2/" $USERDIR"out/" "$wkt" 4
fi fi
#4.- Mask raster with shape file
echo ${RED}"\nApplying Mask..."${NC}
#4.- Delete images with few data
echo "Deleting Images..."
if [ $BYTILE -ne "0" ]; then if [ $BYTILE -ne "0" ]; then
echo "BY TILE" echo "BY TILE"
cd $USERDIR"jp2/" cd $USERDIR"out/"
for d in $(ls -d */) ; do
cd $d
if [ ! -d "mask" ]; then
mkdir mask
fi
for i in $(ls *.tif) ; do
if [ ! -e mask/$i ]; then
gdalwarp -cutline $USERDIR"myshape/wkt.shp" -crop_to_cutline -dstalpha $i mask/$i
else
echo "Passing " $i
fi
done
cd ..
done
else
echo "BY DATE"
cd $USERDIR"out/"
if [ ! -d "mask" ]; then
mkdir mask
fi
for i in $(ls *.tif) ; do
if [ ! -e mask/$i ]; then
gdalwarp -cutline ../myshape/wkt.shp -crop_to_cutline -dstalpha $i mask/$i
else
echo "Passing " $i
fi
done
fi
#5.- Delete images with few data
echo ${RED}"\nDeleting Images..."${NC}
if [ $BYTILE -ne "0" ]; then
echo "BY TILE"
cd $USERDIR"out/"
for i in $(ls -d */) ; do for i in $(ls -d */) ; do
cd ../out/$i cd $i"mask"
echo $i echo $i
if [ ! -d "SCL" ]; then if [ ! -d "SCL" ]; then
mkdir SCL mkdir SCL
fi fi
cp *SCL_60m.tif SCL cp *SCL_60m.tif SCL
cd SCL cd SCL
ls *.tif | parallel --jobs 4 imageMissingData.py {} 50 ls *.tif | parallel --jobs 4 imageMissingData.py {} $USERDIR"findproducts.json" 70
if [ ! -d "../TCI" ]; then if [ ! -d "../TCI" ]; then
mkdir ../TCI mkdir ../TCI
fi fi
for scl in $(ls *.tif) for scl in $(ls *.tif)
do do
fileprefix=$(echo $scl | (cut -d"_" -f1)) fileprefix=$(echo $scl | (cut -d"_" -f1))
cp ../$fileprefix"_TCI_60m.tif" ../TCI/$fileprefix"_TCI_60m.tif" cp ../$fileprefix"_TCI_10m.tif" ../TCI/$fileprefix"_TCI_10m.tif"
done done
cd ../TCI cd $USERDIR"out/"
mogrify -format jpg -quality 100 *.tif
cd $USERDIR"jp2/"
done done
else else
echo "BY DATE" echo "BY DATE"
cd jp2/merge_out/ cd $USERDIR"out/mask/"
if [ ! -d "SCL" ]; then if [ ! -d "SCL" ]; then
mkdir SCL mkdir SCL
fi fi
cp *SCL_60m.tif SCL cp *SCL_60m.tif SCL
cd SCL cd SCL
ls *.tif | parallel --jobs 4 imageMissingData.py {} 50 ls *.tif | parallel --jobs 4 imageMissingData.py {} $USERDIR"findproducts.json" 70
if [ ! -d "../TCI" ]; then if [ ! -d "../TCI" ]; then
mkdir ../TCI mkdir ../TCI
...@@ -102,39 +157,162 @@ else ...@@ -102,39 +157,162 @@ else
do do
fileprefix=$(echo $scl | (cut -d"_" -f1)) fileprefix=$(echo $scl | (cut -d"_" -f1))
cp ../$fileprefix"_TCI_10m.tif" ../TCI/$fileprefix"_TCI_10m.tif" cp ../$fileprefix"_TCI_10m.tif" ../TCI/$fileprefix"_TCI_10m.tif"
done
cd ../TCI cd ../TCI
echo "mogrify -format jpg -quality 100 *.tif" done
fi
#6.-Creating Thumbnails
echo ${RED}"Creating Thumbnails"${NC}
if [ $BYTILE -ne "0" ]; then
echo "BY TILE"
cd $USERDIR"out/"
for i in $(ls -d */) ; do
cd $i"mask"
if [ ! -d "SCL/thumbnails" ]; then
mkdir SCL/thumbnails
fi
if [ ! -d "TCI/thumbnails" ]; then
mkdir TCI/thumbnails
fi
cd TCI
mogrify -format jpg -quality 100 *.tif
mogrify -path ./thumbnails -resize 10% *.jpg
cd $USERDIR"out/"
done
else
echo "BY DATE"
cd $USERDIR"out/mask/"
if [ ! -d "SCL/thumbnails" ]; then
mkdir SCL/thumbnails
fi
if [ ! -d "TCI/thumbnails" ]; then
mkdir TCI/thumbnails
fi
cd TCI
mogrify -format jpg -quality 100 *.tif mogrify -format jpg -quality 100 *.tif
mogrify -path ./thumbnails -resize 10% *.jpg
fi fi
#5.-Extract SCL information
echo "Extracting SCL information..." #7.-Extract SCL information
cd $USERDIR"jp2/" echo ${RED}"Extracting SCL information..."${NC}
if [ $BYTILE -ne "0" ]; then if [ $BYTILE -ne "0" ]; then
echo "BY TILE" echo "BY TILE"
cd $USERDIR"out/"
for i in $(ls -d */) ; do for i in $(ls -d */) ; do
cd ../out/$i"SCL" cd ../out/$i"mask/SCL"
if [ ! -d "json" ]; then if [ ! -d "json" ]; then
mkdir json mkdir json
fi fi
SCLimageListToJSON.sh $USERDIR"out/"$i"SCL/" $USERDIR"findProducts.json" 6 SCLimageListToJSON.sh $USERDIR"out/"$i"mask/SCL/" $USERDIR"findProducts.json" 6
mv *.json json mv *.json json
cd json cd json
tile=$(echo $i | cut -d "/" -f1) tile=$(echo $i | cut -d "/" -f1)
mergeL2ASCL_JSON.py $USERDIR"out/"$i"SCL/json/" > $USERDIR$tile"_sclData.json" mergeL2ASCL_JSON.py $USERDIR"out/"$i"mask/SCL/json/" > $USERDIR$tile"_sclData.json"
cd $USERDIR"jp2/" echo $tile
echo "_sclData.json"
cd $USERDIR"out/"
done done
else else
echo "BY DATE" echo "BY DATE"
cd $USERDIR/jp2/merge_out/SCL cd $USERDIR"out/mask/SCL/"
if [ ! -d "json" ]; then if [ ! -d "json" ]; then
mkdir json mkdir json
fi fi
SCLimageListToJSON.sh $USERDIR"jp2/merge_out/SCL/" $USERDIR"findProducts.json" 6 SCLimageListToJSON.sh $USERDIR"out/mask/SCL/" $USERDIR"findProducts.json" 6
mv *.json json mv *.json json
cd json cd json
mergeL2ASCL_JSON.py $USERDIR"jp2/merge_out/SCL/json/" > $USERDIR"scl_data.json" mergeL2ASCL_JSON.py $USERDIR"out/mask/SCL/json/" > $USERDIR"scl_data.json"
fi
#8.- Split SCL
echo ${RED}"Spliting SCL image..."${NC}
if [ $BYTILE -ne "0" ]; then
echo "BY TILE"
else
echo "BY DATE"
cd $USERDIR"out/mask/SCL/"
if [ ! -d "split" ]; then
mkdir split
fi
for f in $(ls *.tif) ; do
echo $f
splitSCL.py $f split/
done
fi fi
#9.-Create Tiles
echo ${RED}"Creating Tiles..."${NC}
if [ $BYTILE -ne "0" ]; then
echo "BY TILE"
cd $USERDIR"out/"
for d in $(ls -d */) ; do
cd $d"mask/TCI/"
if [ ! -d "tiles" ]; then
mkdir tiles
fi
for i in $(ls *.tif) ;
do
cd tiles
gdal2tiles.py ../$i
cd ..
done
cd $USERDIR"out/"
cd $d"mask/SCL/"
if [ ! -d "tiles" ]; then
mkdir tiles
fi
for i in $(ls *.tif) ;
do
cd tiles
gdal2tiles.py ../$i
cd ..
done
cd $USERDIR"out/"
done
else
echo "BY DATE"
cd $USERDIR"out/mask/TCI/"
if [ ! -d "tiles" ]; then
mkdir tiles
fi
for i in $(ls *.tif) ;
do
cd tiles
gdal2tiles.py --zoom="10-15" ../$i
cd ..
done
cd $USERDIR"out/mask/SCL/split/"
if [ ! -d "tiles" ]; then
mkdir tiles
fi
for i in $(ls *.tif) ;
do
cd tiles
gdal2tiles.py --zoom="10-15" ../$i
cd ..
done
fi
exit
#9.-Create ZIP file
echo "Compresing Files in " $USERDIR "..."
cd $USERDIR
zip -r myzip . --exclude *.json --exclude *.zip --exclude *thumbnails* --exclude *.tif -x ./jp2\* ./L2A\* \*\*/\*/SCL/\*^C
echo "END"
#mv *SCL_60m.tif json #mv *SCL_60m.tif json
#cd json #cd json
#!/bin/sh
#Convert L2A Products to L3A using L2AProductListToL3A Process
DATE=$1
TILE=$2
INDIR=$3
OUTDIR=$4
cd $INDIR
for f in $(ls | egrep $DATE.{21}$TILE) ; do
TMPDIR=MSIL3A_$DATE"_"$TILE"/"
cd $OUTDIR
if [ ! -d $TMPDIR"L2A/TCI/" ]; then
mkdir -p $TMPDIR"L2A/TCI/"
fi
if [ ! -d $TMPDIR"L2A/SCL/" ]; then
mkdir -p $TMPDIR"L2A/SCL/"
fi
JP2DIR=$TMPDIR"GRANULE/L3A_"$TILE"_"$DATE"/IMG_DATA/R60m/"
if [ ! -d $JP2DIR ]; then
mkdir -p $JP2DIR
fi
cd ..
ln -s $INDIR$f $OUTDIR$TMPDIR"L2A/"$f
unzip -n -j $OUTDIR$TMPDIR"L2A/"$f *TCI_60m.jp2 -d $OUTDIR$TMPDIR"L2A/TCI/"
unzip -n -j $OUTDIR$TMPDIR"L2A/"$f *SCL_60m.jp2 -d $OUTDIR$TMPDIR"L2A/SCL/"
echo $f
done
cd $OUTDIR$TMPDIR"L2A/"
L2AtoL3A.py $OUTDIR$TMPDIR"L2A/"
#!/bin/sh #!/bin/sh
INDIR=$1 #Output filename INDIR=$1 #Directory with images to merge
OUTDIR=$2 #Output filename OUTDIR=$2 #Output directory
POLYGON=$3 #Crop Window POLYGON=$3 #Crop Window (minX, maxX, minY, maxY)
DIRNAME=$(echo $INDIR | cut -d"/" -f1) DIRNAME=$(echo $INDIR | cut -d"/" -f1)
echo "INDIR:" $INDIR echo "INDIR:" $INDIR
echo "OUTDIR:" $OUTDIR echo "OUTDIR:" $OUTDIR
...@@ -22,7 +22,7 @@ echo $DIRNAME ...@@ -22,7 +22,7 @@ echo $DIRNAME
for f in $(ls); do echo ${f}; done; for f in $(ls); do echo ${f}; done;
#exit #exit
##------------------------------------------------------------------------------ ##------------------------------------------------------------------------------
MERGEDIMAGE=$DIRNAME"_TCI_60m_merged.tif" MERGEDIMAGE=$DIRNAME"_TCI_10m_merged.tif"
if [ ! -e $MERGEDIMAGE ]; then if [ ! -e $MERGEDIMAGE ]; then
echo "gdal_merge.py -o "$MERGEDIMAGE $(ls *TCI_10m.jp2) echo "gdal_merge.py -o "$MERGEDIMAGE $(ls *TCI_10m.jp2)
gdal_merge.py -o $MERGEDIMAGE $(ls *TCI_10m.jp2) gdal_merge.py -o $MERGEDIMAGE $(ls *TCI_10m.jp2)
...@@ -30,7 +30,7 @@ else ...@@ -30,7 +30,7 @@ else
echo PASSING $MERGEDIMAGE FOUND echo PASSING $MERGEDIMAGE FOUND
fi fi
CROPEDIMAGE=$OUTDIR$DIRNAME"_TCI_60m.tif" CROPEDIMAGE=$OUTDIR$DIRNAME"_TCI_10m.tif"
if [ ! -e $CROPEDIMAGE ]; then if [ ! -e $CROPEDIMAGE ]; then
gdal_translate -projwin $BOX -projwin_srs WGS84 -ot Byte -of GTiff $MERGEDIMAGE $CROPEDIMAGE gdal_translate -projwin $BOX -projwin_srs WGS84 -ot Byte -of GTiff $MERGEDIMAGE $CROPEDIMAGE
else else
......
...@@ -7,6 +7,8 @@ DIRNAME=$(echo $JP2DIR | cut -d"/" -f1) ...@@ -7,6 +7,8 @@ DIRNAME=$(echo $JP2DIR | cut -d"/" -f1)
#echo "POLYGON" #echo "POLYGON"
#echo $POLYGON #echo $POLYGON
BOX=$(polygonToBox.py "$POLYGON") BOX=$(polygonToBox.py "$POLYGON")
echo "Polygon BOX"
echo $BOX
if [ "$#" -le 1 ]; then if [ "$#" -le 1 ]; then
echo "Usage " $0 " <Parent Directory with subdirectories> <POLYGON> [jobs]" echo "Usage " $0 " <Parent Directory with subdirectories> <POLYGON> [jobs]"
exit exit
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment