Update: Gaussian filter (FWHM)

Author: heptaflar

README.md

@@ -91,9 +91,10 @@ will print a detailed usage message to the screen:
                             specify the engine used to write the target product
                             file. (default: None)
       --gaussian-filter {0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5.0}
-                            specify the standard deviation (pixels) of a lateral
-                            Gaussian filter applied to the forecast. If not
-                            specified no filter is applied. (default: None)
+                            specify the full width at half maximum (pixels) of
+                            a lateral Gaussian filter applied to the forecast.
+                            If not specified no filter is applied. (default:
+                            None)
       --horizon {1,2,3,4,5,6,7}
                             specify the forecast horizon (days). (default: None)
       --log-level {debug,info,warning,error,off}

wqf/algorithms/gaussian.py

@@ -22,14 +22,17 @@ class Gaussian(Algorithm):
     """
 
     @override
-    def apply_to(self, cube: da.Array, *, sigma: Number = 1.0) -> da.Array:
+    def apply_to(self, cube: da.Array, *, fwhm: Number = 1.0) -> da.Array:
         """
         Applies a lateral Gaussian filter to a data cube.
 
         :param cube: The data cube.
-        :param sigma: The standard deviation of the Gaussian kernel (pixels).
+        :param fwhm: The full width at half maximum of the Gaussian
+        kernel (pixels).
         :return: The filtered data cube.
         """
+        # noinspection PyTypeChecker
+        sigma = fwhm / 2.35482
         m = da.isnan(cube)
         v = da.where(m, 0.0, cube)
         w = da.where(m, 0.0, 1.0)

wqf/operators/forecastop.py

@@ -67,7 +67,7 @@ def run(self, source: Dataset) -> Dataset:  # noqa: D102
         array = f.apply_to(*feats, names=names)
         if self._args.gaussian_filter is not None:
             g = Gaussian(array.dtype)
-            array = g.apply_to(array, sigma=self._args.gaussian_filter)
+            array = g.apply_to(array, fwhm=self._args.gaussian_filter)
 
         builder = DatasetBuilder()
         for name, value in sorted(

wqf/parser.py

@@ -115,9 +115,9 @@ def _add_options(parser):
         )
         parser.add_argument(
             "--gaussian-filter",
-            help="specify the standard deviation (pixels) of a lateral "
-            "Gaussian filter applied to the forecast. If not specified "
-            "no filter is applied.",
+            help="specify the full width at half maximum (pixels) of "
+            "a lateral Gaussian filter applied to the forecast. "
+            "If not specified no filter is applied.",
             choices=[0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0],
             type=float,
             required=False,

wqf/val/mov.py

@@ -38,7 +38,7 @@ def plot_analysis(chl: DataArray, period: Period):
             continue
         ScenePlot().plot(
             chl[t],
-            title=f"Analysis {time}",
+            title=f"(Re)analysis {time}",
             fn=file,
             cbar_label=r"chlorophyll concentration (mg m$^{-3}$)",
             norm=plc.SymLogNorm(1.0, linscale=0.1, vmin=0.0, vmax=100.0),
@@ -49,14 +49,14 @@ def plot_analysis(chl: DataArray, period: Period):
 
 
 def plot_forecast(
-    chl: DataArray, period: Period, sigma: Number | None, h: int
+    chl: DataArray, period: Period, fwhm: Number | None, h: int
 ):
     chl = period.slice(chl)
 
-    if sigma is not None:
+    if fwhm is not None:
         g = Gaussian(chl.dtype)
         chl = DataArray(
-            data=g.apply_to(chl.data, sigma=sigma),
+            data=g.apply_to(chl.data, fwhm=fwhm),
             coords=chl.coords,
             dims=chl.dims,
             attrs=chl.attrs,
@@ -79,13 +79,13 @@ def plot_forecast(
         )
 
 
-def plot_observed(chl: DataArray, period: Period, sigma: Number | None):
+def plot_observed(chl: DataArray, period: Period, fwhm: Number | None):
     chl = period.slice(chl)
 
-    if sigma is not None:
+    if fwhm is not None:
         g = Gaussian(chl.dtype)
         chl = DataArray(
-            data=g.apply_to(chl.data, sigma=sigma),
+            data=g.apply_to(chl.data, fwhm=fwhm),
             coords=chl.coords,
             dims=chl.dims,
             attrs=chl.attrs,
@@ -114,7 +114,7 @@ def generate_figures(args):
 
     obs = reader.read(args.cube_id, depth_level=3.0)
     ref, pre = XGB(args).predict(obs)
-    plot_observed(ref, period=period, sigma=args.sigma)
+    plot_observed(ref, period=period, fwhm=args.fwhm)
     pre.close()
     ref.close()
     obs.close()
@@ -130,7 +130,7 @@ def generate_figures(args):
     for h in [1, 2, 3, 4, 5, 6, 7]:
         obs = reader.read(args.cube_id, depth_level=3.0)
         ref, pre = XGB(args).predict(obs, h=h)
-        plot_forecast(pre, period=period, sigma=args.sigma, h=h)
+        plot_forecast(pre, period=period, fwhm=args.fwhm, h=h)
         pre.close()
         ref.close()
         obs.close()
@@ -157,29 +157,29 @@ def generate_figures(args):
     )
     parser.add_argument(
         "--analysis",
-        help="plot BGCM analysis",
+        help="plot BGCM (re)analysis",
         action="store_true",
         default=True,
         required=False,
         dest="analysis",
     )
     parser.add_argument(
         "--no-analysis",
-        help="do not plot BGCM analysis",
+        help="do not plot BGCM (re)analysis",
         action="store_false",
         default=False,
         required=False,
         dest="analysis",
     )
     parser.add_argument(
         "--gaussian-filter",
-        help="specify the standard deviation (pixels) of a lateral "
-        "Gaussian filter applied to the forecast. If not specified "
-        "no filter is applied.",
+        help="specify the full width at half maximum (pixels) of "
+        "a lateral Gaussian filter applied to the forecast. "
+        "If not specified no filter is applied.",
         choices=[0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0],
         type=float,
         required=False,
-        dest="sigma",
+        dest="fwhm",
     )
     parser.add_argument(
         "--period-start",