refactor: catalog and forecasts plot methods now uses kwargs instead …

…of plot_args. Kept the `plot_args` for legacy compatibility. docs: Adapted tutorials to the new plot refactoring. Added catalog_plot to Catalog Operations tutorial. ft: Added higher resolutions for basemap autoscale functions, in case extent is lower.
SCECcode · Aug 30, 2024 · c009af4 · c009af4
1 parent 21fcc33
commit c009af4
Show file tree

Hide file tree

Showing 8 changed files with 95 additions and 91 deletions.
diff --git a/csep/core/catalogs.py b/csep/core/catalogs.py
@@ -829,7 +829,8 @@ def b_positive(self):
         """ Implements the b-positive indicator from Nicholas van der Elst """
         pass
 
-    def plot(self, ax=None, show=False, extent=None, set_global=False, plot_args=None):
+    def plot(self, ax=None, show=False, extent=None, set_global=False, plot_args=None,
+             **kwargs):
         """ Plot catalog according to plate-carree projection
 
         Args:
@@ -844,17 +845,7 @@ def plot(self, ax=None, show=False, extent=None, set_global=False, plot_args=Non
 
         # no mutable function arguments
         plot_args_default = {
-             'basemap': 'ESRI_terrain',
-             'markersize': 2,
-             'markercolor': 'red',
-             'alpha': 0.3,
-             'mag_scale': 7,
-             'legend': True,
-             'grid_labels': True,
-             'legend_loc': 3,
-             'figsize': (8, 8),
-             'title': self.name,
-             'mag_ticks': False
+             'basemap': kwargs.pop('basemap', None) or 'ESRI_terrain',
         }
 
         # Plot the region border (if it exists) by default
@@ -869,8 +860,8 @@ def plot(self, ax=None, show=False, extent=None, set_global=False, plot_args=Non
         plot_args_default.update(plot_args)
 
         # this call requires internet connection and basemap
-        ax = plot_catalog(self, ax=ax,show=show, extent=extent,
-                          set_global=set_global, plot_args=plot_args_default)
+        ax = plot_catalog(self, ax=ax, show=show, extent=extent,
+                          set_global=set_global, **plot_args_default, **kwargs)
         return ax
 
 

diff --git a/csep/core/forecasts.py b/csep/core/forecasts.py
@@ -427,7 +427,8 @@ def load_ascii(cls, ascii_fname, start_date=None, end_date=None, name=None, swap
         gds = cls(start_date, end_date, magnitudes=mws, name=name, region=region, data=rates)
         return gds
 
-    def plot(self, ax=None, show=False, log=True, extent=None, set_global=False, plot_args=None):
+    def plot(self, ax=None, show=False, log=True, extent=None, set_global=False, plot_args=None,
+             **kwargs):
         """ Plot gridded forecast according to plate-carree projection
 
         Args:
@@ -446,19 +447,23 @@ def plot(self, ax=None, show=False, log=True, extent=None, set_global=False, plo
             time = f'{round(end-start,3)} years'
 
         plot_args = plot_args or {}
-        plot_args.setdefault('figsize', (10, 10))
-        plot_args.setdefault('title', self.name)
-
+        plot_args.update({
+             'basemap': kwargs.pop('basemap', None) or 'ESRI_terrain',
+             'title': kwargs.pop('title', None) or self.name,
+             'figsize': kwargs.pop('figsize', None) or (8, 8),
+        })
+        plot_args.update(**kwargs)
         # this call requires internet connection and basemap
         if log:
             plot_args.setdefault('clabel', f'log10 M{self.min_magnitude}+ rate per cell per {time}')
             with numpy.errstate(divide='ignore'):
                 ax = plot_gridded_dataset(numpy.log10(self.spatial_counts(cartesian=True)), self.region, ax=ax,
-                                          show=show, extent=extent, set_global=set_global, plot_args=plot_args)
+                                          show=show, extent=extent, set_global=set_global,
+                                          **plot_args)
         else:
             plot_args.setdefault('clabel', f'M{self.min_magnitude}+ rate per cell per {time}')
             ax = plot_gridded_dataset(self.spatial_counts(cartesian=True), self.region, ax=ax, show=show, extent=extent,
-                                      set_global=set_global, plot_args=plot_args)
+                                      set_global=set_global, **plot_args)
         return ax
 
 

diff --git a/csep/utils/plots.py b/csep/utils/plots.py
@@ -559,7 +559,7 @@ def plot_basemap(
     ax = ax or _create_geo_axes(plot_args["figsize"], extent, projection, set_global)
 
     line_autoscaler = cartopy.feature.AdaptiveScaler("110m", (("50m", 50), ("10m", 5)))
-    tile_autoscaler = cartopy.feature.AdaptiveScaler(5, ((6, 50), (7, 15)))
+    tile_autoscaler = cartopy.feature.AdaptiveScaler(5, ((6, 50), (7, 15), (8, 5), (9, 1)))
     tile_depth = (
         4
         if set_global
@@ -1148,9 +1148,11 @@ def plot_comparison_test(
             - **capsize** (`float`): Size of the caps on the error bars.
             - **markersize** (`float`): Size of the markers.
             - **xlabel_fontsize** (`int`): Font size for the X-axis labels.
-            - **ylabel** (`str`): Label for the Y-axis. Defaults to `'Information gain per earthquake'`.
+            - **ylabel** (`str`): Label for the Y-axis. Defaults to
+              `'Information gain per earthquake'`.
             - **ylabel_fontsize** (`int`): Font size for the Y-axis label.
-            - **title** (`str`): Title of the plot. Defaults to the name of the first T-Test result.
+            - **title** (`str`): Title of the plot. Defaults to the name of the first T-Test
+              result.
             - **ylim** (`tuple`): Limits for the Y-axis.
             - **grid** (`bool`): Whether to display grid lines. Defaults to `True`.
             - **hbars** (`bool`): Whether to include horizontal bars for visual separation.
@@ -1682,8 +1684,8 @@ def plot_Molchan_diagram(
     6. Test each ordered and normalized forecasted rate defined in (3) as a threshold value to
        obtain the corresponding contingency table.
     7. Define the "nu" (Miss rate) and "tau" (Fraction of spatial alarmed cells) for each
-       threshold using the information provided by the corresponding contingency table defined in
-       (6).
+       threshold using the information provided by the corresponding contingency table defined
+       in (6).
 
     Note:
     1. The testing catalog and forecast should have exactly the same time-window (duration).
@@ -1696,7 +1698,8 @@ def plot_Molchan_diagram(
         catalog (CSEPCatalog): The evaluation catalog.
         linear (bool): If True, a linear x-axis is used; if False, a logarithmic x-axis is used.
             Defaults to `True`.
-        plot_uniform (bool): If True, include a uniform forecast on the plot. Defaults to `True`.
+        plot_uniform (bool): If True, include a uniform forecast on the plot. Defaults to
+            `True`.
         show (bool): If True, displays the plot. Defaults to `True`.
         ax (matplotlib.axes.Axes): Axes object to plot on. If `None`, creates a new figure.
             Defaults to `None`.

diff --git a/examples/tutorials/catalog_filtering.py b/examples/tutorials/catalog_filtering.py
@@ -92,6 +92,14 @@
 print(catalog)
 
 
+####################################################################################################################################
+# Plot catalog
+# -------------
+#
+# To visualize the catalog, use :meth:`csep.core.catalogs.AbstractBaseCatalog.plot` to plot it spatially
+catalog.plot(show=True)
+
+
 ####################################################################################################################################
 # Write catalog
 # -------------

diff --git a/examples/tutorials/catalog_forecast_evaluation.py b/examples/tutorials/catalog_forecast_evaluation.py
@@ -97,7 +97,7 @@
 print(comcat_catalog)
 
 # Plot the catalog
-comcat_catalog.plot()
+comcat_catalog.plot(show=True)
 
 ####################################################################################################################################
 # Perform number test

diff --git a/examples/tutorials/gridded_forecast_evaluation.py b/examples/tutorials/gridded_forecast_evaluation.py
@@ -103,7 +103,7 @@
 # consistency tests.
 
 ax = plots.plot_consistency_test(spatial_test_result,
-                                 plot_args={'xlabel': 'Spatial likelihood'})
+                                 xlabel='Spatial likelihood')
 plt.show()
 
 ####################################################################################################################################

diff --git a/examples/tutorials/plot_customizations.py b/examples/tutorials/plot_customizations.py
@@ -34,39 +34,35 @@
 forecast = csep.load_gridded_forecast(datasets.hires_ssm_italy_fname,
                                       name='Werner, et al (2010) Italy')
 ####################################################################################################################################
-# **Selecting plotting arguments**
-#
-# Create a dictionary containing the plot arguments
-args_dict = {'title': 'Italy 10 year forecast',
-             'grid_labels': True,
-             'borders': True,
-             'feature_lw': 0.5,
-             'basemap': 'ESRI_imagery',
-             'cmap': 'rainbow',
-             'alpha_exp': 0.8,
-             'projection': cartopy.crs.Mercator()}
+# **Plotting the dataset with fine-tuned arguments**
+
 ####################################################################################################################################
 # These arguments are, in order:
 #
+# * Set an extent
+# * Select ESRI Imagery as a basemap.
 # * Assign a title
 # * Set labels to the geographic axes
 # * Draw country borders
-# * Set a linewidth of 0.5 to country borders
-# * Select ESRI Imagery as a basemap.
-# * Assign ``'rainbow'`` as colormap. Possible values from from ``matplotlib.cm`` library
+# * Set a linewidth of 0.5 to country border
+# * Assign ``'rainbow'`` as colormap. Possible values from ``matplotlib.cm`` library
 # * Defines 0.8 for an exponential transparency function (default is 0 for constant alpha, whereas 1 for linear).
 # * An object cartopy.crs.Projection() is passed as Projection to the map
 #
-# The complete description of plot arguments can be found in :func:`csep.utils.plots.plot_spatial_dataset`
+# The complete description of plot arguments can be found in :func:`csep.utils.plots.plot_gridded_dataset`
 
 ####################################################################################################################################
-# **Plotting the dataset**
-#
-# The map `extent` can be defined. Otherwise, the extent of the data would be used. The dictionary defined must be passed as argument
 
 ax = forecast.plot(extent=[3, 22, 35, 48],
-                   show=True,
-                   plot_args=args_dict)
+                   basemap='ESRI_imagery',
+                   title='Italy 10 year forecast',
+                   grid_labels=True,
+                   borders=True,
+                   borders_linewidth=1.5,
+                   cmap='rainbow',
+                   alpha_exp=0.8,
+                   projection=cartopy.crs.Mercator(),
+                   show=True)
 
 ####################################################################################################################################
 # Example 2: Plot a global forecast and a selected magnitude bin range
@@ -101,27 +97,23 @@
 
 rate_sum = forecast.region.get_cartesian(rate_sum)
 
-####################################################################################################################################
-# **Define plot arguments**
-#
-# We define the arguments and a global projection, centered at $lon=-180$
-
-plot_args = {'figsize': (10,6), 'coastline':True, 'feature_color':'black',
-             'projection': cartopy.crs.Robinson(central_longitude=180.0),
-             'title': forecast.name, 'grid_labels': False,
-             'cmap': 'magma',
-             'clabel': r'$\log_{10}\lambda\left(M_w \in [{%.2f},\,{%.2f}]\right)$ per '
-                       r'${%.1f}^\circ\times {%.1f}^\circ $ per forecast period' %
-                       (low_bound, upper_bound, forecast.region.dh, forecast.region.dh)}
-
 ####################################################################################################################################
 # **Plotting the dataset**
 # To plot a global forecast, we must assign the option ``set_global=True``, which is required by :ref:cartopy to handle
-# internally the extent of the plot
+# internally the extent of the plot. We can further customize the plot using the required arguments from :func:`csep.utils.plots.plot_gridded_dataset`
 
 ax = plots.plot_gridded_dataset(numpy.log10(rate_sum), forecast.region,
-                                show=True, set_global=True,
-                                plot_args=plot_args)
+                                figsize=(10,6),
+                                set_global=True,
+                                coastline_color='black',
+                                projection=cartopy.crs.Robinson(central_longitude=180.0),
+                                title=forecast.name,
+                                grid_labels=False,
+                                colormap='magma',
+                                clabel= r'$\log_{10}\lambda\left(M_w \in [{%.2f},\,{%.2f}]\right)$ per ' \
+                                        r'${%.1f}^\circ\times {%.1f}^\circ $ per forecast period' % \
+                                        (low_bound, upper_bound, forecast.region.dh, forecast.region.dh),
+                                show=True)
 
 ####################################################################################################################################
 
@@ -139,33 +131,37 @@
 min_mag = 4.5
 catalog = csep.query_comcat(start_time, end_time, min_magnitude=min_mag, verbose=False)
 
-# **Define plotting arguments**
-plot_args = {'basemap': csep.datasets.basemap_california,
-             'markersize': 2,
-             'markercolor': 'red',
-             'alpha': 0.3,
-             'mag_scale': 7,
-             'legend': True,
-             'legend_loc': 3,
-             'coastline': False,
-             'mag_ticks': [4.0, 5.0, 6.0, 7.0]}
+# **Define plotting arguments*
 
 ####################################################################################################################################
 # These arguments are, in order:
 #
-# * Assign as basemap the ESRI_terrain webservice
+# * Assign as basemap a TIFF file in the same reference as the plot
 # * Set minimum markersize of 2 with red color
 # * Set a 0.3 transparency
 # * mag_scale is used to exponentially scale the size with respect to magnitude. Recommended 1-8
 # * Set legend True and location in 3 (lower-left corner)
 # * Set a list of Magnitude ticks to display in the legend
 #
 # The complete description of plot arguments can be found in :func:`csep.utils.plots.plot_catalog`
+# The arguments can be stored as a dictionary a priori and then unpacked to the function with `**`.
+
+
+plot_args = {'basemap': csep.datasets.basemap_california,
+             'size': 7,
+             'max_size': 500,
+             'markercolor': 'red',
+             'alpha': 0.3,
+             'mag_scale': 8,
+             'legend': True,
+             'legend_loc': 3,
+             'coastline': False,
+             'mag_ticks': [4.0, 5.0, 6.0, 7.0]}
 
 ####################################################################################################################################
 
 # **Plot the catalog**
-ax = catalog.plot(show=False, plot_args=plot_args)
+ax = catalog.plot(show=True, **plot_args)
 
 
 ####################################################################################################################################
@@ -177,24 +173,25 @@
 # :doc:`gridded_forecast_evaluation` for information on calculating and storing evaluation results)
 
 L_results = [csep.load_evaluation_result(i) for i in datasets.l_test_examples]
-args = {'figsize': (6,5),
-        'title': r'$\mathcal{L}-\mathrm{test}$',
-        'title_fontsize': 18,
-        'xlabel': 'Log-likelihood',
-        'xticks_fontsize': 9,
-        'ylabel_fontsize': 9,
-        'linewidth': 0.8,
-        'capsize': 3,
-        'hbars':True,
-        'tight_layout': True}
+plot_args = {'figsize': (6, 5),
+             'title': r'$\mathcal{L}-\mathrm{test}$',
+             'title_fontsize': 18,
+             'xlabel': 'Log-likelihood',
+             'xticks_fontsize': 9,
+             'ylabel_fontsize': 9,
+             'linewidth': 0.8,
+             'capsize': 3,
+             'hbars':True,
+             'tight_layout': True}
 
 ####################################################################################################################################
 # Description of plot arguments can be found in :func:`plot_poisson_consistency_test`.
 # We set ``one_sided_lower=True`` as usual for an L-test, where the model is rejected if the observed
 # is located within the lower tail of the simulated distribution.
-ax = plots.plot_consistency_test(L_results, one_sided_lower=True, plot_args=args)
+ax = plots.plot_consistency_test(L_results,
+                                 one_sided_lower=True,
+                                 show=True,
+                                 **plot_args)
 
-# Needed to show plots if running as script
-plt.show()
 
 
diff --git a/examples/tutorials/quadtree_gridded_forecast_evaluation.py b/examples/tutorials/quadtree_gridded_forecast_evaluation.py
@@ -193,8 +193,8 @@
 
 stest_result = [spatial_test_single_res_result, spatial_test_multi_res_result]
 ax_spatial = plots.plot_consistency_test(stest_result,
-                                        plot_args={'xlabel': 'Spatial likelihood'})
+                                         xlabel='Spatial likelihood')
 
 ntest_result = [number_test_single_res_result, number_test_multi_res_result]
 ax_number = plots.plot_consistency_test(ntest_result,
-                                        plot_args={'xlabel': 'Number of Earthquakes'})
+                                        xlabel='Number of Earthquakes')