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

Update mass conservation tables to include Ref & Dev in the same file #309

Merged
merged 12 commits into from
May 15, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,10 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
- Added fixed level budget diagnostic to budget operations table
- Function `rename_speciesconc_to_speciesconcvv` in `gcpy/benchmark/modules/benchmark_utils.py`
- Function `copy_file_to_dir` in `gcpy/util.py`. This is a wrapper for `shutil.copyfile`.
- Script `gcpy/benchmark/modules/benchmark_mass_cons_table.py`, with code to create mass conservation tables
- Expanded statistics output in benchmark mass conservation tables
- Function `get_datetimes_from_filenames` in `gcpy/benchmark/modules/benchmark_utils.py`
- Function `replace_whitespace` in `gcpy/util.py`

### Changed
- Bump pip from 23.2.1 to 23.3 (dependabot suggested this)
Expand Down Expand Up @@ -58,6 +62,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
- Environment file `docs/environment_files/environment.yml`
- Environment file `docs/environment_files/requirements.txt`
- Removed `awscli` from the GCPy environment; version 2 is no longer available on conda-forge or PyPi
- Routine `make_benchmark_mass_conservation_table` in `benchmark_funcs.py`; this is now obsolete

## [1.4.2] - 2024-01-26
### Added
Expand Down
166 changes: 0 additions & 166 deletions gcpy/benchmark/modules/benchmark_funcs.py
Original file line number Diff line number Diff line change
Expand Up @@ -4915,172 +4915,6 @@ def make_benchmark_operations_budget(
gc.collect()


def make_benchmark_mass_conservation_table(
datafiles,
runstr,
dst="./benchmark",
overwrite=False,
areapath=None,
spcdb_dir=os.path.dirname(__file__)
):
"""
Creates a text file containing global mass of the PassiveTracer
from Transport Tracer simulations across a series of restart files.

Args:
datafiles: list of str
Path names of restart files.
runstr: str
Name to put in the filename and header of the output file

Keyword Args (optional):
dst: str
A string denoting the destination folder where the file
containing emissions totals will be written.
Default value: "./benchmark"
overwrite: bool
Set this flag to True to overwrite files in the
destination folder (specified by the dst argument).
Default value: False
areapath: str
Path to a restart file containing surface area data.
Default value: None
spcdb_dir: str
Path to the species_database.yml
Default value: points to gcpy/gcpy folder
"""

# ==================================================================
# Initialize
# ==================================================================

# Create the destination folder
util.make_directory(dst, overwrite)

# Load a YAML file containing species properties (such as
# molecular weights), which we will need for unit conversions.
properties = util.read_config_file(
os.path.join(
spcdb_dir,
"species_database.yml"
),
quiet=True
)

# Get the species name
spc_name = 'PassiveTracer'

# Get a list of properties for the given species
species_properties = properties.get(spc_name)

# Specify target units
target_units = "Tg"

dates = []
masses = []

# ==================================================================
# Make sure that surface area data is found
# ==================================================================
with warnings.catch_warnings():
warnings.filterwarnings("ignore", category=xr.SerializationWarning)


# ==================================================================
# Calculate global mass for the tracer at all restart dates
# ==================================================================
for f in datafiles:
ds = xr.open_dataset(f, drop_variables=skip_these_vars)

# Save date in desired format
#datestr = str(pd.to_datetime(ds.time.values[0]))
#dates.append(datestr[:4] + '-' + datestr[5:7] + '-' + datestr[8:10])

# Find the area variable in Dev
if areapath is None:
area = util.get_area_from_dataset(ds)
else:
area = util.get_area_from_dataset(
xr.open_dataset(
areapath,
drop_variables=skip_these_vars
)
)

# Assume typical restart file name format, but avoid using dates
# from within files which may be incorrect for the initial restart
datestr = f.split('/')[-1].split('.')[2][:9]
dates.append(datestr[:4] + '-' + datestr[4:6] + '-' + datestr[6:8])

# Select for GCC or GCHP
delta_p = ds['Met_DELPDRY'] if 'Met_DELPDRY' in list(ds.data_vars) else ds['DELP_DRY']

# ==============================================================
# Convert units of Ref and save to a DataArray
# (or skip if Ref contains NaNs everywhere)
# ==============================================================
# Select for GCC or GCHP
if 'SpeciesRst_PassiveTracer' in list(ds.data_vars):
attrs = ds['SpeciesRst_PassiveTracer'].attrs
da = ds['SpeciesRst_PassiveTracer'].astype(np.float64)
da.attrs = attrs
else:
attrs = ds['SPC_PassiveTracer'].attrs
da = ds['SPC_PassiveTracer'].astype(np.float64)
da.attrs = attrs
da = convert_units(
da,
spc_name,
species_properties,
target_units,
area_m2=area,
delta_p=delta_p
)

# Save total global mass
masses.append(np.sum(da.values))

# Clean up
del ds
del da
gc.collect()

# Calclate max and min mass, absolute diff, percent diff
max_mass = np.max(masses)
min_mass = np.min(masses)
# Convert absdiff to grams
absdiff = (max_mass-min_mass) * 10**12
pctdiff = (max_mass-min_mass)/min_mass * 100

# ==================================================================
# Print masses to file
# ==================================================================
# Create file
outfilename = os.path.join(dst, f"Passive_mass.{runstr}.txt")

with open(outfilename, 'w') as f:
titlestr = ' Global Mass of Passive Tracer in ' + runstr + ' '
#headers
print('%' * (len(titlestr)+4), file=f)
print(titlestr, file=f)
print('%' * (len(titlestr)+4), file=f)
print('', file=f)
print(' Date' + ' ' * 8 + 'Mass [Tg]', file=f)
print(' ' + '-' * 10 + ' ' + '-' * 16, file=f)
#masses
for i in range(len(masses)):
print(f" {dates[i]} {masses[i] : 11.13f}", file=f)
print(' ', file=f)
print(' Summary', file=f)
print(' ' + '-' * 30, file=f)
print(f" Max mass = {max_mass : 2.13f} Tg", file=f)
print(f" Min mass = {min_mass : 2.13f} Tg", file=f)
print(f" Abs diff = {absdiff : >16.3f} g", file=f)
print(f" Pct diff = {pctdiff : >16.10f} %", file=f)

gc.collect()


def get_species_database_dir(config):
"""
Returns the directory in which the species_database.yml file is
Expand Down
Loading