From cb3589da9ecfed22137937d2c39ada4391bbaa05 Mon Sep 17 00:00:00 2001 From: ShixuanZhang Date: Wed, 4 Sep 2024 01:44:08 -0500 Subject: [PATCH 1/7] Code changes to implement PSA1 and PSA2 mode varibility analysis The Pacific South America mode is represented by an eastward-propagating wave train extending from eastern Australia to Argentina, characterized by two invariant empirical orthogonal function (EOF) patterns with associated principal component (PC) time series (PSA1 and PSA2). In this commit, we made code changes to enable the analysis of PSA1 and PSA2 mode varibility with EOF analysis in PCMDI. Specifically, 1. regions.py was revised to include regions used for EOF analysis for PSA1 and PSA2 modes. Note that the region for PSA1 and PSA2 analysis are the same as the SAM mode 2. argparse_functions.py was revised to include PSA1 and PSA2 in the default mode varibility analysis by PCMDI 3. variability_modes_driver.py was revised to include the default setup for EOF analysis of PAS1 or PSA2 in PCMDI 4. plot_map.py was revised to define the setup of mapping plots to generate the graphics of EOF patterns and teleconnections related to PSA1 and PSA2 --- pcmdi_metrics/io/regions.py | 2 ++ .../variability_mode/lib/argparse_functions.py | 2 ++ .../variability_mode/lib/eof_analysis.py | 16 ++++++++++++++++ pcmdi_metrics/variability_mode/lib/plot_map.py | 2 +- .../variability_mode/variability_modes_driver.py | 12 ++++++++++-- 5 files changed, 31 insertions(+), 3 deletions(-) diff --git a/pcmdi_metrics/io/regions.py b/pcmdi_metrics/io/regions.py index 480e1b5b7..b8e4a9f1f 100755 --- a/pcmdi_metrics/io/regions.py +++ b/pcmdi_metrics/io/regions.py @@ -39,6 +39,8 @@ def load_regions_specs() -> dict: "NAM": {"domain": {"latitude": (20.0, 90), "longitude": (-180, 180)}}, "NAO": {"domain": {"latitude": (20.0, 80), "longitude": (-90, 40)}}, "SAM": {"domain": {"latitude": (-20.0, -90), "longitude": (0, 360)}}, + "PSA1": {"domain": {"latitude": (-20.0, -90), "longitude": (0, 360)}}, + "PSA2": {"domain": {"latitude": (-20.0, -90), "longitude": (0, 360)}}, "PNA": {"domain": {"latitude": (20.0, 85), "longitude": (120, 240)}}, "NPO": {"domain": {"latitude": (20.0, 85), "longitude": (120, 240)}}, "PDO": {"domain": {"latitude": (20.0, 70), "longitude": (110, 260)}}, diff --git a/pcmdi_metrics/variability_mode/lib/argparse_functions.py b/pcmdi_metrics/variability_mode/lib/argparse_functions.py index 8be673cec..9ddeadaf9 100644 --- a/pcmdi_metrics/variability_mode/lib/argparse_functions.py +++ b/pcmdi_metrics/variability_mode/lib/argparse_functions.py @@ -224,6 +224,8 @@ def VariabilityModeCheck(mode, P): "NAM", "NAO", "SAM", + "PSA1", + "PSA2", "PNA", "PDO", "NPO", diff --git a/pcmdi_metrics/variability_mode/lib/eof_analysis.py b/pcmdi_metrics/variability_mode/lib/eof_analysis.py index d25a6f712..86dcd7b18 100644 --- a/pcmdi_metrics/variability_mode/lib/eof_analysis.py +++ b/pcmdi_metrics/variability_mode/lib/eof_analysis.py @@ -155,6 +155,22 @@ def arbitrary_checking(mode, eof_Nth): elif mode == "SAM": if eof_Nth.sel({lat_key: slice(-60, -90)}).mean().item() >= 0: reverse_sign = True + elif mode == "PSA1": + if ( + eof_Nth.sel({lat_key: slice(-59.5, -64.5), lon_key: slice(207.5, 212.5)}) + .mean() + .item() + >= 0 + ): + reverse_sign = True + elif mode == "PSA2": + if ( + eof_Nth.sel({lat_key: slice(-57.5, -62.5), lon_key: slice(277.5, 282.5)}) + .mean() + .item() + >= 0 + ): + reverse_sign = True else: # Minimum sign control part was left behind for any future usage.. if not np.isnan(eof_Nth[-1, -1].item()): if eof_Nth[-1, -1].item() >= 0: diff --git a/pcmdi_metrics/variability_mode/lib/plot_map.py b/pcmdi_metrics/variability_mode/lib/plot_map.py index dd9ac3697..d0b8bcc8a 100644 --- a/pcmdi_metrics/variability_mode/lib/plot_map.py +++ b/pcmdi_metrics/variability_mode/lib/plot_map.py @@ -57,7 +57,7 @@ def plot_map( projection = "Lambert" elif mode in ["NAM"]: projection = "Stereo_north" - elif mode in ["SAM"]: + elif mode in ["SAM", "PSA1", "PSA2"]: projection = "Stereo_south" else: sys.exit("Projection for " + mode + "is not defined.") diff --git a/pcmdi_metrics/variability_mode/variability_modes_driver.py b/pcmdi_metrics/variability_mode/variability_modes_driver.py index e7f3ed769..963bd7dea 100755 --- a/pcmdi_metrics/variability_mode/variability_modes_driver.py +++ b/pcmdi_metrics/variability_mode/variability_modes_driver.py @@ -15,6 +15,10 @@ ## EOF2 based variability modes - NPO: North Pacific Oscillation (2nd EOFs of PNA domain) - NPGO: North Pacific Gyre Oscillation (2nd EOFs of PDO domain) +- PSA2: Pacific South America Mode (2nd EOFs of SAM domain) + +## EOF3 based variability modes +- PSA3: Pacific South America Mode (3rd EOFs of SAM domain) ## Reference: Lee, J., K. Sperber, P. Gleckler, C. Bonfils, and K. Taylor, 2019: @@ -178,8 +182,10 @@ if eofn_obs is None: if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: eofn_obs = 1 - elif mode in ["NPGO", "NPO"]: + elif mode in ["NPGO", "NPO", "PSA1"]: eofn_obs = 2 + elif mode in ["PSA2"]: + eofn_obs = 3 else: raise ValueError(f"{eofn_obs} is not given for {mode}") else: @@ -188,8 +194,10 @@ if eofn_mod is None: if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: eofn_mod = 1 - elif mode in ["NPGO", "NPO"]: + elif mode in ["NPGO", "NPO", "PSA2"]: eofn_mod = 2 + elif mode in ["PSA2"]: + eofn_mod = 3 else: raise ValueError(f"{eofn_mod} is not given for {mode}") else: From 186d7b2f227b14903b0c705aa5a7dfb54ac61dbb Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 12:41:23 -0700 Subject: [PATCH 2/7] add description for PSA1 and PSA2 options --- pcmdi_metrics/variability_mode/lib/argparse_functions.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/pcmdi_metrics/variability_mode/lib/argparse_functions.py b/pcmdi_metrics/variability_mode/lib/argparse_functions.py index 9ddeadaf9..99c56c89c 100644 --- a/pcmdi_metrics/variability_mode/lib/argparse_functions.py +++ b/pcmdi_metrics/variability_mode/lib/argparse_functions.py @@ -33,6 +33,8 @@ def AddParserArgument(P): "- NAM: Northern Annular Mode\n" "- NAO: Northern Atlantic Oscillation\n" "- SAM: Southern Annular Mode\n" + "- PSA1: Pacific–South American pattern 1\n" + "- PSA2: Pacific–South American pattern 2\n" "- PNA: Pacific North American Pattern\n" "- PDO: Pacific Decadal Oscillation\n" "- NPO: North Pacific Oscillation\n" From 38e6a422b9d2a1ae2ca60daaea78c5b26a4fd31c Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 13:09:57 -0700 Subject: [PATCH 3/7] remove pre-defined eofn_obs and eofn_mod from demo param to allow flexibility more simply by giving mode name --- doc/jupyter/Demo/basic_mov_param.py.in | 2 -- doc/jupyter/Demo/basic_mov_param_sst.py.in | 2 -- 2 files changed, 4 deletions(-) diff --git a/doc/jupyter/Demo/basic_mov_param.py.in b/doc/jupyter/Demo/basic_mov_param.py.in index b3ebf1203..f46fadf7b 100644 --- a/doc/jupyter/Demo/basic_mov_param.py.in +++ b/doc/jupyter/Demo/basic_mov_param.py.in @@ -23,7 +23,6 @@ varModel = 'psl' ModUnitsAdjust = (True, 'divide', 100.0) # Pa to hPa msyear = 1900 meyear = 2005 -eofn_mod = 1 # OBSERVATIONS SETTINGS reference_data_path = '$INPUT_DIR$/obs4MIPs_PCMDI_monthly/NOAA-ESRL-PSD/20CR/mon/psl/gn/v20210727/psl_mon_20CR_PCMDI_gn_187101-201212.nc' @@ -32,7 +31,6 @@ varOBS = 'psl' ObsUnitsAdjust = (True, 'divide', 100.0) # Pa to hPa; or (False, 0, 0) osyear = 1900 oeyear = 2005 -eofn_obs = 1 # DIRECTORY WHERE TO PUT RESULTS results_dir = os.path.join( diff --git a/doc/jupyter/Demo/basic_mov_param_sst.py.in b/doc/jupyter/Demo/basic_mov_param_sst.py.in index 7e0bb9eae..8a7b0d65b 100644 --- a/doc/jupyter/Demo/basic_mov_param_sst.py.in +++ b/doc/jupyter/Demo/basic_mov_param_sst.py.in @@ -25,7 +25,6 @@ varModel = 'ts' ModUnitsAdjust = (True, "subtract", 273.15) # degK to degC msyear = 1900 meyear = 2005 -eofn_mod = 1 # OBSERVATIONS SETTINGS reference_data_path = '$INPUT_DIR$/obs4MIPs_PCMDI_monthly/MOHC/HadISST-1-1/mon/ts/gn/v20210727/ts_mon_HadISST-1-1_PCMDI_gn_187001-201907.nc' @@ -34,7 +33,6 @@ varOBS = 'ts' ObsUnitsAdjust = (True, "subtract", 273.15) # degK to degC osyear = 1900 oeyear = 2005 -eofn_obs = 1 # DIRECTORY WHERE TO PUT RESULTS results_dir = os.path.join( From 8b49a64c92388b9cf9d21faab30609511ce22bbd Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 16:14:56 -0700 Subject: [PATCH 4/7] update eof number assignment to improve flexibility --- .../variability_modes_driver.py | 43 ++++++++++--------- 1 file changed, 23 insertions(+), 20 deletions(-) diff --git a/pcmdi_metrics/variability_mode/variability_modes_driver.py b/pcmdi_metrics/variability_mode/variability_modes_driver.py index 963bd7dea..7344683b7 100755 --- a/pcmdi_metrics/variability_mode/variability_modes_driver.py +++ b/pcmdi_metrics/variability_mode/variability_modes_driver.py @@ -176,32 +176,35 @@ print("realization: ", realization) # EOF ordinal number -eofn_obs = param.eofn_obs -eofn_mod = param.eofn_mod +eofn_obs = int(param.eofn_obs) +eofn_mod = int(param.eofn_mod) + +if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: + eofn_expected = 1 +elif mode in ["NPGO", "NPO", "PSA1"]: + eofn_expected = 2 +elif mode in ["PSA2"]: + eofn_expected = 3 +else: + raise ValueError( + f"Mode '{mode}' is not defiend with associated expected EOF number" + ) if eofn_obs is None: - if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: - eofn_obs = 1 - elif mode in ["NPGO", "NPO", "PSA1"]: - eofn_obs = 2 - elif mode in ["PSA2"]: - eofn_obs = 3 - else: - raise ValueError(f"{eofn_obs} is not given for {mode}") + eofn_obs = eofn_expected else: - eofn_obs = int(eofn_obs) + if eofn_obs != eofn_expected: + raise ValueError( + f"Observation EOF number ({eofn_obs}) does not match expected EOF number ({eofn_expected}) for mode {mode}" + ) if eofn_mod is None: - if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: - eofn_mod = 1 - elif mode in ["NPGO", "NPO", "PSA2"]: - eofn_mod = 2 - elif mode in ["PSA2"]: - eofn_mod = 3 - else: - raise ValueError(f"{eofn_mod} is not given for {mode}") + eofn_mod = eofn_expected else: - eofn_mod = int(eofn_mod) + if eofn_mod != eofn_expected: + raise ValueError( + f"Model EOF number ({eofn_mod}) does not match expected EOF number ({eofn_expected}) for mode {mode}" + ) print("eofn_obs:", eofn_obs) print("eofn_mod:", eofn_mod) From f241f4a1f0427846728c569202e4b3f138a0c475 Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 16:18:38 -0700 Subject: [PATCH 5/7] minor bug fix --- pcmdi_metrics/variability_mode/variability_modes_driver.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/pcmdi_metrics/variability_mode/variability_modes_driver.py b/pcmdi_metrics/variability_mode/variability_modes_driver.py index 7344683b7..180c812f3 100755 --- a/pcmdi_metrics/variability_mode/variability_modes_driver.py +++ b/pcmdi_metrics/variability_mode/variability_modes_driver.py @@ -176,8 +176,8 @@ print("realization: ", realization) # EOF ordinal number -eofn_obs = int(param.eofn_obs) -eofn_mod = int(param.eofn_mod) +eofn_obs = param.eofn_obs +eofn_mod = param.eofn_mod if mode in ["NAM", "NAO", "SAM", "PNA", "PDO", "AMO"]: eofn_expected = 1 @@ -193,6 +193,7 @@ if eofn_obs is None: eofn_obs = eofn_expected else: + eofn_obs = int(eofn_obs) if eofn_obs != eofn_expected: raise ValueError( f"Observation EOF number ({eofn_obs}) does not match expected EOF number ({eofn_expected}) for mode {mode}" @@ -201,6 +202,7 @@ if eofn_mod is None: eofn_mod = eofn_expected else: + eofn_mod = int(eofn_mod) if eofn_mod != eofn_expected: raise ValueError( f"Model EOF number ({eofn_mod}) does not match expected EOF number ({eofn_expected}) for mode {mode}" From 12eab426b472fc6273a0c4304dc4a6538f1fec3d Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 16:56:16 -0700 Subject: [PATCH 6/7] minor clean up --- pcmdi_metrics/variability_mode/variability_modes_driver.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/pcmdi_metrics/variability_mode/variability_modes_driver.py b/pcmdi_metrics/variability_mode/variability_modes_driver.py index 180c812f3..d7d12f0a1 100755 --- a/pcmdi_metrics/variability_mode/variability_modes_driver.py +++ b/pcmdi_metrics/variability_mode/variability_modes_driver.py @@ -843,7 +843,7 @@ # Conventional EOF approach as supplementary # - - - - - - - - - - - - - - - - - - - - - - - - - if ConvEOF: - eofn_mod_max = 3 + eofn_mod_max = max(3, eofn_mod) # EOF analysis debug_print("conventional EOF analysis start", debug) @@ -989,7 +989,7 @@ debug=debug, ) plot_map( - mode + "_teleconnection", + f"{mode}_teleconnection", f"{mip.upper()} {model} ({run}) - EOF{n + 1}", msyear, meyear, @@ -997,7 +997,7 @@ # eof_lr(longitude=(lon1_global, lon2_global)), eof_lr, frac, - output_img_file + "_teleconnection", + f"{output_img_file}_teleconnection", debug=debug, ) From 74504b13b5997400f3fcfde70c22d76a7062634a Mon Sep 17 00:00:00 2001 From: Jiwoo Lee Date: Mon, 16 Sep 2024 17:00:53 -0700 Subject: [PATCH 7/7] add qc check for monthly time axis --- pcmdi_metrics/variability_mode/lib/lib_variability_mode.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/pcmdi_metrics/variability_mode/lib/lib_variability_mode.py b/pcmdi_metrics/variability_mode/lib/lib_variability_mode.py index 29f224408..e3c203414 100644 --- a/pcmdi_metrics/variability_mode/lib/lib_variability_mode.py +++ b/pcmdi_metrics/variability_mode/lib/lib_variability_mode.py @@ -13,7 +13,7 @@ import pcmdi_metrics from pcmdi_metrics.io import get_time, select_subset, xcdat_open -from pcmdi_metrics.utils import apply_landmask +from pcmdi_metrics.utils import apply_landmask, check_monthly_time_axis def tree(): @@ -67,6 +67,9 @@ def read_data_in( # Open data file ds = xcdat_open(path) + # Data QC check -- time axis check + check_monthly_time_axis(ds) + # Time subset ds_time_subsetted = subset_time(ds, syear, eyear, debug=debug) data_timeseries = ds_time_subsetted[var_in_data]