# coding=utf-8
"""Compute the sea ice extent , area and volume in both hemispheres or a specified region"""
import os
import six
import numpy as np
import iris
import iris.analysis
import iris.coords
import iris.util
from bscearth.utils.log import Log
from earthdiagnostics.constants import Basins
from earthdiagnostics.diagnostic import Diagnostic, DiagnosticBasinListOption, DiagnosticBoolOption
from earthdiagnostics.modelingrealm import ModelingRealms
from earthdiagnostics.utils import Utils, TempFile
# noinspection PyUnresolvedReferences
[docs]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
basin: list of Basin
mask: numpy.array
omit_vol: bool
"""
alias = 'siasiesiv'
"Diagnostic alias for the configuration file"
e1t = None
e2t = None
gphit = None
def __init__(self, data_manager, startdate, member, chunk, masks, var_manager, data_convention, omit_vol):
Diagnostic.__init__(self, data_manager)
self.startdate = startdate
self.member = member
self.chunk = chunk
self.masks = masks
self.generated = {}
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
self.data_convention = data_convention
self.results = {}
for var in ('siarean', 'siareas', 'sivoln', 'sivols', 'siextentn', 'siextents'):
self.results[var] = {}
def __str__(self):
return 'Siasiesiv Startdate: {0.startdate} Member: {0.member} Chunk: {0.chunk} ' \
'Basins: {1} Omit volume: {0.omit_volume}'.format(self,
','.join(str(basin) for basin in self.masks.keys()))
[docs] @classmethod
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
:type options: list[str]
:return:
"""
options_available = (DiagnosticBasinListOption('basins', [Basins().Global]),
DiagnosticBoolOption('omit_volume', False))
options = cls.process_options(options, options_available)
basins = options['basins']
if not basins:
Log.error('Basins not recognized')
return ()
masks = {}
basins.sort()
for basin in basins:
masks[basin] = Utils.get_mask(basin)
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,
options['omit_volume']))
e1t = iris.load_cube('mesh_hgr.nc', 'e1t')
e2t = iris.load_cube('mesh_hgr.nc', 'e2t')
Siasiesiv.area = e1t * e2t
return job_list
[docs] def request_data(self):
"""Request data required by the diagnostic"""
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)
[docs] def declare_data_generated(self):
"""Declare data to be generated by the diagnostic"""
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)
[docs] def compute(self):
"""Run the diagnostic"""
coordinates = ' '.join(('time', 'leadtime', 'time_centered',
self.data_convention.lon_name, self.data_convention.lat_name))
handler = Utils.open_cdf(self.sic.local_file)
handler.variables[self.sic_varname].coordinates = coordinates
handler.close()
sic = iris.load_cube(self.sic.local_file)
if sic.units.origin == '%' and sic.data.max() < 2:
sic.units = '1.0'
extent = sic.copy((sic.data >= 0.15).astype(np.int8)) * Siasiesiv.area.data
sic *= Siasiesiv.area.data
if not self.omit_volume:
handler = Utils.open_cdf(self.sit.local_file)
handler.variables[self.sit_varname].coordinates = coordinates
handler.close()
sit = iris.load_cube(self.sit.local_file)
for basin, mask in six.iteritems(self.masks):
self.results['siarean'][basin] = self.sum(sic, mask, north=True)
self.results['siareas'][basin] = self.sum(sic, mask, north=False)
if not self.omit_volume:
volume = sic * sit
self.results['sivoln'][basin] = self.sum(volume, mask, north=True)
self.results['sivols'][basin] = self.sum(volume, mask, north=False)
self.results['siextentn'][basin] = self.sum(extent, mask, north=True)
self.results['siextents'][basin] = self.sum(extent, mask, north=False)
self.save()
def sum(self, data, mask, north=True):
if north:
condition = data.coord('latitude').points > 0
else:
condition = data.coord('latitude').points < 0
weights = iris.util.broadcast_to_shape(condition, data.shape, data.coord_dims('latitude')) * mask
return data.collapsed(('latitude', 'longitude'), iris.analysis.SUM, weights=weights)
def save(self):
for var, basins in six.iteritems(self.results):
results = iris.cube.CubeList()
for basin, result in six.iteritems(basins):
result.var_name = var
if var.startswith('sivol'):
result.units = 'm^3'
else:
result.units = 'm^2'
result.add_aux_coord(iris.coords.AuxCoord(basin.name, var_name='region'))
results.append(result)
if not results:
continue
self._save_file(results.merge_cube(), var)
def _save_file(self, data, var):
generated_file = self.generated[var]
temp = TempFile.get()
region = data.coord('region').points
data.remove_coord('region')
iris.save(data, temp, zlib=True)
if len(region) > 1:
Utils.rename_variable(temp, 'dim0', 'region', False)
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)