Sd/batch prod

config/clipImage.yaml

@@ -1,21 +1,10 @@
-# clip_image:
-#   working_directory: C:/Users/cmarmy/Documents/STDL/proj-vegroofs/data_test
-#   inputs:
-#     ortho_directory: 01_initial/images/rs_tlm/tiles/
-#     aoi: 01_initial/aoi/STDL_ZH_AOI.shp
-#     epsg: 'epsg:2056'
-#     predicate_sjoin: 'intersects'
-#   outputs: 
-#     clip_ortho_directory: 02_intermediate/images/rs_tlm/tiles
-#     extent_ortho_directory: 01_initial/images/rs_tlm/extent
-
 clip_image:
-  working_directory: C:/Users/cmarmy/Documents/STDL/proj-vegroofs/data
+  working_directory: D:/GitHubProjects/STDL_vegroof_production
   inputs:
-    ortho_directory: 01_initial/images/infer_moitie/tiles/
-    aoi: 01_initial/images/infer_moitie/extent/extent.shp
+    ortho_directory: sources/rasters/ZH
+    aoi: sources/Footprint/ZH/00-AVBodenbedGeb_updated.gpkg
     epsg: 'epsg:2056'
   outputs: 
-    clip_ortho_directory: 02_intermediate/images/infer_moitie/tiles
-    extent_ortho_directory: 01_initial/images/infer_moitie/extent
-
+    result_directory: ML/results
+    clip_ortho_directory: ML/results/infer_moitie/tiles
+    extent_ortho_directory: ML/results/infer_moitie/extent

config/logReg.yaml

@@ -3,25 +3,23 @@ hydra:
     dir: 02_intermediate/th/${now:%Y-%m-%d}/${now:%H-%M-%S}
 
 dev:
-  working_directory: C:/Users/cmarmy/Documents/STDL/proj-vegroofs/data
-  ortho_directory: 02_intermediate/images/infer_moitie/tiles
-  tile_delimitation: 02_intermediate/images/infer_moitie/extent/
-  ndvi_directory: 02_intermediate/images/infer_moitie/ndvi
-  lum_directory: 02_intermediate/images/infer_moitie/lum
-  roofs_file: 02_intermediate/th/2024-08-15/09-12-47/0_500_green_roofs.shp # 02_intermediate/gt/inf_roofs.gpkg # 
+  working_directory: D:/GitHubProjects/STDL_vegroof_production
+  ortho_directory: ML/results/infer_moitie/tiles
+  tile_delimitation: ML/results/infer_moitie/tiles/extent
+  ndvi_directory: ML/results/infer_moitie/ndvi
+  lum_directory: ML/results/infer_moitie/lum
+  roofs_file: sources/Footprint/ZH/00-AVBodenbedGeb_updated.gpkg
   roofs_layer: 
   gt: False
   green_tag: 'veg_new_3'
   green_cls: 'class_3'
-  chm_layer: 02_intermediate/autres/CHM_AOI_inf.gpkg
-  results_directory: 03_results/infer_moitie/
+  chm_layer: sources/CHM/ZH/chm_ZH_total.shp
+  results_directory: ML/results
   egid_train_test: egid_train_test_gt.csv
   th_ndvi: 0 # no thresholding -1
   th_lum: 500 # no thresholding 765 or 210000
   cls_ml: 'binary' # 'binary' 'multi' 'multi_aggreg'
   model_ml: 'LR' # 'LR' 'RF'
-  trained_model_dir: 03_results/scratch_gt/
+  trained_model_dir: ML/models
   epsg: 'epsg:2056'
-
-
-
+  do_overlay: False

scripts/calculate_raster.py

@@ -108,7 +108,7 @@ def calculate_lum(tile, band_nbr_red=1, band_nbr_green=2, band_nbr_blue=3, path=
     tile_list=[]
     tile_list.extend(tile_list_ortho)
 
-    for tile in tqdm(tile_list, 'Processing tiles'):
+    for _, tile in tqdm(enumerate(tile_list), total=len(tile_list), desc='Processing tiles'):
         tile = tile.replace("\\","/") #handle windows path
         ndvi_tile_path=os.path.join(NDVI_DIR, tile.split('/')[-1].replace('.tif', '_NDVI.tif'))
         _ = calculate_ndvi(tile, path=ndvi_tile_path)

scripts/clip_image.py

@@ -2,7 +2,7 @@
 import yaml
 import argparse
 from loguru import logger
-import tqdm as tqdm
+from tqdm import tqdm
 
 import geopandas as gpd
 
@@ -43,28 +43,36 @@
     OUTPUTS=cfg['outputs']
     OUTPUT_DIR=OUTPUTS['clip_ortho_directory']
     TILE_DELIMITATION=OUTPUTS['extent_ortho_directory']
+    RESULT_DIR = OUTPUTS['result_directory']
 
     os.chdir(WORKING_DIR)
     fct_misc.ensure_dir_exists(OUTPUT_DIR)
-
-    fct_misc.generate_extent(ORTHO_DIR, TILE_DELIMITATION, EPSG)
+    if not os.path.isfile(os.path.join(TILE_DELIMITATION,'extent.shp')):
+        fct_misc.generate_extent(ORTHO_DIR, TILE_DELIMITATION, EPSG)
     tiles=gpd.read_file(TILE_DELIMITATION)
 
     logger.info('Reading AOI geometries...')
 
     aoi = gpd.read_file(AOI)
+    # filter out invalid geometries
+    invalid_samples = aoi.loc[~aoi.geometry.is_valid]
+    aoi = aoi.loc[aoi.geometry.is_valid]
+    invalid_samples.to_file(os.path.join(RESULT_DIR, 'invalid_samples.gpkg'), driver='GPKG')
+    aoi.to_file(os.path.join(RESULT_DIR, 'valid_samples.gpkg'), driver='GPKG')
+
     # keep only the geometry column
     aoi = aoi.filter(['geometry'])
     # buffer every geometry by 50 units
-    for index, row in aoi.iterrows():
+    for index, row in tqdm(aoi.iterrows(), total=len(aoi), desc="Buffering geometries"):
         row = row.copy()
-        aoi.loc[index, 'geometry'] = row.geometry.buffer(50,join_style=2)
+        aoi.loc[index, 'geometry'] = row.geometry.buffer(1,join_style=2)
+
 
     aoi_clipped=fct_misc.clip_labels(labels_gdf=aoi, tiles_gdf=tiles, predicate_sjoin='intersects')
     aoi_clipped=aoi_clipped.reset_index(drop=True)
 
     i=1
-    for idx,row in aoi_clipped.iterrows(): 
+    for idx,row in tqdm(aoi_clipped.iterrows(), total=len(aoi_clipped), desc="Clipping rasters"): 
         fct_misc.clip_im(ORTHO_DIR, aoi_clipped.iloc[[idx]], OUTPUT_DIR, i, EPSG)
         i=i+1
     logger.success(f'Successfully clipped {i-1} images.')

scripts/functions/fct_misc.py

@@ -6,7 +6,7 @@
 import numpy as np
 import geopandas as gpd
 import pandas as pd
-from shapely.geometry import mapping, shape
+from shapely.geometry import mapping, shape, MultiPolygon
 from shapely.affinity import scale
 
 import rasterio
@@ -16,6 +16,8 @@
 
 import csv
 import warnings
+
+from tqdm import tqdm
 warnings.filterwarnings('ignore')
 
 def format_logger(logger):
@@ -78,7 +80,6 @@ def clip_labels(labels_gdf: gpd.GeoDataFrame, tiles_gdf: gpd.GeoDataFrame, predi
     labels_tiles_sjoined_gdf = gpd.sjoin(labels_gdf, tiles_gdf, how='inner', predicate=predicate_sjoin)
 
     def clip_row(row, fact=fact):
-
         old_geo = row.geometry
         scaled_tile_geo = scale(row.tile_geometry, xfact=fact, yfact=fact)
         new_geo = old_geo.intersection(scaled_tile_geo)
@@ -152,10 +153,10 @@ def generate_extent(PATH_IN: str, PATH_OUT: str, EPSG: str = 'epsg:2056'):
                 list_name.append(name)
 
     ext_merge=gpd.GeoDataFrame()
-    for _name in list_name:
+    for _, _name in tqdm(enumerate(list_name), total=len(list_name), desc="Computing extent"):
 
         _tif = os.path.join(PATH_IN, _name)
-        logger.info(f'Computing extent of {str(_name)} ...')
+        # logger.info(f'Computing extent of {str(_name)} ...')
 
         with rasterio.open(_tif) as src:
             gdf = gpd.GeoDataFrame.from_features(
@@ -193,7 +194,7 @@ def clip_im(TIFF_FOLDER: str, GPD: str, OUT_FOLDER: str, idx: int, EPSG: str = '
 
     with rasterio.open(os.path.join(TIFF_FOLDER, GPD.iloc[-1]['NAME']+'.tif')) as src:
 
-        logger.info('Clipping ' + GPD.iloc[-1]['NAME'] + '.tif...')
+        # logger.info('Clipping ' + GPD.iloc[-1]['NAME'] + '.tif...')
 
         out_image, out_transform = rasterio.mask.mask(
             src,
@@ -225,4 +226,4 @@ def clip_im(TIFF_FOLDER: str, GPD: str, OUT_FOLDER: str, idx: int, EPSG: str = '
         with rasterio.open(out_path, 'w', **out_meta) as dst:
             dst.write(out_image)
 
-        logger.info('Clipped image ' + GPD.iloc[-1]['NAME']+'_'+str(idx) + ' written...')
+        # logger.info('Clipped image ' + GPD.iloc[-1]['NAME']+'_'+str(idx) + ' written...')

scripts/greenery.py

@@ -10,6 +10,7 @@
 
 import pandas as pd
 import geopandas as gpd
+from time import time
 import fiona
 import rasterio
 from rasterio.features import shapes
@@ -86,6 +87,7 @@ def do_greenery(tile, shapes_roof, roofs):
     TH_NDVI=cfg['th_ndvi']
     TH_LUM=cfg['th_lum']
     EPSG=cfg['epsg']
+    DO_OVERLAY = cfg['do_overlay']
 
     os.chdir(WORKING_DIR)
 
@@ -137,12 +139,13 @@ def do_greenery(tile, shapes_roof, roofs):
     green_roofs_egid['EGID']=green_roofs_egid.index
     green_roofs_egid.index.names = ['Index']
 
-    logger.info('Filtering for overhanging vegetation...')
+    if DO_OVERLAY:
+        logger.info('Filtering for overhanging vegetation...')
+        CHM = os.path.join(WORKING_DIR, CHM_LAYER)
+        CHM_GPD=gpd.read_file(CHM)
+        CHM_GPD['geometry'] = CHM_GPD.buffer(1)
+        green_roofs_egid=gpd.overlay(green_roofs_egid, CHM_GPD, how='difference')
 
-    CHM = os.path.join(WORKING_DIR, CHM_LAYER)
-    CHM_GPD=gpd.read_file(CHM)
-    CHM_GPD['geometry'] = CHM_GPD.buffer(1)
-    green_roofs_egid=gpd.overlay(green_roofs_egid, CHM_GPD, how='difference')
     green_roofs_egid['area_green'] = green_roofs_egid.area