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

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

subsetting to point basic code added

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