Project

General

Profile

Statistics
| Branch: | Tag: | Revision:

pycama / src / pycama / IrradiancePlot.py @ 834:d989df597b80

History | View | Annotate | Download (13.2 KB)

1
#!/usr/bin/env python3
2
# -*- coding: utf-8 -*-
3

    
4
# Copyright 2016-2017 Maarten Sneep, KNMI
5
#
6
# Redistribution and use in source and binary forms, with or without
7
# modification, are permitted provided that the following conditions are met:
8
#
9
# 1. Redistributions of source code must retain the above copyright notice,
10
#    this list of conditions and the following disclaimer.
11
#
12
# 2. Redistributions in binary form must reproduce the above copyright notice,
13
#    this list of conditions and the following disclaimer in the documentation
14
#    and/or other materials provided with the distribution.
15
#
16
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
19
# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
20
# ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
21
# (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
22
# LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
23
# ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24
# (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
25
# SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26

    
27
## @file IrradiancePlot.py
28
#
29
# This file defines a subclass of the pycama.AnalysisAndPlot.AnalysisAndPlot class.
30
# This subclass gathers the scanline independen data.
31
#
32
# @author Maarten Sneep
33

    
34

    
35
import math
36
import logging
37
from collections import OrderedDict
38
import uuid
39
import warnings
40
warnings.filterwarnings("ignore", category=FutureWarning)
41

    
42
import numpy as np
43
import netCDF4
44
import h5py
45

    
46
from .AnalysisAndPlot import AnalysisAndPlot
47

    
48
from .utilities import *
49

    
50
## This class handles irradiance dependent data (i.e. data that has no scanline dimension.).
51
#
52
#  This is a subclass of the pycama.AnalysisAndPlot.AnalysisAndPlot class.
53
#  The irradiance data is gathered once for each unique irradiance input file.
54
#  This requires that the metadata (in the L2 file) is complete. This is handled in
55
#  the pycama.File.File class, both for S5P and OMI files. When not available, all
56
#  variables processed here are treated as being unique.
57
#
58
#  This class relies on the pycama.Reader.Reader class to extract
59
#  the desired data from the level 2 input files. All variables that are determined
60
#  to have no scanline dependence end up here.
61
class IrradiancePlot(AnalysisAndPlot):
62
    ## Define a few variables to collect the irradiance data.
63
    def setup(self, *args, **kwargs):
64
        self.raw_variables = OrderedDict()
65
        self.dimensions = {}
66
        self.has_contents = False
67

    
68
    ## Return True if the variable must be plotted (included) here.
69
    def include_var(self, variable):
70
        return variable.noscanline
71

    
72
    ## Add a variable to the list 'to be processed'
73
    def add_variable(self, var):
74
        self.logger.debug("{0}.add_variable: Adding {1}".format(self.__class__.__name__, var.name))
75
        self.raw_variables[var.name] = var.data
76
        super(self.__class__, self).add_variable(var)
77

    
78
    ## Add the input variables.
79
    #
80
    #  This also ensures uniqueness of the
81
    #  input data (orbit number of the irradiance).
82
    def process(self):
83
        var_names = [k for k in self.raw_variables.keys()]
84
        if len(var_names) == 0:
85
            self.has_contents = False
86
            self.logger.info("No data to process in irradiance action")
87
            return
88
        self.has_contents = True
89

    
90
        orbits = np.asarray([int(k['orbit']) for k in self.raw_variables[var_names[0]].values()], dtype=np.int32)
91
        for k in self.raw_variables[var_names[0]].values():
92
            if k['source'] is None:
93
                k['source'] = "Group_{0}".format(uuid.uuid4()).replace('-', '')
94
        irrad_granules = sorted(list(set([k['source'] for k in self.raw_variables[var_names[0]].values()])))
95
        num_orbits = len(irrad_granules)
96

    
97
        result = OrderedDict()
98
        dimensions = {}
99

    
100
        for k in irrad_granules:
101
            result[name_cleanup(k)] = None
102

    
103
        for orbit in orbits:
104
            key = name_cleanup(self.raw_variables[var_names[0]]["{0:05d}".format(orbit)]['source'])
105
            try:
106
                time = datetime.datetime.strptime(key[11:26], "%Y%m%dT%H%M%S").strftime("%Y-%m-%dT%H:%M:%S")
107
            except:
108
                time = "0001-01-01T00:00:00"
109
            if result[key] is None:
110
                result[key] = {}
111
                for name in var_names:
112
                    data = self.raw_variables[name]["{0:05d}".format(orbit)]
113
                    result[key][name] = {'data': data['data'], 'originator': data['originator'], 'orbit': orbit, 'time': time}
114
                    if data['originator'] == 'KNMI':
115
                        result[key][name]['dimensions'] = ('ground_pixel',)
116
                    else: # DLR
117
                        result[key][name]['dimensions'] = tuple([(k if k != 'number_of_calibrations' else 'ground_pixel') for k in data['dimension_names']])
118
                        for dname, dim in data['dimensions'].items():
119
                            if dname not in dimensions:
120
                                dimensions[dname if (dname != 'number_of_calibrations') else 'ground_pixel'] = dim
121
                        if 'calibration_subwindows_wavelength' in data and 'calibration_subwindows_wavelength' not in result[key]:
122
                            result[key]['calibration_subwindows_wavelength'] = {'data': data['calibration_subwindows_wavelength']}
123
                            result[key]['calibration_subwindows_wavelength']['dimensions'] = ('number_of_subwindows',)
124
        self.variables = result
125
        self.dimensions = dimensions
126

    
127
    ## Write processed data to output file.
128
    #
129
    #  @param fname File to write to
130
    #  @param mode  Writing mode, defaults to append.
131
    #
132
    #  Do nothing if there is no contents.
133
    #
134
    #  Write data (including extraction specific dimensions) to the group with
135
    #  the name given in the `storage_group_name` property ("`irradianceplot_data`" for this class).
136

    
137
    def dump(self, fname, mode='a'):
138
        if not self.has_contents:
139
            return
140

    
141
        compress={'compression':'gzip', 'compression_opts':3, 'shuffle':True, 'fletcher32':True}
142
        with h5py.File(fname, mode=mode) as ref:
143
            time_step = self.time_index_in_output # set in __init__.
144
            try:
145
                grp = ref.create_group(self.storage_group_name)
146
                grp.attrs['comment'] = "irradiance related data"
147
            except:
148
                grp = ref[self.storage_group_name]
149

    
150
            # dimensions
151
            for k, v in self.dimensions.items():
152
                if k == 'ground_pixel':
153
                    continue
154
                try:
155
                    var = grp.create_dataset(k, (len(v),), dtype=np.float32, **compress)
156
                    var[:] = v
157
                except:
158
                    self.logger.warning("Dimension {0} could not be created".format(k))
159
                    pass
160

    
161
            for key, v in self.variables.items():
162
                mod_key = os.path.splitext(key)[0]
163
                try:
164
                    vgrp = grp.create_group(mod_key)
165
                    vgrp.attrs['comment'] = "Data from irradiance file {0}".format(key)
166
                except:
167
                    vgrp = grp[mod_key]
168

    
169
                for name, vv in v.items():
170
                    if 'originator' in vv:
171
                        vgrp.attrs['originator'] = vv['originator']
172
                        vgrp.attrs['orbit'] = vv['orbit']
173
                        vgrp.attrs['time_step'] = time_step
174
                        vgrp.attrs['time'] = vv['time']
175

    
176
                    if name in vgrp:
177
                        continue
178

    
179
                    var = vgrp.create_dataset(name,
180
                                              vv['data'].shape,
181
                                              dtype=np.float32,
182
                                              **compress)
183
                    for i, k in enumerate(vv['dimensions']):
184
                        if k in ref.keys():
185
                            var.dims.create_scale(ref[k], k)
186
                            var.dims[i].attach_scale(ref[k])
187
                        elif k in grp.keys():
188
                            var.dims.create_scale(grp[k], k)
189
                            var.dims[i].attach_scale(grp[k])
190
                        else:
191
                            raise CAMAException("undefined dimension {0}".format(k))
192

    
193
                    var[:] = vv['data']
194
                    try:
195
                        for k, val in self.variables_meta[name].items():
196
                            if val is None:
197
                                continue
198
                            if type(val) == type(True):
199
                                var.attrs[k] = str(val)
200
                            elif k == 'data_range':
201
                                var.attrs[k] = [float(range_value) for range_value in val]
202
                            else:
203
                                var.attrs[k] = val
204
                    except:
205
                        pass
206

    
207
    ## Read processed data from the input file, for specified time index.
208
    #
209
    #  @param fname NetCDF file with input data.
210
    #  @param time_index Time slice to read.
211
    def ingest(self, fname, time_step, exclude=None):
212
        self.time_index_in_output = time_step
213
        with h5py.File(fname, mode='r') as ref:
214
            try:
215
                grp = ref[self.storage_group_name]
216
            except:
217
                self.has_contents = False
218
                return False
219

    
220
            self.has_contents = True
221
            for dim_name in [n for n in grp.keys() if hasattr(grp[n], 'shape') and len(grp[n].shape) == 1]:
222
                try:
223
                    self.dimensions[name] = grp[name][:]
224
                except:
225
                    self.logger.warning("Variable '{0}' not found in '{1}'.".format(name, self.storage_group_name))
226

    
227
            for key, vgrp in [(k,v) for (k,v) in grp.items() if isinstance(v, h5py.Group)]:
228
                if vgrp.attrs['time_step'] == time_step:
229
                    self.variables[key] = {}
230
                    for name, vv in [(kkey,vval) for (kkey, vval) in vgrp.items() if isinstance(vval, h5py.Dataset)]:
231
                        if exclude is not None and name in exclude:
232
                            continue
233
                            
234
                        self.variables[key][name] = {}
235
                        self.variables[key][name]['data'] = vv[:]
236
                        self.variables[key][name]['dimensions'] = [d.keys()[0] for d in vv.dims]
237
                        self.variables[key][name]['originator'] = vgrp.attrs['originator']
238
                        self.variables[key][name]['orbit'] = vgrp.attrs['orbit']
239
                        try:
240
                            self.variables[key][name]['time'] = vgrp.attrs['time']
241
                        except KeyError:
242
                            time = datetime.datetime.strptime(key[11:26], "%Y%m%dT%H%M%S").strftime("%Y-%m-%dT%H:%M:%S")
243
                            self.variables[key][name]['time'] = time
244
                        self.variables_meta[name] = {}
245
                        for atname in vv.attrs.keys():
246
                            self.variables_meta[name][atname] = vv.attrs[atname]
247
                            if isinstance(self.variables_meta[name][atname], bytes):
248
                                self.variables_meta[name][atname] = str(self.variables_meta[name][atname], 'utf-8')
249
        return True
250

    
251
    ## Merge data into a combined dataset.
252
    #
253
    #  @param other The object to be added to self,
254
    #               also an instance of the pycama.IrradiancePlot.IrradiancePlot class.
255
    def __iadd__(self, other):
256
        for key in list(other.variables.keys()):
257
            if key not in self.variables:
258
                self.variables[key] = {}
259
            for name in list(other.variables[key].keys()):
260
                if name not in self.variables[key]:
261
                    self.variables[key][name] = {}
262
                else:
263
                    continue # only collect unique irradiance granules
264
                if name not in self.variables_meta[name]:
265
                    self.variables_meta[name] = {}
266
                for atname in other.variables_meta[name].keys():
267
                    self.variables_meta[name][atname] = other.variables_meta[name][atname]
268
                for part in list(other.variables[key][name].keys()):
269
                    self.variables[key][name][part] = other.variables[key][name][part]
270

    
271

    
272
    ## Make a plot of a specified variable
273
    #
274
    #  @param varname The name of the variable to plot.
275
    #  @param figure  The matplotlib.figure.Figure instance to plot to.
276
    #  @param ax      The matplotlib.axes.Axes object to use. Default is None, in that case `ax = matplotlib.pyplot.gca()` shall be used.
277
    #  @param kwargs  Other keyword parameters (sub-class dependent)
278
    #
279
    #  No plotting is done by this class. We typically have a long list of
280
    #  numbers, not very interesing to plot.
281
    #
282
    #  @return `False` to indicate that no (meaningful) plot was created.
283
    def plot(self, varname, figure, ax=None, **kwargs):
284
        self.logger.warning("Plotting for '{0}' is not implemented.".format(self.__class__.__name__))
285
        return False