Multiple corrections according to new nomenclature and multiple CM versions -...

Multiple corrections according to new nomenclature and multiple CM versions - Introducing Sen2Chain 0.7.0

sen2chain/__init__.py

@@ -32,5 +32,5 @@ from .utils import format_word, grouper, datetime_to_str, str_to_datetime, human
 from .geo_utils import serialise_tiles_index, get_processed_indices_vect, crop_product_by_shp
 from .multi_processing import l2a_multiprocessing, cldidx_multiprocessing, cld_multiprocessing, idx_multiprocessing
 from .tileset import TileSet
-__version__ = "0.6.0"
+__version__ = "0.7.0"
 __author__ = "Jérémy Commins <jebins@openmailbox.org> & Impact <pascal.mouquet@ird.fr>"

sen2chain/library.py

@@ -167,6 +167,45 @@ class Library:
         for tile in self.l2a:
             Tile(tile).update_latest_ql()
 
+    def update_old_cloudmasks(self):
+        """
+        Produce or update the latest quicklook for the L2A library tiles
+        """
+        for tile in self.l2a:
+            Tile(tile).update_old_cloudmasks()
+            
+    def remove_very_old_cloudmasks(self):
+        """
+        Remove very old cloudmasks, matching pattern : *_CLOUD_MASK.tif for the L2A library tiles
+        """
+        for tile in self.l2a:
+            Tile(tile).remove_very_old_cloudmasks()
+
+    def move_old_quicklooks(self):
+        """
+        Move all old quicklooks to QL subfolder for the L2A library tiles
+        """
+        for tile in self.l2a:
+            Tile(tile).move_old_quicklooks()
+
+    def update_old_indices(self):
+        """
+        Rename old indices to match new cloudmask nomenclature
+        for the L2A library tiles
+        """
+        for tile in {f.name for f in list(self._indices_path.glob("*/*/"))}:
+            Tile(tile).update_old_indices()
+            
+    def init_md(self):
+        """
+        Initiate sen2chain metadata for all tile products (l2a, cloudmasks, indices (raw, masked, ql)) 
+        for the L1C, L2A and Indice library tiles
+        """
+        for tile in set([val for sublist in [k for k in [getattr(self,t) for t in self.__dict__]] for val in sublist]):
+           Tile(tile).init_md()
+
+
+
 class TempContainer:
     """Class for managing a downloaded L1C products.
 

sen2chain/products.py

@@ -210,8 +210,10 @@ class L1cProduct(Product):
             # the newly L2A must be moved to his L2A library folder.
             if self._library_product:
                 l2a_prod.archive()
-            L2aProduct(l2a_identifier).setPermissions()
-            L2aProduct(l2a_identifier).update_sen2chain_metadata()
+            l2a_prod.setPermissions()
+            l2a_prod.update_md()
+            #~ L2aProduct(l2a_identifier).setPermissions()
+            #~ L2aProduct(l2a_identifier).update_md()
         return self
 
     def process_ql(self,
@@ -415,7 +417,7 @@ class L2aProduct(Product):
             ql_filename = self.identifier + "_QL.tif"
         
         if out_path is None:
-            ql_folder = self.library_path.parent
+            ql_folder = self.library_path.parent / "QL"
             ql_folder.mkdir(parents=True, exist_ok=True)
             ql_path = ql_folder / ql_filename
         else:
@@ -512,28 +514,29 @@ class L2aProduct(Product):
         
         
     def compute_cloud_mask(self,
-                           version: str = "cm001",
+                           cm_version: str = "cm001",
                            probability: int = 1,
                            iterations: int = 5,
-                           #~ buffering: bool = True,
                            reprocess: bool = False,
                            out_path_mask = None,
                            out_path_mask_b11 = None
                            ) -> "L2aProduct":
         """
         """
-        logger.info("Computing cloudmask version {}: {}".format(version, self.identifier))
+        logger.info("Computing cloudmask version {}: {}".format(cm_version, self.identifier))
 
         cloudmask = NewCloudMaskProduct(l2a_identifier = self.identifier, 
                                         sen2chain_processing_version = self.sen2chain_processing_version,
-                                        version = version, 
+                                        cm_version = cm_version, 
                                         probability = probability, 
                                         iterations = iterations)
         
         if cloudmask.path.exists() and not reprocess:
                 logger.info("{} cloud mask already computed".format(cloudmask.identifier))
         else:
-            if version == "cm001":
+            if not cloudmask.path.parent.exists():
+                cloudmask.path.parent.mkdir(parents=True)
+            if cm_version == "cm001":
                 if cloudmask.path.exists():     # in version 3.8 will be updated using missing_ok = True
                     cloudmask.path.unlink()
                     cloudmask._info_path.unlink()
@@ -542,10 +545,10 @@ class L2aProduct(Product):
                                      dilatation=5,
                                      out_path=cloudmask.path)
                 #~ self.user_cloud_mask = cloudmask.updated_path
-            elif version == "cm002":
+            elif cm_version == "cm002":
                 cloudmask_cm001 = NewCloudMaskProduct(l2a_identifier = self.identifier, 
                                                       sen2chain_processing_version = self.sen2chain_processing_version,
-                                                      version = "cm001")
+                                                      cm_version = "cm001")
                 if cloudmask_cm001.path.exists():
                     if cloudmask.path.exists():    # in version 3.8 will be updated using missing_ok = True
                         cloudmask.path.unlink()
@@ -556,12 +559,12 @@ class L2aProduct(Product):
                                           out_path=Path(str(cloudmask_cm001.path).replace("CM001", "CM002-B11")))
                 else:
                     logger.info("No cloudmask version cm001 found, please compute this one first")
-            elif version == "cm003":
+            elif cm_version == "cm003":
                 toto=12
-            elif version == "cm004":
+            elif cm_version == "cm004":
                 toto=12
             else:
-                logger.info("Wrong cloudmask version {}".format(version))
+                logger.info("Wrong cloudmask version {}".format(cm_version))
 
             cloudmask.init_md()
         

sen2chain/sen2cor.py

@@ -10,7 +10,7 @@ import subprocess
 #~ import re
 from typing import Union
 
-#~ from .config import Config
+from .config import Config
 from .utils import get_current_Sen2Cor_version
 
 

sen2chain/tiles.py

@@ -22,7 +22,7 @@ from itertools import chain
 from .config import Config, SHARED_DATA
 from .utils import str_to_datetime, human_size, getFolderSize
 from .indices import IndicesCollection
-from .products import L1cProduct, L2aProduct, OldCloudMaskProduct, NewCloudMaskProduct
+from .products import L1cProduct, L2aProduct, OldCloudMaskProduct, NewCloudMaskProduct, IndiceProduct
 from .multi_processing import l2a_multiprocessing, cld_multiprocessing, cld_version_probability_iterations_reprocessing_multiprocessing, idx_multiprocessing
 
 logging.basicConfig(level=logging.INFO)
@@ -714,7 +714,7 @@ class Tile:
                     except:
                         logger.error("Can't remove {} from L2A folder".format(f.name))
         # identify 0B cloud masks
-        for f in self._paths["l2a"].glob("*L2A*_CLOUD_MASK.jp2"):
+        for f in self._paths["cloudmasks"].glob("*/*CM*.jp2"):
             if f.stat().st_size == 0:
                 logger.info("Corrupted cloud mask {} in L2A folder".format(f.name))
                 nb_id += 1
@@ -723,15 +723,15 @@ class Tile:
                     f.unlink()
                     nb_rm += 1
         # identify wrong size l2a_QL
-        for f in self._paths["l2a"].glob("*L2A*_QL.tif"):
+        for f in self._paths["l2a"].glob("QL/*_QL.tif"):
             if f.stat().st_size != 3617212:
-                logger.info("Corrupted L2A QL {} in L2A folder".format(f.name))
+                logger.info("Corrupted L2A QL {} in L2A QL folder".format(f.name))
                 nb_id += 1
                 if remove:
                     logger.info("Removing corrupted QL {} from L2A folder".format(f.name))
                     f.unlink()
                     nb_rm += 1
-        # identify 0B or absent indices QL
+        # identify 0B or absent indice QL
         for f in self._paths["indices"]:
             #~ logger.info(f, self._paths["indices"][f])
             for p in self._paths["indices"][f].glob("*_MSIL2A_*/"):
@@ -794,7 +794,7 @@ class Tile:
     def archive_l1c(self, 
                     size_only: bool = False,):
         """
-        Check and move l1c products to archive folder
+        Check and move l1c products to l1c archive folder
         
         """
         
@@ -836,7 +836,7 @@ class Tile:
     def archive_l2a(self,
                     size_only: bool = False,):
         """
-        Check errors and move l2a products to archive folder
+        Check errors and move l2a products to l2a archive folder
         
         """
         
@@ -874,6 +874,9 @@ class Tile:
     
     def archive_all(self, 
                     size_only: bool = False,):
+        """
+        Chain archive_l1c and archive_l2a functions
+        """        
         l1c_size = self.archive_l1c(size_only = size_only)
         l2a_size = self.archive_l2a(size_only = size_only)
         return l1c_size + l2a_size   
@@ -881,12 +884,13 @@ class Tile:
     def update_latest_ql(self):
         """
         Produce or update the latest l2a quicklook for the tile
+        And remove previous ones
         """
         p = self.l2a.last.identifier
         l2a = L2aProduct(p)
 
-        outfullpath = l2a.path.parent / (p[39:44] + '_' + p[0:4] + p[11:19] + '_QL_latest.jpg')
-        old_ql = list(l2a.path.parent.glob('*_QL_latest.jpg*'))
+        outfullpath = l2a.path.parent / "QL" / (l2a.tile + '_' + p[0:4] + Tile._get_date(p).strftime("%Y%m%d") + '_QL_latest.jpg')
+        old_ql = list((l2a.path.parent / "QL").glob('*_QL_latest.jpg*'))
         liste=[a for a in old_ql if str(outfullpath) not in str(a)]
         for f in liste:
             f.unlink()
@@ -895,7 +899,17 @@ class Tile:
             return
         else:
             l2a.process_ql(out_path = outfullpath, out_resolution=(750,750), jpg = True)
-        
+    
+    def move_old_quicklooks(self):
+        """
+        Move all old quicklooks to QL subfolder
+        """
+        logger.info("{}: Moving all quicklooks to QL/ subfolder".format(self.name))
+        (self._paths["l2a"] / "QL").mkdir(exist_ok=True, parents=True)
+        for f in self._paths["l2a"].glob("*_QL*"):
+            if f.is_file():
+                f.replace(f.parent / "QL" / f.name)
+    
     def update_old_cloudmasks(self):
         """
         Move and rename old cloudmasks to new cloudmask folder
@@ -918,6 +932,16 @@ class Tile:
         for f in self._paths["l2a"].glob("*L2A*_CLOUD_MASK*.jp2.aux.xml"):
             f.unlink()
 
+    def remove_very_old_cloudmasks(self):
+        """
+        Remove very old cloudmasks, matching pattern : *_CLOUD_MASK.tif
+        """
+        files = list(self._paths["l2a"].glob("*_CLOUD_MASK.tif"))
+        logger.info("{}: Removing {} very old cloudmasks".format(self.name, len(files)))
+        for f in files:
+            if f.is_file():
+                f.unlink()
+
     def update_old_indices(self):
         """
         Rename old indices to match new cloudmask nomenclature
@@ -926,9 +950,27 @@ class Tile:
         logger.info("Moving and renaming old indices")
         
         for indice, path in self._paths["indices"].items():
-            logger.info("Processing: {}".format(indice.upper()))
+            logger.info("{} - Processing: {}".format(self.name, indice.upper()))
             for f in list(Path(path).glob("*/*MASKED*")) + list(Path(path).glob("*/*QUICKLOOK*")):
                 f_renamed = f.name.replace("MASKED", "CM001").replace("QUICKLOOK", "CM001_QL")
                 f.rename(str(Path(f.parent / f_renamed)))
                 logger.info(f_renamed)
+                IndiceProduct(identifier = f_renamed)
+    
+    def init_md(self):
+        """
+        Initiate sen2chain metadata for all tile products (l2a, cloudmasks, indices (raw, masked, ql)) 
+        """
+        logger.info("{} - Initiating products metadata".format(self.name))
+        
+        for l2a in [product.identifier for product in self.l2a]:
+            L2aProduct(l2a)
+        
+        for cloudmask in [product.identifier for product in self.cloudmasks]:
+            NewCloudMaskProduct(cloudmask)
+        
+        for indice in [val for sublist in [getattr(self, i) for i in [p for p in self.paths["indices"]]] for val in sublist]:
+            IndiceProduct(identifier = indice.identifier)
+        
+      
 

sen2chain/time_series.py

@@ -35,7 +35,7 @@ from .config import Config
 from .tiles import Tile, ProductsList
 from .indices import IndicesCollection
 from .geo_utils import get_tiles_from_file
-from .products import L2aProduct
+from .products import L2aProduct, IndiceProduct
 
 
 logger = logging.getLogger(__name__)
@@ -142,7 +142,7 @@ class TimeSeries:
         :param cover_min: minimum cloud coverage.
         :param cover_max: maximum cloud coverage.
         """
-        products = vars(tile)[indice.lower()].masks
+        products = vars(tile)[indice.lower()].masks.cm001               ##########"" Ajouter ici un truc pour choisir le CM !!!!!
         filt = products.filter_dates(date_min, date_max).filter_clouds(cover_min, cover_max)
         return filt
 
@@ -241,10 +241,14 @@ class TimeSeries:
                                                                self._cover_min, self._cover_max)
 
                         for index1, prod in enumerate(products):
-                            prod_path = tile_indice_path / prod.identifier[:(-12 - len(indice))] / prod.identifier
+                            indice_product = IndiceProduct(prod.identifier)
+                            prod_path = str(indice_product.path)
+                            prod_path_unmasked = prod_path.replace("_CM001", "")
+                            prod_path_cloud_proba = L2aProduct(indice_product.l2a).msk_cldprb_20m
+                            #~ prod_path = tile_indice_path / prod.identifier[:(-12 - len(indice))] / prod.identifier
                             #~ cloud_path = tile_obj.paths["l2a"] / (prod.identifier[:(-12 - len(indice))] + "_CLOUD_MASK.jp2")
-                            prod_path_unmasked = tile_indice_path / prod.identifier[:(-12 - len(indice))] / (prod.identifier[:-11] + '.jp2')
-                            prod_path_cloud_proba = L2aProduct(prod.identifier[:(-12 - len(indice))]).msk_cldprb_20m
+                            #~ prod_path_unmasked = tile_indice_path / prod.identifier[:(-12 - len(indice))] / (prod.identifier[:-11] + '.jp2')
+                            #~ prod_path_cloud_proba = L2aProduct(prod.identifier[:(-12 - len(indice))]).msk_cldprb_20m
                             #~ logger.info(prod_path_cloud_proba)
 
                             
@@ -275,9 +279,14 @@ class TimeSeries:
         logger.info("Execution time: {}".format(timedelta(seconds = end - start)))
 
     def _raster_stats_multi(self, features, shared_list, proc_item):
-        prod_path = proc_item[2] / proc_item[0].identifier[:(-12 - len(proc_item[3]))] / proc_item[0].identifier
-        prod_path_unmasked = proc_item[2] / proc_item[0].identifier[:(-12 - len(proc_item[3]))] / (proc_item[0].identifier[:-11] + '.jp2')
-        prod_path_cloud_proba = L2aProduct(proc_item[0].identifier[:(-12 - len(proc_item[3]))]).msk_cldprb_20m
+        indice_product = IndiceProduct(proc_item[0].identifier)
+        prod_path = str(indice_product.path)
+        prod_path_unmasked = prod_path.replace("_CM001", "")
+        prod_path_cloud_proba = L2aProduct(indice_product.l2a).msk_cldprb_20m
+        
+        #~ prod_path = proc_item[2] / proc_item[0].identifier[:(-12 - len(proc_item[3]))] / proc_item[0].identifier
+        #~ prod_path_unmasked = proc_item[2] / proc_item[0].identifier[:(-12 - len(proc_item[3]))] / (proc_item[0].identifier[:-11] + '.jp2')
+        #~ prod_path_cloud_proba = L2aProduct(proc_item[0].identifier[:(-12 - len(proc_item[3]))]).msk_cldprb_20m
         #~ logger.info(prod_path_cloud_proba)
         
         fid = proc_item[1]
@@ -578,8 +587,9 @@ class TimeSeries:
                     tile_indice_path = tile_obj.paths["indices"][df_name.lower()]
                     tile_l2a_path = tile_obj.paths["l2a"]
                     prod_path = tile_indice_path / prod_id / row['filename']
-                    tci_path = list(Path(str(tile_l2a_path / row['filename'][:(-12 - len(df_name))])+ '.SAFE/')\
-                                             .glob('GRANULE/*/IMG_DATA/R10m/*_TCI_10m.jp2'))
+                    tci_path = L2aProduct(IndiceProduct(row['filename']).l2a).tci_10m
+                    #~ tci_path = list(Path(str(tile_l2a_path / row['filename'][:(-12 - len(df_name))])+ '.SAFE/')\
+                                             #~ .glob('GRANULE/*/IMG_DATA/R10m/*_TCI_10m.jp2'))
                                              
                     crop_extent = gpd.read_file(str(self._vectors_file))
                     raster_tci = rasterio.open(tci_path[0])