source: TI02-CSML/trunk/csml/API/ops_TrajectoryFeature.py @ 2783

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

adding axis interrogation methods to trajectory feature

Line 
1''' ops_TrajectoryFeature  contains operations for TrajectoryFeature'''
2
3import csml, cdms, MV
4import sys
5import copy
6
7def testmethod(self):
8    print 'testmethod for TrajectoryFeature'
9    return 'testmethod - TrajectoryFeature'
10def getAllowedSubsettings(self):
11    return ['subsetToTrajectory'] 
12
13def __getAxis(self, name):
14    '''called by getLongitudeAxis, getTimeAxis and getLatitudeAxis'''   
15    srsname=self.value.trajectoryDomain.srsName
16    cat=csml.csmllibs.csmlcrs.CRSCatalogue()
17    crs=cat.getCRS(srsname,self.value.trajectoryDomain.axisLabels)     
18    axID=None
19    if name == 'lon':
20        axID = crs.lonAxis
21    elif name == 'lat':
22        axID = crs.latAxis
23    elif name == 'time':
24        axID = crs.timeAxis     
25    if axID is not None:
26        return crs.getAxisLabels()[axID]
27    else:
28        return None
29   
30def getLongitudeAxis(self):
31    return self.__getAxis('lon')
32
33def getLatitudeAxis(self):
34    return self.__getAxis('lat')
35
36def getTimeAxis(self):
37    return self.__getAxis('time')
38
39
40def getDomain(self):
41    #returns domain as a dictionary of ordinates {name: [values], ...}
42    self.domain={}
43    self.gridnames={}
44    for gridOrd in self.value.trajectoryDomain.coordTransformTable.gridOrdinates:
45        name=gridOrd.coordAxisLabel.CONTENT
46        self.gridnames[gridOrd.gridAxesSpanned.CONTENT]=name
47        if hasattr(gridOrd.coordAxisValues, 'insertedExtract'):
48            self.domain[name], fill, axisorder, units=gridOrd.coordAxisValues.insertedExtract.getData()
49            #add to cache
50        else:           
51            valList=[]
52            try:
53                vals=gridOrd.coordAxisValues.coordinateList.CONTENT
54                for val in vals.split(): 
55                    valList.append(eval(val))
56            except:
57                vals=gridOrd.coordAxisValues.timePositionList.CONTENT
58                for val in vals.split(): 
59                    valList.append(val)           
60            self.domain[name]=valList   
61    return self.domain   
62
63def _writeNetCDF(self):
64    #Now open a new NetCDF File:   
65    nc=cdms.open(self.pathToSubsetNetCDF,'w')
66   
67    #create the variable
68    var=MV.array(self.data)
69    var.id=self.name.CONTENT
70    var.name=self.name.CONTENT
71    var.units=self.units[0] # hopefully there is just one unit in this list..
72
73   
74    #create the time axis
75    floatTimes=[]
76    tOne=csml.csmllibs.csmltime.getCDtime(self.selectedTimes[0])
77    tbase=csml.csmllibs.csmltime.getBaseUnits(tOne)
78    for time in self.selectedTimes:
79        time=csml.csmllibs.csmltime.getCDtime(time).torel(tbase)
80        floatTimes.append(time.value)
81    timeAx=cdms.createAxis(floatTimes)
82    timeAx.designateTime()
83    timeAx.id='time'
84    timeAx.standard_name='time'
85    timeAx.units=tbase 
86    var.setAxis(0,timeAx)     
87    nc.write(var)
88   
89   
90    #get standard_name table
91    sntable=csml.csmllibs.standardnames.__path__[0]+'/cf-standard-name-table.xml'
92    snt=csml.csmllibs.standardnames.snames.StandardNameTable(tablelocation=sntable)
93   
94    #add auxiliary coordinate variables to hold the locations
95    for ax in self.spatialAxes:
96        var=MV.array(self.spatialAxes[ax])
97        var.id=ax
98        #if id is standard_name, use it as a standard_name       
99        if snt.lookup(ax) is 1:
100            var.standard_name=ax
101        var.setAxis(0,timeAx)
102        nc.write(var)   
103               
104   
105    nc.close()   
106    print 'NetCDF file written to %s'%self.pathToSubsetNetCDF
107
108   
109def subsetToTrajectory(self, outputdir=None, ncname='trajectory.nc' ,time=None, **kwargs):
110    ''' subsetting a trajectory feature by time '''       
111       
112    #work out where to store the resulting files:
113    if outputdir is not None:
114        self.outputdir=outputdir
115    elif csml.API.csmlContainer.globalContainer.OUTPUTDIR is not None:
116        self.outputdir=csml.API.csmlContainer.globalContainer.OUTPUTDIR
117     
118    if self.outputdir is not None:
119        self.pathToSubsetNetCDF=self.outputdir + '/' +ncname
120    else: 
121        self.pathToSubsetNetCDF=ncname
122    if self.outputdir is not None:
123        csml.csmllibs.csmlextra.checkDirExists(self.outputdir)
124 
125 
126    #make sure single time string is converted to a tuple
127    if type(time) is str:
128        selection = tuple(csml.csmllibs.csmlextra.listify(time))
129    else:
130        selection = time 
131     
132    #repurpose the existing trajectoryDomain for the output csml.
133    td=copy.deepcopy(self.value.trajectoryDomain)
134
135     
136    #get the requested times, and also max/min index values, needed to get the corresponding spatial axes 
137    tlist=[]
138    for gridOrd in td.coordTransformTable.gridOrdinates:
139        if hasattr(gridOrd.coordAxisValues, 'timePositionList'):
140            for val in gridOrd.coordAxisValues.timePositionList.CONTENT.split():
141                tlist.append(val)       
142            indices=[]
143            if selection is None:
144                #select all times
145                selection=[tlist[0], tlist[len(tlist)-1]]
146            if  len(selection) ==2:       
147                self.selectedTimes=[]
148                for domaintime in tlist:
149                    if selection[0] == domaintime:
150                        minIndex=tlist.index(domaintime)
151                        self.selectedTimes.append(domaintime)
152                    if domaintime > selection[0]:
153                        if domaintime < selection[1]:
154                            self.selectedTimes.append(domaintime)         
155                    if selection[1] == domaintime:
156                        maxIndex=tlist.index(domaintime)
157                        self.selectedTimes.append(domaintime)         
158                        break                             
159            elif len(selection)==1:   #point only
160                for domaintime in tlist:
161                    if selection[0] == domaintime:               
162                        minIndex=tlist.index(domaintime)
163                        maxIndex=minIndex
164                        self.selectedTimes=selection
165                        break
166            if len(self.selectedTimes) ==1:
167                timestr=self.selectedTimes[0]
168            else:
169                timestr='  '.join(self.selectedTimes) 
170            #now substitute these times into the feature as we will use it later to write out a new feature.
171            gridOrd.coordAxisValues.timePositionList=csml.parser.csString(timestr)
172            gridOrd.coordAxisValues.id = csml.csmllibs.csmlextra.getRandomID()
173            break
174   
175    #change the high limit  (assumes low limit is 0)
176    highlim=''
177    for item in td.limits.high.CONTENT.split():
178        highlim=highlim+ ' ' + str(len(self.selectedTimes)-1)
179    td.limits.high.CONTENT=highlim   
180    #Now get the data from the spatial axes using the indices
181    self.spatialAxes={}
182    for gridOrd in td.coordTransformTable.gridOrdinates:
183        if hasattr(gridOrd.coordAxisValues, 'coordinateList'):
184            if hasattr(gridOrd.coordAxisValues, 'insertedExtract'):
185                data, fillvalue=gridOrd.coordAxisValues.insertedExtract.getDataFromChunks(minIndex, maxIndex)
186            else:
187                data=gridOrd.coordAxisValues.CONTENT
188            if len(data)==1:
189                datastr=data[0]
190            else:       
191                datastr = data.tostring()
192            gridOrd.coordAxisValues.coordinateList=csml.parser.csString(datastr)
193            gridOrd.coordAxisValues.id = csml.csmllibs.csmlextra.getRandomID()
194            self.spatialAxes[gridOrd.coordAxisLabel.CONTENT]=data       
195           
196            #get rid of the old xlink attributes
197            for att in [ 'href', 'role', 'arcrole', 'show','{http://www.w3.org/1999/xlink}href', '{http://www.w3.org/1999/xlink}role', '{http://www.w3.org/1999/xlink}arcrole', '{http://www.w3.org/1999/xlink}show']:
198                try:
199                    delattr(gridOrd.coordAxisValues, att)
200                except:
201                    pass
202   
203    #now get the actual rangeset data
204    if hasattr(self.value.rangeSet.valueArray.valueComponent, 'insertedExtract'):
205        self.data,self.fillvalue=self.value.rangeSet.valueArray.valueComponent.insertedExtract.getDataFromChunks(minIndex, maxIndex)
206    else:
207        print 'warning: inline rangeset not implemented'
208
209   
210    #Now create the file extract
211    fileList=[]
212    try:
213        fileList.append(self.value.rangeSet.valueArray.valueComponent.insertedExtract.fileName.CONTENT)
214        fextract=self.value.rangeSet.valueArray.valueComponent.insertedExtract
215        self.units=self.value.rangeSet.valueArray.valueComponent.uom
216    except:
217        for f in self.value.rangeSet.arrayDescriptor.components.fileList.fileNames.CONTENT.split():
218            fileList.append(f)
219            fextract=self.value.rangeSet.arrayDescriptor.components
220            self.units= self.value.rangeSet.arrayDescriptor.uom.CONTENT
221    #Now write out the CSML Feature:
222    #create a TrajectoryCoverage
223    cvg=csml.parser.TrajectoryCoverage()
224    cvg.id=csml.csmllibs.csmlextra.getRandomID()
225   
226    td.id=csml.csmllibs.csmlextra.getRandomID()
227               
228    sdid=csml.csmllibs.csmlextra.getRandomID() # for the storage descriptor
229    #create a RangeSet
230    rs=csml.parser.RangeSet()
231    va=csml.parser.ValueArray()
232    vc=csml.parser.MeasureOrNullList()
233    vc.href='#%s'%sdid
234    vc.arcrole="http://ndg.nerc.ac.uk/xlinkUsage/insert#QuantityList"
235    vc.role="http://ndg.nerc.ac.uk/fileFormat/csmlStorageDescriptor"
236    vc.show='embed'
237    try:
238        vc.uom=self.value.rangeSet.valueArray.valueComponent.uom
239    except:
240        vc.uom = 'unknown'  #TODO,
241    va.valueComponent=vc
242    va.id=csml.csmllibs.csmlextra.getRandomID()
243    rs.valueArray=va
244   
245    #Add the domain and rangeSet as attributes of the coverage
246    cvg.trajectoryDomain=td
247    cvg.rangeSet=rs
248   
249    totalArraySize=10
250    descriptor=csml.parser.NetCDFExtract(id=sdid,fileName=csml.parser.csString(ncname),variableName=self.name,arraySize=csml.parser.csString(totalArraySize))
251
252   
253    #the parameter of the feature is of type Phenomenon, here href creates "xlink:href=..."
254    param=self.parameter
255   
256   
257    #create a stand alone pointseries feature containing this coverage
258    csmlWrap=csml.csmllibs.csmlfeaturewrap.CSMLWrapper()
259    subsettedFeature=csmlWrap.createTrajectoryFeature(value=cvg,parameter=param,featureID=csml.csmllibs.csmlextra.getRandomID(),description=self.description)
260   
261    self._writeNetCDF()
262   
263
264    return subsettedFeature, self.pathToSubsetNetCDF, descriptor
265
266   
267   
Note: See TracBrowser for help on using the repository browser.