Source code for compass.landice.tests.ismip6_forcing.atmosphere.process_smb_racmo

import os
import shutil

import xarray as xr
from import write_netcdf
from mpas_tools.logging import check_call

from compass.landice.tests.ismip6_forcing.create_mapfile import (
from compass.step import Step

[docs] class ProcessSmbRacmo(Step): """ A step for processing the regional RACMO surface mass balance data """
[docs] def __init__(self, test_case, input_file=None): """ Create the step Parameters ---------- test_case : compass.landice.tests.ismip6_forcing.atmosphere.Atmosphere The test case this step belongs to input_file : file name of ismip6 forcing data processed by this step """ self.input_file = input_file super().__init__(test_case=test_case, name='process_smb_racmo', ntasks=4, min_tasks=1)
[docs] def setup(self): """ Set up this step of the test case """ config = self.config section = config["ismip6_ais"] base_path_mali = section.get("base_path_mali") mali_mesh_name = section.get("mali_mesh_name") mali_mesh_file = section.get("mali_mesh_file") section = config["ismip6_ais_atmosphere"] process_smb_racmo = section.getboolean("process_smb_racmo") self.add_input_file(filename=mali_mesh_file, target=os.path.join(base_path_mali, mali_mesh_file)) if process_smb_racmo: self.add_input_file( target="", database="RACMO2.3p2_ANT27_SMB_yearly_1979_2018") output_file = f"{mali_mesh_name}_RACMO2.3p2_ANT27" \ f"" self.add_output_file(filename=output_file) else: print("\n'Warning: process_smb_racmo' is set to 'False'." " This step will not run unless set 'True' in the" " config file.\n")
[docs] def run(self): """ Run this step of the test case """ logger = self.logger config = self.config section = config["ismip6_ais_atmosphere"] process_smb_racmo = section.getboolean("process_smb_racmo") if not process_smb_racmo: # we don't want to run this step return section = config["ismip6_ais"] mali_mesh_name = section.get("mali_mesh_name") mali_mesh_file = section.get("mali_mesh_file") output_base_path = section.get("output_base_path") section = config["ismip6_ais_atmosphere"] method_remap = section.get("method_remap") racmo_file_temp1 = "" racmo_file_temp2 = "RACMO2.3p2_smb_climatology_1995_2017_" \ "" output_file = f"{mali_mesh_name}_RACMO2.3p2_ANT27" \ f"" output_path = f"{output_base_path}/atmosphere_forcing/" \ f"RACMO_climatology_1995-2017" output_path_final = os.path.join(output_base_path, output_path) if os.path.exists(os.path.join(output_path_final, output_file)):"Processed RACMO SMB data already exists in the " f"output path {output_base_path}. " f"Not processing the source RACMO data...") return input_file = self.inputs[1] # take the time average over the period 1995-2017 args = ["ncra", "-O", "-F", "-d", "time,17,39", input_file, racmo_file_temp1] check_call(args, logger=logger) # interpolate the racmo smb data remapped_file_temp = "" # temporary file name # call the function that reads in, remap and rename the file."Calling the remapping function...") self.remap_source_smb_to_mali(racmo_file_temp1, remapped_file_temp, mali_mesh_name, mali_mesh_file, method_remap) # perform algebraic operation on the source data in unit of kg/m^2 # to be in unit of kg/m^2/s args = ["ncap2", "-O", "-v", "-s", "sfcMassBal=smb/(60*60*24*365)", remapped_file_temp, racmo_file_temp2] check_call(args, logger=logger) # change the unit attribute to kg/m^2/s args = ["ncatted", "-O", "-a", "units,sfcMassBal,m,c,'kg m-2 s-1'", racmo_file_temp2] check_call(args, logger=logger) # call the function that renames the ismip6 variables to MALI variables"Renaming source variables to mali variable names...") self.rename_source_smb_to_mali_vars(racmo_file_temp2, output_file)"Processing done successfully. " "Removing the temporary files...") # remove the temporary remapped and combined files os.remove(remapped_file_temp) os.remove(racmo_file_temp1) os.remove(racmo_file_temp2) # place the output file in appropriate directory output_path = f"{output_base_path}/atmosphere_forcing/" \ f"RACMO_climatology_1995-2017" if not os.path.exists(output_path):"Creating a new directory for the output data") os.makedirs(output_path) src = os.path.join(os.getcwd(), output_file) dst = os.path.join(output_path, output_file) shutil.copy(src, dst)"!---Done processing the file---!")
[docs] def remap_source_smb_to_mali(self, input_file, output_file, mali_mesh_name, mali_mesh_file, method_remap): """ Remap the input ismip6 smb forcing data onto mali mesh Parameters ---------- input_file: str input racmo smb data on its native rotated pole grid output_file : str smb data remapped on mali mesh mali_mesh_name : str name of the mali mesh used to name mapping files mali_mesh_file : str, optional The MALI mesh file if mapping file does not exist method_remap : str, optional Remapping method used in building a mapping file """ mapping_file = f"map_racmo_24km_to_" \ f"{mali_mesh_name}_{method_remap}.nc" # check if a mapfile exists if not os.path.exists(mapping_file): # build a mapping file if it doesn't already exist"Creating a mapping file. " f"Mapping method used: {method_remap}") build_mapping_file(self.config, self.ntasks, self.logger, input_file, mapping_file, scrip_from_latlon=True, mali_mesh_file=mali_mesh_file, method_remap=method_remap) else:"Mapping file exists. " "Remapping the input data...") args = ["ncremap", "-i", input_file, "-o", output_file, "-m", mapping_file, "-v", "smb"] check_call(args, logger=self.logger)
[docs] def rename_source_smb_to_mali_vars(self, remapped_file_temp, output_file): """ Rename variables in the remapped source input data to the ones that MALI uses. Parameters ---------- remapped_file_temp : str temporary ismip6 data remapped on mali mesh output_file : str remapped ismip6 data renamed on mali mesh """ # open dataset in 20 years chunk ds = xr.open_dataset(remapped_file_temp, chunks=dict(time=20), engine="netcdf4") # build dictionary for ismip6 variables that MALI takes in ismip6_to_mali_dims = dict( time="Time", ncol="nCells") ds = ds.rename(ismip6_to_mali_dims) # drop unnecessary variables ds = ds.drop_vars(["height"]) # squeeze unnecessary coordinate variable ds["sfcMassBal"] = ds["sfcMassBal"].squeeze(dim="height") # write to a new netCDF file write_netcdf(ds, output_file) ds.close()
[docs] def correct_smb_anomaly_for_base_smb(self, output_file, mali_mesh_file): """ Apply the MALI base SMB to the ismip6 SMB anomaly field Parameters ---------- output_file : str remapped ismip6 data renamed on mali mesh mali_mesh_file : str initialized MALI mesh file in which the base SMB field exists """ ds = xr.open_dataset(output_file) ds_base = xr.open_dataset(mali_mesh_file) # get the first time index ref_smb = ds_base["sfcMassBal"].isel(Time=0) # correct for the reference smb ds["sfcMassBal"] = ds["sfcMassBal"] + ref_smb # write to a new netCDF file write_netcdf(ds, output_file) ds.close()