From 13103fa920aa7be8d993f2408237d57b74a7f676 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar <132295102+FariborzDaneshvar-NOAA@users.noreply.github.com> Date: Wed, 3 Jul 2024 10:23:02 -0700 Subject: [PATCH 01/40] Update README.md (#56) Add Daneshvar et al (NOAA technical report) citation. --- README.md | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index 8ec0ca3..5430d01 100644 --- a/README.md +++ b/README.md @@ -8,7 +8,10 @@ TO BE COMPLETED... ## References -- Daneshvar, F., et al. Tech Report (TODO) +- Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E. +_Probabilistic Prediction of Tropical Cyclone (TC) Driven Flood Depth_ +_and Extent with a User-Friendly Python Module_ +(NOAA Technical Memorandum NOS CS, under review) - Pringle, W. J., Mani, S., Sargsyan, K., Moghimi, S., Zhang, Y. J., Khazaei, B., Myers, E. (January 2023). _Surrogate-Assisted Bayesian Uncertainty Quantification for From 4604034afc2ef3994448accb88c85b054392ed16 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Thu, 1 Aug 2024 16:19:28 -0500 Subject: [PATCH 02/40] update dask args for Hercules --- stormworkflow/post/analyze_ensemble.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/stormworkflow/post/analyze_ensemble.py b/stormworkflow/post/analyze_ensemble.py index 44a78a3..1228723 100644 --- a/stormworkflow/post/analyze_ensemble.py +++ b/stormworkflow/post/analyze_ensemble.py @@ -415,10 +415,10 @@ def cli(): cluster = SLURMCluster(cores=16, processes=1, memory="500GB", - account="compute", + account="nos-surge", #PW:"compute" ; Hercules:"nos-surge" or "nos-ofs" walltime="08:00:00", header_skip=['--mem'], - interface="eth0") + # interface="eth0") # only for PW cluster.scale(6) client = Client(cluster) print(client) From 8cb64845cdaff793e9a3fa3f41609ea810799530 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Thu, 1 Aug 2024 17:16:19 -0500 Subject: [PATCH 03/40] fix typo in the revised dask form of analyze_ensemble.py --- stormworkflow/post/analyze_ensemble.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stormworkflow/post/analyze_ensemble.py b/stormworkflow/post/analyze_ensemble.py index 1228723..96ec07b 100644 --- a/stormworkflow/post/analyze_ensemble.py +++ b/stormworkflow/post/analyze_ensemble.py @@ -417,7 +417,7 @@ def cli(): memory="500GB", account="nos-surge", #PW:"compute" ; Hercules:"nos-surge" or "nos-ofs" walltime="08:00:00", - header_skip=['--mem'], + header_skip=['--mem'])#, # interface="eth0") # only for PW cluster.scale(6) client = Client(cluster) From 973ec3d8208070f3c1a847a88ca3c3090ad7d901 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Wed, 7 Aug 2024 14:45:38 -0500 Subject: [PATCH 04/40] move dask cluster command from main to cli --- stormworkflow/post/analyze_ensemble.py | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/stormworkflow/post/analyze_ensemble.py b/stormworkflow/post/analyze_ensemble.py index 96ec07b..f12ea8b 100644 --- a/stormworkflow/post/analyze_ensemble.py +++ b/stormworkflow/post/analyze_ensemble.py @@ -408,19 +408,19 @@ def cli(): parser.add_argument('-t', '--tracks-dir', type=Path) parser.add_argument('-s', '--sequential', action='store_true') - main(parser.parse_args()) - - -if __name__ == '__main__': cluster = SLURMCluster(cores=16, processes=1, memory="500GB", - account="nos-surge", #PW:"compute" ; Hercules:"nos-surge" or "nos-ofs" + account="nos-surge", #PW:"compute" ; Hercules:"nos-surge" or "nosofs" walltime="08:00:00", - header_skip=['--mem'])#, - # interface="eth0") # only for PW - cluster.scale(6) - client = Client(cluster) + header_skip=['--mem']) #, +# interface="eth0") # only for PW + cluster.scale(6) + client = Client(cluster) print(client) + main(parser.parse_args()) + + +if __name__ == '__main__': cli() From 098b48e3b387614be937e8118f7a59fee78982cd Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Thu, 8 Aug 2024 16:12:13 -0400 Subject: [PATCH 05/40] Unpin or relax version restrictions due to upstream updates --- environment.yml | 6 +++--- pyproject.toml | 9 ++++----- requirements.txt | 5 ----- 3 files changed, 7 insertions(+), 13 deletions(-) delete mode 100644 requirements.txt diff --git a/environment.yml b/environment.yml index 28cf939..ae94725 100644 --- a/environment.yml +++ b/environment.yml @@ -14,7 +14,7 @@ dependencies: - fiona - gdal - geoalchemy2 - - geopandas>=0.13 + - geopandas - geos - hdf5 - importlib_metadata<8 # Fix issue with esmpy Author import @@ -30,7 +30,7 @@ dependencies: - pyarrow - pygeos - pyproj - - python<3.11 + - python<3.12 - pytz - shapely>=2 - rasterio @@ -42,4 +42,4 @@ dependencies: - tqdm - udunits2 - utm - - xarray==2023.7.0 + - xarray diff --git a/pyproject.toml b/pyproject.toml index 8a8a1e7..e26c3ac 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -21,7 +21,7 @@ description = "A set of scripts to generate probabilistic storm surge results!" license = {file = "LICENSE"} -requires-python = ">= 3.8, < 3.11" +requires-python = ">= 3.8, < 3.12" dependencies = [ "cartopy", @@ -38,8 +38,7 @@ dependencies = [ "ensembleperturbation>=1.1.2", "fiona", "geoalchemy2", - "geopandas==0.14", # to address AttributeError. Should be fixed later in EnsemblePerturbation -# "geopandas", + "geopandas>=0.14", "matplotlib", "mpi4py", "netCDF4", @@ -54,7 +53,7 @@ dependencies = [ "pytz", "pyyaml", "shapely>=2", - "stormevents==2.2.4", + "stormevents>=2.2.4", "rasterio", "requests", "rtree", @@ -63,7 +62,7 @@ dependencies = [ "typing-extensions", "tqdm", "utm", - "xarray==2023.7.0", + "xarray", ] [tool.setuptools_scm] diff --git a/requirements.txt b/requirements.txt deleted file mode 100644 index 5cc4a85..0000000 --- a/requirements.txt +++ /dev/null @@ -1,5 +0,0 @@ -chaospy>=4.2.7 -stormevents==2.2.3 -pyschism>=0.1.15 -coupledmodeldriver>=1.6.6 -ensembleperturbation>=1.1.2 From d3a9d2948102dec2b7c6e6715870d99b8cb1896b Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Thu, 8 Aug 2024 16:33:58 -0400 Subject: [PATCH 06/40] Get countries shape directly from natural earth CDN --- stormworkflow/prep/hurricane_data.py | 11 ++--------- stormworkflow/scripts/workflow.sh | 1 - 2 files changed, 2 insertions(+), 10 deletions(-) diff --git a/stormworkflow/prep/hurricane_data.py b/stormworkflow/prep/hurricane_data.py index 5e3382a..f9977ec 100644 --- a/stormworkflow/prep/hurricane_data.py +++ b/stormworkflow/prep/hurricane_data.py @@ -37,7 +37,7 @@ format='%(asctime)s,%(msecs)d %(levelname)-8s [%(filename)s:%(lineno)d] %(message)s', datefmt='%Y-%m-%d:%H:%M:%S') - +NE_LOW_ADMIN = 'https://naciscdn.org/naturalearth/110m/cultural/ne_110m_admin_0_countries.zip' def trackstart_from_file( leadtime_file: Optional[pathlib.Path], @@ -143,12 +143,11 @@ def main(args): hr_before_landfall = args.hours_before_landfall lead_times = args.lead_times track_dir = args.preprocessed_tracks_dir - countries_shpfile = args.countries_polygon if hr_before_landfall < 0: hr_before_landfall = 48 - ne_low = gpd.read_file(countries_shpfile) + ne_low = gpd.read_file(NE_LOW_ADMIN) shp_US = ne_low[ne_low.NAME_EN.isin(['United States of America', 'Puerto Rico'])].unary_union logger.info("Fetching hurricane info...") @@ -415,12 +414,6 @@ def cli(): help="Existing adjusted track directory", ) - parser.add_argument( - "--countries-polygon", - type=pathlib.Path, - help="Shapefile containing country polygons", - ) - args = parser.parse_args() main(args) diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index e9fab6c..6c88c3e 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -74,7 +74,6 @@ hurricane_data \ --hours-before-landfall "$hr_prelandfall" \ --lead-times "$L_LEADTIMES_DATASET" \ --preprocessed-tracks-dir "$L_TRACK_DIR" \ - --countries-polygon "$L_SHP_DIR/ne_110m_cultural/ne_110m_admin_0_countries.shp" \ $storm $year 2>&1 | tee "${run_dir}/output/head_hurricane_data.out" From 14464d87c1943ff662178a3becbe243a95442a80 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Thu, 18 Jul 2024 15:12:27 +0000 Subject: [PATCH 07/40] Add hotstart ramp fix --- pyproject.toml | 6 +++--- stormworkflow/prep/setup_ensemble.py | 12 ++++++++++++ 2 files changed, 15 insertions(+), 3 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index e26c3ac..a1961ad 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -35,10 +35,10 @@ dependencies = [ "cmocean", "dask", "dask-jobqueue", - "ensembleperturbation>=1.1.2", + "ensembleperturbation==1.2.10", "fiona", "geoalchemy2", - "geopandas>=0.14", + "geopandas", "matplotlib", "mpi4py", "netCDF4", @@ -53,7 +53,7 @@ dependencies = [ "pytz", "pyyaml", "shapely>=2", - "stormevents>=2.2.4", + "stormevents==2.2.5", "rasterio", "requests", "rtree", diff --git a/stormworkflow/prep/setup_ensemble.py b/stormworkflow/prep/setup_ensemble.py index 71d2218..e5cef8e 100644 --- a/stormworkflow/prep/setup_ensemble.py +++ b/stormworkflow/prep/setup_ensemble.py @@ -107,6 +107,17 @@ def _fix_nwm_issues(ensemble_dir, hires_shapefile): _relocate_source_sink(pth, hires_shapefile) +def _fix_hotstart_issue(ensemble_dir): + hotstart_dirs = ensemble_dir.glob('runs/*') + for pth in hotstart_dirs: + nm_list = f90nml.read(pth / 'param.nml') + nm_list['opt']['dramp'] = 0.0 + nm_list['opt']['drampbc'] = 0.0 + nm_list['opt']['dramp_ss'] = 0.0 + nm_list['opt']['drampwind'] = 0.0 + nm_list.write(pth / 'param.nml', force=True) + + def main(args): track_path = args.track_file @@ -267,6 +278,7 @@ def main(args): } ) + _fix_hotstart_issue(workdir) if with_hydrology: _fix_nwm_issues(workdir, hires_reg) if use_wwm: From 7310e710850e98083159316abd86bf1c41738366 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 9 Aug 2024 09:22:27 -0400 Subject: [PATCH 08/40] Update dependency versions --- pyproject.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index a1961ad..a064cd8 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -35,7 +35,7 @@ dependencies = [ "cmocean", "dask", "dask-jobqueue", - "ensembleperturbation==1.2.10", + "ensembleperturbation>=1.3.0", # rmax option "fiona", "geoalchemy2", "geopandas", @@ -53,7 +53,7 @@ dependencies = [ "pytz", "pyyaml", "shapely>=2", - "stormevents==2.2.5", + "stormevents>=2.2.5", "rasterio", "requests", "rtree", From 1bf242889c4ba9e8b1243846b172b5b416fc02bd Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 9 Aug 2024 11:46:31 -0400 Subject: [PATCH 09/40] Add perturb variable argument to input file --- stormworkflow/prep/setup_ensemble.py | 8 ++------ stormworkflow/refs/input.yaml | 6 ++++++ stormworkflow/scripts/workflow.sh | 1 + 3 files changed, 9 insertions(+), 6 deletions(-) diff --git a/stormworkflow/prep/setup_ensemble.py b/stormworkflow/prep/setup_ensemble.py index e5cef8e..29cc81e 100644 --- a/stormworkflow/prep/setup_ensemble.py +++ b/stormworkflow/prep/setup_ensemble.py @@ -189,12 +189,7 @@ def main(args): perturbations=args.num_perturbations, directory=workdir / 'track_files', storm=workdir / 'track_to_perturb.dat', - variables=[ - 'cross_track', - 'along_track', - 'radius_of_maximum_winds', # TODO: add option for persistent - 'max_sustained_wind_speed', - ], + variables=args.variables, sample_from_distribution=args.sample_from_distribution, sample_rule=args.sample_rule, quadrature=args.quadrature, @@ -339,6 +334,7 @@ def parse_arguments(): argument_parser.add_argument('--use-wwm', action='store_true') argument_parser.add_argument('--with-hydrology', action='store_true') argument_parser.add_argument('--pahm-model', choices=['gahm', 'symmetric'], default='gahm') + argument_parser.add_argument('--variables', nargs='+', type=str) argument_parser.add_argument('name', help='name of the storm', type=str) diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index 1e8141e..4b664ba 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -14,6 +14,12 @@ num_perturb: 2 sample_rule: "korobov" spinup_exec: "pschism_PAHM_TVD-VL" hotstart_exec: "pschism_PAHM_TVD-VL" +perturb_vars: + - 'cross_track' + - 'along_track' +# - 'radius_of_maximum_winds' + - 'radius_of_maximum_winds_persistent' + - 'max_sustained_wind_speed' hpc_solver_nnodes: 3 hpc_solver_ntasks: 108 diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index 6c88c3e..3b14816 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -131,6 +131,7 @@ PREP_KWDS+=" --sample-rule $sample_rule" PREP_KWDS+=" --date-range-file $run_dir/setup/dates.csv" PREP_KWDS+=" --nwm-file $L_NWM_DATASET" PREP_KWDS+=" --tpxo-dir $L_TPXO_DATASET" +PREP_KWDS+=" --variables $perturb_vars" if [ $use_wwm == 1 ]; then PREP_KWDS+=" --use-wwm"; fi if [ $hydrology == 1 ]; then PREP_KWDS+=" --with-hydrology"; fi PREP_KWDS+=" --pahm-model $pahm_model" From 9973fbaa785731d887649d47d2fc2d2f86968cf9 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 9 Aug 2024 11:51:13 -0400 Subject: [PATCH 10/40] update cmd dependency version --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index a064cd8..8c7bd32 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -30,7 +30,7 @@ dependencies = [ "cfgrib", "cfunits", "chaospy>=4.2.7", - "coupledmodeldriver>=1.6.6", + "coupledmodeldriver>=1.7.1.post1", "colored-traceback", "cmocean", "dask", From a3a4f2343362e8343c3083055485848fb2a4b766 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 9 Aug 2024 11:57:46 -0400 Subject: [PATCH 11/40] update stormevents to version with rmax forecast stuff --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 8c7bd32..4dd25e5 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -53,7 +53,7 @@ dependencies = [ "pytz", "pyyaml", "shapely>=2", - "stormevents>=2.2.5", + "stormevents>=2.3.2", # rmax forecast "rasterio", "requests", "rtree", From e5ad3e2c935fb6cd79f10b8f8f40e300ae38bddb Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 12 Aug 2024 13:38:10 -0400 Subject: [PATCH 12/40] Bump input version number --- stormworkflow/refs/input.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index 4b664ba..f438167 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -1,5 +1,5 @@ --- -input_version: 0.0.1 +input_version: 0.0.2 storm: "florence" year: 2018 From d66cb721ddd0f0974b21c50a23cf262ed767f05f Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Wed, 14 Aug 2024 15:51:22 -0400 Subject: [PATCH 13/40] Cache geodataset path in the section that runs on headnode --- stormworkflow/prep/hurricane_data.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/stormworkflow/prep/hurricane_data.py b/stormworkflow/prep/hurricane_data.py index f9977ec..8ec20ba 100644 --- a/stormworkflow/prep/hurricane_data.py +++ b/stormworkflow/prep/hurricane_data.py @@ -17,6 +17,7 @@ import pandas as pd import geopandas as gpd +import geodatasets from searvey.coops import COOPS_TidalDatum from searvey.coops import COOPS_TimeZone from searvey.coops import COOPS_Units @@ -147,6 +148,9 @@ def main(args): if hr_before_landfall < 0: hr_before_landfall = 48 + # Caching for next steps! + geodatasets.get_path('naturalearth land') + ne_low = gpd.read_file(NE_LOW_ADMIN) shp_US = ne_low[ne_low.NAME_EN.isin(['United States of America', 'Puerto Rico'])].unary_union From 21ca3e32b4d1598569c5f9cac820a44d45430307 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 11:38:25 -0400 Subject: [PATCH 14/40] Note versions for workflow, coupldemodeldriver and pyschism in output. --- stormworkflow/scripts/workflow.sh | 3 +++ 1 file changed, 3 insertions(+) diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index 3b14816..f1bbf4f 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -48,8 +48,11 @@ function init { done logfile=$run_dir/versions.info + version $logfile stormworkflow version $logfile stormevents version $logfile ensembleperturbation + version $logfile coupledmodeldriver + version $logfile pyschism version $logfile ocsmesh echo "SCHISM: see solver.version each outputs dir" >> $logfile From c2252e04a91bbdb89680da52cf9bc9b79e916e15 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 12:04:28 -0400 Subject: [PATCH 15/40] Add input version check and enforcement --- pyproject.toml | 1 + stormworkflow/main.py | 55 +++++++++++++++++++++++++++++-- stormworkflow/scripts/workflow.sh | 2 +- 3 files changed, 55 insertions(+), 3 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index 4dd25e5..446076b 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -45,6 +45,7 @@ dependencies = [ "numpy", "numba", "ocsmesh==1.5.3", + "packaging", "pandas", "pyarrow", "pygeos", diff --git a/stormworkflow/main.py b/stormworkflow/main.py index e649f6c..1aa5265 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -6,16 +6,62 @@ from argparse import ArgumentParser from pathlib import Path -import stormworkflow import yaml +from packaging.version import Version try: from yaml import CLoader as Loader, CDumper as Dumper except ImportError: from yaml import Loader, Dumper +import stormworkflow _logger = logging.getLogger(__file__) +CUR_INPUT_VER = Version('0.0.2') + + +def _handle_input_v0_0_1_to_v0_0_2(inout_conf): + + # Only update conf if specified version matches the assumed one + if ver != Version('0.0.1'): + return ver + + + _logger.info( + "Adding perturbation variables for persistent RMW perturbation" + ) + inout_conf['perturb_vars'] = [ + 'cross_track', + 'along_track', + 'radius_of_maximum_winds_persistent', + 'max_sustained_wind_speed', + ] + + return Version('0.0.2') + + +def _handle_input_version(inout_conf): + + if 'input_version' not in inout_conf: + ver = CUR_INPUT_VER + _logger.warning( + f"`input_version` is NOT specified in `input.yaml`; assuming {ver}" + ) + + ver = Version(inout_conf['input_version']) + + if ver > CUR_INPUT_VER: + raise ValueError( + f"Input version not supported! Max version supported is {CUR_INPUT_VER}" + ) + + ver = _handle_input_v0_0_1_to_v0_0_2(inout_conf) + + if ver != CUR_INPUT_VER + raise ValueError( + f"Could NOT update input to the latest version! Updated to {ver}" + + def main(): parser = ArgumentParser() @@ -28,12 +74,17 @@ def main(): infile = args.configuration if infile is None: - _logger.warn('No input configuration provided, using reference file!') + _logger.warning( + 'No input configuration provided, using reference file!' + ) infile = refs.joinpath('input.yaml') with open(infile, 'r') as yfile: conf = yaml.load(yfile, Loader=Loader) + _handle_input_version(conf) + # TODO: Write out the updated conf as a yaml file + wf = scripts.joinpath('workflow.sh') run_env = os.environ.copy() diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index f1bbf4f..723f4b0 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -56,7 +56,7 @@ function init { version $logfile ocsmesh echo "SCHISM: see solver.version each outputs dir" >> $logfile - cp $input_file $run_dir/input.yaml + cp $input_file $run_dir/input_asis.yaml echo $run_dir } From bc90ec1d7d5d0a96d40b90a36454785aee1c67fd Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 12:06:33 -0400 Subject: [PATCH 16/40] Syntax error fix --- stormworkflow/main.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 1aa5265..4bdc056 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -57,9 +57,10 @@ def _handle_input_version(inout_conf): ver = _handle_input_v0_0_1_to_v0_0_2(inout_conf) - if ver != CUR_INPUT_VER + if ver != CUR_INPUT_VER: raise ValueError( f"Could NOT update input to the latest version! Updated to {ver}" + ) def main(): From 6382e34bce056e27a48774d4dbb071890ada37ed Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 12:10:18 -0400 Subject: [PATCH 17/40] Update config object version after handling updates --- stormworkflow/main.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 4bdc056..f228960 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -22,7 +22,7 @@ def _handle_input_v0_0_1_to_v0_0_2(inout_conf): - # Only update conf if specified version matches the assumed one + # Only update config if specified version matches the assumed one if ver != Version('0.0.1'): return ver @@ -62,6 +62,8 @@ def _handle_input_version(inout_conf): f"Could NOT update input to the latest version! Updated to {ver}" ) + inout_conf['input_version'] = str(v2) + def main(): @@ -84,7 +86,7 @@ def main(): conf = yaml.load(yfile, Loader=Loader) _handle_input_version(conf) - # TODO: Write out the updated conf as a yaml file + # TODO: Write out the updated config as a yaml file wf = scripts.joinpath('workflow.sh') From e0e95d4b5017769692c8ac7a8f810a9ff8688e0e Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 14:43:48 -0400 Subject: [PATCH 18/40] First iteration of CI tests --- .github/workflows/tests.yml | 36 ++++++++++++++++++++++++ pyproject.toml | 7 ++++- stormworkflow/main.py | 11 +++++--- tests/data/refs/input_v0.0.1.yaml | 40 +++++++++++++++++++++++++++ tests/data/refs/input_v0.0.2.yaml | 46 +++++++++++++++++++++++++++++++ tests/test_input_version.py | 43 +++++++++++++++++++++++++++++ 6 files changed, 178 insertions(+), 5 deletions(-) create mode 100644 .github/workflows/tests.yml create mode 100644 tests/data/refs/input_v0.0.1.yaml create mode 100644 tests/data/refs/input_v0.0.2.yaml create mode 100644 tests/test_input_version.py diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml new file mode 100644 index 0000000..026d464 --- /dev/null +++ b/.github/workflows/tests.yml @@ -0,0 +1,36 @@ +name: tests + +on: + push: + branches: + - main + paths: + - '**.py' + - '.github/workflows/tests.yml' + - 'pyproject.toml' + pull_request: + branches: + - main + +jobs: + test: + name: test + runs-on: ${{ matrix.os }} + strategy: + matrix: + os: [ ubuntu-latest, macos-latest ] + python-version: [ '3.9', '3.10', '3.11' ] + steps: + - name: clone repository + uses: actions/checkout@v4 + - name: conda virtual environment + uses: mamba-org/setup-micromamba@v1 + with: + init-shell: bash + environment-file: environment.yml + - name: install the package + run: pip install ".[dev]" + shell: micromamba-shell {0} + - name: run tests + run: pytest --numprocesses auto + shell: micromamba-shell {0} diff --git a/pyproject.toml b/pyproject.toml index 446076b..2f7cedc 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -21,7 +21,7 @@ description = "A set of scripts to generate probabilistic storm surge results!" license = {file = "LICENSE"} -requires-python = ">= 3.8, < 3.12" +requires-python = ">= 3.9, < 3.12" dependencies = [ "cartopy", @@ -66,6 +66,11 @@ dependencies = [ "xarray", ] +[project.optional-dependencies] +dev = [ + "pytest" +] + [tool.setuptools_scm] version_file = "stormworkflow/_version.py" diff --git a/stormworkflow/main.py b/stormworkflow/main.py index f228960..80fa06f 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -2,6 +2,7 @@ import logging import os import shlex +import warnings from importlib.resources import files from argparse import ArgumentParser from pathlib import Path @@ -40,13 +41,15 @@ def _handle_input_v0_0_1_to_v0_0_2(inout_conf): return Version('0.0.2') -def _handle_input_version(inout_conf): +def handle_input_version(inout_conf): if 'input_version' not in inout_conf: ver = CUR_INPUT_VER - _logger.warning( + warnings.warn( f"`input_version` is NOT specified in `input.yaml`; assuming {ver}" ) + inout_conf['input_version'] = ver + return ver = Version(inout_conf['input_version']) @@ -77,7 +80,7 @@ def main(): infile = args.configuration if infile is None: - _logger.warning( + warnings.warn( 'No input configuration provided, using reference file!' ) infile = refs.joinpath('input.yaml') @@ -85,7 +88,7 @@ def main(): with open(infile, 'r') as yfile: conf = yaml.load(yfile, Loader=Loader) - _handle_input_version(conf) + handle_input_version(conf) # TODO: Write out the updated config as a yaml file wf = scripts.joinpath('workflow.sh') diff --git a/tests/data/refs/input_v0.0.1.yaml b/tests/data/refs/input_v0.0.1.yaml new file mode 100644 index 0000000..1e8141e --- /dev/null +++ b/tests/data/refs/input_v0.0.1.yaml @@ -0,0 +1,40 @@ +--- +input_version: 0.0.1 + +storm: "florence" +year: 2018 +suffix: "" +subset_mesh: 1 +hr_prelandfall: -1 +past_forecast: 1 +hydrology: 0 +use_wwm: 0 +pahm_model: "gahm" +num_perturb: 2 +sample_rule: "korobov" +spinup_exec: "pschism_PAHM_TVD-VL" +hotstart_exec: "pschism_PAHM_TVD-VL" + +hpc_solver_nnodes: 3 +hpc_solver_ntasks: 108 +hpc_account: "" +hpc_partition: "" + +RUN_OUT: "" +L_NWM_DATASET: "" +L_TPXO_DATASET: "" +L_LEADTIMES_DATASET: "" +L_TRACK_DIR: "" +L_DEM_HI: "" +L_DEM_LO: "" +L_MESH_HI: "" +L_MESH_LO: "" +L_SHP_DIR: "" + +TMPDIR: "/tmp" +PATH_APPEND: "" + +L_SOLVE_MODULES: + - "intel/2022.1.2" + - "impi/2022.1.2" + - "netcdf" diff --git a/tests/data/refs/input_v0.0.2.yaml b/tests/data/refs/input_v0.0.2.yaml new file mode 100644 index 0000000..f438167 --- /dev/null +++ b/tests/data/refs/input_v0.0.2.yaml @@ -0,0 +1,46 @@ +--- +input_version: 0.0.2 + +storm: "florence" +year: 2018 +suffix: "" +subset_mesh: 1 +hr_prelandfall: -1 +past_forecast: 1 +hydrology: 0 +use_wwm: 0 +pahm_model: "gahm" +num_perturb: 2 +sample_rule: "korobov" +spinup_exec: "pschism_PAHM_TVD-VL" +hotstart_exec: "pschism_PAHM_TVD-VL" +perturb_vars: + - 'cross_track' + - 'along_track' +# - 'radius_of_maximum_winds' + - 'radius_of_maximum_winds_persistent' + - 'max_sustained_wind_speed' + +hpc_solver_nnodes: 3 +hpc_solver_ntasks: 108 +hpc_account: "" +hpc_partition: "" + +RUN_OUT: "" +L_NWM_DATASET: "" +L_TPXO_DATASET: "" +L_LEADTIMES_DATASET: "" +L_TRACK_DIR: "" +L_DEM_HI: "" +L_DEM_LO: "" +L_MESH_HI: "" +L_MESH_LO: "" +L_SHP_DIR: "" + +TMPDIR: "/tmp" +PATH_APPEND: "" + +L_SOLVE_MODULES: + - "intel/2022.1.2" + - "impi/2022.1.2" + - "netcdf" diff --git a/tests/test_input_version.py b/tests/test_input_version.py new file mode 100644 index 0000000..fb5a0a6 --- /dev/null +++ b/tests/test_input_version.py @@ -0,0 +1,43 @@ +from importlib.resources import files + +import pytest +import yaml +from packaging.version import Version +from yaml import Loader, Dumper + +from stormworkflow.main import handle_input_version, CUR_INPUT_VER + + +refs = files('tests.data.refs') +input_v0_0_1 = refs.joinpath('input_v0.0.1.yaml') +input_v0_0_2 = refs.joinpath('input_v0.0.2.yaml') + + +def read_conf(infile): + with open(infile, 'r') as yfile: + conf = yaml.load(yfile, Loader=Loader) + return conf + + +@pytest.fixture +def conf_v0_0_1(): + return read_conf(input_v0_0_1) + + +@pytest.fixture +def conf_v0_0_2(): + return read_conf(input_v0_0_2) + + +@pytest.fixture +def conf_latest(conf_v0_0_2): + return conf_v0_0_2 + + +def test_no_version_specified(conf_latest): + conf_latest.pop('input_version') + with pytest.warns(UserWarning): + handle_input_version(conf_latest) + + assert conf_latest['input_version'] == CUR_INPUT_VER + From ef073fa807f79eca76cd5fa2703d7ed814155ce6 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 16:25:37 -0400 Subject: [PATCH 19/40] Add more test and bugfix for version check --- stormworkflow/main.py | 6 ++++-- tests/test_input_version.py | 28 ++++++++++++++++++++++++++-- 2 files changed, 30 insertions(+), 4 deletions(-) diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 80fa06f..665dbff 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -23,6 +23,8 @@ def _handle_input_v0_0_1_to_v0_0_2(inout_conf): + ver = Version(inout_conf['input_version']) + # Only update config if specified version matches the assumed one if ver != Version('0.0.1'): return ver @@ -48,7 +50,7 @@ def handle_input_version(inout_conf): warnings.warn( f"`input_version` is NOT specified in `input.yaml`; assuming {ver}" ) - inout_conf['input_version'] = ver + inout_conf['input_version'] = str(ver) return ver = Version(inout_conf['input_version']) @@ -65,7 +67,7 @@ def handle_input_version(inout_conf): f"Could NOT update input to the latest version! Updated to {ver}" ) - inout_conf['input_version'] = str(v2) + inout_conf['input_version'] = str(ver) def main(): diff --git a/tests/test_input_version.py b/tests/test_input_version.py index fb5a0a6..6fa47b3 100644 --- a/tests/test_input_version.py +++ b/tests/test_input_version.py @@ -1,3 +1,4 @@ +from copy import deepcopy from importlib.resources import files import pytest @@ -38,6 +39,29 @@ def test_no_version_specified(conf_latest): conf_latest.pop('input_version') with pytest.warns(UserWarning): handle_input_version(conf_latest) - - assert conf_latest['input_version'] == CUR_INPUT_VER + + assert conf_latest['input_version'] == str(CUR_INPUT_VER) + +def test_invalid_version_specified(conf_latest): + + invalid_1 = deepcopy(conf_latest) + invalid_1['input_version'] = ( + f'{CUR_INPUT_VER.major}.{CUR_INPUT_VER.minor}.{CUR_INPUT_VER.micro + 1}' + ) + with pytest.raises(ValueError) as e: + handle_input_version(invalid_1) + + assert "max" in str(e.value).lower() + + + invalid_2 = deepcopy(conf_latest) + invalid_2['input_version'] = 'a.b.c' + with pytest.raises(ValueError) as e: + handle_input_version(invalid_2) + assert "invalid version" in str(e.value).lower() + + +def test_v0_0_1_to_v0_0_2(conf_v0_0_1, conf_v0_0_2): + handle_input_version(conf_v0_0_1) + assert conf_v0_0_2 == conf_v0_0_1 From f1d82a11aa5e052248678fb3d6d1377f4532d69a Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 16:38:04 -0400 Subject: [PATCH 20/40] Remove macos from tests --- .github/workflows/tests.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 026d464..d03fa52 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -18,7 +18,7 @@ jobs: runs-on: ${{ matrix.os }} strategy: matrix: - os: [ ubuntu-latest, macos-latest ] + os: [ ubuntu-latest ] python-version: [ '3.9', '3.10', '3.11' ] steps: - name: clone repository From edb27e509be6813046303e6b16a927553d61f221 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Mon, 19 Aug 2024 16:46:27 -0400 Subject: [PATCH 21/40] Update pytest arguments --- .github/workflows/tests.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index d03fa52..bcf9649 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -32,5 +32,5 @@ jobs: run: pip install ".[dev]" shell: micromamba-shell {0} - name: run tests - run: pytest --numprocesses auto + run: pytest shell: micromamba-shell {0} From 2e1d955e22ac49206154398130ea06c0e7f838a1 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 23 Aug 2024 08:45:53 -0400 Subject: [PATCH 22/40] Update stormevents version, add rmw input and bump input version --- pyproject.toml | 2 +- stormworkflow/main.py | 29 +++++++++++++++---- stormworkflow/refs/input.yaml | 16 ++++++----- tests/data/refs/input_v0.0.3.yaml | 48 +++++++++++++++++++++++++++++++ tests/test_input_version.py | 19 +++++++++--- 5 files changed, 97 insertions(+), 17 deletions(-) create mode 100644 tests/data/refs/input_v0.0.3.yaml diff --git a/pyproject.toml b/pyproject.toml index 2f7cedc..8646a81 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -54,7 +54,7 @@ dependencies = [ "pytz", "pyyaml", "shapely>=2", - "stormevents>=2.3.2", # rmax forecast + "stormevents>=2.3.4", # rmax option "rasterio", "requests", "rtree", diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 665dbff..7889365 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -18,7 +18,7 @@ _logger = logging.getLogger(__file__) -CUR_INPUT_VER = Version('0.0.2') +CUR_INPUT_VER = Version('0.0.3') def _handle_input_v0_0_1_to_v0_0_2(inout_conf): @@ -43,6 +43,23 @@ def _handle_input_v0_0_1_to_v0_0_2(inout_conf): return Version('0.0.2') +def _handle_input_v0_0_2_to_v0_0_3(inout_conf): + + ver = Version(inout_conf['input_version']) + + # Only update config if specified version matches the assumed one + if ver != Version('0.0.2'): + return ver + + + _logger.info( + "Adding RMW fill method default to persistent" + ) + inout_conf['rmw_fill_method'] = 'persistent' + + return Version('0.0.3') + + def handle_input_version(inout_conf): if 'input_version' not in inout_conf: @@ -60,16 +77,18 @@ def handle_input_version(inout_conf): f"Input version not supported! Max version supported is {CUR_INPUT_VER}" ) - ver = _handle_input_v0_0_1_to_v0_0_2(inout_conf) + for fn in [ + _handle_input_v0_0_1_to_v0_0_2, + _handle_input_v0_0_2_to_v0_0_3, + ]: + ver = fn(inout_conf) + inout_conf['input_version'] = str(ver) if ver != CUR_INPUT_VER: raise ValueError( f"Could NOT update input to the latest version! Updated to {ver}" ) - inout_conf['input_version'] = str(ver) - - def main(): parser = ArgumentParser() diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index f438167..02df434 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -1,5 +1,5 @@ --- -input_version: 0.0.2 +input_version: 0.0.3 storm: "florence" year: 2018 @@ -12,14 +12,16 @@ use_wwm: 0 pahm_model: "gahm" num_perturb: 2 sample_rule: "korobov" +perturb_vars: + - "cross_track" + - "along_track" +# - "radius_of_maximum_winds" + - "radius_of_maximum_winds_persistent" + - "max_sustained_wind_speed" +rmw_fill_method: "persistent" + spinup_exec: "pschism_PAHM_TVD-VL" hotstart_exec: "pschism_PAHM_TVD-VL" -perturb_vars: - - 'cross_track' - - 'along_track' -# - 'radius_of_maximum_winds' - - 'radius_of_maximum_winds_persistent' - - 'max_sustained_wind_speed' hpc_solver_nnodes: 3 hpc_solver_ntasks: 108 diff --git a/tests/data/refs/input_v0.0.3.yaml b/tests/data/refs/input_v0.0.3.yaml new file mode 100644 index 0000000..02df434 --- /dev/null +++ b/tests/data/refs/input_v0.0.3.yaml @@ -0,0 +1,48 @@ +--- +input_version: 0.0.3 + +storm: "florence" +year: 2018 +suffix: "" +subset_mesh: 1 +hr_prelandfall: -1 +past_forecast: 1 +hydrology: 0 +use_wwm: 0 +pahm_model: "gahm" +num_perturb: 2 +sample_rule: "korobov" +perturb_vars: + - "cross_track" + - "along_track" +# - "radius_of_maximum_winds" + - "radius_of_maximum_winds_persistent" + - "max_sustained_wind_speed" +rmw_fill_method: "persistent" + +spinup_exec: "pschism_PAHM_TVD-VL" +hotstart_exec: "pschism_PAHM_TVD-VL" + +hpc_solver_nnodes: 3 +hpc_solver_ntasks: 108 +hpc_account: "" +hpc_partition: "" + +RUN_OUT: "" +L_NWM_DATASET: "" +L_TPXO_DATASET: "" +L_LEADTIMES_DATASET: "" +L_TRACK_DIR: "" +L_DEM_HI: "" +L_DEM_LO: "" +L_MESH_HI: "" +L_MESH_LO: "" +L_SHP_DIR: "" + +TMPDIR: "/tmp" +PATH_APPEND: "" + +L_SOLVE_MODULES: + - "intel/2022.1.2" + - "impi/2022.1.2" + - "netcdf" diff --git a/tests/test_input_version.py b/tests/test_input_version.py index 6fa47b3..5fd615b 100644 --- a/tests/test_input_version.py +++ b/tests/test_input_version.py @@ -12,6 +12,7 @@ refs = files('tests.data.refs') input_v0_0_1 = refs.joinpath('input_v0.0.1.yaml') input_v0_0_2 = refs.joinpath('input_v0.0.2.yaml') +input_v0_0_3 = refs.joinpath('input_v0.0.3.yaml') def read_conf(infile): @@ -31,8 +32,13 @@ def conf_v0_0_2(): @pytest.fixture -def conf_latest(conf_v0_0_2): - return conf_v0_0_2 +def conf_v0_0_3(): + return read_conf(input_v0_0_3) + + +@pytest.fixture +def conf_latest(conf_v0_0_3): + return conf_v0_0_3 def test_no_version_specified(conf_latest): @@ -62,6 +68,11 @@ def test_invalid_version_specified(conf_latest): assert "invalid version" in str(e.value).lower() -def test_v0_0_1_to_v0_0_2(conf_v0_0_1, conf_v0_0_2): +def test_v0_0_1_to_latest(conf_v0_0_1, conf_latest): handle_input_version(conf_v0_0_1) - assert conf_v0_0_2 == conf_v0_0_1 + assert conf_latest == conf_v0_0_1 + + +def test_v0_0_2_to_latest(conf_v0_0_2, conf_latest): + handle_input_version(conf_v0_0_2) + assert conf_latest == conf_v0_0_2 From d827216ceceb2e7049a684ac2d5f557c2efe7ced Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 23 Aug 2024 09:39:35 -0400 Subject: [PATCH 23/40] Actually use the rmw value from the input --- stormworkflow/prep/hurricane_data.py | 45 ++++++++++++++++++++++++---- stormworkflow/prep/setup_ensemble.py | 3 ++ stormworkflow/scripts/workflow.sh | 1 + 3 files changed, 43 insertions(+), 6 deletions(-) diff --git a/stormworkflow/prep/hurricane_data.py b/stormworkflow/prep/hurricane_data.py index 8ec20ba..faf729e 100644 --- a/stormworkflow/prep/hurricane_data.py +++ b/stormworkflow/prep/hurricane_data.py @@ -24,6 +24,7 @@ from shapely.geometry import box, base from stormevents import StormEvent from stormevents.nhc import VortexTrack +from stormevents.nhc.const import RMWFillMethod from stormevents.nhc.track import ( combine_tracks, correct_ofcl_based_on_carq_n_hollandb, @@ -144,6 +145,7 @@ def main(args): hr_before_landfall = args.hours_before_landfall lead_times = args.lead_times track_dir = args.preprocessed_tracks_dir + rmw_fill = RMWFillMethod[args.rmw_fill.lower()] if hr_before_landfall < 0: hr_before_landfall = 48 @@ -183,7 +185,8 @@ def main(args): advisory = 'OFCL' if not local_track_file.is_file(): - # Find and pick a single advisory based on priority + # Find and pick a single advisory based on priority, the + # track is only used to get the available advisories temp_track = event.track(file_deck='a') adv_avail = temp_track.unfiltered_data.advisory.unique() adv_order = ['OFCL', 'HWRF', 'HMON', 'CARQ'] @@ -193,41 +196,65 @@ def main(args): advisory = adv break - # TODO: THIS IS NO LONGER RELEVANT IF WE FAKE RMWP AS OFCL! if advisory == "OFCL" and "CARQ" not in adv_avail: raise ValueError( "OFCL advisory needs CARQ for fixing missing variables!" ) - track = VortexTrack(nhc_code, file_deck='a', advisories=[advisory]) + track = VortexTrack( + nhc_code, + file_deck='a', + advisories=[advisory], + rmw_fill=rmw_fill, + ) else: # read from preprocessed file + advisory = 'OFCL' # If a file exists, use the local file track_raw = pd.read_csv(local_track_file, header=None, dtype=str) - assert len(track_raw[4].unique()) == 1 + # The provided tracks should have a single advisory type, + # e.g. in NHC adjusted track files the value is RMWP + if len(track_raw[4].unique()) != 1: + raise RuntimeError( + "Only single advisory-name track files are supported!" + ) + # Treat the existing advisory as if it's OFCL so that + # stormevents supports reading it track_raw[4] = advisory with tempfile.NamedTemporaryFile() as tmp: track_raw.to_csv(tmp.name, header=False, index=False) + # Track read from file is NOT corrected because it + # does NOT have CARQ advisory unfixed_track = VortexTrack( tmp.name, file_deck='a', advisories=[advisory] ) + # Since we're only getting CARQ, there's no need to + # pass rmw fill method carq_track = event.track(file_deck='a', advisories=['CARQ']) unfix_dict = { **separate_tracks(unfixed_track.data), **separate_tracks(carq_track.data), } - fix_dict = correct_ofcl_based_on_carq_n_hollandb(unfix_dict) + # Fix the file track with the fetched CARQ; if RMW + # is already filled, it fills out other missing values + fix_dict = correct_ofcl_based_on_carq_n_hollandb( + unfix_dict, rmw_fill=rmw_fill + ) fix_track = combine_tracks(fix_dict) + # Create a new VortexTrack object from the datasets. + # Since the values are already filled in, there's + # no need to fill the rmw! track = VortexTrack( fix_track[fix_track.advisory == advisory], file_deck='a', - advisories=[advisory] + advisories=[advisory], + rmw_fill=RMWFillMethod.none, ) @@ -418,6 +445,12 @@ def cli(): help="Existing adjusted track directory", ) + parser.add_argument( + '--rmw-fill', + type=str, + help="Method to use to fill missing RMW data for OFCL track", + ) + args = parser.parse_args() main(args) diff --git a/stormworkflow/prep/setup_ensemble.py b/stormworkflow/prep/setup_ensemble.py index 29cc81e..fe7d889 100644 --- a/stormworkflow/prep/setup_ensemble.py +++ b/stormworkflow/prep/setup_ensemble.py @@ -176,6 +176,9 @@ def main(args): # get has unique forecast time for only the segment we want to # perturb, the preceeding entries are 0-hour forecasts from # previous forecast_times + # + # Here we're working with NA-filled track files, so there's + # no need for rmw fill argument track_to_perturb = VortexTrack.from_file( track_path, start_date=perturb_start, diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index 723f4b0..6956d0e 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -77,6 +77,7 @@ hurricane_data \ --hours-before-landfall "$hr_prelandfall" \ --lead-times "$L_LEADTIMES_DATASET" \ --preprocessed-tracks-dir "$L_TRACK_DIR" \ + --rmw-fill "$rmw_fill_method" \ $storm $year 2>&1 | tee "${run_dir}/output/head_hurricane_data.out" From d6e90c54fd6872a4c7583d32087f987e5bb4c838 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 23 Aug 2024 09:54:48 -0400 Subject: [PATCH 24/40] Add exclusive argument for post jobs --- stormworkflow/slurm/post.sbatch | 1 + 1 file changed, 1 insertion(+) diff --git a/stormworkflow/slurm/post.sbatch b/stormworkflow/slurm/post.sbatch index 4ef59f8..ea61395 100644 --- a/stormworkflow/slurm/post.sbatch +++ b/stormworkflow/slurm/post.sbatch @@ -2,6 +2,7 @@ #SBATCH --parsable #SBATCH --time=05:00:00 #SBATCH --nodes=1 +#SBATCH --exclusive set -ex From 32b9bfbadad741f4c15dafd0675944ad936d7c11 Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Fri, 23 Aug 2024 13:02:13 -0400 Subject: [PATCH 25/40] Cleanup test fixtures --- tests/conftest.py | 39 +++++++++++++++++++++++++++++++++++ tests/test_input_version.py | 41 +++++-------------------------------- 2 files changed, 44 insertions(+), 36 deletions(-) create mode 100644 tests/conftest.py diff --git a/tests/conftest.py b/tests/conftest.py new file mode 100644 index 0000000..edaf46c --- /dev/null +++ b/tests/conftest.py @@ -0,0 +1,39 @@ +from importlib.resources import files + +import pytest +import yaml +from yaml import Loader + + +refs = files('stormworkflow.refs') +test_refs = files('tests.data.refs') +input_v0_0_1 = test_refs.joinpath('input_v0.0.1.yaml') +input_v0_0_2 = test_refs.joinpath('input_v0.0.2.yaml') +input_v0_0_3 = test_refs.joinpath('input_v0.0.3.yaml') +input_latest = refs.joinpath('input.yaml') + + +def read_conf(infile): + with open(infile, 'r') as yfile: + conf = yaml.load(yfile, Loader=Loader) + return conf + + +@pytest.fixture +def conf_v0_0_1(): + return read_conf(input_v0_0_1) + + +@pytest.fixture +def conf_v0_0_2(): + return read_conf(input_v0_0_2) + + +@pytest.fixture +def conf_v0_0_3(): + return read_conf(input_v0_0_3) + + +@pytest.fixture +def conf_latest(): + return read_conf(input_latest) diff --git a/tests/test_input_version.py b/tests/test_input_version.py index 5fd615b..5e2ae38 100644 --- a/tests/test_input_version.py +++ b/tests/test_input_version.py @@ -1,46 +1,10 @@ from copy import deepcopy -from importlib.resources import files import pytest -import yaml -from packaging.version import Version -from yaml import Loader, Dumper from stormworkflow.main import handle_input_version, CUR_INPUT_VER -refs = files('tests.data.refs') -input_v0_0_1 = refs.joinpath('input_v0.0.1.yaml') -input_v0_0_2 = refs.joinpath('input_v0.0.2.yaml') -input_v0_0_3 = refs.joinpath('input_v0.0.3.yaml') - - -def read_conf(infile): - with open(infile, 'r') as yfile: - conf = yaml.load(yfile, Loader=Loader) - return conf - - -@pytest.fixture -def conf_v0_0_1(): - return read_conf(input_v0_0_1) - - -@pytest.fixture -def conf_v0_0_2(): - return read_conf(input_v0_0_2) - - -@pytest.fixture -def conf_v0_0_3(): - return read_conf(input_v0_0_3) - - -@pytest.fixture -def conf_latest(conf_v0_0_3): - return conf_v0_0_3 - - def test_no_version_specified(conf_latest): conf_latest.pop('input_version') with pytest.warns(UserWarning): @@ -76,3 +40,8 @@ def test_v0_0_1_to_latest(conf_v0_0_1, conf_latest): def test_v0_0_2_to_latest(conf_v0_0_2, conf_latest): handle_input_version(conf_v0_0_2) assert conf_latest == conf_v0_0_2 + + +def test_v0_0_3_to_latest(conf_v0_0_3, conf_latest): + handle_input_version(conf_v0_0_3) + assert conf_latest == conf_v0_0_3 From 3bcb0d2f4c066e207c659aefc8ba3208a575519f Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Thu, 29 Aug 2024 12:20:01 -0500 Subject: [PATCH 26/40] call roc_curve in post.sbatch --- stormworkflow/slurm/post.sbatch | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/stormworkflow/slurm/post.sbatch b/stormworkflow/slurm/post.sbatch index ea61395..2e74c23 100644 --- a/stormworkflow/slurm/post.sbatch +++ b/stormworkflow/slurm/post.sbatch @@ -13,3 +13,11 @@ combine_ensemble \ analyze_ensemble \ --ensemble-dir $ENSEMBLE_DIR \ --tracks-dir $ENSEMBLE_DIR/track_files + +storm_roc_curve \ + --storm_name ${storm} \ + --storm_year ${year} \ + --leadtime ${hr_prelandfall} \ + --prob_nc_path $ENSEMBLE_DIR/analyze/linear_k1_p1_n0.025/probabilities.nc \ + --obs_df_path ${NHC_OBS} \ + --save_dir $ENSEMBLE_DIR/analyze/linear_k1_p1_n0.025 From 579ab2cc6d84b2aaad3e16496e2c84cf6bdd61b7 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Fri, 30 Aug 2024 12:39:39 -0500 Subject: [PATCH 27/40] replace ROC_single_run.py with storm_roc_curve.py in [project.scripts] --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 8646a81..d85956d 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -99,4 +99,4 @@ download_data = "stormworkflow.prep.download_data:cli" setup_ensemble = "stormworkflow.prep.setup_ensemble:cli" combine_ensemble = "stormworkflow.post.combine_ensemble:cli" analyze_ensemble = "stormworkflow.post.analyze_ensemble:cli" -storm_roc_curve = "stormworkflow.post.ROC_single_run:cli" +storm_roc_curve = "stormworkflow.post.storm_roc_curve:cli" From dac175329b806c87ea3c7fa24db9ff4cdd9fc34a Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Fri, 30 Aug 2024 12:42:40 -0500 Subject: [PATCH 28/40] add storm_roc_curve with args to the script --- stormworkflow/slurm/post.sbatch | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/stormworkflow/slurm/post.sbatch b/stormworkflow/slurm/post.sbatch index 2e74c23..33ef357 100644 --- a/stormworkflow/slurm/post.sbatch +++ b/stormworkflow/slurm/post.sbatch @@ -15,9 +15,8 @@ analyze_ensemble \ --tracks-dir $ENSEMBLE_DIR/track_files storm_roc_curve \ - --storm_name ${storm} \ - --storm_year ${year} \ + --storm ${storm} \ + --year ${year} \ --leadtime ${hr_prelandfall} \ - --prob_nc_path $ENSEMBLE_DIR/analyze/linear_k1_p1_n0.025/probabilities.nc \ --obs_df_path ${NHC_OBS} \ - --save_dir $ENSEMBLE_DIR/analyze/linear_k1_p1_n0.025 + --ensemble-dir $ENSEMBLE_DIR From da5d501b4838c5cf7611ff003ec5ee7c29a46651 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Fri, 30 Aug 2024 12:49:18 -0500 Subject: [PATCH 29/40] make a copy of ROC_single_run.py and format it to be consistent with the rest of workflow. i.e. variable names and args definition --- stormworkflow/post/storm_roc_curve.py | 293 ++++++++++++++++++++++++++ 1 file changed, 293 insertions(+) create mode 100644 stormworkflow/post/storm_roc_curve.py diff --git a/stormworkflow/post/storm_roc_curve.py b/stormworkflow/post/storm_roc_curve.py new file mode 100644 index 0000000..4758374 --- /dev/null +++ b/stormworkflow/post/storm_roc_curve.py @@ -0,0 +1,293 @@ +import argparse +import logging +import os +import warnings +import numpy as np +import pandas as pd +import xarray as xr +import scipy as sp +import matplotlib.pyplot as plt +from pathlib import Path +from cartopy.feature import NaturalEarthFeature + +os.environ['USE_PYGEOS'] = '0' +import geopandas as gpd + +pd.options.mode.copy_on_write = True + + +def stack_station_coordinates(x, y): + """ + Create numpy.column_stack based on + coordinates of observation points + """ + coord_combined = np.column_stack([x, y]) + return coord_combined + + +def create_search_tree(longitude, latitude): + """ + Create scipy.spatial.CKDTree based on Lat. and Long. + """ + long_lat = np.column_stack((longitude.T.ravel(), latitude.T.ravel())) + tree = sp.spatial.cKDTree(long_lat) + return tree + + +def find_nearby_prediction(ds, variable, indices): + """ + Reads netcdf file, target variable, and indices + Returns max value among corresponding indices for each point + """ + obs_count = indices.shape[0] # total number of search/observation points + max_prediction_index = len(ds.node.values) # total number of nodes + + prediction_prob = np.zeros(obs_count) # assuming all are dry (probability of zero) + + for obs_point in range(obs_count): + idx_arr = np.delete( + indices[obs_point], np.where(indices[obs_point] == max_prediction_index)[0] + ) # len is length of surrogate model array + val_arr = ds[variable].values[idx_arr] + val_arr = np.nan_to_num(val_arr) # replace nan with zero (dry node) + + # # Pick the nearest non-zero probability (option #1) + # for val in val_arr: + # if val > 0.0: + # prediction_prob[obs_point] = round(val,4) #round to 0.1 mm + # break + + # pick the largest value (option #2) + if val_arr.size > 0: + prediction_prob[obs_point] = val_arr.max() + return prediction_prob + + +def plot_probabilities(df, prob_column, gdf_countries, title, save_name): + """ + plot probabilities of exceeding given threshold at obs. points + """ + figure, axis = plt.subplots(1, 1) + figure.set_size_inches(10, 10 / 1.6) + + plt.scatter(x=df.Longitude, y=df.Latitude, vmin=0, vmax=1.0, c=df[prob_column]) + xlim = axis.get_xlim() + ylim = axis.get_ylim() + + gdf_countries.plot(color='lightgrey', ax=axis, zorder=-5) + + axis.set_xlim(xlim) + axis.set_ylim(ylim) + plt.colorbar(shrink=0.75) + plt.title(title) + plt.savefig(save_name) + plt.close() + + +def calculate_hit_miss(df, obs_column, prob_column, threshold, probability): + """ + Reads dataframe with two columns for obs_elev, and probabilities + returns hit/miss/... based on user-defined threshold & probability + """ + hit = len(df[(df[obs_column] >= threshold) & (df[prob_column] >= probability)]) + miss = len(df[(df[obs_column] >= threshold) & (df[prob_column] < probability)]) + false_alarm = len(df[(df[obs_column] < threshold) & (df[prob_column] >= probability)]) + correct_neg = len(df[(df[obs_column] < threshold) & (df[prob_column] < probability)]) + + return hit, miss, false_alarm, correct_neg + + +def calculate_POD_FAR(hit, miss, false_alarm, correct_neg): + """ + Reads hit, miss, false_alarm, and correct_neg + returns POD and FAR + default POD and FAR are np.nan + """ + POD = np.nan + FAR = np.nan + try: + POD = round(hit / (hit + miss), 4) # Probability of Detection + except ZeroDivisionError: + pass + try: + FAR = round(false_alarm / (false_alarm + correct_neg), 4) # False Alarm Rate + except ZeroDivisionError: + pass + return POD, FAR + + +def main(args): + storm = args.storm.capitalize() + year = args.year + leadtime = args.leadtime + obs_df_path = Path(args.obs_df_path) + ensemble_dir = Path(args.ensemble_dir) + + output_directory = ensemble_dir / 'analyze/linear_k1_p1_n0.025' + prob_nc_path = output_directory / 'probabilities.nc' + + if leadtime == -1: + leadtime = 48 + + # *.nc file coordinates + thresholds_ft = [3, 6, 9] # in ft + thresholds_m = [round(i * 0.3048, 4) for i in thresholds_ft] # convert to meter + sources = ['model', 'surrogate'] + probabilities = [0.0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0] + + # attributes of input files + prediction_variable = 'probabilities' + obs_attribute = 'Elev_m_xGEOID20b' + + # search criteria + max_distance = 1000 # [in meters] to set distance_upper_bound + max_neighbors = 10 # to set k + + blank_arr = np.empty((len(thresholds_ft), 1, 1, len(sources), len(probabilities))) + blank_arr[:] = np.nan + + hit_arr = blank_arr.copy() + miss_arr = blank_arr.copy() + false_alarm_arr = blank_arr.copy() + correct_neg_arr = blank_arr.copy() + POD_arr = blank_arr.copy() + FAR_arr = blank_arr.copy() + + # Load obs file, extract storm obs points and coordinates + df_obs = pd.read_csv(obs_df_path) + Event_name = f'{storm}_{year}' + df_obs_storm = df_obs[df_obs.Event == Event_name] + obs_coordinates = stack_station_coordinates( + df_obs_storm.Longitude.values, df_obs_storm.Latitude.values + ) + + # Load probabilities.nc file + ds_prob = xr.open_dataset(prob_nc_path) + + # gdf_countries = gpd.GeoSeries( + # NaturalEarthFeature(category='physical', scale='10m', name='land',).geometries(), + # crs=4326, + # ) + + # Loop through thresholds and sources and find corresponding values from probabilities.nc + threshold_count = -1 + for threshold in thresholds_m: + threshold_count += 1 + source_count = -1 + for source in sources: + source_count += 1 + ds_temp = ds_prob.sel(level=threshold, source=source) + tree = create_search_tree(ds_temp.x.values, ds_temp.y.values) + dist, indices = tree.query( + obs_coordinates, k=max_neighbors, distance_upper_bound=max_distance * 1e-5 + ) # 0.01 is equivalent to 1000 m + prediction_prob = find_nearby_prediction( + ds=ds_temp, variable=prediction_variable, indices=indices + ) + df_obs_storm[f'{source}_prob'] = prediction_prob + + # # Plot probabilities at obs. points + # plot_probabilities( + # df_obs_storm, + # f'{source}_prob', + # gdf_countries, + # f'Probability of {source} exceeding {thresholds_ft[threshold_count]} ft \n {storm}, {year}, {leadtime}-hr leadtime', + # os.path.join( + # output_directory, + # f'prob_{source}_above_{thresholds_ft[threshold_count]}ft_{storm}_{year}_{leadtime}-hr.png', + # ), + # ) + + # Loop through probabilities: calculate hit/miss/... & POD/FAR + prob_count = -1 + for prob in probabilities: + prob_count += 1 + hit, miss, false_alarm, correct_neg = calculate_hit_miss( + df_obs_storm, obs_attribute, f'{source}_prob', threshold, prob + ) + hit_arr[threshold_count, 0, 0, source_count, prob_count] = hit + miss_arr[threshold_count, 0, 0, source_count, prob_count] = miss + false_alarm_arr[threshold_count, 0, 0, source_count, prob_count] = false_alarm + correct_neg_arr[threshold_count, 0, 0, source_count, prob_count] = correct_neg + + pod, far = calculate_POD_FAR(hit, miss, false_alarm, correct_neg) + POD_arr[threshold_count, 0, 0, source_count, prob_count] = pod + FAR_arr[threshold_count, 0, 0, source_count, prob_count] = far + + ds_ROC = xr.Dataset( + coords=dict( + threshold=thresholds_ft, + storm=[storm], + leadtime=[leadtime], + source=sources, + prob=probabilities, + ), + data_vars=dict( + hit=(['threshold', 'storm', 'leadtime', 'source', 'prob'], hit_arr), + miss=(['threshold', 'storm', 'leadtime', 'source', 'prob'], miss_arr), + false_alarm=( + ['threshold', 'storm', 'leadtime', 'source', 'prob'], + false_alarm_arr, + ), + correct_neg=( + ['threshold', 'storm', 'leadtime', 'source', 'prob'], + correct_neg_arr, + ), + POD=(['threshold', 'storm', 'leadtime', 'source', 'prob'], POD_arr), + FAR=(['threshold', 'storm', 'leadtime', 'source', 'prob'], FAR_arr), + ), + ) + ds_ROC.to_netcdf(os.path.join(output_directory, f'{storm}_{year}_{leadtime}hr_POD_FAR.nc')) + + # plot ROC curves + marker_list = ['s', 'x'] + linestyle_list = ['dashed', 'dotted'] + threshold_count = -1 + for threshold in thresholds_ft: + threshold_count += 1 + fig = plt.figure() + ax = fig.add_subplot(111) + plt.axline( + (0.0, 0.0), (1.0, 1.0), linestyle='--', color='grey', label='random prediction' + ) + source_count = -1 + for source in sources: + source_count += 1 + plt.plot( + FAR_arr[threshold_count, 0, 0, source_count, :], + POD_arr[threshold_count, 0, 0, source_count, :], + label=f'{source}', + marker=marker_list[source_count], + linestyle=linestyle_list[source_count], + markersize=5, + ) + plt.legend() + plt.xlabel('False Alarm Rate') + plt.ylabel('Probability of Detection') + + plt.title(f'{storm}_{year}, {leadtime}-hr leadtime, {threshold} ft threshold') + plt.savefig( + os.path.join( + output_directory, + f'ROC_{storm}_{year}_{leadtime}hr_leadtime_{threshold}_ft.png', + ) + ) + plt.close() + + +def cli(): + parser = argparse.ArgumentParser() + + parser.add_argument('--storm', help='name of the storm', type=str) + parser.add_argument('--year', help='year of the storm', type=int) + parser.add_argument('--leadtime', help='OFCL track leadtime hr', type=int) + parser.add_argument('--obs_df_path', help='path to NHC obs data', type=str) + parser.add_argument('--ensemble-dir', help='path to ensemble.dir', type=str) + + main(parser.parse_args()) + + +if __name__ == '__main__': + warnings.filterwarnings('ignore') + # warnings.filterwarnings("ignore", category=DeprecationWarning) + cli() From d3adee26b76c87fc2dac7f457765a789f848a0bd Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Fri, 30 Aug 2024 12:52:30 -0500 Subject: [PATCH 30/40] delete ROC_single_run.py. storm_roc_curve.py is the new format of this file. --- stormworkflow/post/ROC_single_run.py | 300 --------------------------- 1 file changed, 300 deletions(-) delete mode 100644 stormworkflow/post/ROC_single_run.py diff --git a/stormworkflow/post/ROC_single_run.py b/stormworkflow/post/ROC_single_run.py deleted file mode 100644 index 26dbd2c..0000000 --- a/stormworkflow/post/ROC_single_run.py +++ /dev/null @@ -1,300 +0,0 @@ -import argparse -import logging -import os -import warnings -import numpy as np -import pandas as pd -import xarray as xr -import scipy as sp -import matplotlib.pyplot as plt -from pathlib import Path -from cartopy.feature import NaturalEarthFeature - -os.environ['USE_PYGEOS'] = '0' -import geopandas as gpd - -pd.options.mode.copy_on_write = True - - -def stack_station_coordinates(x, y): - """ - Create numpy.column_stack based on - coordinates of observation points - """ - coord_combined = np.column_stack([x, y]) - return coord_combined - - -def create_search_tree(longitude, latitude): - """ - Create scipy.spatial.CKDTree based on Lat. and Long. - """ - long_lat = np.column_stack((longitude.T.ravel(), latitude.T.ravel())) - tree = sp.spatial.cKDTree(long_lat) - return tree - - -def find_nearby_prediction(ds, variable, indices): - """ - Reads netcdf file, target variable, and indices - Returns max value among corresponding indices for each point - """ - obs_count = indices.shape[0] # total number of search/observation points - max_prediction_index = len(ds.node.values) # total number of nodes - - prediction_prob = np.zeros(obs_count) # assuming all are dry (probability of zero) - - for obs_point in range(obs_count): - idx_arr = np.delete( - indices[obs_point], np.where(indices[obs_point] == max_prediction_index)[0] - ) # len is length of surrogate model array - val_arr = ds[variable].values[idx_arr] - val_arr = np.nan_to_num(val_arr) # replace nan with zero (dry node) - - # # Pick the nearest non-zero probability (option #1) - # for val in val_arr: - # if val > 0.0: - # prediction_prob[obs_point] = round(val,4) #round to 0.1 mm - # break - - # pick the largest value (option #2) - if val_arr.size > 0: - prediction_prob[obs_point] = val_arr.max() - return prediction_prob - - -def plot_probabilities(df, prob_column, gdf_countries, title, save_name): - """ - plot probabilities of exceeding given threshold at obs. points - """ - figure, axis = plt.subplots(1, 1) - figure.set_size_inches(10, 10 / 1.6) - - plt.scatter(x=df.Longitude, y=df.Latitude, vmin=0, vmax=1.0, c=df[prob_column]) - xlim = axis.get_xlim() - ylim = axis.get_ylim() - - gdf_countries.plot(color='lightgrey', ax=axis, zorder=-5) - - axis.set_xlim(xlim) - axis.set_ylim(ylim) - plt.colorbar(shrink=0.75) - plt.title(title) - plt.savefig(save_name) - plt.close() - - -def calculate_hit_miss(df, obs_column, prob_column, threshold, probability): - """ - Reads dataframe with two columns for obs_elev, and probabilities - returns hit/miss/... based on user-defined threshold & probability - """ - hit = len(df[(df[obs_column] >= threshold) & (df[prob_column] >= probability)]) - miss = len(df[(df[obs_column] >= threshold) & (df[prob_column] < probability)]) - false_alarm = len(df[(df[obs_column] < threshold) & (df[prob_column] >= probability)]) - correct_neg = len(df[(df[obs_column] < threshold) & (df[prob_column] < probability)]) - - return hit, miss, false_alarm, correct_neg - - -def calculate_POD_FAR(hit, miss, false_alarm, correct_neg): - """ - Reads hit, miss, false_alarm, and correct_neg - returns POD and FAR - default POD and FAR are np.nan - """ - POD = np.nan - FAR = np.nan - try: - POD = round(hit / (hit + miss), 4) # Probability of Detection - except ZeroDivisionError: - pass - try: - FAR = round(false_alarm / (false_alarm + correct_neg), 4) # False Alarm Rate - except ZeroDivisionError: - pass - return POD, FAR - - -def main(args): - storm_name = args.storm_name.capitalize() - storm_year = args.storm_year - leadtime = args.leadtime - prob_nc_path = Path(args.prob_nc_path) - obs_df_path = Path(args.obs_df_path) - save_dir = args.save_dir - - # *.nc file coordinates - thresholds_ft = [3, 6, 9] # in ft - thresholds_m = [round(i * 0.3048, 4) for i in thresholds_ft] # convert to meter - sources = ['model', 'surrogate'] - probabilities = [0.0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0] - - # attributes of input files - prediction_variable = 'probabilities' - obs_attribute = 'Elev_m_xGEOID20b' - - # search criteria - max_distance = 1000 # [in meters] to set distance_upper_bound - max_neighbors = 10 # to set k - - blank_arr = np.empty((len(thresholds_ft), 1, 1, len(sources), len(probabilities))) - blank_arr[:] = np.nan - - hit_arr = blank_arr.copy() - miss_arr = blank_arr.copy() - false_alarm_arr = blank_arr.copy() - correct_neg_arr = blank_arr.copy() - POD_arr = blank_arr.copy() - FAR_arr = blank_arr.copy() - - # Load obs file, extract storm obs points and coordinates - df_obs = pd.read_csv(obs_df_path) - Event_name = f'{storm_name}_{storm_year}' - df_obs_storm = df_obs[df_obs.Event == Event_name] - obs_coordinates = stack_station_coordinates( - df_obs_storm.Longitude.values, df_obs_storm.Latitude.values - ) - - # Load probabilities.nc file - ds_prob = xr.open_dataset(prob_nc_path) - - gdf_countries = gpd.GeoSeries( - NaturalEarthFeature(category='physical', scale='10m', name='land',).geometries(), - crs=4326, - ) - - # Loop through thresholds and sources and find corresponding values from probabilities.nc - threshold_count = -1 - for threshold in thresholds_m: - threshold_count += 1 - source_count = -1 - for source in sources: - source_count += 1 - ds_temp = ds_prob.sel(level=threshold, source=source) - tree = create_search_tree(ds_temp.x.values, ds_temp.y.values) - dist, indices = tree.query( - obs_coordinates, k=max_neighbors, distance_upper_bound=max_distance * 1e-5 - ) # 0.01 is equivalent to 1000 m - prediction_prob = find_nearby_prediction( - ds=ds_temp, variable=prediction_variable, indices=indices - ) - df_obs_storm[f'{source}_prob'] = prediction_prob - - # Plot probabilities at obs. points - plot_probabilities( - df_obs_storm, - f'{source}_prob', - gdf_countries, - f'Probability of {source} exceeding {thresholds_ft[threshold_count]} ft \n {storm_name}, {storm_year}, {leadtime}-hr leadtime', - os.path.join( - save_dir, - f'prob_{source}_above_{thresholds_ft[threshold_count]}ft_{storm_name}_{storm_year}_{leadtime}-hr.png', - ), - ) - - # Loop through probabilities: calculate hit/miss/... & POD/FAR - prob_count = -1 - for prob in probabilities: - prob_count += 1 - hit, miss, false_alarm, correct_neg = calculate_hit_miss( - df_obs_storm, obs_attribute, f'{source}_prob', threshold, prob - ) - hit_arr[threshold_count, 0, 0, source_count, prob_count] = hit - miss_arr[threshold_count, 0, 0, source_count, prob_count] = miss - false_alarm_arr[threshold_count, 0, 0, source_count, prob_count] = false_alarm - correct_neg_arr[threshold_count, 0, 0, source_count, prob_count] = correct_neg - - pod, far = calculate_POD_FAR(hit, miss, false_alarm, correct_neg) - POD_arr[threshold_count, 0, 0, source_count, prob_count] = pod - FAR_arr[threshold_count, 0, 0, source_count, prob_count] = far - - ds_ROC = xr.Dataset( - coords=dict( - threshold=thresholds_ft, - storm=[storm_name], - leadtime=[leadtime], - source=sources, - prob=probabilities, - ), - data_vars=dict( - hit=(['threshold', 'storm', 'leadtime', 'source', 'prob'], hit_arr), - miss=(['threshold', 'storm', 'leadtime', 'source', 'prob'], miss_arr), - false_alarm=( - ['threshold', 'storm', 'leadtime', 'source', 'prob'], - false_alarm_arr, - ), - correct_neg=( - ['threshold', 'storm', 'leadtime', 'source', 'prob'], - correct_neg_arr, - ), - POD=(['threshold', 'storm', 'leadtime', 'source', 'prob'], POD_arr), - FAR=(['threshold', 'storm', 'leadtime', 'source', 'prob'], FAR_arr), - ), - ) - ds_ROC.to_netcdf( - os.path.join(save_dir, f'{storm_name}_{storm_year}_{leadtime}hr_leadtime_POD_FAR.nc') - ) - - # plot ROC curves - marker_list = ['s', 'x'] - linestyle_list = ['dashed', 'dotted'] - threshold_count = -1 - for threshold in thresholds_ft: - threshold_count += 1 - fig = plt.figure() - ax = fig.add_subplot(111) - plt.axline( - (0.0, 0.0), (1.0, 1.0), linestyle='--', color='grey', label='random prediction' - ) - source_count = -1 - for source in sources: - source_count += 1 - plt.plot( - FAR_arr[threshold_count, 0, 0, source_count, :], - POD_arr[threshold_count, 0, 0, source_count, :], - label=f'{source}', - marker=marker_list[source_count], - linestyle=linestyle_list[source_count], - markersize=5, - ) - plt.legend() - plt.xlabel('False Alarm Rate') - plt.ylabel('Probability of Detection') - - plt.title( - f'{storm_name}_{storm_year}, {leadtime}-hr leadtime, {threshold} ft threshold' - ) - plt.savefig( - os.path.join( - save_dir, f'ROC_{storm_name}_{leadtime}hr_leadtime_{threshold}_ft.png' - ) - ) - plt.close() - - -def cli(): - parser = argparse.ArgumentParser() - - parser.add_argument('--storm_name', help='name of the storm', type=str) - - parser.add_argument('--storm_year', help='year of the storm', type=int) - - parser.add_argument('--leadtime', help='OFCL track leadtime hr', type=int) - - parser.add_argument('--prob_nc_path', help='path to probabilities.nc', type=str) - - parser.add_argument('--obs_df_path', help='Path to observations dataframe', type=str) - - # optional - parser.add_argument( - '--save_dir', help='directory for saving analysis', default=os.getcwd(), type=str - ) - - main(parser.parse_args()) - - -if __name__ == '__main__': - warnings.filterwarnings('ignore') - # warnings.filterwarnings("ignore", category=DeprecationWarning) - cli() From d172802c17359a3d1918b2e9bbcaa9d8cdce6292 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Fri, 30 Aug 2024 13:11:04 -0500 Subject: [PATCH 31/40] add tentative path to NHC_OBS for calculation of POD/FAR and plotting ROC ruves --- stormworkflow/refs/input.yaml | 1 + 1 file changed, 1 insertion(+) diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index 02df434..26d3b00 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -38,6 +38,7 @@ L_DEM_LO: "" L_MESH_HI: "" L_MESH_LO: "" L_SHP_DIR: "" +NHC_OBS: "" TMPDIR: "/tmp" PATH_APPEND: "" From 9ae8e0432a04ba40c8c0723a1a8fdaf69f7defe8 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Tue, 3 Sep 2024 09:36:37 -0500 Subject: [PATCH 32/40] increase timelimits of mesh and prep to 1 hour --- stormworkflow/slurm/mesh.sbatch | 2 +- stormworkflow/slurm/prep.sbatch | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/stormworkflow/slurm/mesh.sbatch b/stormworkflow/slurm/mesh.sbatch index 757b9cd..503381a 100644 --- a/stormworkflow/slurm/mesh.sbatch +++ b/stormworkflow/slurm/mesh.sbatch @@ -1,7 +1,7 @@ #!/bin/bash #SBATCH --parsable #SBATCH --exclusive -#SBATCH --time=00:30:00 +#SBATCH --time=01:00:00 #SBATCH --nodes=1 set -ex diff --git a/stormworkflow/slurm/prep.sbatch b/stormworkflow/slurm/prep.sbatch index 8beb298..d5dd455 100644 --- a/stormworkflow/slurm/prep.sbatch +++ b/stormworkflow/slurm/prep.sbatch @@ -1,7 +1,7 @@ #!/bin/bash #SBATCH --parsable #SBATCH --exclusive -#SBATCH --time=00:30:00 +#SBATCH --time=01:00:00 #SBATCH --nodes=1 set -ex From 3111110a04a91e4a50d3b5e7e0cf9a0f99c3518a Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar Date: Wed, 4 Sep 2024 14:05:18 -0500 Subject: [PATCH 33/40] use geodatasets and make spatial plot of probabilities at obs. points --- stormworkflow/post/storm_roc_curve.py | 26 +++++++++++++++----------- 1 file changed, 15 insertions(+), 11 deletions(-) diff --git a/stormworkflow/post/storm_roc_curve.py b/stormworkflow/post/storm_roc_curve.py index 4758374..fe7bde5 100644 --- a/stormworkflow/post/storm_roc_curve.py +++ b/stormworkflow/post/storm_roc_curve.py @@ -10,6 +10,8 @@ from pathlib import Path from cartopy.feature import NaturalEarthFeature +import geodatasets + os.environ['USE_PYGEOS'] = '0' import geopandas as gpd @@ -164,6 +166,8 @@ def main(args): # Load probabilities.nc file ds_prob = xr.open_dataset(prob_nc_path) + gdf_countries = gpd.read_file(geodatasets.get_path('naturalearth land')) + # gdf_countries = gpd.GeoSeries( # NaturalEarthFeature(category='physical', scale='10m', name='land',).geometries(), # crs=4326, @@ -186,17 +190,17 @@ def main(args): ) df_obs_storm[f'{source}_prob'] = prediction_prob - # # Plot probabilities at obs. points - # plot_probabilities( - # df_obs_storm, - # f'{source}_prob', - # gdf_countries, - # f'Probability of {source} exceeding {thresholds_ft[threshold_count]} ft \n {storm}, {year}, {leadtime}-hr leadtime', - # os.path.join( - # output_directory, - # f'prob_{source}_above_{thresholds_ft[threshold_count]}ft_{storm}_{year}_{leadtime}-hr.png', - # ), - # ) + # Plot probabilities at obs. points + plot_probabilities( + df_obs_storm, + f'{source}_prob', + gdf_countries, + f'Probability of {source} exceeding {thresholds_ft[threshold_count]} ft \n {storm}, {year}, {leadtime}-hr leadtime', + os.path.join( + output_directory, + f'prob_{source}_above_{thresholds_ft[threshold_count]}ft_{storm}_{year}_{leadtime}-hr.png', + ), + ) # Loop through probabilities: calculate hit/miss/... & POD/FAR prob_count = -1 From 21e8f9b8731c7be47be14c7db61db07e3d5719bd Mon Sep 17 00:00:00 2001 From: SorooshMani-NOAA Date: Thu, 5 Sep 2024 16:42:25 -0400 Subject: [PATCH 34/40] Update and improve version check --- stormworkflow/main.py | 51 +++++++++++++++++++------------ stormworkflow/refs/input.yaml | 2 +- tests/conftest.py | 6 +++- tests/data/refs/input_v0.0.4.yaml | 49 +++++++++++++++++++++++++++++ tests/test_input_version.py | 5 +++ 5 files changed, 91 insertions(+), 22 deletions(-) create mode 100644 tests/data/refs/input_v0.0.4.yaml diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 7889365..244a3e0 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -18,17 +18,32 @@ _logger = logging.getLogger(__file__) -CUR_INPUT_VER = Version('0.0.3') +CUR_INPUT_VER = Version('0.0.4') +VER_UPDATE_FUNCS = [] -def _handle_input_v0_0_1_to_v0_0_2(inout_conf): +def _input_version(prev, curr): + def decorator(handler): + def wrapper(inout_conf): + ver = Version(inout_conf['input_version']) - ver = Version(inout_conf['input_version']) + # Only update config if specified version matches the + # assumed one + if ver != Version(prev): + return ver - # Only update config if specified version matches the assumed one - if ver != Version('0.0.1'): - return ver + # TODO: Need return values? + handler(inout_conf) + return Version(curr) + global VER_UPDATE_FUNCS + VER_UPDATE_FUNCS.append(wrapper) + return wrapper + return decorator + + +@_input_version('0.0.1', '0.0.2') +def _handle_input_v0_0_1_to_v0_0_2(inout_conf): _logger.info( "Adding perturbation variables for persistent RMW perturbation" @@ -40,24 +55,23 @@ def _handle_input_v0_0_1_to_v0_0_2(inout_conf): 'max_sustained_wind_speed', ] - return Version('0.0.2') - +@_input_version('0.0.2', '0.0.3') def _handle_input_v0_0_2_to_v0_0_3(inout_conf): - ver = Version(inout_conf['input_version']) - - # Only update config if specified version matches the assumed one - if ver != Version('0.0.2'): - return ver - - _logger.info( "Adding RMW fill method default to persistent" ) inout_conf['rmw_fill_method'] = 'persistent' - return Version('0.0.3') + +@_input_version('0.0.3', '0.0.4') +def _handle_input_v0_0_3_to_v0_0_4(inout_conf): + + _logger.info( + "Path to observations" + ) + inout_conf['NHC_OBS'] = '' def handle_input_version(inout_conf): @@ -77,10 +91,7 @@ def handle_input_version(inout_conf): f"Input version not supported! Max version supported is {CUR_INPUT_VER}" ) - for fn in [ - _handle_input_v0_0_1_to_v0_0_2, - _handle_input_v0_0_2_to_v0_0_3, - ]: + for fn in VER_UPDATE_FUNCS: ver = fn(inout_conf) inout_conf['input_version'] = str(ver) diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index 26d3b00..30aefa2 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -1,5 +1,5 @@ --- -input_version: 0.0.3 +input_version: 0.0.4 storm: "florence" year: 2018 diff --git a/tests/conftest.py b/tests/conftest.py index edaf46c..4d3c7e1 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -6,10 +6,11 @@ refs = files('stormworkflow.refs') -test_refs = files('tests.data.refs') +test_refs = files('data.refs') input_v0_0_1 = test_refs.joinpath('input_v0.0.1.yaml') input_v0_0_2 = test_refs.joinpath('input_v0.0.2.yaml') input_v0_0_3 = test_refs.joinpath('input_v0.0.3.yaml') +input_v0_0_4 = test_refs.joinpath('input_v0.0.4.yaml') input_latest = refs.joinpath('input.yaml') @@ -33,6 +34,9 @@ def conf_v0_0_2(): def conf_v0_0_3(): return read_conf(input_v0_0_3) +@pytest.fixture +def conf_v0_0_4(): + return read_conf(input_v0_0_4) @pytest.fixture def conf_latest(): diff --git a/tests/data/refs/input_v0.0.4.yaml b/tests/data/refs/input_v0.0.4.yaml new file mode 100644 index 0000000..30aefa2 --- /dev/null +++ b/tests/data/refs/input_v0.0.4.yaml @@ -0,0 +1,49 @@ +--- +input_version: 0.0.4 + +storm: "florence" +year: 2018 +suffix: "" +subset_mesh: 1 +hr_prelandfall: -1 +past_forecast: 1 +hydrology: 0 +use_wwm: 0 +pahm_model: "gahm" +num_perturb: 2 +sample_rule: "korobov" +perturb_vars: + - "cross_track" + - "along_track" +# - "radius_of_maximum_winds" + - "radius_of_maximum_winds_persistent" + - "max_sustained_wind_speed" +rmw_fill_method: "persistent" + +spinup_exec: "pschism_PAHM_TVD-VL" +hotstart_exec: "pschism_PAHM_TVD-VL" + +hpc_solver_nnodes: 3 +hpc_solver_ntasks: 108 +hpc_account: "" +hpc_partition: "" + +RUN_OUT: "" +L_NWM_DATASET: "" +L_TPXO_DATASET: "" +L_LEADTIMES_DATASET: "" +L_TRACK_DIR: "" +L_DEM_HI: "" +L_DEM_LO: "" +L_MESH_HI: "" +L_MESH_LO: "" +L_SHP_DIR: "" +NHC_OBS: "" + +TMPDIR: "/tmp" +PATH_APPEND: "" + +L_SOLVE_MODULES: + - "intel/2022.1.2" + - "impi/2022.1.2" + - "netcdf" diff --git a/tests/test_input_version.py b/tests/test_input_version.py index 5e2ae38..126b981 100644 --- a/tests/test_input_version.py +++ b/tests/test_input_version.py @@ -45,3 +45,8 @@ def test_v0_0_2_to_latest(conf_v0_0_2, conf_latest): def test_v0_0_3_to_latest(conf_v0_0_3, conf_latest): handle_input_version(conf_v0_0_3) assert conf_latest == conf_v0_0_3 + + +def test_v0_0_4_to_latest(conf_v0_0_4, conf_latest): + handle_input_version(conf_v0_0_4) + assert conf_latest == conf_v0_0_4 From 23c15f430bc08f08a495bfd575278a71f1baf835 Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar <132295102+FariborzDaneshvar-NOAA@users.noreply.github.com> Date: Fri, 27 Sep 2024 10:35:47 -0700 Subject: [PATCH 35/40] Update README.md Updated tech report citation, and added AGU23 presentation. --- README.md | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/README.md b/README.md index 5430d01..fa82219 100644 --- a/README.md +++ b/README.md @@ -8,10 +8,16 @@ TO BE COMPLETED... ## References -- Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E. +- Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E., (2024). _Probabilistic Prediction of Tropical Cyclone (TC) Driven Flood Depth_ -_and Extent with a User-Friendly Python Module_ -(NOAA Technical Memorandum NOS CS, under review) +_and Extent with a User-Friendly Python Module._ +NOAA Technical Memorandum NOS CS 59. + +- Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E., ( 2023). +_Evaluating the effect of inland hydrology and hurricane model on +the probabilistic prediction of tropical cyclone (TC) driven coastal flooding._ +American Geophysical Union (AGU) Fall Meeting, December 2023, San Francisco, CA. + - Pringle, W. J., Mani, S., Sargsyan, K., Moghimi, S., Zhang, Y. J., Khazaei, B., Myers, E. (January 2023). _Surrogate-Assisted Bayesian Uncertainty Quantification for From b166664a86a136898257024aa5bb39d9449f57de Mon Sep 17 00:00:00 2001 From: "Soroosh.Mani" Date: Wed, 9 Oct 2024 16:08:39 -0500 Subject: [PATCH 36/40] Manually add vegetation core param --- stormworkflow/prep/setup_ensemble.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/stormworkflow/prep/setup_ensemble.py b/stormworkflow/prep/setup_ensemble.py index fe7d889..4fa4ac1 100644 --- a/stormworkflow/prep/setup_ensemble.py +++ b/stormworkflow/prep/setup_ensemble.py @@ -117,6 +117,13 @@ def _fix_hotstart_issue(ensemble_dir): nm_list['opt']['drampwind'] = 0.0 nm_list.write(pth / 'param.nml', force=True) +def _fix_veg_parameter_issue(ensemble_dir): + # See https://github.com/schism-dev/pyschism/issues/126 + param_nmls = ensemble_dir.glob('**/param.nml') + for pth in param_nmls: + nm_list = f90nml.read(pth) + nm_list['core']['nbins_veg_vert'] = 2 + nm_list.write(pth, force=True) def main(args): @@ -277,6 +284,7 @@ def main(args): ) _fix_hotstart_issue(workdir) + _fix_veg_parameter_issue(workdir) # For newer SCHISM version if with_hydrology: _fix_nwm_issues(workdir, hires_reg) if use_wwm: From 2107bffd5e0de4dfb4e1d21d7d51314d6c83a39f Mon Sep 17 00:00:00 2001 From: Fariborz Daneshvar <132295102+FariborzDaneshvar-NOAA@users.noreply.github.com> Date: Tue, 15 Oct 2024 11:59:58 -0700 Subject: [PATCH 37/40] Update README.md (#81) Added URL for the technical document --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index fa82219..dcb8066 100644 --- a/README.md +++ b/README.md @@ -11,7 +11,7 @@ TO BE COMPLETED... - Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E., (2024). _Probabilistic Prediction of Tropical Cyclone (TC) Driven Flood Depth_ _and Extent with a User-Friendly Python Module._ -NOAA Technical Memorandum NOS CS 59. +NOAA Technical Memorandum NOS CS 59. https://repository.library.noaa.gov/view/noaa/66123 - Daneshvar, F., Mani, S., Pringle, W. J., Moghimi, S., Myers, E., ( 2023). _Evaluating the effect of inland hydrology and hurricane model on From ebe2eaee29f3ec4ac233ec7279cd415d6a9392b6 Mon Sep 17 00:00:00 2001 From: "Soroosh.Mani" Date: Tue, 22 Oct 2024 10:27:43 -0500 Subject: [PATCH 38/40] Add option to enable disable perturbation features, such as isotach adjustment --- pyproject.toml | 2 +- stormworkflow/main.py | 9 +++++ stormworkflow/prep/setup_ensemble.py | 7 +++- stormworkflow/refs/input.yaml | 4 ++- stormworkflow/scripts/workflow.sh | 6 ++-- tests/conftest.py | 5 +++ tests/data/refs/input_v0.0.5.yaml | 51 ++++++++++++++++++++++++++++ tests/test_input_version.py | 5 +++ 8 files changed, 83 insertions(+), 6 deletions(-) create mode 100644 tests/data/refs/input_v0.0.5.yaml diff --git a/pyproject.toml b/pyproject.toml index d85956d..446219d 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -35,7 +35,7 @@ dependencies = [ "cmocean", "dask", "dask-jobqueue", - "ensembleperturbation>=1.3.0", # rmax option + "ensembleperturbation>=1.3.4", # perturb feature "fiona", "geoalchemy2", "geopandas", diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 244a3e0..4d6f8e0 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -74,6 +74,15 @@ def _handle_input_v0_0_3_to_v0_0_4(inout_conf): inout_conf['NHC_OBS'] = '' +@_input_version('0.0.4', '0.0.5') +def _handle_input_v0_0_4_to_v0_0_5(inout_conf): + + _logger.info("Adding perturbation features") + inout_conf['perturb_features'] = [ + 'isotach_adjustment', + ] + + def handle_input_version(inout_conf): if 'input_version' not in inout_conf: diff --git a/stormworkflow/prep/setup_ensemble.py b/stormworkflow/prep/setup_ensemble.py index 4fa4ac1..4e78cb7 100644 --- a/stormworkflow/prep/setup_ensemble.py +++ b/stormworkflow/prep/setup_ensemble.py @@ -22,7 +22,7 @@ from coupledmodeldriver.generate import SCHISMRunConfiguration from coupledmodeldriver.generate.schism.script import SchismEnsembleGenerationJob from coupledmodeldriver.generate import generate_schism_configuration -from ensembleperturbation.perturbation.atcf import perturb_tracks +from ensembleperturbation.perturbation.atcf import perturb_tracks, PerturberFeatures from pylib_essentials.schism_file import ( read_schism_hgrid_cached, schism_bpfile, @@ -137,6 +137,9 @@ def main(args): use_wwm = args.use_wwm with_hydrology = args.with_hydrology pahm_model = args.pahm_model + setup_features = PerturberFeatures.NONE + for feat in args.perturb_features: + setup_features |= PerturberFeatures[feat.upper()] workdir = out_dir mesh_file = mesh_dir / 'mesh_w_bdry.grd' @@ -208,6 +211,7 @@ def main(args): overwrite=True, file_deck=file_deck, advisories=[advisory], + features=setup_features, ) if perturb_start != model_start_time: @@ -345,6 +349,7 @@ def parse_arguments(): argument_parser.add_argument('--use-wwm', action='store_true') argument_parser.add_argument('--with-hydrology', action='store_true') argument_parser.add_argument('--pahm-model', choices=['gahm', 'symmetric'], default='gahm') + argument_parser.add_argument('--perturb-features', nargs='+', type=str, default=[PerturberFeatures.ISOTACH_ADJUSTMENT.name]) argument_parser.add_argument('--variables', nargs='+', type=str) argument_parser.add_argument('name', help='name of the storm', type=str) diff --git a/stormworkflow/refs/input.yaml b/stormworkflow/refs/input.yaml index 30aefa2..3c629c1 100644 --- a/stormworkflow/refs/input.yaml +++ b/stormworkflow/refs/input.yaml @@ -1,5 +1,5 @@ --- -input_version: 0.0.4 +input_version: 0.0.5 storm: "florence" year: 2018 @@ -19,6 +19,8 @@ perturb_vars: - "radius_of_maximum_winds_persistent" - "max_sustained_wind_speed" rmw_fill_method: "persistent" +perturb_features: + - "isotach_adjustment" spinup_exec: "pschism_PAHM_TVD-VL" hotstart_exec: "pschism_PAHM_TVD-VL" diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index 6956d0e..83d8570 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -139,6 +139,7 @@ PREP_KWDS+=" --variables $perturb_vars" if [ $use_wwm == 1 ]; then PREP_KWDS+=" --use-wwm"; fi if [ $hydrology == 1 ]; then PREP_KWDS+=" --with-hydrology"; fi PREP_KWDS+=" --pahm-model $pahm_model" +PREP_KWDS+=" --perturb-features $perturb_features" export PREP_KWDS # NOTE: We need to wait because run jobs depend on perturbation dirs! setup_id=$(sbatch \ @@ -146,7 +147,7 @@ setup_id=$(sbatch \ --wait \ --job-name=prep_$tag \ --parsable \ - --export=ALL,PREP_KWDS,STORM=$storm,YEAR=$year,IMG="$L_IMG_DIR/prep.sif" \ + --export=ALL,PREP_KWDS,STORM=$storm,YEAR=$year \ $run_dir/slurm/prep.sbatch \ ) @@ -154,7 +155,6 @@ setup_id=$(sbatch \ echo "Launching runs" SCHISM_SHARED_ENV="" SCHISM_SHARED_ENV+="ALL" -SCHISM_SHARED_ENV+=",IMG=$L_IMG_DIR/solve.sif" SCHISM_SHARED_ENV+=",MODULES=$L_SOLVE_MODULES" spinup_id=$(sbatch \ --nodes $hpc_solver_nnodes --ntasks $hpc_solver_ntasks \ @@ -185,5 +185,5 @@ sbatch \ --output "${run_dir}/output/slurm-%j.post.out" \ --job-name=post_$tag \ -d afterok${joblist} \ - --export=ALL,IMG="$L_IMG_DIR/prep.sif",ENSEMBLE_DIR="$run_dir/setup/ensemble.dir/" \ + --export=ALL,ENSEMBLE_DIR="$run_dir/setup/ensemble.dir/" \ $run_dir/slurm/post.sbatch diff --git a/tests/conftest.py b/tests/conftest.py index 4d3c7e1..e62a5af 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -11,6 +11,7 @@ input_v0_0_2 = test_refs.joinpath('input_v0.0.2.yaml') input_v0_0_3 = test_refs.joinpath('input_v0.0.3.yaml') input_v0_0_4 = test_refs.joinpath('input_v0.0.4.yaml') +input_v0_0_5 = test_refs.joinpath('input_v0.0.5.yaml') input_latest = refs.joinpath('input.yaml') @@ -38,6 +39,10 @@ def conf_v0_0_3(): def conf_v0_0_4(): return read_conf(input_v0_0_4) +@pytest.fixture +def conf_v0_0_5(): + return read_conf(input_v0_0_5) + @pytest.fixture def conf_latest(): return read_conf(input_latest) diff --git a/tests/data/refs/input_v0.0.5.yaml b/tests/data/refs/input_v0.0.5.yaml new file mode 100644 index 0000000..3c629c1 --- /dev/null +++ b/tests/data/refs/input_v0.0.5.yaml @@ -0,0 +1,51 @@ +--- +input_version: 0.0.5 + +storm: "florence" +year: 2018 +suffix: "" +subset_mesh: 1 +hr_prelandfall: -1 +past_forecast: 1 +hydrology: 0 +use_wwm: 0 +pahm_model: "gahm" +num_perturb: 2 +sample_rule: "korobov" +perturb_vars: + - "cross_track" + - "along_track" +# - "radius_of_maximum_winds" + - "radius_of_maximum_winds_persistent" + - "max_sustained_wind_speed" +rmw_fill_method: "persistent" +perturb_features: + - "isotach_adjustment" + +spinup_exec: "pschism_PAHM_TVD-VL" +hotstart_exec: "pschism_PAHM_TVD-VL" + +hpc_solver_nnodes: 3 +hpc_solver_ntasks: 108 +hpc_account: "" +hpc_partition: "" + +RUN_OUT: "" +L_NWM_DATASET: "" +L_TPXO_DATASET: "" +L_LEADTIMES_DATASET: "" +L_TRACK_DIR: "" +L_DEM_HI: "" +L_DEM_LO: "" +L_MESH_HI: "" +L_MESH_LO: "" +L_SHP_DIR: "" +NHC_OBS: "" + +TMPDIR: "/tmp" +PATH_APPEND: "" + +L_SOLVE_MODULES: + - "intel/2022.1.2" + - "impi/2022.1.2" + - "netcdf" diff --git a/tests/test_input_version.py b/tests/test_input_version.py index 126b981..9d22429 100644 --- a/tests/test_input_version.py +++ b/tests/test_input_version.py @@ -50,3 +50,8 @@ def test_v0_0_3_to_latest(conf_v0_0_3, conf_latest): def test_v0_0_4_to_latest(conf_v0_0_4, conf_latest): handle_input_version(conf_v0_0_4) assert conf_latest == conf_v0_0_4 + + +def test_v0_0_5_to_latest(conf_v0_0_5, conf_latest): + handle_input_version(conf_v0_0_5) + assert conf_latest == conf_v0_0_5 From dbe85bbee3ee8ea996f0a0371f9c479530ceb13a Mon Sep 17 00:00:00 2001 From: "Soroosh.Mani" Date: Tue, 22 Oct 2024 11:03:06 -0500 Subject: [PATCH 39/40] current version value --- stormworkflow/main.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stormworkflow/main.py b/stormworkflow/main.py index 4d6f8e0..4b4c029 100644 --- a/stormworkflow/main.py +++ b/stormworkflow/main.py @@ -18,7 +18,7 @@ _logger = logging.getLogger(__file__) -CUR_INPUT_VER = Version('0.0.4') +CUR_INPUT_VER = Version('0.0.5') VER_UPDATE_FUNCS = [] From 265392811a386b98e192e6bd8573551c32362399 Mon Sep 17 00:00:00 2001 From: "Soroosh.Mani" Date: Fri, 25 Oct 2024 15:50:04 -0500 Subject: [PATCH 40/40] Put perturb feature (+) arg before pahm model arg (1) --- stormworkflow/scripts/workflow.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stormworkflow/scripts/workflow.sh b/stormworkflow/scripts/workflow.sh index 83d8570..2e22618 100755 --- a/stormworkflow/scripts/workflow.sh +++ b/stormworkflow/scripts/workflow.sh @@ -138,8 +138,8 @@ PREP_KWDS+=" --tpxo-dir $L_TPXO_DATASET" PREP_KWDS+=" --variables $perturb_vars" if [ $use_wwm == 1 ]; then PREP_KWDS+=" --use-wwm"; fi if [ $hydrology == 1 ]; then PREP_KWDS+=" --with-hydrology"; fi -PREP_KWDS+=" --pahm-model $pahm_model" PREP_KWDS+=" --perturb-features $perturb_features" +PREP_KWDS+=" --pahm-model $pahm_model" export PREP_KWDS # NOTE: We need to wait because run jobs depend on perturbation dirs! setup_id=$(sbatch \