esmlab
esmlab copied to clipboard
Function to compute true monthly mean from daily mean.
Thanks for submitting an issue!
Here's a quick checklist in what to include:
-
[x] Include a detailed description of the bug or suggestion
-
[ ]
conda list
of the conda environment you are using -
[x] Minimal, self-contained copy-pastable example that generates the issue if possible. Please be concise with code posted. See guidelines below on how to provide a good bug report:
Bug reports that follow these guidelines are easier to diagnose, and so are often handled much more quickly.
Description
I was trying to fix some bad monthly and daily means from the CESM-CICE history files. @jessluo came up with a script to do a group my year and month. We thought this would be a nice addition to esmlab.
What I Did
The path of the data is:
/gpfs/fs1/p/cesm/pcwg/timeseries-cmip6/b.e21.B1850.f09_g17.CMIP6-piControl.001/ice/proc/tseries/day_1
case = "b.e21.B1850.f09_g17.CMIP6-piControl.001"
year1 = "0061"
year2 = "0099"
year1m = "0001"
year2m = "0099"
filename_month = "../month_1/"+case+".cice.h.aicen."+year1m+"01-"+year2m+"12.nc"
filename_day = case+".cice.h1.aicen_d."+year1+"0101-"+year2+"1231.nc"
filename1 = case+".cice.h1.aicen_d_fix."+year1+"0101-"+year2+"1231.nc"
filename2 = case+".cice.h.aicen_fix."+year1m+"01-"+year2m+"12.nc"
ds = xr.open_dataset(filename1)
ds_day = xr.open_dataset(filename_day)
ds_month = xr.open_dataset(filename_month)
aicen_d = ds['aicen_d']
times = ds_day.time_bounds.values
print(times)
leftbounds_yr = [x[0].timetuple()[0] for x in times]
leftbounds_mo = [x[0].timetuple()[1] for x in times]
yrmo = [cftime.datetime(y, m, 15) for y,m in zip(leftbounds_yr,leftbounds_mo)]
yrmo = np.array(yrmo)
print(yrmo.shape)
aicen_d_val = ds.aicen_d.values
aicen_d = xr.DataArray(aicen_d_val, coords={'time':yrmo, 'NCAT':ds.NCAT, 'TLON':ds.TLON, 'TLAT':ds.TLAT},
dims=('time','nc','nj','ni'))
aicen2 = aicen_d.groupby('time').mean(dim='time')
aicen = ds_month['aicen']
aicen.values[720:][:][:][:] = np.where(aicen2.values[:][:][:][:]>1.,1.0e30,aicen2.values[:][:][:][:])
print(aicen.values[0][4][370][180])
aicen.to_netcdf(filename2)
Paste the command(s) you ran and the output.
If there was a crash, please include the traceback here.