Merge pull request #195 from braden6521/sofast_fringe_phase_unwrap_debug

Added Sofast Fringe debugging plots

opencsp/app/sofast/lib/DebugOpticsGeometry.py

@@ -6,4 +6,6 @@ class DebugOpticsGeometry:
 
     def __init__(self):
         self.debug_active: bool = False
+        """To activate geometry debugging. Default False"""
         self.figures: list = []
+        """List to hold figure objects once created."""

opencsp/app/sofast/lib/ParamsMaskCalculation.py

@@ -11,15 +11,15 @@ class ParamsMaskCalculation(hdf5_tools.HDF5_IO_Abstract):
     """Defines threshold to use when calculating optic mask. Uses a histogram of pixel values
     of the mask difference image (light image - dark image). This is the fraction of the way
     from the first histogram peak (most common dark pixel value) to the the last histogram peak
-    (most common light pixel value)."""
+    (most common light pixel value). (Default 0.5)"""
     filt_width: int = 9
-    """Side length of square kernel used to filter mask image"""
+    """Side length of square kernel used to filter mask image. (Default 9)"""
     filt_thresh: int = 4
-    """Threshold (minimum number of active pixels) to use when removing small active mask areas."""
+    """Threshold (minimum number of active pixels) to use when removing small active mask areas. (Default 4)"""
     thresh_active_pixels: float = 0.05
-    """If number of active mask pixels is below this fraction of total image pixels, thow error."""
+    """If number of active mask pixels is below this fraction of total image pixels, throw error. (Default 0.05)"""
     keep_largest_area: bool = False
-    """Flag to apply processing step that keeps only the largest mask area"""
+    """Flag to apply processing step that keeps only the largest mask area. (Default True)"""
 
     def save_to_hdf(self, file: str, prefix: str = ''):
         """Saves data to given HDF5 file. Data is stored in PREFIX + ParamsMaskCalculation/...

opencsp/app/sofast/lib/ParamsOpticGeometry.py

@@ -9,18 +9,18 @@ class ParamsOpticGeometry(hdf5_tools.HDF5_IO_Abstract):
 
     perimeter_refine_axial_search_dist: float = 50.0
     """The length of the search box (along the search direction) to use when finding optic
-    perimeter. Units pixels. Default 50.0"""
+    perimeter. Units pixels. (Default 50.0)"""
     perimeter_refine_perpendicular_search_dist: float = 50.0
     """The half-width of the search box (perpendicular to the search direction) to use when finding
-    optic perimeter. Units pixels. Default 50.0"""
+    optic perimeter. Units pixels. (Default 50.0)"""
     facet_corns_refine_step_length: float = 10.0
     """The length of the search box (along the search direction) to use when refining facet corner
-    locations (when processing a facet ensemble). Units pixels. Default 10.0"""
+    locations (when processing a facet ensemble). Units pixels. (Default 10.0)"""
     facet_corns_refine_perpendicular_search_dist: float = 10.0
     """The half-width of the search box (perpendicular to the search direction) to use when
-    refining facet corner locations (when processing a facet ensemble). Units pixels. Default 10.0"""
+    refining facet corner locations (when processing a facet ensemble). Units pixels. (Default 10.0)"""
     facet_corns_refine_frac_keep: float = 0.5
-    """The fraction of pixels to consider within search box when finding optic edges. Default 0.5"""
+    """The fraction of pixels to consider within search box when finding optic edges. (Default 0.5)"""
 
     def save_to_hdf(self, file: str, prefix: str = ''):
         """Saves data to given HDF5 file. Data is stored in PREFIX + ParamsOpticGeometry/...

opencsp/app/sofast/lib/ProcessSofastFixed.py

@@ -105,15 +105,6 @@ def _calculate_mask(self) -> ndarray:
         ]
         mask = ip.calc_mask_raw(images, *params)
 
-        if (self.optic_type == 'multi') and self.params.mask.keep_largest_area:
-            lt.warn(
-                '"keep_largest_area" mask processing option cannot be used '
-                'for multifacet ensembles. This will be turned off.'
-            )
-            self.params.mask.keep_largest_area = False
-        elif self.params.mask.keep_largest_area:
-            mask = ip.keep_largest_mask_area(mask)
-
         return mask
 
     def load_measurement_data(self, measurement: MeasurementSofastFixed) -> None:
@@ -266,6 +257,10 @@ def process_single_facet_optic(
         # Calculate mask
         mask_raw = self._calculate_mask()
 
+        # If enabled, fill in holes in mask area
+        if self.params.mask.keep_largest_area:
+            mask_raw = ip.keep_largest_mask_area(mask_raw)
+
         # Process optic geometry (find mask corners, etc.)
         (
             self.data_geometry_general,
@@ -357,6 +352,7 @@ def process_multi_facet_optic(
             self.camera,
             self.measurement.dist_optic_screen,
             self.params.geometry,
+            self.params.mask,
             self.params.debug_geometry,
         )
 

opencsp/app/sofast/lib/ProcessSofastFringe.py

@@ -2,6 +2,8 @@
 to calculate surface slopes.
 """
 
+import matplotlib.pyplot as plt
+from matplotlib import colormaps
 import numpy as np
 
 from opencsp.app.sofast.lib.DefinitionEnsemble import DefinitionEnsemble
@@ -377,13 +379,6 @@ def _process_optic_multifacet_geometry(
         ]
         mask_raw = ip.calc_mask_raw(self.measurement.mask_images, *params)
 
-        if self.params.mask.keep_largest_area:
-            lt.warn(
-                '"keep_largest_area" mask processing option cannot be used '
-                'for multifacet ensembles. This will be turned off.'
-            )
-            self.params.mask.keep_largest_area = False
-
         (
             self.data_geometry_general,
             self.data_image_processing_general,
@@ -399,6 +394,7 @@ def _process_optic_multifacet_geometry(
             self.camera,
             self.measurement.dist_optic_screen,
             self.params.geometry,
+            self.params.mask,
             self.params.debug_geometry,
         )
 
@@ -417,6 +413,17 @@ def _process_display(self) -> None:
         x_periods = self.measurement.fringe_periods_x
         y_periods = self.measurement.fringe_periods_y
 
+        # Prepare for plotting unwrapped phase images
+        if self.params.debug_geometry.debug_active:
+            # X phase RGB image
+            im_phase_x = self.measurement.mask_images[..., 1].copy()
+            im_phase_x = np.stack((im_phase_x,) * 3, 2)
+            im_phase_x = im_phase_x / im_phase_x.max()
+            # Y phase RGB image
+            im_phase_y = self.measurement.mask_images[..., 1].copy()
+            im_phase_y = np.stack((im_phase_y,) * 3, 2)
+            im_phase_y = im_phase_y / im_phase_y.max()
+
         for idx_facet in range(self.num_facets):
             # Get current processed mask layer
             mask_processed = self.data_image_processing_facet[idx_facet].mask_processed
@@ -431,6 +438,25 @@ def _process_display(self) -> None:
             v_screen_points_fractional_screens = Vxy((screen_xs, screen_ys))
             self.data_geometry_facet[idx_facet].v_screen_points_fractional_screens = v_screen_points_fractional_screens
 
+            # Create plot of unwrapped phase (if enabled)
+            if self.params.debug_geometry.debug_active:
+                # Add active pixels as colored pixels
+                cm = colormaps.get_cmap('jet')
+                vals_x_jet = cm(screen_xs)[:, :3]  # remove alpha channel
+                im_phase_x[mask_processed, :] = vals_x_jet
+                vals_y_jet = cm(screen_ys)[:, :3]  # remove alpha channel
+                im_phase_y[mask_processed, :] = vals_y_jet
+                # Plot x image
+                fig = plt.figure(f'ProcessSofastFringe_unwrapped_phase_x_facet_{idx_facet:d}')
+                plt.imshow(im_phase_x)
+                plt.title(f'Unwrapped X Phase Facet {idx_facet:d}')
+                self.params.debug_geometry.figures.append(fig)
+                # Plot y image
+                fig = plt.figure(f'ProcessSofastFringe_unwrapped_phase_y_facet_{idx_facet:d}')
+                plt.imshow(im_phase_y)
+                plt.title(f'Unwrapped Y Phase Facet {idx_facet:d}')
+                self.params.debug_geometry.figures.append(fig)
+
             # Undistort screen points (display coordinates)
             v_screen_points_screen = self.display.interp_func(
                 v_screen_points_fractional_screens
@@ -442,8 +468,9 @@ def _process_display(self) -> None:
             mask_bad_pixels = np.zeros(mask_processed.shape, dtype=bool)
             if np.any(nan_mask):
                 lt.warn(
-                    f'{nan_mask.sum():d} / {nan_mask.size:d} points are NANs in calculated '
-                    'screen points for facet {idx_facet:d}. These data points will be removed.'
+                    'ProcessSofastFringe._process_display(): '
+                    f'{nan_mask.sum():d} / {nan_mask.size:d} screen points are undefined '
+                    f'for facet {idx_facet:d}. These data points will be removed.'
                 )
                 # Make mask of NANs
                 mask_bad_pixels[mask_processed] = nan_mask

opencsp/app/sofast/lib/process_optics_geometry.py

@@ -12,6 +12,7 @@
 from opencsp.app.sofast.lib.DefinitionEnsemble import DefinitionEnsemble
 from opencsp.app.sofast.lib.DefinitionFacet import DefinitionFacet
 from opencsp.app.sofast.lib.ParamsOpticGeometry import ParamsOpticGeometry
+from opencsp.app.sofast.lib.ParamsMaskCalculation import ParamsMaskCalculation
 from opencsp.app.sofast.lib.DebugOpticsGeometry import DebugOpticsGeometry
 import opencsp.app.sofast.lib.image_processing as ip
 from opencsp.app.sofast.lib.SpatialOrientation import SpatialOrientation
@@ -366,7 +367,8 @@ def process_multifacet_geometry(
     orientation: SpatialOrientation,
     camera: Camera,
     dist_optic_screen: float,
-    params: ParamsOpticGeometry = ParamsOpticGeometry(),
+    params_geometry: ParamsOpticGeometry = ParamsOpticGeometry(),
+    params_mask: ParamsMaskCalculation = ParamsMaskCalculation(),
     debug: DebugOpticsGeometry = DebugOpticsGeometry(),
 ) -> tuple[
     cdc.CalculationDataGeometryGeneral,
@@ -393,8 +395,10 @@ def process_multifacet_geometry(
         Camera object
     dist_optic_screen : float
         Optic to screen distance, meters
-    params : ParamsOpticGeometry, optional
+    params_geometry : ParamsOpticGeometry, optional
         ParamsOpticGeometry object, by default ParamsOpticGeometry()
+    params_mask : ParamsMaskCalculation, optional
+        ParamsMaskCalculation object, by default ParamsMaskCalculation()
     debug : DebugOpticsGeometry, optional
         DebugOpticsGeometry object, by default DebugOpticsGeometry()
 
@@ -512,7 +516,10 @@ def process_multifacet_geometry(
         plt.title('Expected Perimeter Points')
 
     # Refine perimeter points
-    args = [params.perimeter_refine_axial_search_dist, params.perimeter_refine_perpendicular_search_dist]
+    args = [
+        params_geometry.perimeter_refine_axial_search_dist,
+        params_geometry.perimeter_refine_perpendicular_search_dist,
+    ]
     loop_ensemble_image_refine = ip.refine_mask_perimeter(loop_ensemble_exp, v_edges_image, *args)
     data_image_processing_general.loop_optic_image_refine = loop_ensemble_image_refine
 
@@ -545,9 +552,9 @@ def process_multifacet_geometry(
 
     # Refine facet corners
     args = [
-        params.facet_corns_refine_step_length,
-        params.facet_corns_refine_perpendicular_search_dist,
-        params.facet_corns_refine_frac_keep,
+        params_geometry.facet_corns_refine_step_length,
+        params_geometry.facet_corns_refine_perpendicular_search_dist,
+        params_geometry.facet_corns_refine_frac_keep,
     ]
     loops_facets_refined: list[LoopXY] = []
     for idx in range(num_facets):
@@ -583,7 +590,11 @@ def process_multifacet_geometry(
     mask_processed *= mask_raw[..., np.newaxis]
     mask_processed = np.logical_and(mask_processed, mask_fitted)
     for idx in range(num_facets):
-        data_image_processing_facet[idx].mask_processed = mask_processed[..., idx]
+        mask = mask_processed[..., idx]
+        # If enabled, keep largest mask area (fill holes) for each individual facet
+        if params_mask.keep_largest_area:
+            mask = ip.keep_largest_mask_area(mask)
+        data_image_processing_facet[idx].mask_processed = mask
 
     # Refine R/T with all refined facet corners
     r_ensemble_cam_refine_2, v_cam_ensemble_cam_refine_2 = sp.calc_rt_from_img_pts(

opencsp/common/lib/deflectometry/SlopeSolverDataDebug.py

@@ -8,9 +8,22 @@ class SlopeSolverDataDebug:
 
     def __init__(self):
         self.debug_active: bool = False
+        """To activate slope solver debugging. (Default False)"""
         self.optic_data: Any = None
+        """Representation of optic (Facet/Mirror) being solved for.
+        The geometry data in this object is used to create visualization plots.
+        This information is updated automatically during SOFAST execution and will
+        overwrite any previously user-given values. (Default None)"""
         self.slope_solver_figures: list = []
+        """List to hold figure objects once created."""
         self.slope_solver_camera_rays_length: float = 0.0
+        """The length (meters) of camera rays to draw when plotting the 3d slope solving scenario plot. (Default 0.0)"""
         self.slope_solver_plot_camera_screen_points: bool = False
+        """To include scatter plot of xyz screen point locations seen by camera in slope solving scenario plot. (Default False)"""
         self.slope_solver_point_downsample: int = 50
+        """The downsample factor (to save computing resources) to apply to screen points
+        Only applicable if plotting screen points is enabled with the
+        `SlopeSolverDataDebug.slope_solver_plot_camera_screen_points` flag). (Default 50)"""
         self.slope_solver_single_plot: bool = False
+        """Flag to plot all iterations of the slope solving algorithm on one plot (True) or create a separate
+        plot for each iteration (False). Default False (new plot for each iteration)"""