# -*- coding: utf-8 -*- """ Sistema de Predicción de Peligro de Incendios Forestales de México Genera el Índice de Verdor Relativo """ # Carpetas de trabajo import glob import os import time import rasterio from rasterio_class import RasterioProcessor import warnings dir_base = 'D:/proyectos/incendios/proceso_diario/DR' dir_NDVI = 'H:/2025-10-03_MADTEC/tif_NDVI_HANTS_LaNOT' # Suprimir warnings temporalmente warnings.filterwarnings('ignore', category=RuntimeWarning) dGlobals = { 'dir_base': dir_base, 'dir_NDVI': dir_NDVI, 'tif_mask': f'{dir_base}/in/mask/maskndvih500.tif', 'tif_NDVI99': f'{dir_base}/in/coefs/NDVI99_V2.tif', 'tif_DIF_NEXT_99': f'{dir_base}/in/coefs/dif_next_99_V2.tif', 'tif_LRMAX': f'{dir_base}/in/coefs/LRmax_V2.tif', 'dir_tmp': f'{dir_base}/tmp', 'dir_out': f'{dir_base}/out' } # Borra contenido de carpeta tmp def borraTMP(ruta): if 'tmp' in ruta: if os.path.isdir(ruta): archivos = glob.glob(ruta + "/*") for archivo in archivos: os.remove(archivo) # # Función con la lógica del proceso general # def proceso_completo(dGlobals, tif_NDVI, tif_mask, tif_NDVI99, tif_DIF_NEXT_99, tif_LRMAX, tif_dr, generar_tif_temporales ): print('-- proceso_completo...') with rasterio.open(tif_NDVI) as src_NDVI, \ rasterio.open(tif_mask) as src_mask, \ rasterio.open(tif_NDVI99) as src_NDVI99, \ rasterio.open(tif_DIF_NEXT_99) as src_DIF_NEXT_99, \ rasterio.open(tif_LRMAX) as src_LRMAX: # Crear instancia del procesador rasterio objRasterioProcessor = RasterioProcessor() # 0. AJUSTAR EXTENT print(' Paso 0/5: Ajustando extent de NDVI...') # Reprojectar el primer dataset para que coincida con el segundo tif_out = os.path.join(dGlobals['dir_tmp'], os.path.basename(tif_NDVI).replace('NDVIR_', 'DR_').replace('.tif', '_extent.tif') ) if generar_tif_temporales else None ds_extent = objRasterioProcessor.reproject_to_match( src_NDVI, src_mask, output_path=tif_out ) # 1. NDVI AJUSTADO X 100 print(' Paso 1/5: Multiplicando NDVI por 100...') tif_out = tif_out.replace('.tif', '_x100.tif') if generar_tif_temporales else None ds_x100 = objRasterioProcessor.multiply_by_scalar( dataset=ds_extent, scalar=100.0, # Escalar flotante output_path=tif_out, preserve_nodata=True ) # Convertiendo en entero print(' Paso 1/5: pasando a entero...') tif_out = tif_out.replace('.tif', '_int.tif') if generar_tif_temporales else None ds_x100int = objRasterioProcessor.float32_to_int16( dataset=ds_x100, output_path=tif_out, nodata_value=128, clamp_values=True ) # 2. DIF NDVI - NDVI99 print(' Paso 2/5: Calculando diferencia NDVI - NDVI99...') tif_out = tif_out.replace('.tif', '_dif.tif') if generar_tif_temporales else None ds_dif = objRasterioProcessor.subtract_datasets( dataset1=ds_x100int, dataset2=src_NDVI99, output_path=tif_out, ignore_nodata=True, # Ignorar áreas sin datos nodata_output=-128 ) # 3. CALCULO DE RG # Calcular RG (EQ3: int(DIF * 100) / DIF_NEXT_99) print(' Paso 3/5: Calculando RG...') # Multiplicación print(' Paso 3/5: Calculando RG (Multiplicación DIF x 100)...') tif_out = tif_out.replace('.tif', '_DIFx100.tif') if generar_tif_temporales else None ds_DIFx100 = objRasterioProcessor.multiply_by_scalar( ds_dif, 100, output_path=tif_out ) # División print(' Paso 3/5: Calculando RG (División (DIF * 100) / DIF_NEXT_99))...') tif_out = tif_out.replace('.tif', '_div.tif') if generar_tif_temporales else None ds_DIFx100divDIF_NEXT99 = objRasterioProcessor.divide_datasets( dataset1=ds_DIFx100, dataset2=src_DIF_NEXT_99, output_path=tif_out, safe_division=True ) # Convertiendo en entero print(' Paso 3/5: Calculando RG (Entero de (DIF * 100) / DIF_NEXT_99))...') tif_out = tif_out.replace('.tif', '_RG.tif') if generar_tif_temporales else None ds_RG = objRasterioProcessor.float32_to_int16( dataset=ds_DIFx100divDIF_NEXT99, output_path=tif_out, nodata_value=-32768, clamp_values=True ) # 4. CALCULO DE LR # Calcular LR (EQ4: int(RG * LRMAX / 100)) print(' Paso 4/5: Calculando LR...') # Multiplicación print(' Paso 4/5: Calculando LR (Multiplicación RG * LRMAX)...') tif_out = tif_out.replace('.tif', '_LRMult.tif') if generar_tif_temporales else None ds_LR_multi = objRasterioProcessor.multiply_datasets( dataset1=ds_RG, dataset2=src_LRMAX, output_path=tif_out, ignore_nodata=True ) # División print(' Paso 4/5: Calculando LR (División (RG * LRMAX / 100))...') tif_out = tif_out.replace('.tif', '_div100.tif') if generar_tif_temporales else None ds_LR_div100 = objRasterioProcessor.divide_by_scalar( dataset=ds_LR_multi, scalar=100, output_path=tif_out, ignore_nodata=True ) # Convertiendo en entero print(' Paso 4/5: Calculando LR (Entero de (RG * LRMAX / 100))...') tif_out = tif_out.replace('.tif', '_int.tif') if generar_tif_temporales else None ds_LR = objRasterioProcessor.float32_to_int16( dataset=ds_LR_div100, output_path=tif_out, nodata_value=-32768, clamp_values=True ) # 5. CALCULO DE DR (en memoria) print(' Paso 5/5: Calculando DR...') # Calcular DR (EQ5: 100 - LR) ds_DR = objRasterioProcessor.subtract_scalar_from_dataset( dataset=ds_LR, scalar=100, output_path=tif_dr, ignore_nodata=True ) #------------------- # Al final se cierran los dataset usados ds_extent.close() ds_x100.close() ds_x100int.close() ds_dif.close() ds_DIFx100.close() ds_DIFx100divDIF_NEXT99.close() ds_RG.close() ds_LR_multi.close() ds_LR_div100.close() ds_LR.close ds_DR.close # # Procesa todos los archivos de disponibles en la carpeta de NDVI # def proceso(dGlobals): # Procesa todos los archivos tif de la carpeta NDVI tif_NDVIs = sorted(glob.glob(f"{dGlobals['dir_NDVI']}/NDVI*20??????.tif"), reverse=True) for tif_NDVI in tif_NDVIs: # Crear nombre base para el archivo de salida DR base_name = os.path.basename(tif_NDVI).split('_')[1] tif_dr = f"{dGlobals['dir_out']}/DR_{base_name}" if not os.path.isfile(tif_dr): print(f'\nProcesando: {tif_NDVI}') # Elimina contenido de carpeta temporal borraTMP(dGlobals['dir_tmp']) # Proceso completo en una sola etapa proceso_completo( dGlobals, tif_NDVI=tif_NDVI, tif_mask=dGlobals['tif_mask'], tif_NDVI99=dGlobals['tif_NDVI99'], tif_DIF_NEXT_99=dGlobals['tif_DIF_NEXT_99'], tif_LRMAX=dGlobals['tif_LRMAX'], tif_dr=tif_dr, generar_tif_temporales=False ) print(f' DR generado: {tif_dr}') ####################################################### # # Proceso principal # ####################################################### fecha = time.strftime("%d/%m/%Y") hora = time.strftime("%H:%M:%S") print (f"Inicio: {fecha} {hora}\n") # Crear directorios de salida si no existen os.makedirs(dGlobals['dir_out'], exist_ok=True) os.makedirs(dGlobals['dir_tmp'], exist_ok=True) proceso(dGlobals) fecha = time.strftime("%d/%m/%Y") hora = time.strftime("%H:%M:%S") print (f"\nProceso terminado:{fecha} {hora}")