source: TI02-CSML/trunk/csml/API/ops_PointSeriesFeature.py @ 2565

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

subset to point feature working

Line 
1''' ops_PointSeriesFeature  contains operations for PointSeriesFeatures'''
2import csml,sys
3import cdms, MV, MA
4
5def testmethod(self):
6    print 'testmethod for pointSeries feature'
7    return 'testmethod pointSeries'
8
9def getAllowedSubsettings(self):
10    return ['subsetToPointSeries']
11   
12def getDomain(self):
13    #returns domain as a dictionary of ordinates {name: [values], ...}
14    self.domain={}
15    self.gridnames={}
16    for gridOrd in self.value.pointSeriesDomain.coordTransformTable.gridOrdinates:
17        name=gridOrd.coordAxisLabel.CONTENT
18        self.gridnames[gridOrd.gridAxesSpanned.CONTENT]=name
19        if hasattr(gridOrd.coordAxisValues, 'insertedExtract'):
20            self.domain[name], fill, axisorder, units=gridOrd.coordAxisValues.insertedExtract.getData()
21        else:           
22            valList=[]
23            try:
24                vals=gridOrd.coordAxisValues.coordinateList.CONTENT
25                for val in vals.split(): 
26                    valList.append(eval(val))
27            except:
28                vals=gridOrd.coordAxisValues.timePositionList.CONTENT
29                for val in vals.split(): 
30                    valList.append(val)           
31            self.domain[name]=valList   
32    return self.domain   
33
34
35
36   
37def subsetToPointSeries(self, csmlpath=None, ncpath=None,times=None):
38    pathToSubsetNetCDF=ncpath
39    if type(times) is str:
40        selection = tuple(csml.csmllibs.csmlextra.listify(times))
41    else:
42        selection = times 
43    tlist=[]   
44    for val in self.value.pointSeriesDomain.timePositionList.CONTENT.split():
45        tlist.append(val)
46   
47    indices=[]
48    if  len(selection) ==2:       
49        selectedTimes=[]
50        for domaintime in tlist:
51            if selection[0] == domaintime:
52                minIndex=tlist.index(domaintime)
53                selectedTimes.append(domaintime)
54            if domaintime > selection[0]:
55                if domaintime < selection[1]:
56                    selectedTimes.append(domaintime)         
57            if selection[1] == domaintime:
58                maxIndex=tlist.index(domaintime)
59                selectedTimes.append(domaintime)         
60                break   
61                   
62    elif len(selection)==1:   #point only
63        for domaintime in tlist:
64            if selection[0] == domaintime:
65               
66                minIndex=tlist.index(domaintime)
67                maxIndex=minIndex
68                selectedTimes=selection
69                break
70   
71   
72   
73    fileList=[]
74    for f in self.value.rangeSet.arrayDescriptor.components.fileList.fileNames.CONTENT.split():
75        fileList.append(f)
76   
77    fextract=self.value.rangeSet.arrayDescriptor.components
78    uom = self.value.rangeSet.arrayDescriptor.uom
79    data, fillvalue=fextract.getDataFromChunks(minIndex, maxIndex)
80   
81    units=[uom.CONTENT]
82   
83    #Now write out the CSML Feature:
84    #create a PointSeriesCoverage
85    cvg=csml.parser.PointSeriesCoverage()
86    cvg.id=csml.csmllibs.csmlextra.getRandomID()
87   
88    #create a TimeSeriesDomain (appropriate domain for a PointSeriesFeature
89    ptsd=csml.parser.TimeSeries()
90    ptsd.id=csml.csmllibs.csmlextra.getRandomID()
91    #create (and populate) a TimePositionList. Using keyword arguements for conciseness.
92    ptsd.timePositionList=csml.parser.TimePositionList(CONTENT=csml.csmllibs.csmlextra.stringify(selectedTimes))
93       
94    sdid=csml.csmllibs.csmlextra.getRandomID() # for the storage descriptor
95    #create a RangeSet
96    rs=csml.parser.RangeSet()
97    va=csml.parser.ValueArray()
98    vc=csml.parser.MeasureOrNullList()
99    vc.href='#%s'%sdid
100    vc.arcrole="http://ndg.nerc.ac.uk/xlinkUsage/insert#QuantityList"
101    vc.role="http://ndg.nerc.ac.uk/fileFormat/csmlStorageDescriptor"
102    vc.show='embed'
103    try:
104        vc.uom=self.value.rangeSet.valueArray.valueComponent.uom
105    except:
106        vc.uom = 'unknown'  #TODO,
107    va.valueComponent=vc
108    va.id=csml.csmllibs.csmlextra.getRandomID()
109    rs.valueArray=va
110   
111   
112    #Add the domain and rangeSet as attributes of the coverage
113    cvg.pointSeriesDomain=ptsd
114    cvg.rangeSet=rs
115   
116    totalArraySize=10
117    descriptor=csml.parser.NetCDFExtract(id=sdid,fileName=csml.parser.csString(pathToSubsetNetCDF),variableName=self.name,arraySize=csml.parser.csString(totalArraySize))
118
119   
120    #the parameter of the feature is of type Phenomenon, here href creates "xlink:href=..."
121    param=self.parameter
122   
123     
124    #create a stand alone pointseries feature containing this coverage
125    csmlWrap=csml.csmllibs.csmlfeaturewrap.CSMLWrapper()
126    subsettedFeature=csmlWrap.createPointSeriesFeature(value=cvg,parameter=param,location=self.location,featureID=self.id,description=self.description)
127   
128     
129   
130    #Now write out the NetCDF File:   
131    nc=cdms.open(pathToSubsetNetCDF,'w')
132    var=MV.array(data)
133    var.id=self.name.CONTENT
134    var.name=self.name.CONTENT
135    var.units=units[0] # hopefully there is just one unit in this list..
136    setattr(var, 'missing_value' ,fillvalue)
137   
138    floatTimes=[]
139   
140    tOne=csml.csmllibs.csmltime.getCDtime(selection[0])
141    tbase=csml.csmllibs.csmltime.getBaseUnits(tOne)
142    for time in selectedTimes:
143        time=csml.csmllibs.csmltime.getCDtime(time).torel(tbase)
144        floatTimes.append(time.value)
145    timeAx=cdms.createAxis(floatTimes)
146    timeAx.designateTime()
147    timeAx.id='time'
148    timeAx.units=tbase 
149    var.setAxis(0,timeAx)     
150    nc.write(var)
151    nc.close()   
152    print 'NetCDF file written to %s'%pathToSubsetNetCDF
153    return subsettedFeature, pathToSubsetNetCDF, descriptor
154   
155   
156def subsetToPoint(self, csmlpath=None, ncpath=None,time=None):
157    pseriesFeature, pathToSubsetNetCDF, descriptor=self.subsetToPointSeries(csmlpath,ncpath,time)
158           
159    #now need to take this point series feature containing one time and rewrite it as a point feature.
160    pointFeature=csml.parser.PointFeature()
161    pointFeature.location=pseriesFeature.location
162    pointFeature.parameter=pseriesFeature.parameter
163    pointFeature.description=pseriesFeature.description
164    cvg=csml.parser.PointCoverage()
165    domain=csml.parser.PointDomain()
166    domain.id=csml.csmllibs.csmlextra.getRandomID()
167    domain.pointMember=csml.parser.csString(pseriesFeature.location.CONTENT)
168    #create (and populate) a TimePositionList. Using keyword arguements for conciseness.
169    cvg.pointDomain=domain
170    pointFeature.value=cvg
171    pointFeature.time=csml.parser.csString(time)
172    return pointFeature, pathToSubsetNetCDF, descriptor
Note: See TracBrowser for help on using the repository browser.