source: TI02-CSML/trunk/csml/API/ops_GridSeriesFeature.py @ 3704

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

Added Ag's file cache code to DataInterface?. The Data Interface is read only

Line 
1''' ops_GridSeriesFeature  contains operations for GridSeriesFeatures. These methods are attached to GridFeature instances at runtime.'''
2
3import csml.parser
4import csml.csmllibs.csmltime
5import csml.csmllibs.csmldocument
6import csml.API.ops_AbstractFeature
7import csml.API.genSubset
8import csml.csmllibs.netCDFWriter
9import csml.csmllibs.csmlcrs
10import csmlutils
11import os
12
13def testmethod(self):
14    ''' testmethod for gridseries feature'''
15    return 'testmethod - gridseries'
16
17def getAllowedSubsettings(self):
18    '''get the allowed operations
19    @return:    list of operations'''
20    return ['subsetToGridSeries', 'subsetToProfileSeries']  #other operations
21
22
23def getBoundingBox(self):
24    ''' get the simple boundingBox '''
25    try:
26        lower1=eval(bbox.boundedBy.lowerCorner.CONTENT.split()[0])
27        lower2=eval(bbox.boundedBy.lowerCorner.CONTENT.split()[1])
28        upper1=eval(bbox.boundedBy.upperCorner.CONTENT.split()[0])
29        upper2=eval(bbox.boundedBy.upperCorner.CONTENT.split()[1])
30        return [lower1,lower2,upper1,upper2]
31    except:
32        return None
33   
34   
35def getDomain(self):
36    ''' get domain of the feature
37    @return:   domain as a dictionary of ordinates {name: [values], ...}'''
38    self.domain={}
39    self.gridnames={}
40    for gridOrd in self.value.gridSeriesDomain.coordTransformTable.gridOrdinates:
41        name=gridOrd.coordAxisLabel.CONTENT
42        self.gridnames[gridOrd.gridAxesSpanned.CONTENT]=name
43        if hasattr(gridOrd.coordAxisValues, 'insertedExtract'):
44            self.domain[name], fill, axisorder, units=gridOrd.coordAxisValues.insertedExtract.getData()
45            #add to cache
46        else:           
47            valList=[]
48            try:
49                vals=gridOrd.coordAxisValues.coordinateList.CONTENT
50                for val in vals.split(): 
51                    valList.append(eval(val))
52            except:
53                vals=gridOrd.coordAxisValues.timePositionList.CONTENT
54                for val in vals.split(): 
55                    valList.append(val)           
56            self.domain[name]=valList   
57    return self.domain
58
59def getNativeCRS(self):
60    ''' get the the coordinate reference system of the domain of the feature
61    return:   a string, the srsName'''
62    return self.value.gridSeriesDomain.srsName
63
64def getUom(self):
65    ''' get the uom of the phenomenon.
66    @return: string - uom attribute'''
67    uom=None
68    try: 
69        uom=self.value.rangeSet.valueArray.valueComponent.quantityList.insertedExtract.uom.CONTENT
70    except AttributeError:
71        uom =None
72    return uom
73   
74def getDomainUnits(self):
75    ''' get the units used in the domain
76    @return:     list of units'''   
77    #srsname=self.value.gridSeriesDomain.srsName
78    #cat=csml.csmllibs.csmlcrs.CRSCatalogue()
79    #crs=cat.getCRS(srsname)
80    #return crs.units
81    return self.value.gridSeriesDomain.uomLabels.split()
82   
83def __getAxis(self, name):
84    '''called by getLongitudeAxis, getTimeAxis and getLatitudeAxis'''   
85    srsname=self.value.gridSeriesDomain.srsName
86    cat=csml.csmllibs.csmlcrs.CRSCatalogue()
87    crs=cat.getCRS(srsname) 
88    axID=None
89    if name == 'lon':
90        axID = crs.lonAxis
91    elif name == 'lat':
92        axID = crs.latAxis
93    elif name == 'time':
94        axID = crs.timeAxis     
95    if axID is not None:
96        return crs.getAxisLabels()[axID]
97    else:
98        return None
99   
100def getLongitudeAxis(self):
101    ''' get the index of the longitude axis'''
102    return self.__getAxis('lon')
103
104def getLatitudeAxis(self):
105    ''' get the index of the latitude axis'''
106    return self.__getAxis('lat')
107
108def getTimeAxis(self):
109    ''' get the index of the time axis'''
110    return self.__getAxis('time')
111
112def getAxisLabels(self):
113    '''get the names of the axes'''
114    labels=[]
115    for label in self.value.gridSeriesDomain.axisLabels.split():
116        labels.append(label)
117    return labels
118
119def _subsetGrid(self, **kwargs):
120    '''this takes a selection from a gridseries object (be it a profile, a grid, whatever, it is still just a selection).  It is called internally  by the subsetToGridSeries method.
121    @return: a load of things which can then be used to rebuild a new csml feature: pathToSubsetNetCDF, crs, frame, timeName, times, strTimes, caltype, axisorder,units, fulldata, fillvalue, kwargs '''
122           
123    #set self.domain:
124    self.getDomain()     
125   
126    #get the CRS from a  the  catalogue
127    cat=csml.csmllibs.csmlcrs.CRSCatalogue()
128    crs=cat.getCRS(self.value.gridSeriesDomain.srsName) 
129   
130    #non-feature specific setup code, mainly handles the time dimension/calendar
131   
132    pathToSubsetNetCDF, kwargs, timeAxis, timeName, caltype, times=csml.API.genSubset.genericSubset(self, self.outputdir, self.ncname, self.domain, kwargs)
133    ##Get names of variables in file and relate them to the subset selection
134    selection={}
135    frame=''
136   
137    for gridOrd in self.value.gridSeriesDomain.coordTransformTable.gridOrdinates:
138        try:
139            selection[gridOrd.gridAxesSpanned.CONTENT]=kwargs[gridOrd.coordAxisLabel.CONTENT]           
140        except KeyError:
141            allValues=tuple(self.domain[gridOrd.coordAxisLabel.CONTENT])
142        if hasattr(gridOrd.coordAxisValues, 'timePositionList'):
143            if hasattr(gridOrd.coordAxisValues, 'frame'):
144                frame=gridOrd.coordAxisValues.frame
145    try:
146        del selection[timeAxis] #NOTE: Haven't resolved the single CRS issue so cannot subset by time within an individual file for now
147    except KeyError:
148        pass
149   
150    #get slice indices
151    selection=self.getSliceIndices(selection)                 #commented out while fixing subsetting
152
153    strTimes, axisorder, units, fulldata, fillvalue =csml.API.genSubset.getTheData(self, selection, times, timeName, timeAxis)
154    return pathToSubsetNetCDF, crs, frame, timeName, times, strTimes, caltype, axisorder,units, fulldata, fillvalue, kwargs
155   
156def subsetToGridSeries(self, outputdir=None, ncname='gridseries.nc',**kwargs):
157    '''Subset a GridSeriesFeature to another smaller  GridSeriesFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs
158    @param outputdir:    name of directory to write data files to
159    @param ncname:     name of netcdf file to write out
160    @param kwargs:      subset selection by axis name.
161    @return:     subsettedFeature (GridSeriesFeature instance)  pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)
162   
163    '''
164    if outputdir is not None:
165        self.outputdir=outputdir
166    elif csml.API.csmlContainer.globalContainer.OUTPUTDIR is not None:
167        self.outputdir=csml.API.csmlContainer.globalContainer.OUTPUTDIR
168    csml.csmllibs.csmlextra.checkDirExists(self.outputdir)
169    self.ncname=ncname
170    pathToSubsetNetCDF, crs, frame, timeName, times, strTimes, caltype, axisorder,units, fulldata, fillvalue, kwargs=self._subsetGrid(**kwargs) 
171    #Okay, got the data now. Need to write CSML feature and NetCDF files.
172    #Writing out the CSML feature   
173   
174    # define domain/coverage  to use in 'value' attribute   
175    newdomain=csml.parser.GridSeriesDomain()
176    domainSubset, totalArraySize=csml.API.genSubset.subsetDomain(timeName,strTimes,self.domain, **kwargs)
177    cTT=csml.API.genSubset.getCoordTransformTable(domainSubset, crs, frame)
178    newdomain.id=csml.csmllibs.csmlextra.getRandomID()
179    newdomain.coordTransformTable=cTT
180    newdomain.srsName=self.value.gridSeriesDomain.srsName 
181    newdomain.axisLabels=self.value.gridSeriesDomain.axisLabels
182    newdomain.srsDimension=self.value.gridSeriesDomain.srsDimension
183    newdomain.dimension=self.value.gridSeriesDomain.dimension
184    env=csml.parser.GridEnvelope()
185    env.low=csml.parser.csString('0 0 0') #TODO
186    env.high=csml.parser.csString('0 0 0')
187    newdomain.limits=env
188    newdomain.aLabels=self.value.gridSeriesDomain.aLabels
189    rs=csml.parser.RangeSet()
190    sdid=csml.csmllibs.csmlextra.getRandomID()
191    descriptor=csml.parser.NetCDFExtract(id=sdid,fileName=csml.parser.csString(self.ncname),variableName=self.name,arraySize=csml.parser.csString(totalArraySize))
192    rs=csml.parser.RangeSet()
193    va=csml.parser.ValueArray()
194    vc=csml.parser.MeasureOrNullList()
195    vc.href='#%s'%sdid
196    vc.arcrole="http://ndg.nerc.ac.uk/xlinkUsage/insert#QuantityList"
197    vc.role="http://ndg.nerc.ac.uk/fileFormat/csmlStorageDescriptor"
198    vc.show='embed'
199    try:
200        vc.uom=self.value.rangeSet.valueArray.valueComponent.uom
201    except:
202        vc.uom = 'unknown'  #TODO,
203    va.valueComponent=vc
204    va.id=csml.csmllibs.csmlextra.getRandomID()
205    rs.valueArray=va
206    #gridseries coverage
207    cvg=csml.parser.GridSeriesCoverage()
208    cvg.id=csml.csmllibs.csmlextra.getRandomID()
209    cvg.rangeSet=rs
210    cvg.gridSeriesDomain=newdomain   
211   
212    #parameter, as before subsetting.
213    param = self.parameter
214       
215    #create a stand alone gridseries feature containing this coverage
216    csmlWrap=csml.csmllibs.csmlfeaturewrap.CSMLWrapper()
217    subsettedFeature=csmlWrap.createGridSeriesFeature(value=cvg,parameter=param,featureID=csml.csmllibs.csmlextra.getRandomID(),name=self.name,description=self.description)
218 
219    ### write netcdf using NCWriter class (wraps cdms) ###
220    nc=csml.csmllibs.netCDFWriter.NCwriter(pathToSubsetNetCDF)
221    ords=cTT.gridOrdinates
222    axislist=[]
223    for a in axisorder:
224        axislist.append(self.gridnames[a])
225    stdname=param.getStandardName()
226    nc.genWriteVar(self.name.CONTENT,ords, times, caltype, axislist, units, stdname, fulldata, fillvalue)
227    nc.closeFinishedFile()
228    print 'NetCDF file written to %s'%pathToSubsetNetCDF
229    return subsettedFeature, pathToSubsetNetCDF, descriptor
230   
231
232def subsetToProfileSeries(self, outputdir=None, ncname='profileseries.nc',**kwargs):
233    '''Subset a GridSeriesFeature to a ProfileSeriesFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs
234    @param outputdir:    name of directory to write data files to
235    @param ncname:     name of netcdf file to write out
236    @param kwargs:      subset selection by axis name.
237    @return:     subsettedFeature (ProfileSeriesFeature instance)  pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)'''
238   
239   
240    #TODO   !!!!!!!!! Need to perform some sort of testing on the kwargs to check it is a profileseries request.
241    if outputdir is not None:
242        self.outputdir=outputdir
243    elif csml.API.csmlContainer.globalContainer.OUTPUTDIR is not None:
244        self.outputdir=csml.API.csmlContainer.globalContainer.OUTPUTDIR
245    csml.csmllibs.csmlextra.checkDirExists(self.outputdir)
246    self.ncname=ncname   
247   
248    #perform the subset (note this included nearest neighbour searching, so may return a different set of kwargs
249    pathToSubsetNetCDF, crs,frame, timeName, times, strTimes, caltype, axisorder,units, fulldata, fillvalue, kwargs=self._subsetGrid(**kwargs) 
250
251    latName=crs.axes[crs.latAxis]
252    lonName=crs.axes[crs.lonAxis]
253    #Okay, got the data now. Need to write CSML feature and NetCDF files.
254    #Writing out the CSML feature
255   
256    # define domain/coverage  to use in 'value' attribute   
257    newdomain=csml.parser.ProfileSeriesDomain()
258   
259   
260    domainSubset, totalArraySize=csml.API.genSubset.subsetDomain(timeName,strTimes,self.domain, **kwargs)
261    cTT=csml.API.genSubset.getCoordTransformTable(domainSubset, crs, frame) 
262   
263    #Need to remove location (lat/lon) axes from the cTT
264    newGridOrds=[]
265    for gridOrd in cTT.gridOrdinates:
266        if gridOrd.coordAxisLabel.CONTENT not in  [latName,lonName]:
267            newGridOrds.append(gridOrd)
268    cTT.gridOrdinates=newGridOrds
269    axes=[]
270    for ord in newGridOrds:
271        axes.append(ord.coordAxisLabel.CONTENT)
272    cat = csml.csmllibs.csmlcrs.CRSCatalogue()
273    crsInfo=cat.determineCRS(knownCRSAxes=axes)   
274    crs=cat.getCRS(crsInfo[0].srsName)
275   
276    newdomain.coordTransformTable=cTT   
277    newdomain.id=csml.csmllibs.csmlextra.getRandomID()
278    #TODO need to remove 'location' from this srs
279    newdomain.srsName=crsInfo[0].srsName
280    axisLabels=csml.csmllibs.csmlextra.stringify(crsInfo[1].keys())
281    newdomain.axisLabels=axisLabels
282    newdomain.srsDimension=2
283    newdomain.dimension=2
284    env=csml.parser.GridEnvelope()
285    env.low=csml.parser.csString('0 0 0') #TODO
286    env.high=csml.parser.csString('0 0 0')
287    newdomain.limits=env
288    newdomain.aLabels=self.value.gridSeriesDomain.aLabels #todo 
289    rangeSet=csml.parser.RangeSet()
290    descid=csml.csmllibs.csmlextra.getRandomID()
291    descriptor=csml.parser.NetCDFExtract(id=descid,fileName=csml.parser.csString(ncname),variableName=self.name,arraySize=csml.parser.csString(totalArraySize))
292    rs=csml.parser.RangeSet()
293    va=csml.parser.ValueArray()
294    vc=csml.parser.MeasureOrNullList()
295    vc.href='#%s'%descid
296    vc.arcrole="http://ndg.nerc.ac.uk/xlinkUsage/insert#QuantityList"
297    vc.role="http://ndg.nerc.ac.uk/fileFormat/csmlStorageDescriptor"
298    vc.show='embed'
299    vc.insertedExtract=descriptor
300    try:
301        vc.uom=self.value.rangeSet.valueArray.valueComponent.uom
302    except:
303        vc.uom = 'unknown'  #TODO,
304   
305    va.valueComponent=vc
306    va.id=csml.csmllibs.csmlextra.getRandomID()
307    rs.valueArray=va
308    #profileseries coverage
309    cvg=csml.parser.ProfileSeriesCoverage()
310    cvg.id=csml.csmllibs.csmlextra.getRandomID()
311    cvg.rangeSet=rs
312    cvg.profileSeriesDomain=newdomain   
313    csmlWrap=csml.csmllibs.csmlfeaturewrap.CSMLWrapper()
314   
315    #parameter, as before subsetting.
316    param = self.parameter
317   
318    #create 'location' attribute:
319    if (latName, lonName) is not (None, None):
320        locStr ='%s %s'%(kwargs[latName], kwargs[lonName])
321        loc=csml.parser.csString(locStr)
322    else:
323        loc = csml.parser.csString('unknown location')
324       
325   
326    #create a stand alone gridseries feature containing this coverage
327    subsettedFeature=csmlWrap.createProfileSeriesFeature(value=cvg,parameter=param, location=loc, featureID=csml.csmllibs.csmlextra.getRandomID(),name=self.name,description=self.description)
328 
329   
330    ### write netcdf using NCWriter class (wraps cdms) ###
331    nc=csml.csmllibs.netCDFWriter.NCwriter(pathToSubsetNetCDF)
332    ords=cTT.gridOrdinates
333    axislist=[]
334       
335    for a in axisorder:
336        axislist.append(self.gridnames[a])
337    stdname=param.getStandardName()
338    nc.genWriteVar(self.name.CONTENT,ords, times, caltype, axislist, units, stdname, fulldata, fillvalue, latitude=kwargs[latName], longitude=kwargs[lonName])
339    nc.closeFinishedFile()
340   
341    return subsettedFeature, pathToSubsetNetCDF, descriptor
342
343def subsetToProfile(self, outputdir = None, ncname='profile.nc',**kwargs):
344    '''Subset a GridSeriesFeature to a ProfileFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs.
345       
346        Is a two step process - subsets GridSeries to ProfileSeries, then ProfileSeries to Profile
347    @param outputdir:    name of directory to write data files to
348    @param ncname:     name of netcdf file to write out
349    @param kwargs:      subset selection by axis name.
350    @return:     subsettedFeature (ProfileFeature instance)  pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)'''
351
352    odir=outputdir
353    tmpfile= csml.csmllibs.csmlextra.getRandomID()+'.nc'
354    profileSeries, pSfile, descriptor=self.subsetToProfileSeries(odir, tmpfile, **kwargs)   
355    del kwargs['latitude']    #TODO - need to remove excess kwargs based on the domain of the temporary profileSeries feature
356    del kwargs['longitude']
357    subsettedFeature, pathToSubsetNetCDF, descriptor =profileSeries.subsetToProfile(odir,ncname, **kwargs)   
358    tmpfile=odir + '/' + tmpfile
359    os.remove(tmpfile)
360    return subsettedFeature, pathToSubsetNetCDF, descriptor
361
362def subsetToPointSeries(self, outputdir =None, ncname='pointseries.nc',**kwargs):
363    '''Subset a GridSeriesFeature to a PointSeriesFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs.
364       
365        Is a two step process - subsets GridSeries to ProfileSeries, then ProfileSeries to PointSeries
366    @param outputdir:    name of directory to write data files to
367    @param ncname:     name of netcdf file to write out
368    @param kwargs:      subset selection by axis name.
369    @return:     subsettedFeature (PointSeriesFeature instance)  pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)'''
370    tmpfile= csml.csmllibs.csmlextra.getRandomID()+'.nc'
371    profileSeries, pSfile, descriptor =self.subsetToProfileSeries(outputdir, tmpfile, **kwargs)   
372    del kwargs['latitude']    #TODO - need to remove excess kwargs based on the domain of the temporary profileSeries feature
373    del kwargs['longitude']
374    subsettedFeature, pathToSubsetNetCDF, descriptor=profileSeries.subsetToPointSeries(outputdir, ncname, **kwargs)   
375    tmpfile=outputdir + '/' + tmpfile
376    os.remove(tmpfile)
377    return subsettedFeature, pathToSubsetNetCDF, descriptor
Note: See TracBrowser for help on using the repository browser.