Newer
Older
import netCDF4
import os
from earthdiagnostics.constants import Basins
from earthdiagnostics.diagnostic import Diagnostic
from earthdiagnostics.utils import Utils, TempFile
import cdftoolspython
Javier Vegas-Regidor
committed
class Siasiesiv(Diagnostic):
Compute the sea ice extent , area and volume in both hemispheres or a specified region.
:original author: Virginie Guemas <virginie.guemas@bsc.es>
:contributor: Neven Fuckar <neven.fuckar@bsc.es>
:contributor: Ruben Cruz <ruben.cruzgarcia@bsc.es>
:contributor: Javier Vegas-Regidor <javier.vegas@bsc.es>
:created: April 2012
:last modified: June 2016
alias = 'siasiesiv'
"Diagnostic alias for the configuration file"
e1t = None
e2t = None
gphit = None
def __init__(self, data_manager, basin, startdate, member, chunk, mask):
"""
:param data_manager: data management object
:type data_manager: DataManager
:param startdate: startdate
:type startdate: str
:param member: member number
:type member: int
:param chunk: chunk's number
:type chunk: int
:param mask: mask to use
:type mask: numpy.array
"""
Javier Vegas-Regidor
committed
Diagnostic.__init__(self, data_manager)
self.basin = basin
self.startdate = startdate
self.member = member
self.chunk = chunk
Javier Vegas-Regidor
committed
self.required_vars = ['sit', 'sic']
self.generated_vars = ['siextents', 'sivols', 'siareas', 'siextentn', 'sivoln', 'siarean']
def __str__(self):
return 'Siasiesiv Startdate: {0} Member: {1} Chunk: {2} Basin: {3}'.format(self.startdate, self.member,
self.chunk, self.basin.fullname)
def generate_jobs(cls, diags, options):
"""
Creates a job for each chunk to compute the diagnostic
:param diags: Diagnostics manager class
:type diags: Diags
:param options: basin
:return:
"""
if len(options) != 2:
raise Exception('You must specify the basin for the siasiesiv diagnostic (and nothing else)')
Javier Vegas-Regidor
committed
basin = Basins.parse(options[1])
mask = Utils.get_mask(basin)
Javier Vegas-Regidor
committed
job_list = list()
Javier Vegas-Regidor
committed
for startdate, member, chunk in diags.exp_manager.get_chunk_list():
job_list.append(Siasiesiv(diags.data_manager, basin, startdate, member, chunk, mask))
mesh_handler = Utils.openCdf('mesh_hgr.nc')
Siasiesiv.e1t = np.asfortranarray(mesh_handler.variables['e1t'][0, :])
Siasiesiv.e2t = np.asfortranarray(mesh_handler.variables['e2t'][0, :])
Siasiesiv.gphit = np.asfortranarray(mesh_handler.variables['gphit'][0, :])
mesh_handler.close()
Javier Vegas-Regidor
committed
return job_list
def compute(self):
sit_file = self.data_manager.get_file('seaIce', 'sit', self.startdate, self.member, self.chunk)
sit_handler = Utils.openCdf(sit_file)
sit = np.asfortranarray(sit_handler.variables['sit'][:])
timesteps = sit_handler.dimensions['time'].size
sit_handler.close()
sic_file = self.data_manager.get_file('seaIce', 'sic', self.startdate, self.member, self.chunk)
sic = np.asfortranarray(sic_handler.variables['sic'][:])
sic_handler.close()
result = np.empty((8, timesteps))
for t in range(0, timesteps):
result[:, t] = cdftoolspython.icediag.icediags(Siasiesiv.e1t, Siasiesiv.e2t, self.mask,
Siasiesiv.gphit, sit[t, :], sic[t, :])
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[4, :], 'sivols', '10^9 m3'),
'seaIce', 'sivols', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[5, :], 'siareas', '10^9 m2'),
'seaIce', 'siareas', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[7, :], 'siextents', '10^9 m2'),
'seaIce', 'siextents', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
Javier Vegas-Regidor
committed
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[0, :], 'sivoln', '10^9 m3'),
'seaIce', 'sivoln', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[1, :], 'siarean', '10^9 m2'),
'seaIce', 'siarean', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
Javier Vegas-Regidor
committed
self.data_manager.send_file(self._extract_variable_and_rename(sit_file, result[3, :], 'siextentn', '10^9 m2'),
'seaIce', 'siextentn', self.startdate, self.member, self.chunk,
region=self.basin.fullname)
@staticmethod
def _extract_variable_and_rename(reference_file, values, cmor_name, units):
temp = TempFile.get()
os.remove(temp)
handler = netCDF4.Dataset(temp, 'w')
# Create dimensions
handler.createDimension('bnds', 2)
# Copy time variable
Utils.copy_variable(reference_handler, handler, 'time')
Utils.copy_variable(reference_handler, handler, 'time_bnds')
Utils.copy_variable(reference_handler, handler, 'leadtime')
reference_handler.close()
Javier Vegas-Regidor
committed
new_var = handler.createVariable(cmor_name, float, 'time', fill_value=1.0e20)
new_var.units = units
new_var.short_name = cmor_name
Javier Vegas-Regidor
committed
new_var[:] = values
handler.close()
Javier Vegas-Regidor
committed
return temp