diff --git a/app/Archive.py b/app/Archive.py
index cb9596b57f1b5e483b315302723fe7537284b750..058bba0467453724b4272566d3df484b1dd857e6 100644
--- a/app/Archive.py
+++ b/app/Archive.py
@@ -15,6 +15,7 @@ from uuid import uuid4
from collections import defaultdict, UserDict
import app.Satellites as Satellites
import app.Outils as Outils
+import app.Constantes as Constantes
class Archive():
"""
@@ -39,14 +40,23 @@ class Archive():
:type annee_fin: Entier
"""
- def __init__(self, capteur, extensions, niveau, emprise, sortie, annee_debut, annee_fin):
+ def __init__(self, capteur, bandes, niveau, emprise, zone_etude, sortie, annee_debut, annee_fin, seuil):
"""
Créé une instance de la classe 'Archive'.
"""
self._capteur = capteur
self.niveau = niveau
- self.extent_img = extensions
+
+ if bandes == "RGB" :
+ self.extent_img = Satellites.SATELLITE[self._capteur][niveau][:-1]
+ else :
+ self.extent_img = Satellites.SATELLITE[self._capteur][niveau]
+
+ if self.niveau == "LEVEL2A" :
+ self.nuage = Satellites.SATELLITE[self._capteur]["NUAGE"]
+ else :
+ self.nuage = None
self.liste_annees = []
@@ -56,9 +66,10 @@ class Archive():
self.liste_annees.append(annee_fin)
self.emprise = emprise
+ self.zone_etude = zone_etude
self.dossier_sortie = sortie
- self.emprise = emprise
-
+ self.seuil_nuage = seuil
+
self.liste_archive = []
self.logger = Outils.Log("log", "Archive")
@@ -392,6 +403,43 @@ class Archive():
self.logger.info("All images have been downloaded !")
+ def calcul_ennuagement(self, tuiles_nuage):
+
+ option_warp = gdal.WarpOptions(format='Mem', cropToCutline=True, outputType=gdal.GDT_Float32,\
+ cutlineDSName=self.emprise, dstNodata="NaN")
+ mosaic_nuage = gdal.Warp("", tuiles_nuage, options=option_warp)
+
+ rows = mosaic_nuage.RasterYSize # Rows number
+ cols = mosaic_nuage.RasterXSize # Columns number
+
+ data = mosaic_nuage.GetRasterBand(1).ReadAsArray(0, 0, cols, rows)
+
+ nombre_pixel_zone = 0
+ nombre_pixel_null = 0
+ nombre_pixel_non_nuage = 0
+
+ for donnees in data :
+ d = donnees.tolist()
+ nombre_pixel_zone += len(d)
+ nombre_pixel_null += d.count(-10000) + np.isnan(d).sum()
+ nombre_pixel_non_nuage += d.count(0)
+
+ nombre_pixel_non_null = nombre_pixel_zone - nombre_pixel_null
+
+ # self.logger.debug("nombre_pixel_non_null : {}".format(nombre_pixel_non_null))
+ # self.logger.debug("nombre_pixel_null : {}".format(nombre_pixel_null))
+ # self.logger.debug("nombre_pixel_non_nuage : {}".format(nombre_pixel_non_nuage))
+
+ del data
+
+ ennuagement = (float(nombre_pixel_non_nuage)/float(nombre_pixel_non_null)) if nombre_pixel_non_null != 0 else 0.0
+ self.logger.debug("Ennuagement : {}%".format(100 - round(ennuagement*100, 2)))
+
+ # Computer cloud's percentage with dist (sum of cloud) by sum of the image's extent
+ # return mosaic_nuage, 1.0 - (float(dist)/(np.sum(mask_spec)) if np.sum(mask_spec) != 0 else 0)
+ return 1.0 - ennuagement
+
+
def calcul_aire(self, tuiles_image):
"""
Méthode pour fusionner un ensemble d'images et calculer l'aire de l'image obtenue.
@@ -409,10 +457,10 @@ class Archive():
cols = mosaic_image.RasterXSize # Columns number
data = mosaic_image.GetRasterBand(1).ReadAsArray(0, 0, cols, rows).astype(np.float16)
- mask_spec = np.isin(data, [-10000, math.isnan], invert=True)
+ mask_spec = ~(np.isnan(data) | np.isin(data, [-10000]))
area = float((np.sum(mask_spec) * mosaic_image.GetGeoTransform()[1] * abs(mosaic_image.GetGeoTransform()[-1]) )/10000)
- self.logger.info("Aire = {0} ha".format(area))
+ # self.logger.info("Aire = {0} ha".format(area))
return mosaic_image, area
@@ -458,10 +506,11 @@ class Archive():
self.logger.info("Date : {0} -> {1} image(s)".format(date, len(liste_content)))
tuiles_image = []
+ tuiles_nuage = []
#Â Options Gdal pour la fusion des bandes
- options_vrt = gdal.BuildVRTOptions(separate=True)
- options_translate = gdal.TranslateOptions(format="Mem")
+ options_vrt = gdal.ParseCommandLine('-resolution highest -srcnodata -10000 -vrtnodata NaN -separate')
+ options_translate = gdal.ParseCommandLine('-of Mem -ot Float32 -a_nodata NaN')
self.logger.info("Extraction des images")
@@ -493,12 +542,32 @@ class Archive():
# on découpe l'image selon l'emprise
# et on conserve l'image obtenue
vrt = gdal.BuildVRT("", liste_bandes, options=options_vrt)
- tuiles_image.append(Outils.clip(gdal.Translate("", vrt, options=options_translate), self.emprise))
- # On libère la mémoire
+ image = Outils.clip(vrt, self.zone_etude)
+ _, aire = self.calcul_aire(image)
+
+ if aire > 0.0 :
+ tuiles_image.append(vrt)
+
+ # On libère la mémoire
+ for mmap_name in liste_mem :
+ gdal.Unlink(mmap_name)
+
+ liste_mem = []
+
+ if self.nuage is not None :
+ image_nuage = [f for f in zip_img if self.nuage[0] in f][0]
+ mmap_name = "/vsimem/"+uuid4().hex
+
+ gdal.FileFromMemBuffer(mmap_name, tzip.read(image_nuage))
+ liste_mem.append(mmap_name)
+ tuiles_nuage.append(Outils.clip(gdal.Open(mmap_name), self.zone_etude))
+
for mmap_name in liste_mem :
gdal.Unlink(mmap_name)
+ liste_mem = []
+
liste_bandes = None
del tzip
@@ -506,21 +575,22 @@ class Archive():
del liste_content
- # On effectue une mosaïque des images qu'on découpe selon l'emprise,
- # puis on calcule l'aire de l'image ainsi obtenue
- self.logger.info("Calcul de l'aire")
- image, aire = self.calcul_aire(tuiles_image)
- del tuiles_image
+ if len(tuiles_image) > 0 :
+
+ if len(tuiles_nuage) == 0 or self.calcul_ennuagement(tuiles_nuage) <= self.seuil_nuage:
+
+ # On effectue une mosaïque des images qu'on découpe selon l'emprise.
+ image = Outils.clip(tuiles_image, self.emprise)
- # Si l'aire est positive on enregistre l'image
- if aire > 0.0 :
+ self.logger.info("Sauvegarde des images")
- self.logger.info("Sauvegarde des images")
+ dossier = "{0}/{1}/Images".format(self.dossier_sortie, date[:4])
- dossier = "{0}/{1}/Images".format(self.dossier_sortie, date[:4])
+ self.logger.debug("Dossier image : {0}".format(dossier))
- self.logger.debug("Dossier image : {0}".format(dossier))
+ self.ecriture_geotiff(image, "{0}/{1}.tif".format(dossier,date))
- self.ecriture_geotiff(image, "{0}/{1}.tif".format(dossier,date))
+ del tuiles_nuage
- del image
\ No newline at end of file
+ del tuiles_image
+ del image
\ No newline at end of file
diff --git a/app/Constantes.py b/app/Constantes.py
index e7b1d076377aa234c53b8b073febdc3564fd1e10..7307edc7106f46dd1d856445a11ee8f22a3681d6 100644
--- a/app/Constantes.py
+++ b/app/Constantes.py
@@ -4,13 +4,4 @@
import logging
NIVEAU_DEFAUT = logging.DEBUG
-
-CLOUD_THRESHOLD = 0.0
-
EPSG_PHYMOBAT = 2154
-
-# Extensions des bandes bleue, verte et rouge
-BANDES_RGB = ['_FRC_B2.tif', '_FRC_B3.tif', '_FRC_B4.tif']
-
-# Extensions des bandes bleue, verte et rouge + proche infra-rouge
-BANDES_RGBI = ['_FRC_B2.tif', '_FRC_B3.tif', '_FRC_B4.tif', '_FRC_B8.tif']
diff --git a/app/Outils.py b/app/Outils.py
index bda1c573b7e9792e7d050b92cbc34c3ddcae6241..efec4bc9f938a84ed13cd076876ce555958f6e33 100644
--- a/app/Outils.py
+++ b/app/Outils.py
@@ -97,7 +97,7 @@ def limitation_temporelle(secondes):
# except TimeoutException:
# pass
-def clip(image, cut, form="Mem", dst=""):
+def clip(image, cut, form="Mem", dst="", type=gdal.GDT_Float32, nodata="NaN"):
"""
Découpe une image selon un shapefile donné. Si le type de l'image en sortie
n'est pas spécifié, l'image sera retournée en pointeur mémoire et ne sera pas sauvegardée.
@@ -116,6 +116,6 @@ def clip(image, cut, form="Mem", dst=""):
"""
option_clip = gdal.WarpOptions(cropToCutline=True,\
- cutlineDSName=cut, outputType=gdal.GDT_Float32 , format=form, dstNodata='NaN')
+ cutlineDSName=cut, outputType=type , format=form, dstNodata=nodata)
return gdal.Warp(dst, image, options=option_clip)
\ No newline at end of file
diff --git a/app/Processing.py b/app/Processing.py
index 04010fb4475d85f76870bc08917b27874600551d..3c72971bd47528e2db3e62a850446640cfac3ca3 100644
--- a/app/Processing.py
+++ b/app/Processing.py
@@ -27,7 +27,8 @@ class Processing(object):
if not self.annee_fin :
self.annee_fin = datetime.datetime.now().year
- self.check_download = Archive.Archive(self.capteur, self.extensions, self.niveau, self.emprise, self.resultats, self.annee_debut, int(self.annee_fin))
+ self.check_download = Archive.Archive(self.capteur, self.bandes, self.niveau, self.emprise, self.zone_etude,\
+ self.resultats, self.annee_debut, int(self.annee_fin), self.seuil_nuage)
self.check_download.listing()
self.check_download.download_auto(self.id, self.mdp, self.proxy)
@@ -99,5 +100,5 @@ class Processing(object):
"""
Calul le ndvi, fusionnne en une seule image puis lance le module OTBPhenology
"""
- self.calcul_ndvi()
- self.otbPhenologie()
\ No newline at end of file
+ # self.calcul_ndvi()
+ # self.otbPhenologie()
\ No newline at end of file
diff --git a/app/Satellites.py b/app/Satellites.py
index 2f3eb969221b7dc85e44fd96ca5970d4d7067b58..568009590cc76730b74d2c970a03d61b69435e39 100644
--- a/app/Satellites.py
+++ b/app/Satellites.py
@@ -14,6 +14,9 @@ SATELLITE["SENTINEL2"]["resto"] = "resto2"
SATELLITE["SENTINEL2"]["token_type"] = "text"
SATELLITE["SENTINEL2"]["R"] = 2
SATELLITE["SENTINEL2"]["PIR"] = 3
+SATELLITE["SENTINEL2"]["LEVEL2A"] = ['_SRE_B2.tif', '_SRE_B3.tif', '_SRE_B4.tif', '_SRE_B8.tif']
+SATELLITE["SENTINEL2"]["NUAGE"] = ['_CLM_R1.tif']
+SATELLITE["SENTINEL2"]["LEVEL3A"] = ['_FRC_B2.tif', '_FRC_B3.tif', '_FRC_B4.tif', '_FRC_B8.tif']
########################## LANDSAT #################################
diff --git a/main.py b/main.py
index 8d616da9b10024769665ff1cd46e2d253067f61c..055fe537d3a8613070831911991f9671fc0d7ca1 100644
--- a/main.py
+++ b/main.py
@@ -34,11 +34,7 @@ class Telechargement(Processing):
# Capteur utilisé
self.capteur = "{}".format(configfile["satellite"]["capteur"])
self.niveau = "{}".format(configfile["satellite"]["processingLevel"])
-
- if configfile["satellite"]["bandes"] == "RGB" :
- self.extensions = Constantes.BANDES_RGB
- else :
- self.extensions = Constantes.BANDES_RGBI
+ self.bandes = "{}".format(configfile["satellite"]["bandes"])
# Date de début et de fin de la recherche
try:
@@ -47,9 +43,14 @@ class Telechargement(Processing):
raise "L'année de départ est requise."
self.annee_fin = "{}".format(configfile["donnees"]["annee_fin"])
+ self.seuil_nuage = float("{}".format(configfile["donnees"]["seuil_nuage"]))/100.0
- # Emprise de l'étude
- self.emprise = "{}".format(configfile["donnees"]["chemin_emprise"])
+ # Emprise et zone de l'étude
+ self.emprise = "{}".format(configfile["donnees"]["chemin_emprise"])
+ self.zone_etude = "{}".format(configfile["donnees"]["chemin_zone_etude"])
+
+ if not self.zone_etude :
+ self.zone_etude = self.emprise
# Identifiant, mot de passe et proxy pour le téléchargement des images Théia
self.id = "{}".format(configfile["theia"]["identifiant"])
@@ -60,9 +61,8 @@ class Telechargement(Processing):
def run(self):
"""
Fonction pour lancer le programme
- """
-
- # Début du processus
+ """
+ Début du processus
debut = time.time()
# Recherche de nouvelles images non traitées et téléchargement de celles-ci le cas échéant
@@ -77,7 +77,16 @@ class Telechargement(Processing):
self.logger.info('Programme terminé en {} secondes'.format(fin - debut))
nb_jours = int(time.strftime('%d', time.gmtime(fin - debut))) - 1
self.logger.info("Cela représente {} jour(s) {}".format(nb_jours, time.strftime('%Hh %Mmin%S', time.gmtime(fin - debut))))
-
+
if __name__ == "__main__":
+
app = Telechargement()
- sys.exit(app.run())
\ No newline at end of file
+
+ SUCCESS = False
+
+ while not SUCCESS :
+ try:
+ sys.exit(app.run())
+ SUCCESS = True
+ except Exception as e:
+ pass
\ No newline at end of file