Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add plotting script for thermal forcing and surface mass balance #566

Open
wants to merge 7 commits into
base: master
Choose a base branch
from

Conversation

trhille
Copy link
Collaborator

@trhille trhille commented Apr 23, 2024

This merge adds a script that plots thermal forcing and surface mass balance time series. Currently the thermal forcing-related fields are hard-coded to ISMIP6 variable names, but that could be changed in the future. The user must specify files from which to plot thermal forcing and SMB; a mesh file containing coordinates, ice geometry, etc; and optionally a regions file that contains the regionCellMasks field.
There are many (maybe too many) user-defined options for plotting:

  1. The option to plot at (i) a specified depth (or average over a depth range) using the -d argument, (ii) to interpolate the thermal forcing field to the seafloor (--seafloor), ice-shelf base (--shelf_base), or (iii) average over the full ice-shelf depth range (--shelf_range).
  2. The option to specify a csv file containing x and y coordinates for thermal forcing using the -c flag, or to specify coordinates in the command using -x and -y as comma-separated lists.
  3. The option to set start and end time indices using --start_time and --end_time.
  4. The option to take a random sample of point from the thermal forcing coordinates using --n_samples. This is just for prototyping, and the random samples are not currently reproducible.
  5. The option to average thermal forcing over all cells using the --average flag. If this flag is not set, a separate time series will be plotted for every cell, which will get messy when plotting more than a few cells. Note that averaging is not currently weighted by cell area.
  6. The option to save .png and .pdf files using the --save flag.
  7. The option to plot SMB for a specified region using the -n flag, which corresponds to the region number in regionCellMasks.

Plot ISMIP6 thermal forcing time series for Amery/Lambert UQ paper.
Add option to plot SMB across a region or the entire domain.
Update linestyles based on ESM, save a pdf version, use lowercase
axis labels, etc.
Support plotting average over depth range, for example using '-d -200,-500'
in the command.
Add option to plot over whole shelf depth range, using the
--shelf_range flag. This may not be very accurate because it will
sample deep temperatures far from the grounding line, where they
are unlikely to ever be felt by the ice shelf during a simulation,
but it is much faster than interpolating to a specific depth or
to the ice-shelf base or seafloor.
@trhille
Copy link
Collaborator Author

trhille commented Apr 23, 2024

Example plot, created using the following command:

python plot_forcing_time_series.py -t /global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE02/AIS_4to20km_r01_20220907_TF_CCSM4_RCP85_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE03/AIS_4to20km_r01_20220907_TF_HadGEM2-ES_RCP85_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE04/AIS_4to20km_r01_20220907_TF_CESM2-WACCM_SSP585_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE05/AIS_4to20km_r01_20220907_TF_UKESM1-0-LL_SSP585_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier2_experiments/expAE10/AIS_4to20km_TF_UKESM1-0-LL_SSP126_2300.nc \
-s /global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE02/AIS_4to20km_r01_20220907_SMB_CCSM4_RCP85_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE03/AIS_4to20km_r01_20220907_SMB_HadGEM2-ES_RCP85_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE04/AIS_4to20km_r01_20220907_SMB_CESM2-WACCM_SSP585_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier1_experiments/expAE05/AIS_4to20km_r01_20220907_smb_UKESM1-0-LL_SSP585_2300.nc,/global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/forcing/ais_mesh_4to20km_res/tier2_experiments/expAE10/AIS_4to20km_SMB_UKESM1-0-LL_SSP126_2300.nc -m /global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/initial_conditions/AIS_4to20km_20230105/relaxation_0TGmelt_10yr/relaxed_10yrs_4km.nc \
-r /global/cfs/cdirs/fanssie/MALI_projects/ISMIP6-2300/initial_conditions/AIS_4to20km_20230105/AIS_4to20km_r01_20220907.regionMask_ismip6.nc\
-c /global/cfs/cdirs/fanssie/MALI_projects/Amery_UQ/amery.csv -n 2 --shelf_base \
--average --start_time 5 --end_time -1 --save /pscratch/sd/t/trhille/Amery_TF_SMB

Amery_TF_SMB

Copy link
Member

@matthewhoffman matthewhoffman left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@trhille , thanks for making a PR to capture this script. I have some cleanup suggestions and a couple questions about the functionality.

landice/output_processing_li/plot_forcing_time_series.py Outdated Show resolved Hide resolved
landice/output_processing_li/plot_forcing_time_series.py Outdated Show resolved Hide resolved
landice/output_processing_li/plot_forcing_time_series.py Outdated Show resolved Hide resolved
landice/output_processing_li/plot_forcing_time_series.py Outdated Show resolved Hide resolved
Comment on lines +152 to +153
# Assume this is floating ice
depth = -1. * rhoi / rhosw * thk[0, nearest_cell]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you want to do this? Because TF is defined everywhere, it would be easy to inadvertently include values inside grounded ice. Maybe that's desired, though generally probably not?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Because this only plots TF at user-specified locations and the concept of the shelf base is meaningless for grounded ice, I don't think this is a problem. But I could add an assertion that depth >= bed.

Comment on lines +34 to +41
parser.add_option('-c', dest='coords_file', default=None,
help='CSV file containing x in first column, y in second. No header.')
parser.add_option('-x', dest='x_coords', default=None,
help='List of x coordinates of transect if not \
using a csv file. Comma-separated, no spaces.')
parser.add_option('-y', dest='y_coords', default=None,
help='List of y coordinates of transect if not \
using a csv file. Comma-separated, no spaces.')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

From what I understand of the code below, it looks like it is necessary to specify x,y locations in one of these two ways. Is that right? Is there not a way to just use all data in a region, say? I'm not asking for you to generalize that behavior if it doesn't already exist, but clarify, probably in the main doc strings, how this works. Also that all the depth options apply over just the horizontal locations as defined by these options, i.e. plot_shelf_base_thermal_forcing is not the entire shelf base but the shelf base at the x,y locations determined by these horizontal options. (if I follow correctly)

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There's not currently a way to use all data in a region, but that would be good to add. I'll update the doc strings.

Comment on lines +48 to +49
parser.add_option('-n', dest='region_number', default=None,
help='Region number to plot. If None, use entire domain.')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It looks like the region only applies to the SMB plot - is that intended?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's correct. I do want to generalize it, but it wasn't necessary for the Lambert paper figure (although it might have made things easier if I had started with that). I can add a comment for the region is just SMB for now.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants