diff --git a/.github/workflows/python-package.yml b/.github/workflows/python-package.yml index 5cea8fb..69f53f3 100644 --- a/.github/workflows/python-package.yml +++ b/.github/workflows/python-package.yml @@ -16,18 +16,28 @@ jobs: strategy: fail-fast: false matrix: - python-version: ["3.10", "3.11", "3.12", "3.13",] + include: + - python-version: "3.10" + pixi-env: py310dev + - python-version: "3.11" + pixi-env: py311dev + - python-version: "3.12" + pixi-env: py312dev + - python-version: "3.13" + pixi-env: py313dev + - python-version: "3.14" + pixi-env: py314dev steps: - - uses: actions/checkout@v4 - - name: Set up Python ${{ matrix.python-version }} - uses: actions/setup-python@v5 + - uses: actions/checkout@v6 + - name: Set up Pixi ${{ matrix.python-version }} + uses: prefix-dev/setup-pixi@v0.9.4 with: - python-version: ${{ matrix.python-version }} - - name: Install dependencies + pixi-version: v0.62.2 + # Project environment caching uses pixi.lock as the cache key input. + # If you don't commit a lockfile, disable caching to avoid ENOENT. + # cache: ${{ hashFiles('pixi.lock') != '' }} + environments: ${{ matrix.pixi-env }} + - name: Run tests run: | - python3 -m pip install --upgrade pip - python3 -m pip install -e .\[dev\] - - name: Test with pytest - run: | - pytest -vvv + pixi run -e ${{ matrix.pixi-env }} pytest -vvv diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 528d937..c67f382 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -10,6 +10,7 @@ repos: - id: debug-statements - id: destroyed-symlinks - id: detect-private-key + - repo: https://github.com/astral-sh/ruff-pre-commit # Ruff version. rev: v0.13.0 @@ -18,7 +19,3 @@ repos: - id: ruff-check # Run the formatter. - id: ruff-format - - repo: https://github.com/pycqa/isort - rev: 6.0.1 - hooks: - - id: isort diff --git a/README.md b/README.md index c20d467..152030d 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ ![](https://raw.githubusercontent.com/wiki/LIVVkit/LIVVkit/imgs/livvkit.png) -# LIVVkit Extensions (LEX) +# LIVVkit Extensions (LivvEX) This repository holds a collection of extensions to [LIVVkit](https://livvkit.github.io/Docs/index.html) for validation and @@ -19,41 +19,43 @@ machines, including Perlmutter at NERSC, and Chrysalis at ANL\'s LCRC. The Python package itself is described in `pyproject.toml`, which is used by `pip` to install this package -Currently, LEX is designed to run on Perlmutter, but future work is -planned to support other machines. +Currently, LEX is designed to run on NERSC's Perlmutter, and ANL-LCRC's Chrysalis, +but future work is planned to support other machines where E3SM runs. ## Environment setup For setting up an environment to which lex and dependencies will be -installed, conda and Python virtualenv are documented here. **NB** this -will only currently work on NERSC's Perlmutter, the environment should be +installed, Pixi and conda are documented here. **NB** this +will only currently work on Perlmutter and Chrysalis, the environment should be created there. -### Conda environment +### Pixi environment +[pixi](https://pixi.sh/latest/) is a package management tool, and the primary +environment management tool for LEX. +First, it must be [installed](https://pixi.prefix.dev/latest/installation/) locally. +Then an enviornment for lex development can be created: ```bash $ git clone https://github.com/LIVVkit/lex.git $ cd lex -$ {conda, mamba} create -n lex_env python --file requirements.txt -$ {conda, mamba} activate lex_env -$ pip install -e . # Installs the lex module as an editable Python package to the lex_env environment. +$ pixi install # This will install the default environment +$ pixi shell -e default # To activate the default (runtime) environment ``` -### Python virtualenv +If development tools are needed (e.g. `pytest`), install and activate the `dev` environment: ```bash -$ git clone https://github.com/LIVVkit/lex.git -$ cd lex python -m venv .env -$ source .env/bin/activate -$ pip install --upgrade pip # Needed if the system pip version < 21.3 -$ pip install -e . +$ pixi install -e dev # Installs the environment +$ pixi shell -e dev # Activates the dev environment (for code testing and checks) ``` +**N.B.** it is not recommended to do `pixi install --all`, as this will install the environment versions used for testing. -This will create a virtual environment at `lex/.env`, and install the -LEX package as editable with all its Python requirements to run. - -### Other available environment management solutions -Not documented here, but also available for environment management -- [uv](https://docs.astral.sh/uv/) -- [pixi](https://pixi.sh/latest/) +### Conda environment +```bash +$ git clone https://github.com/LIVVkit/lex.git +$ cd lex +$ {conda, mamba} env create -n lex_env python --file env.yml +$ {conda, mamba} activate lex_env +$ pip install -e . # Installs the lex module as an editable Python package to the lex_env environment. +``` ## Basic usage @@ -72,7 +74,9 @@ website you\'d run this command: ```bash $ livv -V config/example/example.yml -o vv_test -s ``` -This will create a directory in the current directory called `vv_test` (~7.5 MB), and spawn an HTTP server, which should only be used for testing purposes. (This works best if the output is in the current directory) +This will create a directory in the current directory called `vv_test` (~7.5 MB), and +spawn an HTTP server, which should only be used for testing purposes. +(This works best if the output is in the current directory) *Note:* All the extension configurations files assume you are working from the top level `lex` directory. You *can* run any of these @@ -103,20 +107,24 @@ $ sbatch run_lex_pm-cpu.sbatch ## Running new cases on PM-CPU ### Generate a single timeseries file from ELM h0 outputs -- `ncrcat -v topo,landfrac,QSNOFRZ,FSRND,FSRVD,FSDSVD,FSDSND,EFLX_LH_TOT,FIRA,FLDS,FSA,FSDS,FSH,QICE,QRUNOFF,QSNOMELT,QSOIL,RAIN,SNOW,TSA elm*h0*.nc -o ${CASE}.nc` +- `CASE="The case name"` +- `ncrcat -v topo,landfrac,QSNOFRZ,FSRND,FSRVD,FSDSVD,FSDSND,EFLX_LH_TOT,FIRA,FLDS,FSA,FSDS,FSH,QICE,QRUNOFF,QSNOMELT,QSOIL,RAIN,SNOW,TSA,SNOWICE,SNOWLIQ,H2OSNO ${CASE}.elm.h0*.nc -o ${CASE}.nc` ### Perform post-processing on a single time series ELM h0 output - Edit the `lex/lex/postproc/e3sm/postproc.sbatch` batch file to mach the new run Key variables: - `INDIR`: Path which contains single output time series file - - `OUTCASE`: Name of the new case which is the name of the netCDF file without extension (e.g. `v2.1.r025.IGERA5ELM_MLI-deep_firn_1980_2020`) + - `OUTCASE`: Name of the new case which is the name of the netCDF file without + extension (e.g. `v2.1.r025.IGERA5ELM_MLI-deep_firn_1980_2020`) - `RES`: ELM output resolution (currently accepts `R05` and `R025`) - - `OUTDIR`: Scratch directory into which climatology files will be written, defaults to `${SCRATCH}/lex/data/e3sm/${OUTCASE}` + - `OUTDIR`: Scratch directory into which climatology files will be written, + defaults to `${SCRATCH}/lex/data/e3sm/${OUTCASE}` - Run the post-processing script: - `cd lex/lex/postproc/e3sm; sbatch postproc.sbatch` -**NB**: the `postproc.sbatch` script will create the configuration for your case (based on `OUTCASE` and `OUTDIR`), then run LIVVkit on it with `lex/run_livv.sh` +**NB**: the `postproc.sbatch` script will create the configuration for your case +(based on `OUTCASE` and `OUTDIR`), then run LIVVkit on it with `lex/run_livv.sh` ## Developing a custom extension diff --git a/config/template_r025/all.jinja b/config/template_r025/all.jinja new file mode 100644 index 0000000..93fb066 --- /dev/null +++ b/config/template_r025/all.jinja @@ -0,0 +1,737 @@ +common: &common + meta: &meta + Case ID: [{{ case_id }}] + Climatology years: [1980-2020] + Model: [E3SM-ELM] + climo: {{ case_out_dir }}/{{ case_id }}.{clim}_mean.nc + latlon: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + elevation: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + lnd_climo: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + topo: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + glc_surf: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + latv: lat + lonv: lon + topov: topo + landfracv: landfrac + img_height: 300 + +common_smb: &common_smb + smbv: QICE + maskv: gis_mask2 + smbscale: 31536000 + cmap: BrBG + cmap_diff: BrBG + clim_even: 1 + units: mm w. e. yr^-1 + References: + - {{ livvproj_dir }}/data/smb/smb_icecores.bib + - {{ livvproj_dir }}/data/livvkit.bib + + +common_energy: &common_energy + desc: Surface energy balance component {component}{comment} from {data_var_names} + references: + - {{ livvproj_dir }}/data/racmo/Noel2015.bib + - {{ livvproj_dir }}/data/cism/glissade/cism-glissade.bib + - {{ livvproj_dir }}/data/e3sm/Evans2019.bib + - {{ livvproj_dir }}/data/livvkit.bib + cmap: plasma + cmap_diff: RdBu_r + +common_racmo: &common_racmo + clim_years: + dset_a: { year_s: 1980, year_e: 2020 } + dataset_names: { model: ELM, dset_a: RACMO 2.4 } + in_dirs: { dset_a: "{{ racmo_root_dir }}/clm/{_var}" } + file_patterns: + dset_a: "{_var}_{icesheet}_{season}_{sea_s}_{sea_e}_climo.nc" + timeseries_dirs: + model: "{{ model_ts_dir }}" + dset_a: "{{ racmo_root_dir }}/ts" + ts_file_patterns: + model: "{{ case_id }}.nc" + dset_a: "{_var}_{icesheet}_{year_s}01_{year_e}12.nc" + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_rcm.nc + model: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_rcm.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r025.nc + +common_era5: &common_era5 + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_era5.nc + dataset_names: {model_remap: 'ELM [ERA5 Grid]', model_native: ELM, dset_a: ERA5} + aavg_sort: [ELM, 'ELM [ERA5 Grid]', ERA5] + scales: {model: 1, dset_a: 1} + clim_years: + dset_a: { year_s: 1979, year_e: 2019} + in_dirs: {dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/ERA5} + file_patterns: {dset_a: 'ERA5_{season}_{sea_s}_{sea_e}_climo.nc'} + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_era5.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_era5.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r025.nc + +common_merra: &common_merra + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_merra2.nc + in_dirs: {dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/MERRA2} + file_patterns: {dset_a: 'MERRA2_{season}_{sea_s}_{sea_e}_climo.nc'} + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_merra2.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_merra2.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r025.nc + dataset_names: {model_remap: 'ELM [MERRA2 Grid]', model_native: ELM, dset_a: MERRA2} + aavg_sort: [ELM, 'ELM [MERRA2 Grid]', MERRA2] + scales: {model: 1, dset_a: 1} + clim_years: + dset_a: { year_s: 1980, year_e: 2016} + +common_ceres: &common_ceres + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_cmip6.nc + dataset_names: + { model_native: ELM, model_remap: "ELM [CERES Grid]", dset_a: CERES } + aavg_sort: [ELM, "ELM [CERES Grid]", CERES] + scales: { model: 1, dset_a: 1 } + clim_years: + dset_a: { year_s: 2001, year_e: 2018 } + in_dirs: + dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/ceres_ebaf_surface_v4.1 + file_patterns: + dset_a: ceres_ebaf_surface_v4.1_{season}_{sea_s}_{sea_e}_climo.nc + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_cmip6.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_cmip6.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r025.nc + +common_racmo_gis: &common_racmo_gis + core_year_s: 1980 + core_year_e: 2021 + ib_year_e: 1987-1976 + ib_year_s: 2014-2004 + preprocess: [] + preproc_dir: {{ livvproj_dir }}/data/smb/processed + zwally_file: model_zwally_basins_elm_r025.csv + ib_file: IceBridge_modelXY_elm{}_r025.csv + smb_cf_file: SMB_CoreFirnEstimates_elm{}_r025.csv + smb_mo_file: SMB_Obs_Model_elm{}_r025.csv + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_rcmgis.nc + icesheet: gis + +common_racmo_ais: &common_racmo_ais + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_rcmais.nc + icesheet: ais + +common_racmo_smb_vars: &common_racmo_smb_vars + # Surface mass balance variables for dset_a: RACMO 2.4, model: ELM + # Each field must have: + # - title: Human readable field name + # - dset_a: definition for the field in dset_a dataset (i.e. RACMO) + # - model: definition for the field in model dataset (i.e. ELM) + # - units: Units of this field + # - ac_contrib_sign: For each dataset, multiply the annual cycle by +/- 1 to achieve + # correct contribution to annual cycle of SMB + # - aavg: parameters for ice sheet-area-averaging + # - scale: Further scale the area-average, (this typically converts mm w.e./yr to GT/yr) + # - units: Units of area average (only applicable if different than units for the field) + # - sum: Perform area-weighted sum if true (area weighted mean if false) + + # Additional parameters which may be defined are: + # - mask_weight: mask variable acts as a weight (e.g. partial mask at edge of the ice sheet) + # - primary_var: True for the primary field of interest for an annual cycle plot + # - cmap: Override for colourmap for the fields + # - cmap_diff: Override for colourmap for the difference plot + # - cmin, cmax, cmin_d, cmax_d: Override for min / max values for fields and differences respectively + # - comment: added to figure caption on output (usually used to indicate signedness of fluxes) + + # Surface mass balance + - title: Surface Mass Balance + dset_a: smbgl + model: QICE + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + mask_weight: true + primary_var: true + + # Total precipitation + - title: Total precip + dset_a: prgl + model: [+, SNOW, RAIN] + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + + # Snowfall + - title: Snowfall + dset_a: sf + model: SNOW + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Rainfall + - title: Rain + dset_a: [+, crrate, lsrrate] + model: RAIN + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Re-freeze + - title: Re-freeze + dset_a: rfrzgl + model: QSNOFRZ + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + cmin: 0 + cmap: YlGnBu + + # Runoff + - title: Runoff + dset_a: totrunoff + model: QRUNOFF + ac_contrib_sign: { model: -1, dset_a: -1 } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + # cmin: -20.0 + cmin: 0 + # cmax: 20.0 + cmap: YlGnBu + + # Snow & ice melt + - title: Snow + ice melt + dset_a: mltgl + model: ["+", QSNOMELT, QICE_MELT] + ac_contrib_sign: { model: -1, dset_a: -1 } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + # cmin: -20.0 + # cmax: 20.0 + cmin: 0 + # cmax: 20.0 + cmap: YlGnBu + + # Sublimation + - title: Sublimation + dset_a: sublgl + model: QSOIL + comment: " (Positive to atmosphere)" + ac_contrib_sign: { model: -1, dset_a: 1 } + scales: { model: 365 * 24 * 3600, dset_a: -365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + +common_racmo_cmb_vars: &common_racmo_cmb_vars + # Climatic mass balance variables for dset_a: RACMO 2.4, model: ELM + + # Climatic mass balance + - title: Climatic Mass Balance + dset_a: ["-", "prgl", ["-", "totrunoff", "sublgl"]] + model: ["-", ["+", "SNOW", "RAIN"], ["+", "QRUNOFF", "QSOIL"]] + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + mask_weight: true + primary_var: true + + # Snowfall + - title: Snowfall + dset_a: sf + model: SNOW + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Rainfall + - title: Rain + dset_a: [+, crrate, lsrrate] + model: RAIN + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Runoff + - title: Runoff + dset_a: totrunoff + model: QRUNOFF + ac_contrib_sign: { model: -1, dset_a: -1, } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + cmin: 0 + cmap: YlGnBu + + # Sublimation + - title: Sublimation + dset_a: sublgl + model: QSOIL + comment: " (Positive to atmosphere)" + ac_contrib_sign: { model: -1, dset_a: 1, } + scales: + { model: 365 * 24 * 3600, dset_a: -365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + +common_racmo_energy_vars: &common_racmo_energy_vars + - title: Surface temperature + dset_a: tas + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsusgl, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 0.9 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsusgl] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave net + dset_a: strgl + model: FIRA + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Sensible heat + dset_a: hfssgl + model: FSH + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hflsgl + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - + + - ['-', rsds, rsusgl] + - - + + - strgl + - [+, hfssgl, hflsgl] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_era_energy_vars: &common_era_energy_vars + - title: Surface temperature + dset_a: ts + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 1.0 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsus] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - model: FIRA + title: Longwave net + dset_a: ['-', rlus, rlds] + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - model: FSH + dset_a: hfss + title: Sensible heat + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hfls + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - '-' + - ['-', rsds, rsus] + - - + + - ['-', rlus, rlds] + - [+, hfss, hfls] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_merra_energy_vars: &common_merra_energy_vars + - title: Surface temperature + dset_a: tas + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 0.83 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsus] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave net + dset_a: ['-', rlus, rlds] + model: FIRA + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Sensible heat + dset_a: hfss + dset_b: FSH + model: FSH + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hfls + dset_b: EFLX_LH_TOT + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - '-' + - ['-', rsds, rsus] + - - + + - ['-', rlus, rlds] + - [+, hfss, hfls] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_ceres_energy_vars: &common_ceres_energy_vars + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: { dset_a: 1, model: 1 } + units: unitless + cmin: 0.5 + cmax: 0.83 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: sfc_net_sw_all_mon + model: FSA + comment: " (Positive to surface)" + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + + - title: Longwave net + dset_a: sfc_net_lw_all_mon + model: FIRA + comment: " (Positive to atmosphere)" + sign: 1 + scales: { dset_a: -1, model: 1 } + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + +{% if run_gis %} +# Greenland +{% if set_cmb %} +Climatic_Mass_Balance_GIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_cmb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_gis] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: Climatic Mass Balance + desc: "{component} component of CMB from {data_var_names}" +{% endif %} + +{% if set_smb %} +Surface_Mass_Balance_GIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_smb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_gis] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: smbgl + desc: "{component} component of SMB from {data_var_names}" +{% endif %} + +{% if set_energy_racmo %} +Energy_Balance_RACMO_GIS: + module: lex/energy/energy.py + <<: [*common, *common_energy, *common_racmo, *common_racmo_gis] + scales: {model: 1, dset_a: 1} + data_vars: *common_racmo_energy_vars +{% endif %} + +{% if set_energy_era5 %} +Energy_Balance_ERA5_GIS: + module: lex/energy/energy.py + icesheet: gis + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + <<: [*common, *common_energy, *common_era5] + data_vars: *common_era_energy_vars +{% endif %} + +{% if set_energy_merra2 %} +Energy_Balance_MERRA2_GIS: + module: lex/energy/energy.py + icesheet: gis + <<: [*common, *common_energy, *common_merra] + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + data_vars: *common_merra_energy_vars +{% endif %} + +{% if set_energy_ceres %} +Energy_Balance_CERES_GIS: + module: lex/energy/energy.py + icesheet: gis + <<: [*common, *common_energy, *common_ceres] + mask_ocean: { model_native: false, model_remap: true, dset_a: true } + data_vars: *common_ceres_energy_vars +{% endif %} +{% endif %} + +{% if run_ais %} +# Antarctica +{% if set_cmb %} +Climatic_Mass_Balance_AIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_cmb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_ais] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: Climatic Mass Balance + desc: "{component} component of CMB from {data_var_names}" +{% endif %} + +{% if set_smb %} +Surface_Mass_Balance_AIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_smb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_ais] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: smbgl + desc: "{component} component of SMB from {data_var_names}" +{% endif %} + +{% if set_energy_racmo %} +Energy_Balance_RACMO_AIS: + module: lex/energy/energy.py + <<: [*common, *common_energy, *common_racmo, *common_racmo_ais] + scales: {model: 1, dset_a: 1} + data_vars: *common_racmo_energy_vars +{% endif %} + +{% if set_energy_era5 %} +Energy_Balance_ERA5_AIS: + module: lex/energy/energy.py + icesheet: ais + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + <<: [*common, *common_energy, *common_era5] + data_vars: *common_era_energy_vars +{% endif %} + +{% if set_energy_merra2 %} +Energy_Balance_MERRA2_AIS: + module: lex/energy/energy.py + icesheet: ais + <<: [*common, *common_energy, *common_merra] + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + data_vars: *common_merra_energy_vars +{% endif %} + +{% if set_energy_ceres %} +Energy_Balance_CERES_AIS: + module: lex/energy/energy.py + icesheet: ais + <<: [*common, *common_energy, *common_ceres] + mask_ocean: { model_native: false, model_remap: true, dset_a: true } + data_vars: *common_ceres_energy_vars +{% endif %} +{% endif %} diff --git a/config/template_r05/all.jinja b/config/template_r05/all.jinja new file mode 100644 index 0000000..06b4a7d --- /dev/null +++ b/config/template_r05/all.jinja @@ -0,0 +1,737 @@ +common: &common + meta: &meta + Case ID: [{{ case_id }}] + Climatology years: [1980-2020] + Model: [E3SM-ELM] + climo: {{ case_out_dir }}/{{ case_id }}.{clim}_mean.nc + latlon: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + elevation: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + lnd_climo: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + topo: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + glc_surf: {{ case_out_dir }}/{{ case_id }}.ANN_mean.nc + latv: lat + lonv: lon + topov: topo + landfracv: landfrac + img_height: 300 + +common_smb: &common_smb + smbv: QICE + maskv: gis_mask2 + smbscale: 31536000 + cmap: BrBG + cmap_diff: BrBG + clim_even: 1 + units: mm w. e. yr^-1 + References: + - {{ livvproj_dir }}/data/smb/smb_icecores.bib + - {{ livvproj_dir }}/data/livvkit.bib + + +common_energy: &common_energy + desc: Surface energy balance component {component}{comment} from {data_var_names} + references: + - {{ livvproj_dir }}/data/racmo/Noel2015.bib + - {{ livvproj_dir }}/data/cism/glissade/cism-glissade.bib + - {{ livvproj_dir }}/data/e3sm/Evans2019.bib + - {{ livvproj_dir }}/data/livvkit.bib + cmap: plasma + cmap_diff: RdBu_r + +common_racmo: &common_racmo + clim_years: + dset_a: { year_s: 1980, year_e: 2020 } + dataset_names: { model: ELM, dset_a: RACMO 2.4 } + in_dirs: { dset_a: "{{ racmo_root_dir }}/clm/{_var}" } + file_patterns: + dset_a: "{_var}_{icesheet}_{season}_{sea_s}_{sea_e}_climo.nc" + timeseries_dirs: + model: "{{ model_ts_dir }}" + dset_a: "{{ racmo_root_dir }}/ts" + ts_file_patterns: + model: "{{ case_id }}.nc" + dset_a: "{_var}_{icesheet}_{year_s}01_{year_e}12.nc" + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_rcm.nc + model: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_rcm.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r05.nc + +common_era5: &common_era5 + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_era5.nc + dataset_names: {model_remap: 'ELM [ERA5 Grid]', model_native: ELM, dset_a: ERA5} + aavg_sort: [ELM, 'ELM [ERA5 Grid]', ERA5] + scales: {model: 1, dset_a: 1} + clim_years: + dset_a: { year_s: 1979, year_e: 2019} + in_dirs: {dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/ERA5} + file_patterns: {dset_a: 'ERA5_{season}_{sea_s}_{sea_e}_climo.nc'} + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_era5.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_era5.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r05.nc + +common_merra: &common_merra + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_merra2.nc + in_dirs: {dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/MERRA2} + file_patterns: {dset_a: 'MERRA2_{season}_{sea_s}_{sea_e}_climo.nc'} + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_merra2.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_merra2.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r05.nc + dataset_names: {model_remap: 'ELM [MERRA2 Grid]', model_native: ELM, dset_a: MERRA2} + aavg_sort: [ELM, 'ELM [MERRA2 Grid]', MERRA2] + scales: {model: 1, dset_a: 1} + clim_years: + dset_a: { year_s: 1980, year_e: 2016} + +common_ceres: &common_ceres + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_cmip6.nc + dataset_names: + { model_native: ELM, model_remap: "ELM [CERES Grid]", dset_a: CERES } + aavg_sort: [ELM, "ELM [CERES Grid]", CERES] + scales: { model: 1, dset_a: 1 } + clim_years: + dset_a: { year_s: 2001, year_e: 2018 } + in_dirs: + dset_a: {{ e3sm_diags_data_dir }}/observations/Atm/climatology/ceres_ebaf_surface_v4.1 + file_patterns: + dset_a: ceres_ebaf_surface_v4.1_{season}_{sea_s}_{sea_e}_climo.nc + masks: + dset_a: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_cmip6.traave.nc + model_remap: {{ livvproj_dir }}/grids/msk_{icesheet}_r025_remap_cmip6.traave.nc + model_native: {{ livvproj_dir }}/grids/msk_{icesheet}_rcm_r05.nc + +common_racmo_gis: &common_racmo_gis + core_year_s: 1980 + core_year_e: 2021 + ib_year_e: 1987-1976 + ib_year_s: 2014-2004 + preprocess: [] + preproc_dir: {{ livvproj_dir }}/data/smb/processed + zwally_file: model_zwally_basins_elm_r05.csv + ib_file: IceBridge_modelXY_elm{}_r05.csv + smb_cf_file: SMB_CoreFirnEstimates_elm{}_r05.csv + smb_mo_file: SMB_Obs_Model_elm{}_r05.csv + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_rcmgis.nc + icesheet: gis + +common_racmo_ais: &common_racmo_ais + climo_remap: {{ case_out_dir }}/remap/{{ case_id }}.{clim}_mean_rcmais.nc + icesheet: ais + +common_racmo_smb_vars: &common_racmo_smb_vars + # Surface mass balance variables for dset_a: RACMO 2.4, model: ELM + # Each field must have: + # - title: Human readable field name + # - dset_a: definition for the field in dset_a dataset (i.e. RACMO) + # - model: definition for the field in model dataset (i.e. ELM) + # - units: Units of this field + # - ac_contrib_sign: For each dataset, multiply the annual cycle by +/- 1 to achieve + # correct contribution to annual cycle of SMB + # - aavg: parameters for ice sheet-area-averaging + # - scale: Further scale the area-average, (this typically converts mm w.e./yr to GT/yr) + # - units: Units of area average (only applicable if different than units for the field) + # - sum: Perform area-weighted sum if true (area weighted mean if false) + + # Additional parameters which may be defined are: + # - mask_weight: mask variable acts as a weight (e.g. partial mask at edge of the ice sheet) + # - primary_var: True for the primary field of interest for an annual cycle plot + # - cmap: Override for colourmap for the fields + # - cmap_diff: Override for colourmap for the difference plot + # - cmin, cmax, cmin_d, cmax_d: Override for min / max values for fields and differences respectively + # - comment: added to figure caption on output (usually used to indicate signedness of fluxes) + + # Surface mass balance + - title: Surface Mass Balance + dset_a: smbgl + model: QICE + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + mask_weight: true + primary_var: true + + # Total precipitation + - title: Total precip + dset_a: prgl + model: [+, SNOW, RAIN] + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + + # Snowfall + - title: Snowfall + dset_a: sf + model: SNOW + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Rainfall + - title: Rain + dset_a: [+, crrate, lsrrate] + model: RAIN + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Re-freeze + - title: Re-freeze + dset_a: rfrzgl + model: QSNOFRZ + ac_contrib_sign: { model: 1, dset_a: 1 } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + cmin: 0 + cmap: YlGnBu + + # Runoff + - title: Runoff + dset_a: totrunoff + model: QRUNOFF + ac_contrib_sign: { model: -1, dset_a: -1 } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + # cmin: -20.0 + cmin: 0 + # cmax: 20.0 + cmap: YlGnBu + + # Snow & ice melt + - title: Snow + ice melt + dset_a: mltgl + model: ["+", QSNOMELT, QICE_MELT] + ac_contrib_sign: { model: -1, dset_a: -1 } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + # cmin: -20.0 + # cmax: 20.0 + cmin: 0 + # cmax: 20.0 + cmap: YlGnBu + + # Sublimation + - title: Sublimation + dset_a: sublgl + model: QSOIL + comment: " (Positive to atmosphere)" + ac_contrib_sign: { model: -1, dset_a: 1 } + scales: { model: 365 * 24 * 3600, dset_a: -365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + +common_racmo_cmb_vars: &common_racmo_cmb_vars + # Climatic mass balance variables for dset_a: RACMO 2.4, model: ELM + + # Climatic mass balance + - title: Climatic Mass Balance + dset_a: ["-", "prgl", ["-", "totrunoff", "sublgl"]] + model: ["-", ["+", "SNOW", "RAIN"], ["+", "QRUNOFF", "QSOIL"]] + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + mask_weight: true + primary_var: true + + # Snowfall + - title: Snowfall + dset_a: sf + model: SNOW + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Rainfall + - title: Rain + dset_a: [+, crrate, lsrrate] + model: RAIN + ac_contrib_sign: { model: 1, dset_a: 1, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + units: mm w. e. yr^-1 + cmap: YlGnBu + cmin: 0 + mask_weight: true + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + + # Runoff + - title: Runoff + dset_a: totrunoff + model: QRUNOFF + ac_contrib_sign: { model: -1, dset_a: -1, } + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + cmin: 0 + cmap: YlGnBu + + # Sublimation + - title: Sublimation + dset_a: sublgl + model: QSOIL + comment: " (Positive to atmosphere)" + ac_contrib_sign: { model: -1, dset_a: 1, } + scales: + { model: 365 * 24 * 3600, dset_a: -365 * 24 * 3600, } + aavg: { scale: 1e-06, units: GT yr^-1, sum: true } + mask_weight: true + units: mm w. e. yr^-1 + +common_racmo_energy_vars: &common_racmo_energy_vars + - title: Surface temperature + dset_a: tas + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsusgl, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 0.9 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsusgl] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave net + dset_a: strgl + model: FIRA + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Sensible heat + dset_a: hfssgl + model: FSH + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hflsgl + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: -1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - + + - ['-', rsds, rsusgl] + - - + + - strgl + - [+, hfssgl, hflsgl] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_era_energy_vars: &common_era_energy_vars + - title: Surface temperature + dset_a: ts + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 1.0 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsus] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - model: FIRA + title: Longwave net + dset_a: ['-', rlus, rlds] + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - model: FSH + dset_a: hfss + title: Sensible heat + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hfls + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - '-' + - ['-', rsds, rsus] + - - + + - ['-', rlus, rlds] + - [+, hfss, hfls] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_merra_energy_vars: &common_merra_energy_vars + - title: Surface temperature + dset_a: tas + model: TSA + sign: 1 + scales: {dset_a: 1, model: 1} + units: K + + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: {dset_a: 1, model: 1} + units: unitless + cmin: 0.5 + cmax: 0.83 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: ['-', rsds, rsus] + model: FSA + comment: ' (Positive to surface)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave net + dset_a: ['-', rlus, rlds] + model: FIRA + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Sensible heat + dset_a: hfss + dset_b: FSH + model: FSH + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Latent Heat Flux + dset_a: hfls + dset_b: EFLX_LH_TOT + model: EFLX_LH_TOT + comment: ' (Positive to atmosphere)' + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + + - title: Net energy balance + dset_a: + - '-' + - ['-', rsds, rsus] + - - + + - ['-', rlus, rlds] + - [+, hfss, hfls] + model: + - '-' + - FSA + - - + + - FIRA + - [+, FSH, EFLX_LH_TOT] + sign: 1 + scales: {dset_a: 1, model: 1} + units: W m^-2 + cmin: -10 + cmax: 10 + cmap: RdBu_r + +common_ceres_energy_vars: &common_ceres_energy_vars + - title: Albedo + dset_a: [/, rsus, rsds] + model: + - / + - [+, FSRND, FSRVD] + - [+, FSDSVD, FSDSND] + sign: 1 + scales: { dset_a: 1, model: 1 } + units: unitless + cmin: 0.5 + cmax: 0.83 + cmin_d: -0.2 + cmax_d: 0.2 + + - title: Shortwave net + dset_a: sfc_net_sw_all_mon + model: FSA + comment: " (Positive to surface)" + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + + - title: Shortwave down + dset_a: rsds + model: FSDS + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + + - title: Longwave net + dset_a: sfc_net_lw_all_mon + model: FIRA + comment: " (Positive to atmosphere)" + sign: 1 + scales: { dset_a: -1, model: 1 } + units: W m^-2 + + - title: Longwave down + dset_a: rlds + model: FLDS + sign: 1 + scales: { dset_a: 1, model: 1 } + units: W m^-2 + +{% if run_gis %} +# Greenland +{% if set_cmb %} +Climatic_Mass_Balance_GIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_cmb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_gis] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: Climatic Mass Balance + desc: "{component} component of CMB from {data_var_names}" +{% endif %} + +{% if set_smb %} +Surface_Mass_Balance_GIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_smb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_gis] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: smbgl + desc: "{component} component of SMB from {data_var_names}" +{% endif %} + +{% if set_energy_racmo %} +Energy_Balance_RACMO_GIS: + module: lex/energy/energy.py + <<: [*common, *common_energy, *common_racmo, *common_racmo_gis] + scales: {model: 1, dset_a: 1} + data_vars: *common_racmo_energy_vars +{% endif %} + +{% if set_energy_era5 %} +Energy_Balance_ERA5_GIS: + module: lex/energy/energy.py + icesheet: gis + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + <<: [*common, *common_energy, *common_era5] + data_vars: *common_era_energy_vars +{% endif %} + +{% if set_energy_merra2 %} +Energy_Balance_MERRA2_GIS: + module: lex/energy/energy.py + icesheet: gis + <<: [*common, *common_energy, *common_merra] + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + data_vars: *common_merra_energy_vars +{% endif %} + +{% if set_energy_ceres %} +Energy_Balance_CERES_GIS: + module: lex/energy/energy.py + icesheet: gis + <<: [*common, *common_energy, *common_ceres] + mask_ocean: { model_native: false, model_remap: true, dset_a: true } + data_vars: *common_ceres_energy_vars +{% endif %} +{% endif %} + +{% if run_ais %} +# Antarctica +{% if set_cmb %} +Climatic_Mass_Balance_AIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_cmb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_ais] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: Climatic Mass Balance + desc: "{component} component of CMB from {data_var_names}" +{% endif %} + +{% if set_smb %} +Surface_Mass_Balance_AIS: + module: lex/smb/smb_icecores.py + data_vars: *common_racmo_smb_vars + <<: [*common, *common_smb, *common_racmo, *common_racmo_ais] + scales: + { model: 365 * 24 * 3600, dset_a: 365 * 24 * 3600 } + primary_var: smbgl + desc: "{component} component of SMB from {data_var_names}" +{% endif %} + +{% if set_energy_racmo %} +Energy_Balance_RACMO_AIS: + module: lex/energy/energy.py + <<: [*common, *common_energy, *common_racmo, *common_racmo_ais] + scales: {model: 1, dset_a: 1} + data_vars: *common_racmo_energy_vars +{% endif %} + +{% if set_energy_era5 %} +Energy_Balance_ERA5_AIS: + module: lex/energy/energy.py + icesheet: ais + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + <<: [*common, *common_energy, *common_era5] + data_vars: *common_era_energy_vars +{% endif %} + +{% if set_energy_merra2 %} +Energy_Balance_MERRA2_AIS: + module: lex/energy/energy.py + icesheet: ais + <<: [*common, *common_energy, *common_merra] + mask_ocean: {model_native: false, model_remap: true, dset_a: true} + data_vars: *common_merra_energy_vars +{% endif %} + +{% if set_energy_ceres %} +Energy_Balance_CERES_AIS: + module: lex/energy/energy.py + icesheet: ais + <<: [*common, *common_energy, *common_ceres] + mask_ocean: { model_native: false, model_remap: true, dset_a: true } + data_vars: *common_ceres_energy_vars +{% endif %} +{% endif %} diff --git a/env.yml b/env.yml new file mode 100644 index 0000000..93ab51f --- /dev/null +++ b/env.yml @@ -0,0 +1,32 @@ +# +# This file is autogenerated by pyproject2conda +# with the following command: +# +# $ pyproject2conda yaml -o env.yml +# +# You should not manually edit this file. +# Instead edit the corresponding pyproject.toml file. +# +channels: + - conda-forge +dependencies: + - cartopy>=0.25.0 + - cftime>=1.6.5 + - dask>=2026.1.2 + - f90nml>=1.5 + - jinja2>=3.1.6 + - livvkit>=3.3.0 + - loguru>=0.7.3 + - mache>=2.1.0 + - matplotlib>=3.10.8 + - nc-time-axis>=1.4.1 + - netcdf4>=1.7.4 + - nco>=5.3.6 + - numpy>=2.2.6 + - pandas>=2.2.3 + - pybtex>=0.25.1 + - ruamel.yaml>=0.19.1 + - scikit-learn>=1.7.2 + - scipy>=1.15.3 + - seaborn>=0.13.2 + - xarray>=2025.6.1 diff --git a/lex/annual_cycle.py b/lex/annual_cycle.py index 5a34470..15e8358 100644 --- a/lex/annual_cycle.py +++ b/lex/annual_cycle.py @@ -6,6 +6,7 @@ import pandas as pd import xarray as xr from livvkit import elements as el +from loguru import logger import lex.common as lxc @@ -18,8 +19,9 @@ def main(args, config): + _files = [lxc.proc_climo_file(config, "climo_remap", mon) for mon in range(1, 13)] model_data = xr.open_mfdataset( - [config["climo_remap"].format(clim=f"{mon:02d}") for mon in range(1, 13)], + _files, combine="nested", concat_dim="time", ) @@ -40,6 +42,7 @@ def main(args, config): model_data_out = {} for idx, data_var in enumerate(config["data_vars"]): + logger.info(f"WORKING ON {data_var['title']}") _obs_in = {} aavg_config = data_var.get("aavg", None) @@ -64,7 +67,7 @@ def main(args, config): / 365 ) except KeyError: - print(f"MODEL DATA NOT FOUND FOR {data_var['model']}") + logger.error(f"MODEL DATA NOT FOUND FOR {data_var['model']}") continue obs_aavg[data_var["title"]] = {} @@ -77,7 +80,9 @@ def main(args, config): lxc.area_avg( _obs_in[_vers], {}, - area_file=config["masks"][_vers], + area_file=config["masks"][_vers].format( + icesheet=config["icesheet"] + ), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -86,7 +91,7 @@ def main(args, config): model_aavg[data_var["title"]], _, _, _ = lxc.area_avg( _model_plt, {}, - area_file=config["masks"]["model"], + area_file=config["masks"]["model"].format(icesheet=config["icesheet"]), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -125,6 +130,7 @@ def main(args, config): marker=".", lw=lw, ) + logger.info(f"DONE - WORKING ON {data_var['title']}") obs_data_out["month"] = np.arange(1, 12 + 1) model_data_out["month"] = np.arange(1, 12 + 1) diff --git a/lex/common.py b/lex/common.py index 3c28c00..5d1cc63 100644 --- a/lex/common.py +++ b/lex/common.py @@ -18,6 +18,7 @@ from cartopy import feature as cfeature from livvkit import elements as el from livvkit.util.LIVVDict import LIVVDict +from loguru import logger from numpy import ma import lex.utils as lxu @@ -81,15 +82,18 @@ def check_longitude(data, lon_coord="lon"): def get_season_bounds(season, year_s, year_e, sep="_"): """Determine season bounds for climatology files.""" - _seasons = {"DJF": (1, 12), "MAM": (3, 5), "JJA": (6, 8), "SON": (9, 11)} + _seasons = { + "DJF": (1, 12), + "MAM": (3, 5), + "JJA": (6, 8), + "SON": (9, 11), + "ANN": (1, 12), + } _annual = (1, 12) if season in _seasons: _lb, _ub = _seasons[season] bound_l = f"{year_s:04d}{_lb:02d}" bound_u = f"{year_e:04d}{_ub:02d}" - elif season.upper() == "ANN": - bound_l = f"{year_s:04d}01" - bound_u = f"{year_e:04d}12" else: # Months if isinstance(season, str): @@ -101,6 +105,42 @@ def get_season_bounds(season, year_s, year_e, sep="_"): return bound_l, bound_u +def proc_climo_file(config, file_tag, sea): + """ + Process the climatology file to maintain backward compatibility with standalone LEX. + + Parameters + ---------- + config : dict + LIVVkit /LEX configuration dict + file_tag : str + Configuration item which points to climatology filename to be formatted, usually + `climo` or `climo_remap` + sea : str + Season identifier + + Returns + ------- + climo_file : str + Formatted name of climatology file + + """ + _filename = config[file_tag] + if "{sea_s}" in _filename: + sea_s, sea_e = get_season_bounds( + sea, config.get("year_s", None), config.get("year_e", None) + ) + if isinstance(sea, int): + sea = f"{sea:02d}" + climo_file = _filename.format(clim=sea, sea_s=sea_s, sea_e=sea_e) + else: + if isinstance(sea, int): + sea = f"{sea:02d}" + climo_file = _filename.format(clim=sea) + + return climo_file + + def get_cycle(sea): """Get the name of which type of averaging is used (annual, seasonal, monthly).""" if sea == "ANN": @@ -343,7 +383,10 @@ def _fcn_filt(_var): clim_years = config.get("clim_years", None) - clim_years_default = {"year_s": None, "year_e": None} + clim_years_default = { + "year_s": config.get("year_s", None), + "year_e": config.get("year_e", None), + } if clim_years: clim_years = clim_years.get(overs, clim_years_default) else: @@ -377,6 +420,7 @@ def _fcn_filt(_var): return var_files +@logger.catch def load_timeseries_data(config): """Load data for timeseries.""" files = {} @@ -391,17 +435,28 @@ def load_timeseries_data(config): if len(set(files[overs])) == 1: files[overs] = files[overs][0] - + _nfiles = 1 + else: + _nfiles = len(set(files[overs])) + _dsname = config["dataset_names"].get( + overs, config["dataset_names"].get("model_native") + ) + logger.info(f"LOAD TIMESERIES DATA FOR {overs}: {_dsname} NFILES: {_nfiles}") try: - obs_data[overs] = xr.open_mfdataset(files[overs]).squeeze() + obs_data[overs] = xr.open_mfdataset(files[overs]).squeeze().load() except (xr.MergeError, ValueError): if isinstance(files[overs], Path): - obs_data[overs] = xr.open_dataset(files[overs]).squeeze() + obs_data[overs] = xr.open_dataset(files[overs]).squeeze().load() else: - obs_data[overs] = xr.open_mfdataset( - files[overs], - combine="nested", - ).squeeze() + obs_data[overs] = ( + xr.open_mfdataset( + files[overs], + combine="nested", + ) + .squeeze() + .load() + ) + logger.info(f"DONE - LOAD TIMESERIES DATA FOR {overs}: {_dsname}") return obs_data @@ -468,7 +523,12 @@ def area_avg( Input `data` masked by `isheet_mask` """ - area_data = xr.open_dataset(area_file) + try: + area_data = xr.open_dataset(area_file) + except ValueError: + logger.error(f"INCOMPATABLE FILE {area_file}") + raise + area_data = check_longitude(area_data) if mask_file is None: mask_data = area_data diff --git a/lex/compare_gridded.py b/lex/compare_gridded.py index c8349fa..e0640dc 100644 --- a/lex/compare_gridded.py +++ b/lex/compare_gridded.py @@ -48,6 +48,11 @@ def add_colorbar(color_field, fig, axis, cblabel, cbar_span=False, ndsets=3): _ = fig.colorbar( color_field, cax=cbar_ax, orientation="horizontal", label=cblabel ) + elif ndsets == 1: + _ = fig.colorbar( + color_field, location="right", shrink=0.85, pad=0.01, label=cblabel + ) + else: raise NotImplementedError( f"THIS NUMBER OF PLOTS ({ndsets}) NOT IMPLEMENTED (yet)" @@ -74,11 +79,22 @@ def annotate_plot( xlabel=None, cblabel=None, cnrtxt=None, + gridline_args=None, icesheet="gis", ): """Add land / ocean, gridlines, colourbar.""" axis.coastlines(linewidth=0.1) - axis.gridlines(linestyle="--", linewidth=0.1) + gridline_args_def = {"linestyle": "--", "linewidth": 0.1} + + if gridline_args is None: + gridline_args = gridline_args_def + else: + for _arg, _value in gridline_args_def.items(): + if _arg not in gridline_args: + gridline_args[_arg] = _value + + axis.gridlines(**gridline_args) + if icesheet.lower() == "gis": axis.set_extent([-60, -27, 59, 84], ccrs.PlateCarree()) elif icesheet.lower() == "ais": @@ -116,12 +132,22 @@ def annotate_plot( ) -def get_figure(n_dsets, proj, icesheet="gis"): +def get_figure(n_dsets, proj=None, icesheet="gis"): """Set up figure based on number of datasets to be plotted.""" fig_size = { - "gis": {3: (10, 10), 2: (10, 8)}, - "ais": {3: (10, 10), 2: (16, 7)}, + "gis": {3: (10, 10), 2: (10, 8), 1: (7, 10)}, + "ais": {3: (10, 10), 2: (16, 7), 1: (8, 8)}, } + if proj is None: + if icesheet == "gis": + lon_0 = -45 + lat_0 = 70 + proj = ccrs.LambertConformal( + central_longitude=lon_0, central_latitude=lat_0 + ) + elif icesheet == "ais": + proj = ccrs.SouthPolarStereo(central_longitude=0) + fig = plt.figure(figsize=fig_size[icesheet][n_dsets], dpi=90) if n_dsets == 3: @@ -131,8 +157,10 @@ def get_figure(n_dsets, proj, icesheet="gis"): ) elif n_dsets == 2: axes = [fig.add_subplot(1, 3, i + 1, projection=proj) for i in range(3)] + elif n_dsets == 1: + axes = [fig.add_subplot(1, 1, 1, projection=proj)] - return fig, axes + return fig, axes, proj def main(args, config, sea="ANN"): @@ -159,13 +187,15 @@ def main(args, config, sea="ANN"): if "model_native" in config["dataset_names"]: all_data = { - "model_remap": xr.open_dataset(config["climo_remap"].format(clim=sea)), - "model_native": xr.open_dataset(config["climo"].format(clim=sea)), + "model_remap": xr.open_dataset( + lxc.proc_climo_file(config, "climo_remap", sea) + ), + "model_native": xr.open_dataset(lxc.proc_climo_file(config, "climo", sea)), **lxc.load_obs(config, sea, mode=mode), } else: all_data = { - "model": xr.open_dataset(config["climo_remap"].format(clim=sea)), + "model": xr.open_dataset(lxc.proc_climo_file(config, "climo_remap", sea)), **lxc.load_obs(config, sea, mode=mode), } for _vers in all_data: @@ -258,7 +288,7 @@ def main(args, config, sea="ANN"): all_aavg[_vers], mask_r[_vers], area_r[_vers], _ = lxc.area_avg( _plt_data[_vers], config, - area_file=config["masks"][_vers], + area_file=config["masks"][_vers].format(icesheet=config["icesheet"]), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -274,7 +304,7 @@ def main(args, config, sea="ANN"): diffs_aavg[_diffname], _, _, _ = lxc.area_avg( diffs[_diffname], config, - area_file=config["masks"][_ds2], + area_file=config["masks"][_ds2].format(icesheet=config["icesheet"]), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -303,7 +333,7 @@ def main(args, config, sea="ANN"): else: raise NotImplementedError(f"ICESHEET {icesheet} NOT FOUND USE ais / gis") - fig, axes = get_figure(n_dsets_to_plot, proj, icesheet=icesheet) + fig, axes, _ = get_figure(n_dsets_to_plot, proj, icesheet=icesheet) for _vers in _plt_data: try: diff --git a/lex/energy/energy.py b/lex/energy/energy.py index 8dc5cd9..11ef8e2 100644 --- a/lex/energy/energy.py +++ b/lex/energy/energy.py @@ -38,6 +38,7 @@ import livvkit import pandas as pd from livvkit import elements as el +from loguru import logger from lex import compare_gridded, time_series_plot, utils from lex.common import SEASON_NAME @@ -64,6 +65,7 @@ def run(name, config): """ img_dir = Path(livvkit.output_dir, "validation", "imgs", name) + logger.info(f"Starting ENERGY BALANCE WITH OUTPUT TO {img_dir}") if not img_dir.exists(): img_dir.mkdir(parents=True) @@ -86,6 +88,9 @@ def run(name, config): tables = {} for season in ["ANN", "DJF", "MAM", "JJA", "SON"]: + logger.info( + f"PLOTTING COMPARE GRIDDED FOR {config.get('icesheet', '')} {season}" + ) _plots, aavgs = compare_gridded.main(args, config, sea=season) images[season] = _plots @@ -99,10 +104,16 @@ def run(name, config): transpose=True, ) tables[season] = table_el + logger.info( + "FINISHED PLOTTING COMPARE GRIDDED FOR " + f"{config.get('icesheet', '')} {season}" + ) timeseries_img = [] if "timeseries_dirs" in config: + logger.info(f"PLOTTING TIMESERIES FOR {config.get('icesheet', '')}") timeseries_img.extend(time_series_plot.main(args, config)) + logger.info(f"FINISHED PLOTTING TIMESERIES FOR {config.get('icesheet', '')}") tabs = {} @@ -132,6 +143,7 @@ def run(name, config): tabs["References"] = [ref_ele] elements = [run_summary, el.Tabs(tabs)] + logger.info(f"FINISHED ENERGY BALANCE WITH OUTPUT TO {img_dir}") return el.Page(name, PAGE_DOCS[config.get("icesheet", "gis")], elements) diff --git a/lex/generate_cfg.py b/lex/generate_cfg.py new file mode 100644 index 0000000..9e1fa7d --- /dev/null +++ b/lex/generate_cfg.py @@ -0,0 +1,175 @@ +"""Generate a LIVVkit Extensions (LEX) config based on template information and defaults.""" + +import argparse +from pathlib import Path + +import jinja2 +import mache + +import lex + +ALL_SHEETS = "gis,ais" +ALL_SETS = "cmb,smb,energy_racmo,energy_era5,energy_merra2,energy_ceres" + + +def args(): + parser = argparse.ArgumentParser( + description=__doc__, formatter_class=argparse.ArgumentDefaultsHelpFormatter + ) + + parser.add_argument( + "--template", + "-t", + type=Path, + help="Path to configuration template file.", + ) + + parser.add_argument( + "--case", + "-c", + type=str, + help="Case ID", + required=True, + ) + + parser.add_argument( + "--casedir", + "-d", + type=str, + help="Output directory where climatology files for `case` are stored", + required=True, + ) + + parser.add_argument( + "--cfg_out", + "-o", + type=Path, + help="Output directory for config file.", + default=Path("./").resolve(), + ) + + parser.add_argument( + "--sets", + "-s", + type=str, + help=( + "Analysis sets to run: cmb, smb, energy_racmo, energy_era5, " + "energy_merra2, energy_ceres, or all to run all available" + ), + ) + + parser.add_argument( + "--icesheets", + "-i", + type=str, + help=( + "Comma separated icesheets to analyse (ais for Antarctica," + " gis for Greenland), or all to run both" + ), + ) + + parser.add_argument( + "-z", + "--from-zppy", + action="store_true", + default=False, + help="Flag to be called from zppy, makes grid parsed", + ) + + return parser.parse_args() + + +def gen_cfg(cfg_template, params, cfg_out): + jenv = jinja2.Environment( + loader=jinja2.FileSystemLoader(cfg_template.resolve().parent) + ) + template = jenv.get_template(cfg_template.name) + cfg = template.render(**params) + + if not Path(cfg_out.parent).exists(): + print(f"CREATE {cfg_out.parent}") + Path(cfg_out.parent).mkdir(parents=True) + print(f"WRITE: {cfg_out}") + with open(cfg_out, "w", encoding="utf-8") as _cfgout: + _cfgout.write(cfg) + return cfg_out + + +def parse_sets(sheets, sets): + """Parse comma separated strings of sets / icesheets to analyse.""" + + params = {} + if sheets.lower() == "all": + sheets = ALL_SHEETS + if sets.lower() == "all": + sets = ALL_SETS + + _sheets = sheets.lower().split(",") + _sets = sets.lower().split(",") + + sheets = [_sheet.strip() for _sheet in _sheets] + sets = [_set.strip() for _set in _sets] + + for _sheet in sheets: + params[_sheet] = True + + for _set in sets: + params[_set] = True + return params + + +def main(): + cl_args = args() + mach = mache.discover_machine() + mach_info = mache.MachineInfo() + + defaults = { + "chrysalis": { + "livvproj_dir": Path("/lcrc/group/e3sm/livvkit"), + "model_ts_dir": Path("/lcrc/group/e3sm/ac.zender/scratch/livvkit"), + "grid_dir": Path("/lcrc/group/e3sm/zender/grids"), + }, + "pm-cpu": { + "livvproj_dir": Path("/global/cfs/cdirs/e3sm/livvkit"), + "model_ts_dir": Path("/global/cfs/projectdirs/e3sm/zender/livvkit"), + "grid_dir": Path("/global/cfs/cdirs/e3sm/zender/grids"), + "racmo_root_dir": Path("/global/cfs/cdirs/fanssie/racmo/2.4.1"), + }, + } + climo_dirs = {} + ts_dirs = {} + + if cl_args.from_zppy: + for _grid in ["native", "racmo_ais", "racmo_gis", "era5", "ceres", "merra2"]: + if _grid == "ceres": + _ceres_grid = "180x360_traave" + _grid_dir = cl_args.casedir.replace("DATA_GRID", _ceres_grid) + else: + _grid_dir = cl_args.casedir.replace("DATA_GRID", _grid) + + _ts_dir = _grid_dir.replace("/clim/", "/ts/monthly/") + climo_dirs[f"dir_{_grid}"] = Path(_grid_dir) + ts_dirs[f"dir_ts_{_grid}"] = Path(_ts_dir) + + case_dir = Path(cl_args.casedir) + + _mach_defaults = defaults[mach] + _mach_defaults["e3sm_diags_data_dir"] = Path( + mach_info.config.get("diagnostics", "base_path") + ) + params = { + **_mach_defaults, + "case_id": cl_args.case, + "case_out_dir": case_dir, + "lex_root": Path(lex.__file__).parent, + **parse_sets(cl_args.icesheets, cl_args.sets), + **climo_dirs, + **ts_dirs, + } + out_cfg = Path(cl_args.cfg_out, "livvkit.yml") + gen_cfg(cl_args.template, params, out_cfg) + print(f"CONFIGURATION WRITTEN TO: {out_cfg}") + + +if __name__ == "__main__": + main() diff --git a/lex/postproc/e3sm/postproc.sbatch b/lex/postproc/e3sm/postproc.sbatch index b4f8fbb..dc6553d 100644 --- a/lex/postproc/e3sm/postproc.sbatch +++ b/lex/postproc/e3sm/postproc.sbatch @@ -8,8 +8,22 @@ conda_activate_script=/global/common/software/e3sm/anaconda_envs/base/etc/profile.d/conda.sh echo $conda_activate_script source $conda_activate_script -conda activate lex +if [ -n ${LEX_ENV} ]; then + conda activate ${LEX_ENV} +else if [ -d ${HOME}/.conda/envs/lex_env ]; then + conda activate ${HOME}/.conda/envs/lex_env +else if [ -d ${HOME}/anaconda/envs/lex_env ]; then + conda activate ${HOME}/anaconda/envs/lex_env +else + echo "LEX ENV NOT FOUND AT EITHER $HOME/.conda or $HOME/anaconda " + echo "SET LEX_ENV variable to point to \$CONDA_PREFIX for lex_env" + exit 1 +fi +fi +fi + +LEXDIR=$HOME/lex INDIR=/global/cfs/projectdirs/e3sm/zender/livvkit LIVVPROJ=/global/cfs/cdirs/e3sm/livvkit @@ -86,30 +100,34 @@ ANN_MEAN_FILE=${OUTDIR}/${OUTCASE}.ANN_mean.nc ncra $INFILE ${ANN_MEAN_FILE} for var in topo landfrac do - ncks -m -v ${var} ${INFILE} > /dev/null || ncks -A -C -v ${var} ${INFILE_REF} ${ANN_MEAN_FILE} + ncks -m -v ${var} ${INFILE} > /dev/null 2>&1 || ncks -A -C -v ${var} ${INFILE_REF} ${ANN_MEAN_FILE} done -srun parallel --jobs 12 ./split_files.sh ::: {0..11} +srun parallel --jobs 12 ${LEXDIR}/lex/postproc/e3sm/split_files.sh ::: {0..11} ncra \ + -O \ ${OUTDIR}/${OUTCASE}.01.nc \ ${OUTDIR}/${OUTCASE}.02.nc \ ${OUTDIR}/${OUTCASE}.12.nc \ ${OUTDIR}/${OUTCASE}.DJF_mean.nc & ncra \ + -O \ ${OUTDIR}/${OUTCASE}.03.nc \ ${OUTDIR}/${OUTCASE}.04.nc \ ${OUTDIR}/${OUTCASE}.05.nc \ ${OUTDIR}/${OUTCASE}.MAM_mean.nc & ncra \ + -O \ ${OUTDIR}/${OUTCASE}.06.nc \ ${OUTDIR}/${OUTCASE}.07.nc \ ${OUTDIR}/${OUTCASE}.08.nc \ ${OUTDIR}/${OUTCASE}.JJA_mean.nc & ncra \ + -O \ ${OUTDIR}/${OUTCASE}.09.nc \ ${OUTDIR}/${OUTCASE}.10.nc \ ${OUTDIR}/${OUTCASE}.11.nc \ @@ -122,7 +140,7 @@ do do # Check if the mean output file has ${var} (topography or landfrac) if it # doesn't, then get that field from a reference file (set by the resolution) - ncks -m -v ${var} ${MEAN_FILE} > /dev/null || ncks -A -C -v ${var} ${INFILE_REF} ${MEAN_FILE} + ncks -m -v ${var} ${MEAN_FILE} > /dev/null 2>&1 || ncks -A -C -v ${var} ${INFILE_REF} ${MEAN_FILE} done ncremap \ @@ -163,7 +181,7 @@ echo "#########################" echo "POST PROCESSING COMPLETE" echo "#########################" -pushd $HOME/lex +pushd ${LEXDIR} if [ ! -d config/${OUTCASE} ]; then echo "CREATING CONFIG DIRECTORY FOR ${OUTCASE}" cp -R config/template_${RES} config/${OUTCASE} diff --git a/lex/postproc/e3sm/split_files.sh b/lex/postproc/e3sm/split_files.sh index 6d70e57..954eb4d 100755 --- a/lex/postproc/e3sm/split_files.sh +++ b/lex/postproc/e3sm/split_files.sh @@ -7,14 +7,14 @@ MON=$(printf "%02d" $(($MON_IDX + 1))) # Monthly split MONFILE=${OUTCASE}.${YEAR_STR}${MON} -ncks -t 16 -d time,$MON_IDX,$MAXTIME,12 $INFILE ${OUTDIR}/$MONFILE.nc +ncks -O -t 16 -d time,$MON_IDX,$MAXTIME,12 $INFILE ${OUTDIR}/$MONFILE.nc # echo "MEAN FOR ${MON}: ${MON_IDX}-${MAXTIME}" ncra ${OUTDIR}/${MONFILE}.nc ${OUTDIR}/${MONFILE}_mean.nc for var in topo landfrac do - ncks -m -v ${var} ${OUTDIR}/${MONFILE}_mean.nc > /dev/null || ncks -A -C -v ${var} ${INFILE_REF} ${OUTDIR}/${MONFILE}_mean.nc + ncks -m -v ${var} ${OUTDIR}/${MONFILE}_mean.nc > /dev/null 2>&1 || ncks -A -C -v ${var} ${INFILE_REF} ${OUTDIR}/${MONFILE}_mean.nc done ncremap \ diff --git a/lex/smb/smb/plot_spatial.py b/lex/smb/smb/plot_spatial.py index 6dca4df..d31a808 100644 --- a/lex/smb/smb/plot_spatial.py +++ b/lex/smb/smb/plot_spatial.py @@ -7,17 +7,21 @@ import os from pathlib import Path +import cartopy.crs as ccrs import matplotlib.pyplot as plt import numpy as np import numpy.ma as ma import pandas as pd import smb.preproc as preproc from livvkit import elements as el +from loguru import logger from matplotlib import colors as c -from mpl_toolkits.basemap import Basemap from netCDF4 import Dataset from scipy.interpolate import griddata +import lex.common as lxc +import lex.compare_gridded as lxcg + DESCRIBE_CORE = """ Filled contour of modeled annual surface mass balance of the Greenland ice sheet, with firn and core field estimates overlaid as filled circles. Data were @@ -57,6 +61,13 @@ IMG_GROUP = "Spatial" +GRIDLINE_ARGS = { + "draw_labels": ["top", "bottom", "left"], + "x_inline": False, + "y_inline": False, + "linewidth": 1.0, +} + def mali_to_contour(lon_cell, lat_cell, data_cell): """Convert MALI unstructured to gridded data.""" @@ -75,27 +86,35 @@ def mali_to_contour(lon_cell, lat_cell, data_cell): return x_grid, y_grid, z_grid -def load_model_data(args, config, regrid=True): +def load_model_data(config, regrid=True): """Load Model data.""" - nc1 = Dataset(config["climo"].format(m_s=1, m_e=12, clim="ANN"), mode="r") - nc2 = Dataset(args.latlon, mode="r") - nc3 = Dataset(args.elevation, mode="r") - - lats_model = nc1.variables[config["latv"]][:] - lons_model = nc1.variables[config["lonv"]][:] - if nc2.variables[args.lonv].units == "radians": + sea_s, sea_e = lxc.get_season_bounds( + "ANN", config.get("year_s", 0), config.get("year_e", 1) + ) + _climfile = config["climo"].format(sea_s=sea_s, sea_e=sea_e, clim="ANN") + _elevfile = config["elevation"].format(sea_s=sea_s, sea_e=sea_e, clim="ANN") + _gridfile = config["latlon"].format(sea_s=sea_s, sea_e=sea_e, clim="ANN") + logger.info(f"LOADING CLIMO FILE: {_climfile}") + clim_nc = Dataset(_climfile, mode="r") + grid_nc = Dataset(_gridfile, mode="r") + elev_nc = Dataset(_elevfile, mode="r") + + lats_model = clim_nc.variables[config["latv"]][:] + lons_model = clim_nc.variables[config["lonv"]][:] + if grid_nc.variables[config["lonv"]].units == "radians": lats_model = np.degrees(lats_model) lons_model = np.degrees(lons_model) - smb_model = nc1.variables[args.smbv][:] + smb_model = clim_nc.variables[config["smbv"]][:] + + thk_model = clim_nc.variables[config["topov"]][:] + if config["landfracv"] in clim_nc.variables: + smb_model *= clim_nc.variables[config["landfracv"]][:] + smb_model *= config["smbscale"] - thk_model = nc1.variables[args.topov][:] - if args.landfracv in nc1.variables: - smb_model *= nc1.variables[args.landfracv][:] - smb_model *= args.smbscale - nc1.close() - nc2.close() - nc3.close() + clim_nc.close() + grid_nc.close() + elev_nc.close() mask = thk_model.flatten() < 0.0001 smb_flat = smb_model.flatten() @@ -112,110 +131,58 @@ def load_model_data(args, config, regrid=True): return lons_model, lats_model, msmb -def gen_map(lon_0, lat_0): - """ - Create Basemap on which to plot. - - Parameters - ---------- - lon_0, lat_0 : float - Center location of map to be drawn - - Returns - ------- - pmap : basemap.Basemap - Lambert conformal conic `Basemap` centered on Greenland - (future: do a stereographic for Antarctica) - - """ - pmap = Basemap( - width=2000000, - height=3000000, - rsphere=(6378137.00, 6356752.3142), - resolution="l", - projection="lcc", - lat_0=lat_0, - lon_0=lon_0, - ) - return pmap - - -def annotate_map(pmap, axis=None): - """ - Add common details to basemap. - - Draw grid and coastlines, and fill the continents with colour. - - Parameters - ---------- - pmap : basemap.Basemap - Basemap on which to draw annotations - - axis : matplotlib.pyplot.Axis, optional - Axis where Basemap is drawn, default is None - - """ - pmap.drawcoastlines(ax=axis) - pmap.drawcountries(ax=axis) - pmap.fillcontinents(color="gainsboro", lake_color="aqua", zorder=1) - pmap.drawparallels( - np.arange(-80.0, 81.0, 5.0), - labels=[1, 0, 0, 0], - fontsize=10, - color="lightgrey", - zorder=1, - ax=axis, - ) - pmap.drawmeridians( - np.arange(-180.0, 181.0, 10.0), - labels=[0, 0, 0, 1], - fontsize=10, - color="lightgrey", - zorder=1, - ax=axis, - ) - - def plot_core(args, config): """Plot Ice Core data on map with model data.""" img_list = [] _, _, _, core_file = preproc.core_file(config) smb_avg = pd.read_csv(core_file) - plt.figure(figsize=(12, 12)) - - # lon_0 = lons_model.mean() - # lat_0 = lats_model.mean() - lon_0 = -40.591 - lat_0 = 71.308 + tform = ccrs.PlateCarree() + lons_model, lats_model, msmb = load_model_data(config) + fig, axes, proj = lxcg.get_figure(1, icesheet="gis") - lons_model, lats_model, msmb = load_model_data(args, config) - pmap = gen_map(lon_0, lat_0) - annotate_map(pmap) - xi, yi = pmap(lons_model.squeeze(), lats_model.squeeze()) vabsmax = 2000 - # smb = pmap.scatter(xi, yi, c=msmb.squeeze(), vmin=-1500, vmax=1500, - # cmap='Spectral', zorder=2) - smb = pmap.pcolormesh( - xi, yi, msmb.squeeze(), vmin=-vabsmax, vmax=vabsmax, cmap="Spectral", zorder=2 + cf_smb_model = axes[0].pcolormesh( + lons_model, + lats_model, + msmb.squeeze(), + vmin=-vabsmax, + vmax=vabsmax, + cmap="Spectral", + zorder=2, + transform=tform, + ) + + lxcg.annotate_plot( + axes[0], + icesheet="gis", + gridline_args=GRIDLINE_ARGS, + ) + + lxcg.add_colorbar( + cf_smb_model, + fig, + axes[0], + "Surface mass balance (kg m$^{-2}$ a$^{-1}$)", + cbar_span=False, + ndsets=1, ) - cbar = pmap.colorbar(smb, location="bottom", pad="5%") - cbar.set_label("Surface mass balance (kg m$^{-2}$ a$^{-1}$)") lat_obs = smb_avg["Y"].values lon_obs = smb_avg["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) smbobs = smb_avg.b - _ = pmap.scatter( - xobs, - yobs, + + _ = axes[0].scatter( + lon_obs, + lat_obs, c=smbobs, - vmin=-1500, - vmax=1500, + vmin=-vabsmax, + vmax=vabsmax, marker="o", edgecolor="black", cmap="Spectral", zorder=3, + transform=tform, ) plt.tight_layout() @@ -244,34 +211,36 @@ def plot_ib_spatial(args, config): img_list = [] _, _, ib_file = preproc.ib_outfile(config) ice_bridge = pd.read_csv(ib_file) - # lons_model, lats_model, msmb = load_model_data(args, config) - fig = plt.figure(figsize=(12, 12)) - axis = fig.add_subplot(1, 1, 1) - # lon_0 = lons_model.mean() - # lat_0 = lats_model.mean() - lon_0 = -40.591 - lat_0 = 71.30 - - # print(f"LON0: {lon_0} LAT0: {lat_0}") - pmap = gen_map(lon_0, lat_0) - annotate_map(pmap, axis) + # Plot the IceBridge data + tform = ccrs.PlateCarree() + _, _, ib_file = preproc.ib_outfile(config) + ice_bridge = pd.read_csv(ib_file) + fig, axes, proj = lxcg.get_figure(1, icesheet="gis") lat_obs = ice_bridge["Y"].values lon_obs = ice_bridge["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) - # xmod, ymod = m(lons_model, lats_model) smbobs = ice_bridge["b"].values - obs = pmap.scatter( - xobs, yobs, c=smbobs, marker="o", cmap="viridis", zorder=3, lw=0, ax=axis + cf_obs = axes[0].scatter( + lon_obs, + lat_obs, + c=smbobs, + marker="o", + cmap="Spectral", + zorder=3, + lw=0, + transform=tform, + ) + lxcg.annotate_plot(axes[0], icesheet="gis", gridline_args=GRIDLINE_ARGS) + lxcg.add_colorbar( + cf_obs, + fig, + axes[0], + "Surface mass balance (kg m$^{-2}$ a$^{-1}$)", + cbar_span=False, + ndsets=1, ) - # _map = pmap.contour(xmod, ymod, climo_data[0], colors="k") - # _map = pmap.contour(xmod, ymod, climo_data, color="k") - - cbar = pmap.colorbar(obs, location="bottom", pad="5%", ax=axis) - cbar.set_label("Surface mass balance (kg m$^-2$ a$^{-1}$)") - plt.tight_layout() img_file = os.path.join(args.out, "IB_spatial.png") plt.savefig(img_file) @@ -290,30 +259,33 @@ def plot_ib_spatial(args, config): ) img_list.append(img_elem) - fig = plt.figure(figsize=(12, 12)) - axis = fig.add_subplot(1, 1, 1) - annotate_map(pmap, axis) - - lat_obs = ice_bridge["Y"].values - lon_obs = ice_bridge["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) - smbobs = ice_bridge["mod_b"].values - ice_bridge["b"].values - obs = pmap.scatter( - xobs, - yobs, - c=smbobs, + # Plot IB - Model difference + fig, axes, proj = lxcg.get_figure(1, icesheet="gis") + smbobs_diff = ice_bridge["mod_b"].values - ice_bridge["b"].values + cf_diff = axes[0].scatter( + lon_obs, + lat_obs, + c=smbobs_diff, marker="o", cmap="RdBu", zorder=3, vmin=-200, vmax=200, lw=0, - ax=axis, + transform=tform, ) - - cbar = pmap.colorbar(obs, location="bottom", pad="5%") - cbar.set_label( - "Difference in surface mass balance \n (Model - IceBridge) (kg m$^-2$ a$^{-1}$)" + lxcg.annotate_plot( + axes[0], + icesheet="gis", + gridline_args=GRIDLINE_ARGS, + ) + lxcg.add_colorbar( + cf_diff, + fig, + axes[0], + "Surface mass balance difference (kg m$^{-2}$ a$^{-1}$)", + cbar_span=False, + ndsets=1, ) plt.tight_layout() @@ -340,6 +312,8 @@ def plot_metadata(args, config): """ Create plot of basin locations, ice bridge transects, and core locations. """ + tform = ccrs.PlateCarree() + img_list = [] _, _, ib_file = preproc.ib_outfile(config) _, _, _, core_file = preproc.core_file(config) @@ -348,11 +322,11 @@ def plot_metadata(args, config): except pd.errors.EmptyDataError: print(ib_file) raise + smb_avg = pd.read_csv(core_file) zwally_data = pd.read_csv(Path(config["preproc_dir"], config["zwally_file"])) - lons_model, lats_model, smb_model = load_model_data(args, config, regrid=False) - plt.figure(figsize=(12, 12)) + lons_model, lats_model, smb_model = load_model_data(config, regrid=False) forcolors = c.ListedColormap( [ @@ -368,13 +342,7 @@ def plot_metadata(args, config): ] ) - # This will be the center of our map - # lon_0 = lons_model.mean() - # lat_0 = lats_model.mean() - lon_0 = -40.591 - lat_0 = 71.308 - pmap = gen_map(lon_0, lat_0) - annotate_map(pmap) + fig, axes, proj = lxcg.get_figure(1, icesheet="gis") # Read in the zwally basins and mask out model cells that are missing a basin designation basins = np.floor(zwally_data.zwally_basin.values) @@ -384,15 +352,15 @@ def plot_metadata(args, config): mbasins = ma.masked_invalid(basins) # Plotting the basins pseudocolor - if lons_model.ndim == 1 and lons_model.shape[0] < 1441: + if lons_model.ndim == 1 and "elm" in config["meta"]["Model"][0].lower(): lons_model, lats_model = np.meshgrid(lons_model, lats_model) - elif lons_model.shape[0] >= 1441: + elif "mali" in config["meta"]["Model"][0].lower(): lons_model, lats_model, mbasins = mali_to_contour( lons_model, lats_model, mbasins ) - - xi, yi = pmap(lons_model.squeeze(), lats_model.squeeze()) - _ = pmap.pcolormesh(xi, yi, mbasins, cmap=forcolors, zorder=2) + _ = axes[0].pcolormesh( + lons_model, lats_model, mbasins, cmap=forcolors, zorder=2, transform=tform + ) basin_labels = { "1": (-48, 80.5), @@ -405,43 +373,76 @@ def plot_metadata(args, config): "8": (-55, 75), } for lbl, loc in basin_labels.items(): - xi, yi = pmap(loc[0], loc[1]) - plt.text(xi, yi, lbl, fontsize=16, fontweight="bold", color="#FF7900") + # xi, yi = pmap(loc[0], loc[1]) + plt.text( + loc[0], + loc[1], + lbl, + fontsize=16, + fontweight="bold", + color="#FF7900", + transform=tform, + ) # Adding in firn/core measurement locations. Size by the number of years in the record lat_obs = smb_avg["Y"].values lon_obs = smb_avg["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) - smb_avg.loc[smb_avg["nyears"] > 50] = 50 - smb_avg.loc[smb_avg["nyears"] < 5] = 5 - _ = pmap.scatter( - xobs, yobs, marker="o", edgecolor="black", s=4 * smb_avg["nyears"], zorder=4 + # xobs, yobs = pmap(lon_obs, lat_obs) + smb_avg[smb_avg["nyears"] > 50].loc[:, "nyears"] = 50 + smb_avg[smb_avg["nyears"] < 5].loc[:, "nyears"] = 5 + + _ = axes[0].scatter( + lon_obs, + lat_obs, + marker="o", + edgecolor="black", + s=4 * smb_avg["nyears"], + zorder=4, + transform=tform, ) # Color the ablation zone (PROMICE) measurements yellow smb_promice = smb_avg[smb_avg.source == "promice"].copy() lat_obs = smb_promice["Y"].values lon_obs = smb_promice["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) - smb_promice.loc[smb_promice["nyears"] > 50] = 50 - smb_promice.loc[smb_promice["nyears"] < 5] = 5 - _ = pmap.scatter( - xobs, - yobs, + + smb_promice[smb_promice["nyears"] > 50].loc[:, "nyears"] = 50 + smb_promice[smb_promice["nyears"] < 5].loc[:, "nyears"] = 5 + + _ = axes[0].scatter( + lon_obs, + lat_obs, marker="^", color="yellow", edgecolor="black", s=4 * smb_promice["nyears"], zorder=4, + transform=tform, ) # Add IceBridge transects + _, _, ib_file = preproc.ib_outfile(config) + ice_bridge = pd.read_csv(ib_file) lat_obs = ice_bridge["Y"].values lon_obs = ice_bridge["X"].values - xobs, yobs = pmap(lon_obs, lat_obs) - _ = pmap.scatter(xobs, yobs, marker="o", lw=0, zorder=3, s=3, color="white") + _ = axes[0].scatter( + lon_obs, + lat_obs, + marker="o", + lw=0, + zorder=3, + s=2, + color="white", + transform=tform, + ) + lxcg.annotate_plot( + axes[0], + icesheet="gis", + gridline_args=GRIDLINE_ARGS, + ) plt.tight_layout() + img_file = os.path.join(args.out, "plot_meta_old.png") plt.savefig(img_file) plt.close() diff --git a/lex/smb/smb_icecores.py b/lex/smb/smb_icecores.py index 8cb88ad..95756d6 100644 --- a/lex/smb/smb_icecores.py +++ b/lex/smb/smb_icecores.py @@ -47,6 +47,8 @@ import smb.preproc as preproc import smb.utils as utils +from loguru import logger + PAGE_DOCS = { "gis": """Validation of the Greenland Ice Sheet (GrIS) surface mass balance by comparing modeled surface mass balance to estimates from in situ measurements @@ -84,8 +86,8 @@ def run(name, config): Returns: A LIVVkit page element containing the LIVVkit elements to display on a webpage """ - img_dir = os.path.join(livvkit.output_dir, "validation", "imgs", name) + logger.info(f"Starting SMB_ICECORES OUTPUT TO {img_dir}") fn.mkdir_p(img_dir) config_arg_list = [] for key, val in config.items(): @@ -93,31 +95,41 @@ def run(name, config): config_arg_list.extend(["--out", img_dir]) args = parse_args(config_arg_list) - if config["preprocess"]: + if config.get("preprocess", False): preproc.main(args, config) spatial_img = [] statistic_img = [] timeseries_img = [] if "smb_cf_file" in config and "smb_mo_file" in config and "ib_file" in config: + logger.info("PLOT SPATIAL METADATA") spatial_img.extend(plt_spatial.plot_metadata(args, config)) + logger.info("DONE - PLOT SPATIAL METADATA") if "smb_cf_file" in config and "smb_mo_file" in config: + logger.info("PLOT SPATIAL CORE DATA") spatial_img.extend(plt_spatial.plot_core(args, config)) transects = c_transects.main(args, config) statistic_img.extend(transects[3:]) statistic_img.extend(IB_scatter.main(args, config)) statistic_img.extend(c_hists.main(args, config)) + logger.info("DONE - PLOT SPATIAL CORE DATA") if "ib_file" in config: + logger.info("PLOT SPATIAL IB DATA") spatial_img.extend(plt_spatial.plot_ib_spatial(args, config)) statistic_img.extend(IB_hist.main(args, config)) + logger.info("DONE - PLOT SPATIAL IB DATA") if "smb_cf_file" in config and "smb_mo_file" in config: + logger.info("PLOT STATSTICAL DATA") statistic_img.extend(transects[:3]) + logger.info("DONE - PLOT STATSTICAL DATA") if "timeseries_dirs" in config: + logger.info("PLOT TIMESERIES DATA") timeseries_img.extend(time_series_plot.main(args, config)) + logger.info("DONE - PLOT TIMESERIES DATA") seasons = ["ANN", "DJF", "MAM", "JJA", "SON"] seasonal_components = {} @@ -129,12 +141,16 @@ def _format_table(x): return pd.Series([f"{xi:.2f}" for xi in x], index=x.index) for season in seasons: + logger.info(f"COMPARE GRIDDED {season} DATA") _img, _aavg = compare_gridded.main(args, config, sea=season) + logger.info(f"DONE - COMPARE GRIDDED {season} DATA") seasonal_components[season] = [] if season == "ANN": + logger.info("PLOT ANNUAL CYCLE DATA") seasonal_components[season].extend(annual_cycle.main(args, config)) + logger.info("DONE - PLOT ANNUAL CYCLE DATA") seasonal_components[season].extend(_img) seasonal_tables[season] = el.Table( @@ -178,6 +194,7 @@ def _format_table(x): tabs["References"] = [refs] + logger.info(f"FINISHED SMB_ICECORES WITH OUTPUT TO {img_dir}") return el.Page( name, PAGE_DOCS[config.get("icesheet", "gis")], diff --git a/lex/time_series_plot.py b/lex/time_series_plot.py index b21e559..2f85c41 100644 --- a/lex/time_series_plot.py +++ b/lex/time_series_plot.py @@ -12,6 +12,7 @@ import numpy as np import xarray as xr from livvkit import elements as el +from loguru import logger import lex.common as lxc import lex.utils as lxu @@ -39,7 +40,7 @@ def assemble_outdata(args, config, dataset, aavg_data, ts_data, aavg_units): else: _aavg = aavg_data.get(data_var["title"]) else: - print(f"DATA NOT FOUND FOR {data_var['title']} in {dataset}") + logger.error(f"DATA NOT FOUND FOR {data_var['title']} in {dataset}") continue aavg_out[data_var["title"].replace(" ", "_")] = xr.DataArray( @@ -79,6 +80,7 @@ def main(args, config): config_names[_dset] = _dset img_elem = [] for idx, data_var in enumerate(config["data_vars"]): + logger.info(f" PLOTTING {config.get('icesheet', '')} TS: {data_var['title']}") _obs_in = {} aavg_config = data_var.get("aavg", None) @@ -114,7 +116,9 @@ def main(args, config): lxc.area_avg( _obs_in[_vers], {}, - area_file=config["masks"][_vers], + area_file=config["masks"][_vers].format( + icesheet=config["icesheet"] + ), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -124,7 +128,9 @@ def main(args, config): model_aavg[data_var["title"]], _, _, _ = lxc.area_avg( _model_plt, {}, - area_file=config["masks"]["model_native"], + area_file=config["masks"]["model_native"].format( + icesheet=config["icesheet"] + ), area_var="area", mask_var="Icemask", sum_out=_do_sum, @@ -259,7 +265,9 @@ def main(args, config): relative_to="", ) img_elem.append(_img_elem) - + logger.info( + f" DONE - PLOTTING {config.get('icesheet', '')} TS: {data_var['title']}" + ) # assemble_outdata(args, config, "model", model_aavg, ts_data, _aavg_units) # assemble_outdata(args, config, "dset_a", obs_aavg, ts_data, _aavg_units) diff --git a/pyproject.toml b/pyproject.toml index 790e972..c12a67b 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -37,11 +37,32 @@ authors = [ { name = "Joseph H. Kennedy", email = "kennedyjh@ornl.gov" }, ] name = "lex" -requires-python = ">=3.9" +requires-python = ">=3.10" description = "Extensions for the land ice verification and validation toolkit (Livvkit EXtensions)" readme = { file = "README.rst", content-type = "text/x-rst" } license = { text = "BSD-3-Clause" } +dependencies = [ + "numpy>=2.2.6", + "scipy>=1.15.3", + "xarray>=2025.6.1", + "dask>=2026.1.2", + "seaborn>=0.13.2", + "cftime>=1.6.5", + "nc-time-axis>=1.4.1", + "scikit-learn>=1.7.2", + "ruamel.yaml>=0.19.1", + "jinja2>=3.1.6", + "matplotlib>=3.10.8", + "netcdf4>=1.7.4", + "pandas>=2.2.3", + "cartopy>=0.25.0", + "pybtex>=0.25.1", + "f90nml>=1.5", + "loguru>=0.7.3", + "livvkit>=3.3.0", +] + classifiers = [ "Development Status :: 5 - Production/Stable", "Intended Audience :: Science/Research", @@ -51,25 +72,67 @@ classifiers = [ "Natural Language :: English", "Operating System :: OS Independent", "Programming Language :: Python :: 3", - "Programming Language :: Python :: 3.9", "Programming Language :: Python :: 3.10", "Programming Language :: Python :: 3.11", "Programming Language :: Python :: 3.12", "Programming Language :: Python :: 3.13", + "Programming Language :: Python :: 3.14", ] - -dynamic = ["version", "dependencies", "optional-dependencies"] +dynamic = ["version"] [tool.setuptools] packages = ["lex", "lex.smb", "lex.energy"] [tool.setuptools.dynamic] -dependencies = { file = ["requirements.txt"] } version = { attr = "lex.__version__" } -[tool.setuptools.dynamic.optional-dependencies] -dev = { file = ["requirements-dev.txt"] } +[project.urls] +Homepage = "https://github.com/LIVVkit/lex" + +[project.scripts] +lex-cfg = "lex.generate_cfg:main" +[tool.pixi.workspace] +channels = ["conda-forge"] +platforms = ["linux-64"] -[project.urls] -Homepage = "https://code.ornl.gov/LIVVkit/lex" +[tool.pixi.pypi-dependencies] +lex = { path = ".", editable = true } + +[tool.pixi.tasks] + +[tool.pixi.dependencies] +mache = ">=2.1.0" +nco = ">=5.3.6" + +# Development tools +[tool.pixi.feature.dev.dependencies] +isort = "*" +pytest = "*" +ruff = "*" +pre-commit = "*" + +# Different Python versions for CI testing +[tool.pixi.feature.py310dev.dependencies] +python = "3.10.*" + +[tool.pixi.feature.py311dev.dependencies] +python = "3.11.*" + +[tool.pixi.feature.py312dev.dependencies] +python = "3.12.*" + +[tool.pixi.feature.py313dev.dependencies] +python = "3.13.*" + +[tool.pixi.feature.py314dev.dependencies] +python = "3.14.*" + +[tool.pixi.environments] +default = { solve-group = "default" } +dev = ["dev"] +py310dev = ["py310dev", "dev"] +py311dev = ["py311dev", "dev"] +py312dev = ["py312dev", "dev"] +py313dev = ["py313dev", "dev"] +py314dev = ["py314dev", "dev"] diff --git a/requirements-dev.txt b/requirements-dev.txt deleted file mode 100644 index db616d0..0000000 --- a/requirements-dev.txt +++ /dev/null @@ -1,3 +0,0 @@ -isort -pytest -ruff diff --git a/requirements.txt b/requirements.txt deleted file mode 100644 index 7acf04b..0000000 --- a/requirements.txt +++ /dev/null @@ -1,19 +0,0 @@ - numpy - scipy - xarray - dask - seaborn - cftime - nc-time-axis - scikit-learn - ruamel.yaml - pylint - matplotlib - basemap - netcdf4 - pandas - nco - cartopy - pybtex - f90nml - livvkit>=3.2.0 diff --git a/run_lex_pm-cpu.sbatch b/run_lex_pm-cpu.sbatch index 0b1054a..18dce3d 100644 --- a/run_lex_pm-cpu.sbatch +++ b/run_lex_pm-cpu.sbatch @@ -4,14 +4,9 @@ #SBATCH --nodes=1 #SBATCH --constraint=cpu -conda_activate_script=/global/common/software/e3sm/anaconda_envs/base/etc/profile.d/conda.sh -echo $conda_activate_script -source $conda_activate_script - # NERSC Portal project location export WEBDIR=/global/cfs/projectdirs/e3sm/www/${USER} -conda activate lex CASES=( "v2.1.r025.IGERA5ELM_MLI-deep_firn_1980_2020" "v3.LR.piControl-deepfirn_0001_0100" diff --git a/run_livv.sh b/run_livv.sh index 32367a6..5c6f9d9 100755 --- a/run_livv.sh +++ b/run_livv.sh @@ -1,28 +1,63 @@ #!/bin/bash CASE=$1 +if [ $(command -v pixi) ];then + echo "ACTIVATING PIXI SHELL" + eval $(pixi shell-hook -e default) +else if [ $(command -v conda) ]; then + conda_activate_script=/global/common/software/e3sm/anaconda_envs/base/etc/profile.d/conda.sh + echo "ACTIVATE: ${conda_activate_script}" + source $conda_activate_script + if [ ! -z ${LEX_ENV+x} ]; then + conda activate ${LEX_ENV} + else if [ -d ${HOME}/.conda/envs/lex_env ]; then + conda activate ${HOME}/.conda/envs/lex_env + else if [ -d ${HOME}/anaconda/envs/lex_env ]; then + conda activate ${HOME}/anaconda/envs/lex_env + else + echo "LEX ENV NOT FOUND AT EITHER $HOME/.conda or $HOME/anaconda " + echo "SET LEX_ENV variable to point to \$CONDA_PREFIX for lex_env" + exit 1 + fi + fi + fi +fi +fi + +livv_cmd='livv' +livv_exe=`which ${livv_cmd}` +if [ -z "${livv_exe}" ]; then + echo "ERROR: Unable to find LIVV binary executable for command \"${livv_cmd}\"" + exit 1 +else + echo "Running LIVV from ${livv_exe}" +fi + echo "LEX ON ${CASE}" +if [[ ${CASE} == *"r05"* || ${CASE} == *".LR."* ]]; then + template=config/template_r05/all.jinja +else if [[ ${CASE} == *"r025"* || ${CASE} == *".HR."* ]]; then + template=config/template_r025/all.jinja +fi +fi # Allow for standalone (outside of batch script) by setting WEBDIR if it's not already set WEBDIR="${WEBDIR:-/global/cfs/projectdirs/e3sm/www/${USER}}" +# (Re-)generate the config file for this run +lex-cfg \ + --template ${template} \ + --casedir $PSCRATCH/lex/data/e3sm/${CASE} \ + --case ${CASE} \ + --cfg_out ./config \ + --icesheets run_gis,run_ais \ + --sets set_all + # Run LIVVkit with all the configs we want for this case # writes the output website to the user's scratch directory mkdir -p ${SCRATCH}/lex - -livv -V \ - config/${CASE}/cmb_gis.yml \ - config/${CASE}/smb_gis.yml \ - config/${CASE}/energy_e3sm_racmo_gis.yml \ - config/${CASE}/energy_e3sm_era5_gis.yml \ - config/${CASE}/energy_e3sm_merra2_merra_grid_gis.yml \ - config/${CASE}/energy_e3sm_ceres_gis.yml \ - config/${CASE}/cmb_ais.yml \ - config/${CASE}/smb_ais.yml \ - config/${CASE}/energy_e3sm_racmo_ais.yml \ - config/${CASE}/energy_e3sm_era5_ais.yml \ - config/${CASE}/energy_e3sm_merra2_merra_grid_ais.yml \ - config/${CASE}/energy_e3sm_ceres_ais.yml \ - -o $SCRATCH/lex/${CASE} &> livv_log_${CASE}.log +livv \ + --validate config/${CASE}/livvkit.yml \ + --out-dir $SCRATCH/lex/${CASE} >> livv_stdoe_${CASE}.log 2>&1 # Backup the existing published version of this analysis mv ${WEBDIR}/${CASE} ${WEBDIR}/${CASE}_bkd_$(date +'%Y%m%dT%H%M') @@ -42,3 +77,5 @@ for htmlfile in validation/*.html do sed -i "s/\(<\!-- GROUP LINK-->\)/\\nCurrent\ runs\n\<\/a\>/g" ${htmlfile} done +echo "LIVVkit results availalble at:" +echo "https://portal.nersc.gov/project/e3sm/${USER}/${CASE}/index.html" diff --git a/tests/simple_test.yml b/tests/simple_test.yml index bfd88d0..94336fe 100644 --- a/tests/simple_test.yml +++ b/tests/simple_test.yml @@ -1,4 +1,9 @@ --- +common: + meta: &meta + Case ID: TestCaseId + Clim Years: 1890-1906 ExtensionTest: module: tests/extension_simple.py description: A minimal LIVVkit extensions test. + meta: *meta diff --git a/tests/test_common.py b/tests/test_common.py new file mode 100644 index 0000000..83b56a1 --- /dev/null +++ b/tests/test_common.py @@ -0,0 +1,54 @@ +import numpy as np +import xarray as xr + +import lex.common as lxc + + +def test_img_file_prefix(): + config_1 = {"module": "energy.py"} + config_2 = {"EMPTY": "empty.py"} + + assert "energy" == lxc.img_file_prefix(config_1), "PREFIX DOES NOT MATCH" + try: + _mod = lxc.img_file_prefix(config_2) + except KeyError as _err: + assert "module" in _err.args, "UN-EXPECTED ERROR" + + +def test_check_longitude(): + lons = np.arange(-180, 180, 45.0) + data = {"lon": lons} + b_coord = "longitude" + data_b = {b_coord: lons} + assert lxc.check_longitude(data) == {"lon": lons} + assert lxc.check_longitude(data_b, lon_coord=b_coord) == {b_coord: lons} + + lons_2 = np.arange(0, 360, 45.0) + _shift_lons = np.array([0.0, 45.0, 90.0, 135.0, -180.0, -135.0, -90.0, -45.0]) + _check = lxc.check_longitude({"lon": lons_2}) + assert np.all(_check["lon"] == _shift_lons), "SHIFTED LONGITUDE DOES NOT MATCH" + + test_data = [ + 0, + 1, + 2, + 3, + 4, + 5, + 6, + 7, + ] + test_data_shift = [4, 5, 6, 7, 0, 1, 2, 3] + _ds1 = _ds1 = xr.Dataset( + data_vars={"x": (["lon"], test_data)}, coords={"lon": lons_2} + ) + _check = lxc.check_longitude(_ds1) + assert isinstance(_check, xr.Dataset), "CHECK_LON DID NOT RETURN XARRAY.DATASET" + assert (_check["lon"] == lons).all(), "LONGITUDE NOT SHIFTED PROPERLY" + assert (_check["x"] == test_data_shift).all(), "DATA NOT ROLLED PROPERLY" + + _da1 = xr.DataArray(test_data, coords={"lon": lons_2}, dims=("lon",)) + _check = lxc.check_longitude(_da1) + assert isinstance(_check, xr.DataArray), "CHECK_LON DID NOT RETURN XARRAY.DATAARRAY" + assert (_check["lon"] == lons).all() + assert (_check.values == test_data_shift).all()