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

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

handling standard name in phenomenon

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