Merge pull request #19 from luhipi/enhancement/colormaps

Add scientific colormaps.

HISTORY.rst

@@ -2,6 +2,11 @@
 History
 =======
 
+1.1.0 (2024-11-06)
+------------------
+
+* Use Scientific colour maps from Crameri.
+
 1.0.0 (2024-08-12) Strawberry Pie
 ---------------------------------
 

sarvey/objects.py

@@ -38,6 +38,7 @@
 from pyproj.aoi import AreaOfInterest
 from pyproj.database import query_utm_crs_info
 from logging import Logger
+from cmcrameri import cm
 
 from miaplpy.objects.slcStack import slcStack
 from mintpy.utils import readfile
@@ -130,7 +131,7 @@ def plot(self, *, ax: plt.Axes = None, logger: Logger):
         if ax is None:
             fig = plt.figure(figsize=(15, 5))
             ax = fig.add_subplot()
-        ax.imshow(self.background_map, cmap=plt.cm.get_cmap("gray"))
+        ax.imshow(self.background_map, cmap=cm.grayC)
         meta = {"ORBIT_DIRECTION": self.orbit_direction}
         auto_flip_direction(meta, ax=ax, print_msg=False)
 

sarvey/preparation.py

@@ -198,7 +198,7 @@ def selectPixels(*, path: str, selection_method: str, thrsh: float,
         cand_mask = quality >= thrsh
         grid_min_val = False
         unit = "Temporal\nCoherence [ ]"
-        cmap = "autumn"
+        cmap = "lajolla"
 
     if selection_method == "miaplpy":
         raise NotImplementedError("This part is not developed yet. MiaplPy data is read in another way.")
@@ -208,7 +208,7 @@ def selectPixels(*, path: str, selection_method: str, thrsh: float,
         # quality = pl_coherence
         # grid_min_val = False
         # unit = "Phase-Linking\nCoherence [ ]"
-        # cmap = "autumn"
+        # cmap = "lajolla"
 
     if grid_size is not None:  # -> sparse pixel selection
         coord_utm_obj = CoordinatesUTM(file_path=join(path, "coordinates_utm.h5"), logger=logger)

sarvey/processing.py

@@ -30,9 +30,9 @@
 """Processing module for SARvey."""
 from os.path import join, exists
 import matplotlib.pyplot as plt
-from matplotlib import colormaps
 import numpy as np
 from logging import Logger
+import cmcrameri as cmc
 
 from miaplpy.objects.slcStack import slcStack
 from mintpy.utils import readfile
@@ -207,7 +207,7 @@ def runPreparation(self):
 
         fig = plt.figure(figsize=(15, 5))
         ax = fig.add_subplot()
-        im = ax.imshow(temp_coh, cmap=colormaps["gray"], vmin=0, vmax=1)
+        im = ax.imshow(temp_coh, cmap=cmc.cm.cmaps["grayC"], vmin=0, vmax=1)
         auto_flip_direction(slc_stack_obj.metadata, ax=ax, print_msg=True)
         ax.set_xlabel("Range")
         ax.set_ylabel("Azimuth")
@@ -243,11 +243,11 @@ def runConsistencyCheck(self):
 
         fig = plt.figure(figsize=(15, 5))
         ax = fig.add_subplot()
-        ax.imshow(mask_valid_area, cmap=plt.cm.get_cmap("gray"), alpha=0.5, zorder=10, vmin=0, vmax=1)
+        ax.imshow(mask_valid_area, cmap=cmc.cm.cmaps["grayC"], alpha=0.5, zorder=10, vmin=0, vmax=1)
         bmap_obj.plot(ax=ax, logger=self.logger)
         coord_xy = np.array(np.where(cand_mask1)).transpose()
         val = np.ones_like(cand_mask1)
-        sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=val[cand_mask1], s=0.5, cmap=plt.get_cmap("autumn_r"),
+        sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=val[cand_mask1], s=0.5, cmap=cmc.cm.cmaps["lajolla_r"],
                         vmin=1, vmax=2)  # set min, max to ensure that points are yellow
         cbar = plt.colorbar(sc, pad=0.03, shrink=0.5)
         cbar.ax.set_visible(False)  # make size of axis consistent with all others
@@ -322,7 +322,7 @@ def runConsistencyCheck(self):
             ax, cbar = viewer.plotColoredPointNetwork(x=point_obj.coord_xy[:, 1], y=point_obj.coord_xy[:, 0],
                                                       arcs=net_par_obj.arcs[arc_mask, :],
                                                       val=net_par_obj.gamma[arc_mask],
-                                                      ax=ax, linewidth=1, cmap_name="autumn", clim=(0, 1))
+                                                      ax=ax, linewidth=1, cmap="lajolla", clim=(0, 1))
             ax.set_title("Coherence from temporal unwrapping\n"
                          r"(only arcs with $\gamma \leq$ {} "
                          "shown)\nBefore outlier removal".format(thrsh_visualisation))
@@ -347,7 +347,7 @@ def runConsistencyCheck(self):
             ax, cbar = viewer.plotColoredPointNetwork(x=coord_xy[:, 1], y=coord_xy[:, 0],
                                                       arcs=net_par_obj.arcs[arc_mask, :],
                                                       val=net_par_obj.gamma[arc_mask],
-                                                      ax=ax, linewidth=1, cmap_name="autumn", clim=(0, 1))
+                                                      ax=ax, linewidth=1, cmap="lajolla", clim=(0, 1))
             ax.set_title("Coherence from temporal unwrapping\n"
                          r"(only arcs with $\gamma \leq$ {} "
                          "shown)\nAfter outlier removal".format(thrsh_visualisation))
@@ -403,7 +403,8 @@ def runUnwrappingTimeAndSpace(self):
         #                                               max_rm_fraction=0.001)
         fig = viewer.plotScatter(value=-demerr, coord=point_obj.coord_xy,
                                  ttl="Parameter integration: DEM correction in [m]",
-                                 bmap_obj=bmap_obj, s=3.5, cmap="jet_r", symmetric=True, logger=self.logger)[0]
+                                 bmap_obj=bmap_obj, s=3.5, cmap="vanimo", symmetric=True,
+                                 logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_2_estimation_dem_correction.png"), dpi=300)
         plt.close(fig)
 
@@ -422,7 +423,8 @@ def runUnwrappingTimeAndSpace(self):
         #                                            max_rm_fraction=0.001)
         fig = viewer.plotScatter(value=-vel, coord=point_obj.coord_xy,
                                  ttl="Parameter integration: mean velocity in [m / year]",
-                                 bmap_obj=bmap_obj, s=3.5, cmap="jet_r", symmetric=True, logger=self.logger)[0]
+                                 bmap_obj=bmap_obj, s=3.5, cmap="roma", symmetric=True,
+                                 logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_2_estimation_velocity.png"), dpi=300)
         plt.close(fig)
 
@@ -509,7 +511,7 @@ def runUnwrappingSpace(self):
                                                   y=point_obj.coord_xy[:, 0],
                                                   arcs=arcs,
                                                   val=np.zeros(arcs.shape[0], dtype=np.float32),
-                                                  ax=ax, linewidth=0.5, cmap_name="hot", clim=(0, 1))
+                                                  ax=ax, linewidth=0.5, cmap="lajolla", clim=(0, 1))
         cbar.ax.set_visible(False)
         ax.set_xlabel("Range")
         ax.set_ylabel("Azimuth")
@@ -584,7 +586,7 @@ def runFiltering(self):
         auto_corr_img[~mask] = np.inf
 
         fig = viewer.plotScatter(value=auto_corr, coord=point1_obj.coord_xy, bmap_obj=bmap_obj,
-                                 ttl="Temporal autocorrelation", unit="[ ]", s=3.5, cmap="autumn_r",
+                                 ttl="Temporal autocorrelation", unit="[ ]", s=3.5, cmap="lajolla",
                                  vmin=0, vmax=1, logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_3_temporal_autocorrelation.png"), dpi=300)
         plt.close(fig)
@@ -618,7 +620,7 @@ def runFiltering(self):
         # store plot for quality control during processing
         fig, ax = viewer.plotScatter(value=auto_corr_img[cand_mask_sparse], coord=point1_obj.coord_xy,
                                      bmap_obj=bmap_obj, ttl="Selected pixels for APS estimation",
-                                     unit="Auto-correlation\n[ ]", s=5, cmap="autumn_r", vmin=0, vmax=1,
+                                     unit="Auto-correlation\n[ ]", s=5, cmap="lajolla", vmin=0, vmax=1,
                                      logger=self.logger)[:2]
         viewer.plotGridFromBoxList(box_list=box_list, ax=ax, edgecolor="k", linewidth=0.2)
         fig.savefig(join(self.path, "pic", "step_3_stable_points.png"), dpi=300)
@@ -669,12 +671,12 @@ def runFiltering(self):
 
                 fig = plt.figure(figsize=(15, 5))
                 ax = fig.add_subplot()
-                ax.imshow(mask_pl_aoi, cmap=plt.cm.get_cmap("gray"), alpha=0.5, zorder=10, vmin=0, vmax=1)
+                ax.imshow(mask_pl_aoi, cmap=cmc.cm.cmaps["grayC"], alpha=0.5, zorder=10, vmin=0, vmax=1)
                 bmap_obj.plot(ax=ax, logger=self.logger)
                 coord_xy = np.array(np.where(cand_mask_pl)).transpose()
                 val = np.ones_like(cand_mask_pl)
                 sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=val[cand_mask_pl], s=0.5,
-                                cmap=plt.get_cmap("autumn_r"),
+                                cmap=cmc.cm.cmaps["lajolla_r"],
                                 vmin=1, vmax=2)  # set min, max to ensure that points are yellow
                 cbar = plt.colorbar(sc, pad=0.03, shrink=0.5)
                 cbar.ax.set_visible(False)  # make size of axis consistent with all others
@@ -707,12 +709,12 @@ def runFiltering(self):
 
         fig = plt.figure(figsize=(15, 5))
         ax = fig.add_subplot()
-        ax.imshow(mask_valid_area, cmap=plt.cm.get_cmap("gray"), alpha=0.5, zorder=10, vmin=0, vmax=1)
+        ax.imshow(mask_valid_area, cmap=cmc.cm.cmaps["grayC"], alpha=0.5, zorder=10, vmin=0, vmax=1)
         bmap_obj.plot(ax=ax, logger=self.logger)
         coord_xy = np.array(np.where(cand_mask2)).transpose()
         val = np.ones_like(cand_mask2)
-        sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=val[cand_mask2], s=0.5, cmap=plt.get_cmap("autumn_r"),
-                        vmin=1, vmax=2)  # set min, max to ensure that points are yellow
+        sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=val[cand_mask2], s=0.5, cmap=cmc.cm.cmaps["lajolla_r"],
+                        vmin=0, vmax=10)  # set min, max to ensure that points are yellow
         cbar = plt.colorbar(sc, pad=0.03, shrink=0.5)
         cbar.ax.set_visible(False)  # make size of axis consistent with all others
         plt.tight_layout()
@@ -951,8 +953,8 @@ def runDensificationTimeAndSpace(self):
 
         bmap_obj = AmplitudeImage(file_path=join(self.path, "background_map.h5"))
         fig = viewer.plotScatter(value=gamma, coord=point2_obj.coord_xy, bmap_obj=bmap_obj,
-                                 ttl="Coherence from temporal unwrapping\nBefore outlier removal", s=3.5, cmap="autumn",
-                                 vmin=0, vmax=1, logger=self.logger)[0]
+                                 ttl="Coherence from temporal unwrapping\nBefore outlier removal", s=3.5,
+                                 cmap="lajolla", vmin=0, vmax=1, logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_4_temporal_unwrapping_p2_coh{}.png".format(coh_value)), dpi=300)
         plt.close(fig)
 
@@ -976,20 +978,22 @@ def runDensificationTimeAndSpace(self):
         plt.close(fig)
 
         fig = viewer.plotScatter(value=gamma[mask_gamma], coord=point2_obj.coord_xy, bmap_obj=bmap_obj,
-                                 ttl="Coherence from temporal unwrapping\nAfter outlier removal", s=3.5, cmap="autumn",
-                                 vmin=0, vmax=1, logger=self.logger)[0]
+                                 ttl="Coherence from temporal unwrapping\nAfter outlier removal", s=3.5,
+                                 cmap="lajolla", vmin=0, vmax=1, logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_4_temporal_unwrapping_p2_coh{}_reduced.png".format(coh_value)),
                     dpi=300)
         plt.close(fig)
 
         fig = viewer.plotScatter(value=-vel[mask_gamma], coord=point2_obj.coord_xy,
                                  ttl="Mean velocity in [m / year]",
-                                 bmap_obj=bmap_obj, s=3.5, cmap="jet_r", symmetric=True, logger=self.logger)[0]
+                                 bmap_obj=bmap_obj, s=3.5, cmap="roma", symmetric=True,
+                                 logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_4_estimation_velocity_p2_coh{}.png".format(coh_value)), dpi=300)
         plt.close(fig)
 
-        fig = viewer.plotScatter(value=-demerr[mask_gamma], coord=point2_obj.coord_xy, ttl="DEM error in [m]",
-                                 bmap_obj=bmap_obj, s=3.5, cmap="jet_r", symmetric=True, logger=self.logger)[0]
+        fig = viewer.plotScatter(value=-demerr[mask_gamma], coord=point2_obj.coord_xy, ttl="DEM correction in [m]",
+                                 bmap_obj=bmap_obj, s=3.5, cmap="vanimo", symmetric=True,
+                                 logger=self.logger)[0]
         fig.savefig(join(self.path, "pic", "step_4_estimation_dem_correction_p2_coh{}.png".format(coh_value)), dpi=300)
         plt.close(fig)
 

sarvey/sarvey_plot.py

@@ -34,14 +34,13 @@
 from os.path import join, basename, dirname
 import matplotlib
 import matplotlib.pyplot as plt
-from matplotlib import colormaps
 import numpy as np
 import logging
 from logging import Logger
 import sys
+import cmcrameri as cmc
 
 from mintpy.utils import ptime
-from mintpy.objects.colors import ColormapExt
 from mintpy.utils.plot import auto_flip_direction
 
 from sarvey.ifg_network import IfgNetwork
@@ -98,7 +97,7 @@ def plotMap(*, obj_name: str, save_path: str, interactive: bool = False, input_p
 
     ax = bmap_obj.plot(logger=logger)
     sc = ax.scatter(point_obj.coord_xy[:, 1], point_obj.coord_xy[:, 0], c=demerr, s=scatter_size,
-                    cmap=colormaps["jet_r"])
+                    cmap=cmc.cm.cmaps["vanimo"])
     plt.colorbar(sc, label="[m]", pad=0.03, shrink=0.5)
     plt.title("DEM correction")
     plt.ylabel('Azimuth')
@@ -110,25 +109,11 @@ def plotMap(*, obj_name: str, save_path: str, interactive: bool = False, input_p
     else:
         plt.close(plt.gcf())
 
-    ax = bmap_obj.plot(logger=logger)
-    sc = ax.scatter(point_obj.coord_xy[:, 1], point_obj.coord_xy[:, 0], c=omega, s=scatter_size,
-                    cmap=colormaps["autumn_r"])
-    plt.colorbar(sc, label="", pad=0.03, shrink=0.5)
-    plt.title("Squared sum of residuals")
-    plt.ylabel('Azimuth')
-    plt.xlabel('Range')
-    plt.tight_layout()
-    plt.gcf().savefig(join(save_path, "map_squared_sum_of_residuals.png"), dpi=300)
-    if interactive:
-        plt.show()
-    else:
-        plt.close(plt.gcf())
-
     v_range = np.max(np.abs(vel * 100))
 
     ax = bmap_obj.plot(logger=logger)
     sc = ax.scatter(point_obj.coord_xy[:, 1], point_obj.coord_xy[:, 0], c=vel * 100, s=scatter_size,
-                    cmap=colormaps["jet_r"],
+                    cmap=cmc.cm.cmaps["roma"],
                     vmin=-v_range, vmax=v_range)
     plt.colorbar(sc, label="[cm / year]", pad=0.03, shrink=0.5)
     plt.title("Mean Velocity")
@@ -143,7 +128,7 @@ def plotMap(*, obj_name: str, save_path: str, interactive: bool = False, input_p
 
     ax = bmap_obj.plot(logger=logger)
     sc = ax.scatter(point_obj.coord_xy[:, 1], point_obj.coord_xy[:, 0], c=coherence, vmin=0, vmax=1, s=scatter_size,
-                    cmap=colormaps["autumn"])
+                    cmap=cmc.cm.cmaps["lajolla"])
     plt.colorbar(sc, label="[-]", pad=0.03, shrink=0.5)
     plt.title("Temporal coherence")
     plt.ylabel('Azimuth')
@@ -161,7 +146,7 @@ def plotMap(*, obj_name: str, save_path: str, interactive: bool = False, input_p
 
     ax = bmap_obj.plot(logger=logger)
     sc = ax.scatter(point_obj.coord_xy[:, 1], point_obj.coord_xy[:, 0], c=stc * 100, s=scatter_size,
-                    cmap=colormaps["autumn_r"])
+                    cmap=cmc.cm.cmaps["lajolla"])
     plt.colorbar(sc, label="[cm]", pad=0.03, shrink=0.5)
     plt.title("Spatiotemporal consistency")
     plt.ylabel('Azimuth')
@@ -319,11 +304,11 @@ def plotAllIfgs(*, obj_name: str, save_path: str, interactive: bool = False, log
     start_time = time.time()
     logger.info(msg="plot and save figures of ifgs.")
     for i in range(num_ifgs):
-        fig = plt.figure(figsize=[15, 5])
+        fig = plt.figure(figsize=(15, 5))
         ax = fig.add_subplot()
         ifg = np.angle(ifgs[:, :, i])
         ifg[ifg == 0] = np.nan
-        im = plt.imshow(ifg, cmap=ColormapExt('cmy').colormap, interpolation='nearest', vmin=-np.pi, vmax=np.pi)
+        im = plt.imshow(ifg, cmap=cmc.cm.cmaps["romaO"], interpolation='nearest', vmin=-np.pi, vmax=np.pi)
         auto_flip_direction(ifg_stack_obj.metadata, ax=ax, print_msg=False)
         ax.set_xlabel("Range")
         ax.set_ylabel("Azimuth")
@@ -441,7 +426,7 @@ def main(iargs=None):
 
     config = loadConfiguration(path=config_file_path)
 
-    folder_name = "p1" if "p1" in basename(args.input_file) else basename(args.input_file)[:5]
+    folder_name = "p1" if "p1" in basename(args.input_file) else basename(args.input_file)[:8]
     folder_name = "ifgs" if "ifg_stack" in basename(args.input_file) else folder_name
 
     save_path = join(dirname(args.input_file), "pic", folder_name)

sarvey/triangulation.py

@@ -137,7 +137,7 @@ def triangulateKnn(self, *, k: int):
             idx = tree.query(self.coord_xy[p1, :], k)[1]
             self.adj_mat[p1, idx] = True
             count += 1
-            prog_bar.update(value=count + 1, every=np.int16(num_points / 250),
+            prog_bar.update(value=count + 1, every=np.int16(num_points / (num_points / 5)),
                             suffix='{}/{} points triangulated'.format(count + 1, num_points + 1))
         prog_bar.close()
         m, s = divmod(time.time() - start_time, 60)

sarvey/unwrapping.py

@@ -41,6 +41,7 @@
 from scipy.sparse.linalg import lsqr
 from scipy.optimize import minimize
 from logging import Logger
+import cmcrameri as cmc
 
 from mintpy.utils import ptime
 
@@ -1001,7 +1002,7 @@ def parameterBasedNoisyPointRemoval(*, net_par_obj: NetworkParameter, point_id:
             # vel
             ax = bmap_obj.plot(logger=logger)
             sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=rmse_vel * 1000, s=3.5,
-                            cmap=plt.cm.get_cmap("autumn_r"), vmin=0, vmax=rmse_thrsh * 1000)
+                            cmap=cmc.cm.cmaps["lajolla"], vmin=0, vmax=rmse_thrsh * 1000)
             plt.colorbar(sc, pad=0.03, shrink=0.5)
             ax.set_title("{}. iteration\nmean velocity - RMSE per point in [mm / year]".format(it_count))
             fig = ax.get_figure()
@@ -1013,7 +1014,7 @@ def parameterBasedNoisyPointRemoval(*, net_par_obj: NetworkParameter, point_id:
             # demerr
             ax = bmap_obj.plot(logger=logger)
             sc = ax.scatter(coord_xy[:, 1], coord_xy[:, 0], c=rmse_demerr, s=3.5,
-                            cmap=plt.cm.get_cmap("autumn_r"))
+                            cmap=cmc.cm.cmaps["lajolla"])
             plt.colorbar(sc, pad=0.03, shrink=0.5)
             ax.set_title("{}. iteration\nDEM correction - RMSE per point in [m]".format(it_count))
             fig = ax.get_figure()

sarvey/version.py

@@ -29,6 +29,6 @@
 
 """Version module for SARvey."""
 
-__version__ = '1.0.0'
-__versiondate__ = '2024-08-12_01'
+__version__ = '1.1.0'
+__versiondate__ = '2024-11-06_01'
 __versionalias__ = 'Strawberry Pie'