Commit 9731a821 authored by Clément Haëck's avatar Clément Haëck
Browse files

Minor refactorization

parent 0df763c4
......@@ -23,8 +23,8 @@ from scipy import stats
import lib
import lib.data.hi
import lib.data.hists
import lib.data.SN_separation
import lib.zones
from lib import m_start, m_end, m_progress
REDUCE = False
......@@ -33,24 +33,32 @@ HI_TARGET = 9.5
def main(args):
m_start("Getting data")
ds = lib.data.hi.get_data(args)
m_end()
if REDUCE:
print("DATASET REDUCED FOR TESTING")
ds = ds.isel(time=slice(None, None, 50))
m_start("Getting static masks")
# Remove land and select standard study region.
args['grid'] = '4km_EPSG32662'
zone = lib.zones.get_data(args)['total']
land = lib.zones.get_land(args)['land_large']
ds = ds.where(zone*~land)
# Compute V, S, B coefs
coef = coef_components(ds)
m_end()
# If coef argument is 0, recompute normalisation coefs for V, S, B.
if args['coef'] == 0:
coef = get_components_coef(ds)
else:
# Else just use those already computed
coef = lib.data.hi.get_coef(args)
print('components:', coef)
# Compute HI coef
coef_HI = coef_hi(ds, args, coef)
coef_HI = coef_hi(ds, coef)
print('HI:', coef_HI)
coef.update(coef_HI)
......@@ -68,20 +76,25 @@ def main(args):
return s
def coef_components(ds):
def get_components_coef(ds):
m_start("Computing SVB coef")
coef = {}
for variable in 'SVB':
m_progress()
data = ds[variable]
if variable == 'S':
data = np.fabs(data)
std = data.std().values
coef[variable] = 1./std
m_end()
return coef
def coef_hi(ds, args, coef):
def coef_hi(ds, coef):
m_start("Getting HI coef")
coef['HI'] = 1.
ds['HI'] = lib.data.hi.apply_coef(ds, coef)
......@@ -92,6 +105,7 @@ def coef_hi(ds, args, coef):
hi_x = h.ppf(Q_TARGET)
coef_hi = HI_TARGET/hi_x
m_end()
return {'HI': coef_hi}
......
......@@ -22,14 +22,14 @@ import xarray as xr
import lib
import lib.data.globcolour
import lib.data.hi
import lib.data.hists
import lib.data.hists as lh
import lib.data.ostia
import lib.data.p_frt_mask
import lib.data.SN_separation
import lib.zones
from lib import m_start, m_end, m_next
from lib import m_start, m_next, m_end, m_progress
# Width of SST bins
# Width of SST discretization step
SST_STEP = 1e-2
# Number of Chl bins
CHL_NBINS = 500
......@@ -90,8 +90,7 @@ def main(args):
up[mask] = up[mask] * op(ds.sst, ds.threshold)
zones['GS3_'+zone] = up
m_end()
print("Setting up histogram computations", end='', flush=True)
m_next("Setting up histogram computations")
hists_var = []
for var in VARS:
hists = []
......@@ -99,7 +98,7 @@ def main(args):
bins_name = 'bins_' + var
for zone, zone_ds in zones.items():
for m, mask in zip(MASKS, MASKS_VAR):
print('.', end='', flush=True)
m_progress()
h = zone_ds[var].where(zone_ds[mask]).groupby('time').map(
hist_grp, shortcut=True, args=[bins, bins_name]
)
......@@ -108,7 +107,7 @@ def main(args):
hist = xr.combine_by_coords(hists)
# Get a DataArray with the correct name
hist = hist[var].rename('hist_' + var)
hist = hist[var].rename(lh.var_name('hist', var))
hist = hist.assign_coords({bins_name: bins.edges[:-1]})
hist[bins_name].attrs['right_edge'] = VAR_RANGE[var][1]
hist[bins_name].attrs['VAR'] = var
......@@ -124,8 +123,9 @@ def main(args):
# We have at most every pixel of an image in one bin,
# so approx 1000*1000 = 1e6 values. Uint32 stores up to 2**32~4e9.
encoding = {v: {'dtype': 'uint32', '_FillValue': 2**30-1}
for v in ['hist_' + v for v in VARS]}
encoding = {lh.var_name('hist', v): {'dtype': 'uint32',
'_FillValue': 2**30-1}
for v in VARS}
m_next("Executing computations / Writing to disk")
ofile = lib.data.hists.get_filename(args)
......
......@@ -22,7 +22,7 @@ import xarray as xr
import lib
import lib.data.globcolour
import lib.data.hi
import lib.data.hists
import lib.data.hists as lh
import lib.data.ostia
import lib.data.p_frt_mask
import lib.data.SN_separation
......@@ -96,7 +96,7 @@ def main(args):
for var in VARS:
hists = []
bins = get_bins(var)
bins_name = 'bins_' + var
bins_name = lh.var_name('bins', var)
for zone, zone_ds in zones.items():
for m, mask in zip(MASKS, MASKS_VAR):
print('.', end='', flush=True)
......@@ -108,7 +108,7 @@ def main(args):
hist = xr.combine_by_coords(hists)
# Get a DataArray with the correct name
hist = hist[var].rename('hist_' + var)
hist = hist[var].rename(lh.var_name('hist', var))
hist = hist.assign_coords({bins_name: bins.edges[:-1]})
hist[bins_name].attrs['right_edge'] = VAR_RANGE[var][1]
hist[bins_name].attrs['VAR'] = var
......
......@@ -57,7 +57,7 @@ lib.data.name.get_filename(Y=2007, m=3, d=8)
#+end_src
Any matcher that was fixed in the ~fixes~ argument does not need to be specified again here.
*** ~get_data~
This returns a dataset created with [[https://xarray.pydata.org/en/stable/generated/xarray.open_mfdataset.html][~xarray.open_mfdataset()~]]. The files are supplied by the Finder object. Other arguments to the function are supplied by the ~lib.data.OPEN_MF_KW_DEF~ dictionnary, completed by the argument ~open_mf_kw~ of the ~lib.data.create_dataset~ function (see below).
This returns a dataset created with [[https://xarray.pydata.org/en/stable/generated/xarray.open_mfdataset.html][~xarray.open_mfdataset()~]]. The files are supplied by the Finder object. Other arguments to the function are supplied by the ~lib.data.OPEN_MF_KW_DEF~ dictionnary, completed by the argument ~open_mf_kw~ of the ~lib.data.create_dataset~ and ~lib.data.<dataset>.get_data~ functions (see below).
** ~lib.data.create_dataset()~
This function automatically populates the module with the functions above. It also returns a dictionnary containing those functions; this is useful if one desires to have a more specific implementation of one of the functions.
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment