Newer
Older
"""Compute the sea ice extent , area and volume in both hemispheres or a specified region"""
Javier Vegas-Regidor
committed
import iris
import iris.analysis
import iris.coords
Javier Vegas-Regidor
committed
import iris.util
from earthdiagnostics.constants import Basins
from earthdiagnostics.diagnostic import Diagnostic, \
DiagnosticBasinListOption, DiagnosticBoolOption
from earthdiagnostics.modelingrealm import ModelingRealms
from earthdiagnostics.utils import Utils, TempFile
Javier Vegas-Regidor
committed
import diagonals.siasie as siasie
from diagonals.mesh_helpers.nemo import Nemo
# noinspection PyUnresolvedReferences
Javier Vegas-Regidor
committed
class Siasiesiv(Diagnostic):
Compute the sea ice extent , area and volume in both hemispheres or a
specified region.
Parameters
----------
data_manager: DataManager
startdate: str
member: int
chunk: init
domain: ModellingRealm
variable: str
Javier Vegas-Regidor
committed
basin: list of Basin
alias = 'siasiesiv'
"Diagnostic alias for the configuration file"
def __init__(self, data_manager, startdate, member, chunk, masks,
var_manager, data_convention, omit_vol):
Javier Vegas-Regidor
committed
Diagnostic.__init__(self, data_manager)
self.startdate = startdate
self.member = member
self.chunk = chunk
Javier Vegas-Regidor
committed
self.masks = masks
self.var_manager = var_manager
self.omit_volume = omit_vol
self.sic_varname = self.var_manager.get_variable('sic').short_name
self.sit_varname = self.var_manager.get_variable('sit').short_name
Javier Vegas-Regidor
committed
Javier Vegas-Regidor
committed
self.results = {}
for var in ('siarean', 'siareas', 'siextentn', 'siextents'):
Javier Vegas-Regidor
committed
self.results[var] = {}
return 'Siasiesiv Startdate: {0.startdate} Member: {0.member} Chunk: {0.chunk} ' \
Javier Vegas-Regidor
committed
'Basins: {1} Omit volume: {0.omit_volume}'.format(self,
','.join(str(basin) for basin in self.masks.keys()))
def generate_jobs(cls, diags, options):
Create a job for each chunk to compute the diagnostic
:param diags: Diagnostics manager class
:type diags: Diags
:param options: basin
options_available = (DiagnosticBasinListOption('basins',
[Basins().Global]),
options = cls.process_options(options, options_available)
basins = options['basins']
if not basins:
Javier Vegas-Regidor
committed
Log.error('Basins not recognized')
return ()
Javier Vegas-Regidor
committed
masks = {}
basins.sort()
for basin in basins:
Javier Vegas-Regidor
committed
masks[basin] = Utils.get_mask(basin)
Javier Vegas-Regidor
committed
job_list = list()
for startdate, member, chunk in diags.config.experiment.get_chunk_list():
job_list.append(Siasiesiv(diags.data_manager, startdate, member,
chunk, masks, diags.config.var_manager,
diags.config.data_convention,
Javier Vegas-Regidor
committed
Javier Vegas-Regidor
committed
return job_list
if not self.omit_volume:
self.sit = self.request_chunk(ModelingRealms.seaIce,
self.sit_varname,
self.startdate,
self.member,
self.chunk)
self.sic = self.request_chunk(ModelingRealms.seaIce, self.sic_varname,
self.startdate, self.member, self.chunk)
if not self.omit_volume:
self._declare_var('sivols')
self._declare_var('sivoln')
self._declare_var('siareas')
self._declare_var('siextents')
self._declare_var('siarean')
self._declare_var('siextentn')
def _declare_var(self, var_name):
self.generated[var_name] = self.declare_chunk(ModelingRealms.seaIce,
var_name, self.startdate,
self.member, self.chunk)
Javier Vegas-Regidor
committed
def compute(self):
self._fix_coordinates_attribute(
self.sic.local_file, self.sic_varname
)
Javier Vegas-Regidor
committed
sic = iris.load_cube(self.sic.local_file)
if sic.units.origin == '%' and sic.data.max() < 2:
Javier Vegas-Regidor
committed
sic.units = '1.0'
sic_slices = []
for sic_data in sic.slices_over('time'):
sic_data.data = np.ma.filled(sic_data.data, 0.0).astype(np.float32)
sic_slices.append(sic_data)
mesh = Nemo('mesh_hgr.nc', 'mask_regions.nc')
areacello = mesh.get_areacello(cell_point='T')
gphit = mesh.get_grid_latitude(cell_point='T')
self.results['siextentn'], self.results['siextents'], self.results['siarean'], self.results['siareas'] = siasie.compute(gphit, areacello, sic_slices, self.masks)
Javier Vegas-Regidor
committed
self.save()
def _fix_coordinates_attribute(self, filepath, var_name):
add_coordinates = {
'time', 'leadtime', 'time_centered',
self.data_convention.lon_name, self.data_convention.lat_name
}
handler = Utils.open_cdf(filepath)
coordinates = handler.variables[var_name].coordinates.split()
handler.variables[var_name].coordinates = \
' '.join(set(coordinates) | add_coordinates)
handler.close()
Javier Vegas-Regidor
committed
def save(self):
Javier Vegas-Regidor
committed
results = iris.cube.CubeList()
for time, basin in six.iteritems(time):
for basin, result in six.iteritems(basin):
cube = iris.cube.Cube(result)
cube.var_name = var
cube.units = 'm^2'
cube.add_aux_coord(iris.coords.AuxCoord(basin.name,
var_name='region'))
cube.add_aux_coord(iris.coords.AuxCoord(time,
var_name='time'))
self._save_file(results.merge_cube(), var)
Javier Vegas-Regidor
committed
def _save_file(self, data, var):
generated_file = self.generated[var]
temp = TempFile.get()
Javier Vegas-Regidor
committed
region = data.coord('region').points
data.remove_coord('region')
iris.save(data, temp, zlib=True)
handler = Utils.open_cdf(temp)
var = handler.createVariable('region2', str, ('region',))
var[...] = region
handler.close()
Utils.rename_variable(temp, 'region2', 'region', True)
else:
handler = Utils.open_cdf(temp)
if 'region' not in handler.dimensions:
new_file = TempFile.get()
new_handler = Utils.open_cdf(new_file, 'w')
new_handler.createDimension('region', 1)
for dimension in handler.dimensions:
Utils.copy_dimension(handler, new_handler, dimension)
for variable in handler.variables.keys():
if variable in (var, 'region'):
continue
Utils.copy_variable(handler, new_handler, variable)
old_var = handler.variables[var]
new_var = new_handler.createVariable(var, old_var.dtype,
('region',) + old_var.dimensions,
zlib=True, fill_value=1.0e20)
Utils.copy_attributes(new_var, old_var)
new_var[0, :] = old_var[:]
new_var = new_handler.createVariable('region', str,
('region',))
new_var[0] = region[0]
new_handler.close()
os.remove(temp)
temp = new_file
handler.close()
generated_file.set_local_file(temp)