Merge pull request #603 from E3SM-Project/diurnal_cycle

Fix and Enhancement for Diurnal Cycle of Precipitation

e3sm_diags/driver/utils/diurnal_cycle.py

@@ -46,12 +46,21 @@ def composite_diurnal_cycle(var, season, fft=True):
     #    tbounds = var_time.getBounds()
     #    var_time[:] = 0.5*(tbounds[:,0]+tbounds[:,1]) #time bounds for h1-h4 are problematic
     var_time_absolute = var_time.asComponentTime()
-    time_freq = int(
-        24 / (var_time_absolute[1].hour - var_time_absolute[0].hour)
-    )  # This only valid for time interval >= 1hour
-    start_time = var_time_absolute[0].hour
-    logger.info(f"start_time {var_time_absolute[0]} {var_time_absolute[0].hour}")
-    logger.info("var_time_freq={}".format(time_freq))
+    # i.e. var_time_absolute[0] = "2000-1-1 1:30:0.0"
+    time_0 = (
+        var_time_absolute[0].hour
+        + var_time_absolute[0].minute / 60
+        + var_time_absolute[0].second / 3600
+    )
+    time_1 = (
+        var_time_absolute[1].hour
+        + var_time_absolute[1].minute / 60
+        + var_time_absolute[1].second / 3600
+    )
+    time_freq = int(24 / (time_1 - time_0))
+    start_time = time_0
+    logger.info(f"start_time {var_time_absolute[0]} {start_time}")
+    logger.info(f"var_time_freq={time_freq}")
 
     # Convert to masked array
     v = var.asma()

e3sm_diags/parameter/diurnal_cycle_parameter.py

@@ -5,6 +5,7 @@ class DiurnalCycleParameter(CoreParameter):
     def __init__(self):
         super(DiurnalCycleParameter, self).__init__()
         self.print_statements = False
-        self.normalize_test_amp = True
+        self.normalize_test_amp = False
         self.ref_timeseries_input = False
         self.test_timeseries_input = False
+        self.normalize_amp_int = 0

e3sm_diags/parser/diurnal_cycle_parser.py

@@ -50,5 +50,15 @@ def load_default_args(self, files=[]):
             "--normalize_test_amp",
             dest="normalize_test_amp",
             help="Normalize test data by maximum diurnal cycle amplitude from reference data",
+            action="store_const",
+            const=True,
+            required=False,
+        )
+
+        self.add_argument(
+            "--normalize_amp_int",
+            dest="normalize_amp_int",
+            type=int,
+            help="Normalize both test and obs by a specified integer",
             required=False,
         )

e3sm_diags/plot/cartopy/diurnal_cycle_plot.py

@@ -67,26 +67,32 @@ def determine_tick_step(degrees_covered):
 def plot_panel(n, fig, proj, var, amp, amp_ref, title, parameter):
 
     normalize_test_amp = parameter.normalize_test_amp
+    specified_max_amp = parameter.normalize_amp_int
+
     lat = var.getLatitude()
     var = ma.squeeze(var.asma())
-    max_amp = amp.max()
-    max_amp_ref = amp_ref.max()
+    max_amp = round(amp.max())
+    max_amp_ref = round(amp_ref.max())
     amp = ma.squeeze(amp.asma())
     amp_ref = ma.squeeze(amp_ref.asma())
 
     if normalize_test_amp:
         img = np.dstack(
-            (
-                (var / 24 - 0.5) % 1,
-                (amp / max_amp * (max_amp / max_amp_ref)) ** 0.5,
-                np.ones_like(amp),
-            )
+            ((var / 24 - 0.5) % 1, (amp / max_amp_ref) ** 0.5, np.ones_like(amp))
         )
         max_amp = max_amp_ref
+        logger.info(
+            f"Scale test diurnal cycle amplitude to max of reference ({max_amp_ref}) mm/day"
+        )
     else:
+        if specified_max_amp != 0:
+            max_amp = specified_max_amp
+
         img = np.dstack(
             ((var / 24 - 0.5) % 1, (amp / max_amp) ** 0.5, np.ones_like(amp))
         )
+        logger.info(f"Scale test diurnal cycle amplitude to specified {max_amp} mm/day")
+    # Note: hsv_to_rgb would clipping input data to the valid range for imshow with RGB data ([0..1]
     img = hsv_to_rgb(img)
 
     # imshow plot
@@ -189,8 +195,8 @@ def plot_panel(n, fig, proj, var, amp, amp_ref, title, parameter):
     bar_ax.set_theta_direction(-1)
     bar_ax.set_theta_offset(np.pi / 2)
     bar_ax.set_xticklabels(["0h", "3h", "6h", "9h", "12h", "15h", "18h", "21h"])
-    bar_ax.set_yticklabels(["", "", "{:.2f}".format(max_amp)])
-    bar_ax.set_rlabel_position(340)
+    bar_ax.set_yticklabels(["", "", f"{int(max_amp)}"])
+    bar_ax.set_rlabel_position(350)
     bar_ax.get_yticklabels()[-2].set_weight("bold")
     # We change the fontsize of minor ticks label
     bar_ax.tick_params(axis="both", labelsize=7, pad=0, length=0)