#!/usr/bin/env python

## Just run code in recipes, use it carefully
## example: use netcdf4 to make cache file. Since python is faster than ncl.
#;DIAG_NORCPM; RUNTHESECODES: print('    no codes run here.')

import xarray as xa
import numpy as np
import sys,glob,os

reg = 'GLandEast'
latb  = 66.
late  = 75.
lonb  = -25.
lone  = -12.
var = 'aice'
fns = glob.glob(f'./{var}_??.nc')
fns.sort()


if lonb < 0: lonb += 360
if lone < 0: lone += 360

gridds = xa.open_dataset(fns[0])
wgt = gridds.tarea
lat = gridds.TLAT
lon = gridds.TLON
lon = lon.where(lon>0, lon+360) ## 0-360

wgt = wgt.where(latb<lat ,0)
wgt = wgt.where(late>lat ,0)
wgt = wgt.where(lonb<lon ,0)
wgt = wgt.where(lone>lon ,0)
wgt.name = 'weighting'
wgt.to_netcdf(f'{reg}_weighting.nc')  ## save for plot

tsfns = list()
for fn in fns:
    n = fn[-5:-3]
    ofn = f'{var}_{n}_{reg}_ts.nc'
    tsfns.append(ofn)
    if os.path.exists(ofn): continue
    ds = xa.open_dataset(fn)
    dd = ds[var].weighted(wgt)
    ovar  = dd.mean(('nj','ni'),skipna=True)
    ovar.to_netcdf(ofn)

tsfs = xa.open_mfdataset(tsfns,combine='nested',concat_dim='member')
tsmean = tsfs[var].mean(dim='member')
tsstd  = tsfs[var].std(dim='member')
outds = xa.Dataset({
    'avg': tsmean
    ,'std': tsstd
    })
outds.to_netcdf(f'{var}_{reg}_ens_ts.nc')


