fix log1 duplicate call and introduce product_ananalysis_mask for NAS…

rapida/ntl/outage.py

@@ -1,4 +1,4 @@
-
+from matplotlib import pyplot as plt
 from datetime import datetime, timedelta
 import numbers
 import logging
@@ -148,25 +148,31 @@ async def detect_outage(
                 level = product.split('_')[0][-2:] if 'NRT' in deliverable else product[-2:]
                 sub_dataset_name = nasa_const.SUB_DATASETS[level]
                 daily_data = extract_bb(image_files=local_image_files, bbox=bbox,sds_name=sub_dataset_name,progress=progress)
+                if len(set(daily_data.ravel())) ==  1:
+                    logger.info(f'Invalid data for product {product} for {timestamp}. Skipping...')
+                    continue
                 positive = daily_data>=0
                 log_daily_data = np.zeros_like(daily_data)
                 log_daily_data[positive] = np.log1p(daily_data[positive])
-                log_daily_data = np.log1p(daily_data)
 
                 daily_data_label = f'{product}_{timestamp}'
                 arrays[daily_data_label] = log_daily_data
 
                 if mask_clouds:
+                    print('MC', product)
                     qf_array = extract_bb(image_files=local_image_files, sds_name='QF_Cloud_Mask',
                                           bbox=bbox, progress=progress).astype('u2')
                     cloud_confidence = (qf_array >> 6) & 0b11
                     is_cloudy = cloud_confidence == 3
                     arrays['CLOUD_MASK'] = is_cloudy
-                    analysis_mask |= is_cloudy
-                arrays[f'{daily_data_label}_MASK'] = analysis_mask
+                    product_analysis_mask =  (analysis_mask | is_cloudy).astype(bool)
+
+                else:
+                    product_analysis_mask = analysis_mask
+                arrays[f'{daily_data_label}_MASK'] = product_analysis_mask
                 log_difference, zscore, outage = utils.logdiff_outage(
                     log_monthly_data=log_monthly_data,log_daily_data=log_daily_data,
-                    analysis_mask=analysis_mask,percentage_drop=percentage_drop
+                    analysis_mask=product_analysis_mask,percentage_drop=percentage_drop
 
                 )
                 arrays[f'{daily_data_label}_LOGDIFF'] = log_difference
@@ -175,9 +181,11 @@ async def detect_outage(
 
     # --- 5. UNIFIED DISPLAY & EXPORT ---
     #file_name = utils.get_custom_bbox_label(bbox)
+
     file_name = get_best_semantic_label(bbox=bbox)
     outage_tif_path = os.path.join(dst_dir, f'{deliverable}_{file_name}.tif')
     outage_gpkg_path = os.path.join(dst_dir, f'{deliverable}_{file_name}.gpkg')
+    logger.info(f'Writing outage results to {outage_tif_path}')
     write_outage_tif(src_arrays=arrays, gt=gt, dst_path=outage_tif_path)
 
     if pop_vars:

rapida/ntl/vis.py

@@ -155,35 +155,51 @@ def display2(data=dict(), interpolation='nearest', title='', max_discrete_vals=5
         ax.set_aspect('equal')
 
         # --- THE NEW DYNAMIC LOGIC ---
-        # Filter out NaNs to safely check unique values
-        valid_vals = a[~np.isnan(a)]
+        # 1. Safely filter NaNs
+        if np.issubdtype(a.dtype, np.floating):
+            valid_vals = a[~np.isnan(a)]
+        else:
+            valid_vals = a
+
         unique_vals = np.unique(valid_vals)
 
-        # Check if this is a discrete/classification map (e.g., grid_health or masks)
-        if len(unique_vals) <= max_discrete_vals:
-            # It's discrete! Use specific colors and a custom legend.
-            # (Matches your Black/Yellow/Red layout for 3 values)
+        # 2. Handle empty (all NaNs) or uniform (all 1s) arrays natively
+        if len(unique_vals) < 2:
+            cmap = 'viridis' if 'Mask' in iname else 'magma'
+            # Lock the scale to 0-1 for masks so 'all 1s' doesn't default to purple
+            vmin, vmax = (0, 1) if ('Mask' in iname or a.dtype == bool) else (None, None)
+
+            im = ax.imshow(a, interpolation='none', cmap=cmap, vmin=vmin, vmax=vmax)
+
+            # Only draw a colorbar if there is at least one valid value to label
+            if len(unique_vals) == 1:
+                divider = make_axes_locatable(ax)
+                cax = divider.append_axes("right", size="5%", pad=0.1)
+                plt.colorbar(im, cax=cax, ticks=unique_vals)
+
+        # 3. Check if this is a discrete/classification map (2 to max_discrete_vals)
+        elif len(unique_vals) <= max_discrete_vals:
+            display_a = np.zeros_like(a, dtype=int)
+            for idx, val in enumerate(unique_vals):
+                display_a[a == val] = idx
+
             discrete_palette = ['black', 'yellow', 'red', 'cyan', 'magenta']
             colors = discrete_palette[:len(unique_vals)]
             cmap = ListedColormap(colors)
 
-            # BoundaryNorm ensures the colorbar is split perfectly by the number of classes
             bounds = np.arange(len(unique_vals) + 1) - 0.5
             norm = BoundaryNorm(bounds, cmap.N)
 
-            # Force interpolation='none' so classes don't blur at the edges
-            im = ax.imshow(a, interpolation='none', cmap=cmap, norm=norm)
+            im = ax.imshow(display_a, interpolation='none', cmap=cmap, norm=norm)
 
             divider = make_axes_locatable(ax)
             cax = divider.append_axes("right", size="5%", pad=0.1)
 
-            # Build a discrete colorbar with centered ticks
             cbar = plt.colorbar(im, cax=cax, ticks=np.arange(len(unique_vals)))
-            # Label the ticks with their actual array values
             cbar.ax.set_yticklabels([f'Val: {v}' for v in unique_vals])
 
         else:
-            # It's continuous! (e.g., raw radiance or log_diff)
+            # 4. Continuous maps
             cmap = 'viridis' if 'Mask' in iname else 'magma'
             im = ax.imshow(a, interpolation=interpolation, cmap=cmap)