source: TI02-CSML/trunk/parser/API/ops_GridSeriesFeature.py @ 1086

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/TI02-CSML/trunk/parser/API/ops_GridSeriesFeature.py@1086
Revision 1086, 6.2 KB checked in by domlowe, 14 years ago (diff)

time and calendar working in netcdf output. Time subsetting by min/max or explicit values. Ready for attempt at data extractor integration.

Line 
1''' ops_GridSeriesFeature  contains operations for GridSeriesFeatures'''
2from API import *
3from CSMLDocument import *
4from Scientific.IO.NetCDF import *  #use this instead of cdms for now for it's simple write interface..
5from Numeric import *
6
7
8def testmethod(self):
9    print 'testmethod for gridseries feature'
10    return 'testmethod - gridseries'
11
12def getAllowedSubsettings(self):
13    return ['subsetToGridSeries']
14
15def getDomainReference(self):
16    #Inspects a time position list for the domain reference.
17    #TODO, does not handle a file extract in place of a list.
18    if isinstance(self.domain.domainReference,Parser.TimePositionList):
19        time = {}
20        time['t'] = self.domain.domainReference.timePositions
21        domainref  = time
22    return domainref
23   
24def getDomainComplement(self):
25    #This will return a list containing one or more ordinates:
26    #currently in form [Name, values]
27    #assumes ordinate.axisValues is a file extract
28    #TODO axisValues may be inline -see also domainReference for similar problem
29    domaincomp ={}
30    dc = self.domain.domainComplement
31    #dc should be a grid!
32    if isinstance(dc, Parser.Grid):
33        for ordinate in dc.ordinates:
34            domaincomp[ordinate.definesAxis]=ordinate.axisValues.getData()
35    return domaincomp
36
37def getDomain(self):
38    #returns both the domain reference axes and domain compliment axes in a single domain dictionary
39    #axes are in no particular order
40    domain = {}
41    dr=getDomainReference(self)
42    dc=getDomainComplement(self)
43    for key in dc.keys():
44        domain[key]=dc[key]
45    for key in dr.keys():
46        domain[key]=dr[key]
47    return domain
48
49
50       
51
52def subsetToGridSeries(self, timeSubset,  **kwargs):
53    pathToSubsetCSML = 'temp.xml'
54    pathToSubsetNetCDF='temp.nc'
55    domainref = getDomainReference(self) 
56    self.times=timeSubset
57    self.files=[]
58    strTimes=''
59    fulldata=[]
60    if len(self.times) == 2:
61        tone=ops_AbstractFeature.__getCDtime(self.times[0])
62        ttwo=ops_AbstractFeature.__getCDtime(self.times[1])
63        dr=getDomainReference(self)
64        self.times=[]
65        for time in dr['t'].split():
66            timeok=ops_AbstractFeature.__compareTimes(tone,time,ttwo)
67            if timeok ==1:
68                self.times.append(time)
69    for time in self.times:
70        listPosition=domainref['t'].split().index(time)
71        strTimes= strTimes + ' ' + time
72        for comp in self.rangeSet.aggregatedArray.component:
73            data=comp.getData(fileposition=listPosition, **kwargs)
74            self.files.append(comp.fileName.split()[listPosition])
75            if fulldata ==[]:
76                fulldata = data.tolist()
77            else:
78                for item in data.tolist():
79                    fulldata.append(item)
80     #get the calendar type
81    try:
82        caltype, calunits = ops_AbstractFileExtract.__calendar(self.rangeSet.aggregatedArray.component[0].fileName.split()[0], 't')    #TODO should accept any time dim!!
83        csmltime.setcdtimeCalendar(caltype)
84    except:
85        caltype=cdtime.DefaultCalendar   
86    ### define domain and rangeSet to use for feature in csml document####
87    domain=Parser.GridSeriesDomain()
88    domain.domainReference=Parser.TimePositionList(timePositions=strTimes) 
89    grid=Parser.Grid()
90    dc = self.getDomainComplement()
91    ordinates= []
92    i=0
93    valueStore=[]  # use the values again later to generate netcdf
94    arraySize=0
95    totalArraySize=1
96    for key in dc.keys():
97        arraySize=0
98        i=i+1
99        god=Parser.GridOrdinateDescription()
100        god.gridAxesSpanned='dim%s'%i
101        god.sequenceRule='+x+y+z'
102        god.definesAxis=key
103        straxisValues=''
104        if key in kwargs:
105            for val in dc[key]:
106                if val >= kwargs[key][0]:
107                    if val <= kwargs[key] [1]:
108                        arraySize=arraySize+1
109                        straxisValues=straxisValues+ str(val) + ', '
110        else: # this dimension has not been subsetted
111            for val in dc[key]:
112                arraySize=arraySize+1
113                straxisValues=straxisValues+ str(val) + ', '
114        totalArraySize=totalArraySize*arraySize
115        god.axisValues=straxisValues[:-2]
116        ordinates.append(god)
117    totalArraySize=totalArraySize*len(self.times)
118    grid.ordinates=ordinates
119    domain.domainComplement=grid
120    rangeSet=Parser.RangeSet()
121    rangeSet.arrayDescriptor=Parser.NetCDFExtract(id=self.id,fileName='temp.nc',variableName=self.id,arraySize=[arraySize])
122
123    #### write csml document ##### -move this to the csmldocument module?
124    subsetCSML=CSMLDocument()
125    subsetCSML=subsetCSML.makeGridSeries(domain,rangeSet)
126    output=open(pathToSubsetCSML,'w')
127    output.write(subsetCSML)
128    output.close()
129
130    ### create and write netcdf - uses scientific python####
131    ncfile=NetCDFFile(pathToSubsetNetCDF,'w')
132    # create the dimensions       
133    ncfile.createDimension ( 'time', len(self.times))
134    time_var = ncfile.createVariable ( 'time', Float, ('time',) )
135    time_var.longname = 'time'
136    floatTimes=[]
137    print len(fulldata[0])
138    for time in self.times:
139        time=ops_AbstractFeature.__getCDtime(time).torel(calunits)
140        floatTimes.append(time.value)
141    time_var[:] =floatTimes[:]
142    time_var.units=calunits
143    time_var.calendar=caltype
144
145    for ordinate in ordinates:
146        ncfile.createDimension(ordinate.definesAxis, len(ordinate.axisValues.split())) 
147        item_var = ncfile.createVariable (ordinate.definesAxis, Float, (ordinate.definesAxis,) )
148        #convert to list
149        vals=[]
150        for val in ordinate.axisValues.split(','):
151            vals.append(float(val))
152        ordinate.axisValues=vals
153        item_var[:]=vals[:]
154        print ordinate.definesAxis
155    #this needs reconsidering - do the shapes always match up??
156    if len(ordinates)==3:
157        feature_var = ncfile.createVariable (self.id, Float, ('time',ordinates[0].definesAxis,ordinates[1].definesAxis,ordinates[2].definesAxis))
158    elif len(grid.ordinates)==2:
159        feature_var = ncfile.createVariable (self.id, Float, ('time',ordinates[0].definesAxis,ordinates[1].definesAxis))
160    print shape(feature_var)
161    print shape(fulldata)
162    feature_var[:]=fulldata[:]
163    ncfile.close()
164
165    return pathToSubsetCSML, pathToSubsetNetCDF, totalArraySize
Note: See TracBrowser for help on using the repository browser.