1 | # Adapted for numpy/ma/cdms2 by convertcdms.py |
---|
2 | ''' ops_PointSeriesFeature contains operations for PointSeriesFeatures. These methods are attached to PointSeriesFeature instances at runtime.''' |
---|
3 | import csml,sys |
---|
4 | import cdms2 as cdms, MV2 as MV, numpy.oldnumeric.ma as MA |
---|
5 | |
---|
6 | def testmethod(self): |
---|
7 | ''' testmethod for gridseries feature''' |
---|
8 | print 'testmethod for pointSeries feature' |
---|
9 | return 'testmethod pointSeries' |
---|
10 | |
---|
11 | def getAllowedSubsettings(self): |
---|
12 | '''get the allowed operations |
---|
13 | @return: list of operations''' |
---|
14 | return ['subsetToPointSeries'] |
---|
15 | |
---|
16 | def getDomain(self): |
---|
17 | ''' returns domain as a dictionary of ordinates {name: [values], ...} ''' |
---|
18 | self.domain={} |
---|
19 | valList=[] |
---|
20 | for val in self.value.pointSeriesDomain.timePositionList.CONTENT.split(): |
---|
21 | valList.append(val) |
---|
22 | self.domain['time']=valList |
---|
23 | [lat,lon]=self.location.CONTENT.split() |
---|
24 | self.domain['longitude']=lon |
---|
25 | self.domain['latitude']=lat |
---|
26 | |
---|
27 | return self.domain |
---|
28 | |
---|
29 | def getLongitudeAxis(self): |
---|
30 | ''' get the index of the longitude axis''' |
---|
31 | return 'longitude' #dummy value not used, subsetting is by time only |
---|
32 | |
---|
33 | def getLatitudeAxis(self): |
---|
34 | ''' get the index of the latitude axis''' |
---|
35 | return 'latitude' #dummy value not used, subsetting is by time only |
---|
36 | |
---|
37 | def getTimeAxis(self): |
---|
38 | ''' get the index of the time axis''' |
---|
39 | return 'times' |
---|
40 | |
---|
41 | def subsetToPointSeries(self, outputdir=None, ncname='pointseries.nc' ,times=None): |
---|
42 | '''Subset a PointSeriesFeature to a PointSeriesFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs. |
---|
43 | |
---|
44 | @param outputdir: name of directory to write data files to |
---|
45 | @param ncname: name of netcdf file to write out |
---|
46 | @param times: selection by times |
---|
47 | @return: subsettedFeature (PointSeriesFeature instance) pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)''' |
---|
48 | if outputdir is not None: |
---|
49 | self.outputdir=outputdir |
---|
50 | elif csml.API.csmlContainer.globalContainer.OUTPUTDIR is not None: |
---|
51 | self.outputdir=csml.API.csmlContainer.globalContainer.OUTPUTDIR |
---|
52 | |
---|
53 | if self.outputdir is not None: |
---|
54 | pathToSubsetNetCDF=self.outputdir + '/' +ncname |
---|
55 | else: |
---|
56 | pathToSubsetNetCDF=ncname |
---|
57 | if self.outputdir is not None: |
---|
58 | csml.csmllibs.csmlextra.checkDirExists(self.outputdir) |
---|
59 | |
---|
60 | if type(times) is str: |
---|
61 | selection = tuple(csml.csmllibs.csmlextra.listify(times)) |
---|
62 | else: |
---|
63 | selection = times |
---|
64 | tlist=[] |
---|
65 | for val in self.value.pointSeriesDomain.timePositionList.CONTENT.split(): |
---|
66 | tlist.append(val) |
---|
67 | |
---|
68 | indices=[] |
---|
69 | if len(selection) ==2: |
---|
70 | selectedTimes=[] |
---|
71 | for domaintime in tlist: |
---|
72 | if selection[0] == domaintime: |
---|
73 | minIndex=tlist.index(domaintime) |
---|
74 | selectedTimes.append(domaintime) |
---|
75 | if domaintime > selection[0]: |
---|
76 | if domaintime < selection[1]: |
---|
77 | selectedTimes.append(domaintime) |
---|
78 | if selection[1] == domaintime: |
---|
79 | maxIndex=tlist.index(domaintime) |
---|
80 | selectedTimes.append(domaintime) |
---|
81 | break |
---|
82 | |
---|
83 | elif len(selection)==1: #point only |
---|
84 | for domaintime in tlist: |
---|
85 | if selection[0] == domaintime: |
---|
86 | |
---|
87 | minIndex=tlist.index(domaintime) |
---|
88 | maxIndex=minIndex |
---|
89 | selectedTimes=selection |
---|
90 | break |
---|
91 | |
---|
92 | |
---|
93 | |
---|
94 | #fileList=[] |
---|
95 | |
---|
96 | if hasattr(self.value.rangeSet, 'quantityList'): |
---|
97 | #inline |
---|
98 | data=self.value.rangeSet.quantityList.CONTENT.split()[minIndex:maxIndex+1] |
---|
99 | fillvalue=None |
---|
100 | units=self.value.rangeSet.quantityList.uom |
---|
101 | else: |
---|
102 | #file extract |
---|
103 | if hasattr(self.value.rangeSet, 'valueArray'): |
---|
104 | fextract=self.value.rangeSet.valueArray.valueComponent.insertedExtract.components |
---|
105 | uom = self.value.rangeSet.valueArray.valueComponent.insertedExtract.uom |
---|
106 | else: |
---|
107 | fextract=self.value.rangeSet.arrayDescriptor.components |
---|
108 | uom = self.value.rangeSet.arrayDescriptor.uom |
---|
109 | #for f in fextract.fileList.fileNames.CONTENT.split(): |
---|
110 | #fileList.append(f) |
---|
111 | |
---|
112 | data, fillvalue=fextract.getDataFromChunks(minIndex, maxIndex) |
---|
113 | units=[uom.CONTENT] |
---|
114 | |
---|
115 | #Now write out the CSML Feature: |
---|
116 | #create a PointSeriesCoverage |
---|
117 | cvg=csml.parser.PointSeriesCoverage() |
---|
118 | cvg.id=csml.csmllibs.csmlextra.getRandomID() |
---|
119 | |
---|
120 | #create a TimeSeriesDomain (appropriate domain for a PointSeriesFeature |
---|
121 | ptsd=csml.parser.TimeSeries() |
---|
122 | ptsd.id=csml.csmllibs.csmlextra.getRandomID() |
---|
123 | #create (and populate) a TimePositionList. Using keyword arguements for conciseness. |
---|
124 | ptsd.timePositionList=csml.parser.TimePositionList(CONTENT=csml.csmllibs.csmlextra.stringify(selectedTimes)) |
---|
125 | |
---|
126 | sdid=csml.csmllibs.csmlextra.getRandomID() # for the storage descriptor |
---|
127 | #create a RangeSet |
---|
128 | rs=csml.parser.RangeSet() |
---|
129 | va=csml.parser.ValueArray() |
---|
130 | vc=csml.parser.MeasureOrNullList() |
---|
131 | vc.href='#%s'%sdid |
---|
132 | vc.arcrole="http://ndg.nerc.ac.uk/xlinkUsage/insert#QuantityList" |
---|
133 | vc.role="http://ndg.nerc.ac.uk/fileFormat/csmlStorageDescriptor" |
---|
134 | vc.show='embed' |
---|
135 | try: |
---|
136 | vc.uom=self.value.rangeSet.valueArray.valueComponent.uom |
---|
137 | except: |
---|
138 | vc.uom = 'unknown' #TODO, |
---|
139 | va.valueComponent=vc |
---|
140 | va.id=csml.csmllibs.csmlextra.getRandomID() |
---|
141 | rs.valueArray=va |
---|
142 | |
---|
143 | |
---|
144 | #Add the domain and rangeSet as attributes of the coverage |
---|
145 | cvg.pointSeriesDomain=ptsd |
---|
146 | cvg.rangeSet=rs |
---|
147 | |
---|
148 | totalArraySize=10 |
---|
149 | if hasattr(self, 'name'): |
---|
150 | varname=self.name.CONTENT |
---|
151 | else: |
---|
152 | varname=self.id |
---|
153 | descriptor=csml.parser.NetCDFExtract(id=sdid,fileName=csml.parser.csString(ncname),variableName=csml.parser.csString(varname),arraySize=csml.parser.csString(totalArraySize)) |
---|
154 | |
---|
155 | #the parameter of the feature is of type Phenomenon, here href creates "xlink:href=..." |
---|
156 | param=self.parameter |
---|
157 | |
---|
158 | |
---|
159 | #create a stand alone pointseries feature containing this coverage |
---|
160 | csmlWrap=csml.csmllibs.csmlfeaturewrap.CSMLWrapper() |
---|
161 | subsettedFeature=csmlWrap.createPointSeriesFeature(value=cvg,parameter=param,location=self.location,featureID=self.id,description=self.description) |
---|
162 | |
---|
163 | |
---|
164 | #Now write out the NetCDF File: |
---|
165 | nc=cdms.open(pathToSubsetNetCDF,'w') |
---|
166 | var=MV.array(data) |
---|
167 | var.id=varname |
---|
168 | var.name=varname |
---|
169 | var.units=units[0] # hopefully there is just one unit in this list.. |
---|
170 | setattr(var, 'missing_value' ,fillvalue) |
---|
171 | |
---|
172 | floatTimes=[] |
---|
173 | |
---|
174 | tOne=csml.csmllibs.csmltime.getCDtime(selection[0]) |
---|
175 | tbase=csml.csmllibs.csmltime.getBaseUnits(tOne) |
---|
176 | for time in selectedTimes: |
---|
177 | time=csml.csmllibs.csmltime.getCDtime(time).torel(tbase) |
---|
178 | floatTimes.append(time.value) |
---|
179 | timeAx=cdms.createAxis(floatTimes) |
---|
180 | timeAx.designateTime() |
---|
181 | timeAx.id='time' |
---|
182 | timeAx.units=tbase |
---|
183 | var.setAxis(0,timeAx) |
---|
184 | nc.write(var) |
---|
185 | nc.close() |
---|
186 | print 'NetCDF file written to %s'%pathToSubsetNetCDF |
---|
187 | return subsettedFeature, pathToSubsetNetCDF, descriptor |
---|
188 | |
---|
189 | |
---|
190 | def subsetToPoint(self, outputdir=None, ncname='point.nc',time=None): |
---|
191 | '''Subset a PointSeriesFeature to a PointFeature. Performs the subset (note this includes nearest neighbour searching, so may return a different set of kwargs. |
---|
192 | |
---|
193 | @param outputdir: name of directory to write data files to |
---|
194 | @param ncname: name of netcdf file to write out |
---|
195 | @param time: selection by single time |
---|
196 | @return: subsettedFeature (PointFeature instance) pathToSubsetNetCDF (filepath), descriptor (array descriptor instance)''' |
---|
197 | if outputdir is not None: |
---|
198 | self.outputdir=outputdir |
---|
199 | elif csml.API.csmlContainer.globalContainer.OUTPUTDIR is not None: |
---|
200 | self.outputdir=csml.API.csmlContainer.globalContainer.OUTPUTDIR |
---|
201 | pseriesFeature, pathToSubsetNetCDF, descriptor=self.subsetToPointSeries(self.outputdir, ncname,time) |
---|
202 | #now need to take this point series feature containing one time and rewrite it as a point feature. |
---|
203 | pointFeature=csml.parser.PointFeature() |
---|
204 | pointFeature.location=pseriesFeature.location |
---|
205 | pointFeature.parameter=pseriesFeature.parameter |
---|
206 | pointFeature.description=pseriesFeature.description |
---|
207 | cvg=csml.parser.PointCoverage() |
---|
208 | domain=csml.parser.PointDomain() |
---|
209 | domain.id=csml.csmllibs.csmlextra.getRandomID() |
---|
210 | domain.pointMember=csml.parser.csString(pseriesFeature.location.CONTENT) |
---|
211 | #create (and populate) a TimePositionList. Using keyword arguements for conciseness. |
---|
212 | cvg.pointDomain=domain |
---|
213 | pointFeature.value=cvg |
---|
214 | pointFeature.time=csml.parser.csString(time) |
---|
215 | return pointFeature, pathToSubsetNetCDF, descriptor |
---|