source: DPPP/kml/csml2kml/python/csml2kmlpylon/csml2kmlpylon/controllers/csmlGrapher.py @ 3636

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/DPPP/kml/csml2kml/python/csml2kmlpylon/csml2kmlpylon/controllers/csmlGrapher.py@3636
Revision 3636, 14.6 KB checked in by mkochan, 11 years ago (diff)

Made csmlGrapher use mapping of times to values, to avoid duplicates.

Line 
1# Pylons-specific imports
2import logging
3from csml2kmlpylon.lib.base import *
4from paste.httpexceptions import HTTPBadRequest, HTTPNotFound
5
6# Other imports
7import Image
8import pylab
9from pylab import *
10from datetime import datetime, timedelta
11import re
12import urllib
13from cStringIO import StringIO
14from tempfile import NamedTemporaryFile
15from cElementTree import Element, SubElement, ElementTree, XML
16from urllib import quote
17
18import csml
19import csml2kml.Station
20from csml2kml.utils import wget
21
22log = logging.getLogger(__name__)
23
24class CsmlgrapherController(BaseController):
25
26    def __call__(self, environ, start_response):
27        '''
28        Initialise the web service by overriding a method that always get called upon controller construction.
29        More specifically, read in the config file.
30        '''
31
32        # The name of the config file is set in the server's "development.ini" file.
33        configFileName = config['app_conf']['csmlGrapher.configfile']
34       
35        # Load the configuration XML element (but only the part pertaining to the web service)
36        grapherConfig = (ElementTree().parse(configFileName)).find('GrapherWebService')
37
38        # Set the configurable variables
39        self.displayIntervalStart = dates.dateutil.parser.parse(grapherConfig.find('DisplayIntervalStart').text)
40        self.displayIntervalEnd = dates.dateutil.parser.parse(grapherConfig.find('DisplayIntervalEnd').text)
41        self.geoServerUrl = grapherConfig.find('GeoServerURL').text
42
43        # Call the __call__ method of the parent class
44        return BaseController.__call__(self, environ, start_response)
45
46    def _datetimeToGeoServerDate(self, datetime):
47        monthCodes = ['JAN', 'FEB', 'MAR', 'APR', 'MAY', 'JUN', 'JUL', 'AUG', 'SEP', 'OCT', 'NOV', 'DEC']
48        return repr(datetime.day) + '-' + monthCodes[datetime.month-1] + '-' + repr(datetime.year)
49
50    def _retrieveCsmlPointSeriesFeature(self, station_name, feature_id):
51        '''
52        @return: A C{PointSeriesFeature} object representing the single CSML feature if the feature was found;
53        C{None} if a feature collection has been returned, but contains no CSML features (this happens when
54        there are no measured time points in the used time interval).
55        @throws: A C{LookupError} exception with a message if there is no response from the server.
56        '''
57       
58        # Prepare the request for the GeoServer.
59        geoServerRequestUrl = urllib.quote(self.geoServerUrl + '?request=getfeature&service=wfs&version=1.1.0&typename=csml:PointSeriesFeature&filter=<ogc:Filter xmlns:ogc="http://www.opengis.net/ogc" xmlns:gml="http://www.opengis.net/gml"><ogc:And><ogc:PropertyIsEqualTo><ogc:PropertyName>gml:description</ogc:PropertyName><ogc:Literal>'+ station_name +'</ogc:Literal></ogc:PropertyIsEqualTo><ogc:PropertyIsEqualTo><ogc:PropertyName>csml:parameter/swe:Phenomenon/gml:name</ogc:PropertyName><ogc:Literal>'+ feature_id +'</ogc:Literal></ogc:PropertyIsEqualTo><ogc:PropertyIsBetween><ogc:PropertyName>csml:value/csml:PointSeriesCoverage/csml:pointSeriesDomain/csml:TimeSeries/csml:timePositionList</ogc:PropertyName><ogc:LowerBoundary><ogc:Literal>'+ self._datetimeToGeoServerDate(self.displayIntervalStart) +'</ogc:Literal></ogc:LowerBoundary><ogc:UpperBoundary><ogc:Literal>'+ self._datetimeToGeoServerDate(self.displayIntervalEnd) +'</ogc:Literal></ogc:UpperBoundary></ogc:PropertyIsBetween></ogc:And></ogc:Filter>', '/._:?&=')
60
61        # Acquire the GeoServer response and parse it into a CSMLFeatureCollection object
62        geoServerResponse = wget(geoServerRequestUrl)
63        if not geoServerResponse:
64            raise LookupError('Cannot acquire response from server (wrong URL or server down)')
65        csmlFeatureCollection = csml.parser.CSMLFeatureCollection()
66        csmlFeatureCollection.fromXML( XML(geoServerResponse) )
67
68        # If csmlFeatureCollection contains a single CSML feature, the isolate it. If it contains no feature,
69        # probably no data points are available for that interval.
70        try:
71            csmlFeature = csmlFeatureCollection.featureMembers
72        except AttributeError:
73            raise LookupError(
74                'No data points available for the given time interval. ' + 
75                'Perhaps the interval is shorter than the frequency of data points, ' +
76                'or collection of data for this feature was interrupted at the time.'
77                )
78
79        # Make sure that the feature is a PointSeriesFeature
80        if not isinstance(csmlFeature, csml.parser.PointSeriesFeature):
81            raise LookupError('CSML feature not a PointSeriesFeature')
82        csmlPointSeriesFeature = csmlFeature
83
84        return csmlPointSeriesFeature
85
86    def plot(self):
87        '''
88        Handler for plotting a specific PointSeries feature given by URL parameter "feature_id".
89        Retrieves the data from a web service, which should preferably run on the same machine.
90
91        Request parameters:
92                               feature_id      Unique identifier of the requested CSML PointSeries feature.
93                               station_name    Name of the requested station.
94                               
95        Response:              An image/png of the time series for the CSML PointSeries feature.
96        '''
97        def _plot_feature(feature):
98            '''
99            Plot a PointSeries feature using matplotlib, into a temporary file
100            '''
101           
102            # Make sure that time is considered to be UTC time. This is because all CSML data is supposed to be in UTC,
103            # but actually pylab is "time zone aware" and will interpret the time as local time zone time.
104            # Not quite sure if this would be a problem, but rather enforce it is UTC.
105            def _enforce_UTC_timezone(datestr):
106                if datestr[-1] == 'Z':
107                    return datestr
108                else:
109                    return datestr + 'Z'
110           
111            # Isolate time points
112            times=feature.value.pointSeriesDomain.timePositionList.CONTENT.split()
113            times = map(_enforce_UTC_timezone, times)
114
115            # Make a list of times as float numbers being days passed since start of epoch (here 01-01-0001)
116            elapsed_times = map(datestr2num, times)
117
118            # Get values of the measured quantity
119            ql = feature.value.rangeSet.quantityList
120            vals = map(float, ql.CONTENT.split())
121
122            print '---times (' + str(len(times)) + ' of them):\n' + str(times)
123            print '---vals (' + str(len(times)) + ' of them):\n' + str(vals)
124
125            # Make sure the times are sorted in ascending order (and rearrange the values accordingly),
126            # this is because we get them potentially unsorted from GeoServer. Also have to handle the sad fact,
127            # that sometimes we are getting duplicate time-value samples from the MIDAS dataset!
128            sortTimesExplicitly = True
129            if sortTimesExplicitly:
130                reorder = {}
131                for i in range(len(elapsed_times)):
132                    reorder[elapsed_times[i]] = i
133                elapsed_times = reorder.keys()
134                elapsed_times.sort()
135                times2 = []
136                vals2 = []
137                for et in elapsed_times:
138                    times2.append(times[reorder[et]])
139                    vals2.append(vals[reorder[et]])
140                times = times2
141                vals = vals2
142
143            print '---times (' + str(len(times)) + ' of them):\n' + str(times)
144            print '---vals (' + str(len(times)) + ' of them):\n' + str(vals)
145
146            # Note the start and end date (they are actually date+time), and time span between them
147
148            start_date = dates.num2date(elapsed_times[0])
149            end_date = dates.num2date(elapsed_times[-1])
150            span = elapsed_times[-1] - elapsed_times[0]
151
152            # Generate tick locator and formatter; these determine the time axis: locator determines tick positions
153            # and formatter format of time
154            print '---start_date:' + str(start_date)
155            print '---end_date:' + str(end_date)
156            print '---span: ' + str(span)
157            (tickLocator, tickFormatter) = dates.date_ticker_factory(span, numticks=8)
158
159            # Determine the units of measurement
160            uom=ql.uom.title()
161            if ql.uom.islower():
162                uom = uom.lower()
163            if ql.uom.isupper():
164                uom = uom.upper()
165
166
167            # Prepare the plot the figure (actual plotting actions are postponed until save)
168            fig = figure()
169            plot_date(elapsed_times, vals, 'b-', xdate=True, lw=1)
170            plot_date(elapsed_times, vals, 'go', markeredgecolor = 'g', xdate=True, lw=2)
171            ax = gca()
172            ax.xaxis.set_major_locator(tickLocator)
173            ax.xaxis.set_major_formatter(tickFormatter)
174            fig.autofmt_xdate()  # show times without overlaps
175            time_format = '%d-%b-%Y %H:%M:%S UTC'
176            xlabel('Times between %s and %s ' % (start_date.strftime(time_format), end_date.strftime(time_format)))
177            ylabel('Values [%s]' % uom)
178            title('"%s" (%s)' % (feature.id, feature.description.CONTENT))
179            grid(True)
180
181            # Save the figure to a temporary file
182            tempFile = NamedTemporaryFile(suffix='.png')
183            savefig(tempFile.name)
184           
185            return tempFile
186
187        def _set_response(tempFile):
188            '''
189            Set the WSGI response to an image, containing image read from a temporary location.
190            '''
191            img = Image.open(tempFile.name)
192            buf = StringIO()
193            img.save(buf, 'PNG')
194            response.content_type = 'image/png'
195            response.content = buf.getvalue()
196           
197        #----------------------
198
199        # Get parameters from the request object
200        try:
201            feature_id = str(request.params['feature_id'])     # convert back from Unicode
202            station_name = str(request.params['station_name']) # convert back from Unicode
203        except KeyError:
204            raise HTTPBadRequest('Parameters "feature_id" and "station_name" must be supplied.')
205
206        # Try to retrieve the csmlPointSeriesFeature object
207        try:
208            csmlPointSeriesFeature = self._retrieveCsmlPointSeriesFeature(station_name, feature_id)
209        except LookupError, e:
210            errorMessage = e.args[0]
211            httpNotFound = HTTPNotFound()
212            httpNotFound.explanation = ''
213            httpNotFound.detail = errorMessage
214            raise httpNotFound
215
216        # Try to plot the feature into a temporary file, and put the contents of that file into the response
217        try:
218            tempFile = _plot_feature(csmlPointSeriesFeature)    # plot the feature into a temporary file
219            _set_response(tempFile)                             # set the response as an image containing the plot
220        finally:
221            try:
222                tempFile.close()
223            except NameError:
224                pass  # tempFile undefined -- no need to close
225
226    def list(self):
227        '''
228        Handler for generating a list of CSML features per station, with dynamic links to the plot service.
229        Request parameters:
230                               station_name    Name of the requested station.
231                               
232        Response:              A text/html with the page containing links to the plot service on the same server.
233        '''
234
235        # Get parameters from the request object
236        try:
237            station_name = str(request.params['station_name']) # convert back from Unicode
238        except KeyError:
239            raise HTTPBadRequest('Parameter "station_name" must be supplied.')
240
241        # Prepare the request for the GeoServer -- to return a collection of np:Station, containing a single np:Station,
242        # which contains a list of CSML features.
243        geoServerRequestUrl = urllib.quote(self.geoServerUrl + '?request=getFeature&service=wfs&version=1.1.0&typename=np:Station&filter=<ogc:Filter xmlns:ogc="http://www.opengis.net/ogc" xmlns:gml="http://www.opengis.net/gml"><ogc:PropertyIsEqualTo><ogc:PropertyName>stationName</ogc:PropertyName><ogc:Literal>' + station_name + '</ogc:Literal></ogc:PropertyIsEqualTo></ogc:Filter>', '/._:?&=')
244
245        # Acquire the GeoServer response and parse it into a WFSStationCollection object
246        geoServerResponse = wget(geoServerRequestUrl)
247        if not geoServerResponse:
248            raise LookupError('Cannot acquire response from server (wrong URL or server down)')
249        wfsStationCollection = csml2kml.Station.WFSStationCollection()
250        wfsStationCollection.parseString(geoServerResponse)
251        if len(wfsStationCollection.stations) != 1:
252            raise ValueError('Multiple stations match OGC selection filter (only one must match)')
253        wfsStation = wfsStationCollection.stations[0]
254
255        # Generate HTML with a list of CSML features contained in the station.
256        # Currently, this is done simply by directly generating HTML (i.e. no template was used).
257        htmlElement = Element('html')
258        SubElement(htmlElement, 'title').text = 'List of CSML features for station ' + station_name
259        bodyElement = SubElement(htmlElement, 'body')
260        SubElement(bodyElement, 'h2').text = 'List of CSML features for station ' + station_name
261        tableElement = SubElement(bodyElement, 'table')
262        tableElement.set('border', '1')
263        headingRowElement = SubElement(tableElement, 'tr')
264        SubElement(headingRowElement, 'th').text = 'Feature name'
265        SubElement(headingRowElement, 'th').text = 'Data collected since*'
266        SubElement(headingRowElement, 'th').text = 'Data collected until* '
267        for stationFeature in wfsStation.stationFeatures:
268            if self.displayIntervalStart >= stationFeature.collectBeginDate and self.displayIntervalEnd <= stationFeature.collectEndDate:
269                rowElement = SubElement(tableElement, 'tr')
270                featureNameElement = SubElement(rowElement, 'td')
271                anchorElement = SubElement(featureNameElement, 'a')
272                linkToGrapher = 'http://bond.badc.rl.ac.uk:8089/csmlGrapher/plot?station_name=' + station_name + '&feature_id=' + stationFeature.featureId
273                anchorElement.set('href', linkToGrapher)
274                anchorElement.text = stationFeature.featureId
275                SubElement(rowElement, 'td').text = self._datetimeToGeoServerDate(stationFeature.collectBeginDate)
276                SubElement(rowElement, 'td').text = self._datetimeToGeoServerDate(stationFeature.collectEndDate)
277
278        htmlStringIO = StringIO()
279        ElementTree(htmlElement).write(htmlStringIO)
280
281        response.content_type = 'text/html'
282        response.content = htmlStringIO.getvalue()
Note: See TracBrowser for help on using the repository browser.