Skip to content
Snippets Groups Projects
Commit 6a1513f8 authored by Peters, Wouter's avatar Peters, Wouter
Browse files

TransCom files are now created based on a 1x1 to TC conversion, note that... 

TransCom files are now created based on a 1x1 to TC conversion, note that there is no covaraince yet
parent 691ac1b0
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
da/analysis/expand_fluxes.py
+ 81
− 64
View file @ 6a1513f8
...@@ -6,7 +6,8 @@ import getopt ...@@ -6,7 +6,8 @@ import getopt
from pylab import array, arange, transpose, date2num from pylab import array, arange, transpose, date2num
from runinfo import * from runinfo import *
from datetime import datetime from datetime import datetime
from ct_netcdf import CT_CDF, std_savedict, GetVariable, CreateDirs from da.tools.io4 import CT_CDF, std_savedict, GetVariable
from da.tools.general import CreateDirs
""" """
Author: Wouter Peters (Wouter.Peters@noaa.gov) Author: Wouter Peters (Wouter.Peters@noaa.gov)
...@@ -43,7 +44,7 @@ def SaveWeeklyAvg1x1Data(rundat): ...@@ -43,7 +44,7 @@ def SaveWeeklyAvg1x1Data(rundat):
# #
dirname = 'data_flux1x1_weekly' dirname = 'data_flux1x1_weekly'
# #
dirname = CreateDirs(rundat,dirname) dirname = CreateDirs(os.path.join(rundat.outputdir,dirname))
# #
# Some help variables # Some help variables
# #
...@@ -172,16 +173,15 @@ def SaveWeeklyAvgEcoData(rundat): ...@@ -172,16 +173,15 @@ def SaveWeeklyAvgEcoData(rundat):
dirname='data_eco_weekly' dirname='data_eco_weekly'
dirname=CreateDirs(rundat,dirname) dirname = CreateDirs(os.path.join(rundat.outputdir,dirname))
# Create NetCDF output file # Create NetCDF output file
# #
saveas=os.path.join(rundat.outputdir,dirname,'ecofluxes.nc') saveas = os.path.join(dirname,'ecofluxes.nc')
ncf=CT_CDF(saveas,'write') ncf = CT_CDF(saveas,'write')
dimdate=ncf.AddDateDim() dimdate = ncf.AddDateDim()
dimidateformat=ncf.AddDateDimFormat() dimidateformat = ncf.AddDateDimFormat()
dimregs=ncf.AddRegionDim(rundat,type='eco') dimregs = ncf.AddRegionDim(type='eco')
dummy=ncf.AddRegionDefs(rundat,type='olson')
# #
dectime0=date2num(datetime(2000,1,1)) dectime0=date2num(datetime(2000,1,1))
...@@ -334,23 +334,22 @@ def SaveWeeklyAvgTCData(rundat): ...@@ -334,23 +334,22 @@ def SaveWeeklyAvgTCData(rundat):
*** Inputs *** *** Inputs ***
rundat : a RunInfo object rundat : a RunInfo object
""" """
from da.analysis.tools_region import globarea from da.analysis.tools_regions import globarea, StateToGrid
from da.analysis.tools_transcom import StateToTranscom, StateCovToTranscom from da.analysis.tools_transcom import StateToTranscom, StateCovToTranscom, transcommask
import da.tools.io4 as io
from numpy import dot from numpy import dot
import pycdf as CDF
dirname='data_tc_weekly' dirname='data_tc_weekly'
dirname=CreateDirs(rundat,dirname) dirname = CreateDirs(os.path.join(rundat.outputdir,dirname))
# Create NetCDF output file # Create NetCDF output file
# #
saveas=os.path.join(rundat.outputdir,dirname,'tcfluxes.nc') saveas = os.path.join(dirname,'tcfluxes.nc')
ncf=CT_CDF(saveas,'create') ncf = io.CT_CDF(saveas,'create')
dimdate=ncf.AddDateDim() dimdate = ncf.AddDateDim()
dimidateformat=ncf.AddDateDimFormat() dimidateformat = ncf.AddDateDimFormat()
dimregs=ncf.AddRegionDim(rundat,type='tc') dimregs = ncf.AddRegionDim(type='tc')
dummy=ncf.AddRegionDefs(rundat,type='tc')
# #
dectime0=date2num(datetime(2000,1,1)) dectime0=date2num(datetime(2000,1,1))
...@@ -358,33 +357,45 @@ def SaveWeeklyAvgTCData(rundat): ...@@ -358,33 +357,45 @@ def SaveWeeklyAvgTCData(rundat):
area=globarea() area=globarea()
dt=rundat.dt dt=rundat.dt
infile=os.path.join(rundat.outputdir,'data_eco_weekly','ecofluxes.nc') infile=os.path.join(rundat.outputdir,'data_flux1x1_weekly','flux_1x1.nc')
if not os.path.exists(infile): if not os.path.exists(infile):
print "Needed input file (%s) does not exist yet, please create weekly eco flux files first, returning..."%infile print "Needed input file (%s) does not exist yet, please create weekly eco flux files first, returning..."%infile
return None return None
ncf_in=CDF.CDF(infile) ncf_in = io.CT_CDF(infile)
vardict=ncf_in.variables() vardict = ncf_in.variables
for vname, vprop in vardict.iteritems(): for vname, vprop in vardict.iteritems():
data=array(ncf_in.var(vname)[:])
atts=ncf_in.var(vname).attributes() data = ncf_in.GetVariable(vname)[:]
if vname=='latitude': continue
elif vname=='longitude': continue
elif vname=='date': dims=dimdate
elif vname=='idate': dims=dimdate+dimidateformat
if vname not in ['date','idate']: if vname not in ['date','idate']:
if 'cov' in vname: if 'cov' in vname:
data=Pto_tc(data.shape[-1],data) alldata=[]
for dd in data:
dd=StateCovToGrid(dd*area,transcommask,reverse=True)
alldata.append(dd)
data=array(alldata)
dims=dimdate+dimregs+dimregs dims=dimdate+dimregs+dimregs
else: else:
data=xto_tc(data.shape[-1],data) alldata=[]
for dd in data:
dd=StateToGrid(dd*area,transcommask,reverse=True)
alldata.append(dd)
data=array(alldata)
dims=dimdate+dimregs dims=dimdate+dimregs
elif vname=='date': dims=dimdate
elif vname=='idate': dims=dimdate+dimidateformat
else:
print 'Dataset with unknown dimensions encountered in file: %s'%vname
savedict=ncf.StandardVar(varname=vname) print vname,data.shape
savedict['values']=data.tolist() savedict = ncf.StandardVar(varname=vname)
savedict['dims']=dims savedict['values'] = data.tolist()
savedict['units']=atts['units'] savedict['dims'] = dims
savedict['count']=0 savedict['units'] = 'mol/region/s'
savedict['count'] = 0
ncf.AddData(savedict,nsets=data.shape[0]) ncf.AddData(savedict,nsets=data.shape[0])
ncf.close() ncf.close()
...@@ -402,8 +413,8 @@ def SaveTCDataExt(rundat): ...@@ -402,8 +413,8 @@ def SaveTCDataExt(rundat):
*** Example *** *** Example ***
./expand_savestate project=enkf_release sd=20000101 ed=20010101 """ ./expand_savestate project=enkf_release sd=20000101 ed=20010101 """
from da.analysis.tools_transcom import ExtendedRCRegions from da.analysis.tools_transcom import *
import pycdf as CDF import da.tools.io4 as io
infile=os.path.join(rundat.outputdir,'data_tc_weekly','tcfluxes.nc') infile=os.path.join(rundat.outputdir,'data_tc_weekly','tcfluxes.nc')
if not os.path.exists(infile): if not os.path.exists(infile):
...@@ -412,36 +423,43 @@ def SaveTCDataExt(rundat): ...@@ -412,36 +423,43 @@ def SaveTCDataExt(rundat):
# Create NetCDF output file # Create NetCDF output file
# #
saveas=os.path.join(rundat.outputdir,'data_tc_weekly','tc_extfluxes.nc') saveas = os.path.join(rundat.outputdir,'data_tc_weekly','tc_extfluxes.nc')
ncf=CT_CDF(saveas,'create') ncf = io.CT_CDF(saveas,'create')
dimdate=ncf.AddDateDim() dimdate = ncf.AddDateDim()
dimidateformat=ncf.AddDateDimFormat() dimidateformat = ncf.AddDateDimFormat()
dimregs=ncf.AddRegionDim(rundat,type='tc_ext') dimregs = ncf.AddRegionDim(type='tc_ext')
dummy=ncf.AddRegionDefs(rundat,type='tc_ext')
ncf_in=CDF.CDF(infile) ncf_in=io.CT_CDF(infile,'read')
vardict=ncf_in.variables() vardict = ncf_in.variables
for vname, vprop in vardict.iteritems(): for vname, vprop in vardict.iteritems():
dims=()
data=array(ncf_in.var(vname)[:]) data = ncf_in.GetVariable(vname)[:]
atts=ncf_in.var(vname).attributes()
if vname=='latitude': continue
elif vname=='longitude': continue
elif vname=='date': dims=dimdate
elif vname=='idate': dims=dimdate+dimidateformat
if vname not in ['date','idate']: if vname not in ['date','idate']:
if 'cov' in vname: if 'cov' in vname:
data=xtc(data[:,0:23,0:23],cov=True) alldata=[]
for dd in data:
dd=StateCovToGrid(dd*area,transcommask,reverse=True)
alldata.append(dd)
data=array(alldata)
dims=dimdate+dimregs+dimregs dims=dimdate+dimregs+dimregs
else: else:
data=xtc(data[:,0:23]) data=ExtendedTCRegions(data)
dims=dimdate+dimregs dims=dimdate+dimregs
elif vname=='date': dims=dimdate
elif vname=='idate': dims=dimdate+dimidateformat
else:
print 'Dataset with unknown dimensions encountered in file: %s'%vname
savedict=ncf.StandardVar(varname=vname) print vname,data.shape
savedict['units']=atts['units'] savedict = ncf.StandardVar(varname=vname)
savedict['values']=data.tolist() savedict['values'] = data.tolist()
savedict['dims']=dims savedict['dims'] = dims
savedict['units'] = 'mol/region/s'
savedict['count'] = 0
ncf.AddData(savedict,nsets=data.shape[0]) ncf.AddData(savedict,nsets=data.shape[0])
ncf.close() ncf.close()
...@@ -474,8 +492,7 @@ def SaveEcoDataExt(rundat): ...@@ -474,8 +492,7 @@ def SaveEcoDataExt(rundat):
ncf=CT_CDF(saveas,'create') ncf=CT_CDF(saveas,'create')
dimdate=ncf.AddDateDim() dimdate=ncf.AddDateDim()
dimidateformat=ncf.AddDateDimFormat() dimidateformat=ncf.AddDateDimFormat()
dimregs=ncf.AddRegionDim(rundat,type='eco_ext') dimregs = ncf.AddRegionDim(type='tc_ext')
dummy=ncf.AddRegionDefs(rundat,type='eco_ext')
ncf_in=CDF.CDF(infile) ncf_in=CDF.CDF(infile)
...@@ -522,7 +539,7 @@ def SaveTimeAvgData(rundat,infile,avg='monthly'): ...@@ -522,7 +539,7 @@ def SaveTimeAvgData(rundat,infile,avg='monthly'):
import da.tools.io4 as io import da.tools.io4 as io
dirname,filename = os.path.split(infile) dirname,filename = os.path.split(infile)
dirname = CreateDirs(rundat,dirname.replace('weekly',avg) ) dirname = CreateDirs(os.path.join(rundat.outputdir,dirname.replace('weekly',avg) ) )
dectime0 = date2num(datetime(2000,1,1)) dectime0 = date2num(datetime(2000,1,1))
...@@ -582,7 +599,7 @@ def SaveTimeAvgData(rundat,infile,avg='monthly'): ...@@ -582,7 +599,7 @@ def SaveTimeAvgData(rundat,infile,avg='monthly'):
dims += (dimgrid[1],) dims += (dimgrid[1],)
for type in ['eco','eco_ext','tc','tc_ext','olson']: for type in ['eco','eco_ext','tc','tc_ext','olson']:
if 'regions_%s' % type == d: if 'regions_%s' % type == d:
dimregs = ncf.AddRegionDim(rundat,type) dimregs = ncf.AddRegionDim(type)
dims += dimregs dims += dimregs
# #
# If the variable name from the infile is not defined in the standard, simply copy attributes from the old to the new sds # If the variable name from the infile is not defined in the standard, simply copy attributes from the old to the new sds
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment