Project

General

Profile

Statistics
| Branch: | Tag: | Revision:

pycama / src / pycama / transform / Filters.py @ 840:a913bf0bc65f

History | View | Annotate | Download (40.1 KB)

1 318:50d30ced7cb7 Maarten
#!/usr/bin/env python3
2
# -*- coding: utf-8 -*-
3
4
# Copyright 2016-2017 Maarten Sneep, KNMI
5 366:87a959e11e51 Maarten
#
6
# Redistribution and use in source and binary forms, with or without
7 318:50d30ced7cb7 Maarten
# modification, are permitted provided that the following conditions are met:
8 366:87a959e11e51 Maarten
#
9
# 1. Redistributions of source code must retain the above copyright notice,
10 318:50d30ced7cb7 Maarten
#    this list of conditions and the following disclaimer.
11
#
12 366:87a959e11e51 Maarten
# 2. Redistributions in binary form must reproduce the above copyright notice,
13
#    this list of conditions and the following disclaimer in the documentation
14 318:50d30ced7cb7 Maarten
#    and/or other materials provided with the distribution.
15
#
16 366:87a959e11e51 Maarten
# 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 318:50d30ced7cb7 Maarten
# SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26
27
## \file Filters.py
28
#  Filter functions for PyCAMA.
29
#
30
#  These subclasses of AbstractTransformer are filters available to PyCAMA by default.
31
#  User filters can be placed in a direcory pointed to by the `PYCAMA_PLUGIN_PATH`
32 366:87a959e11e51 Maarten
#  environment variable (can be a list of directories). These can be accessed
33
#  through the same name space as the standard filters here, i.e. by specifying
34
#  <tt>transform._Classname_</tt> in the configuration.
35 318:50d30ced7cb7 Maarten
# @author Maarten Sneep
36
37 5:3bfdb29d0ed0 maarten
import numpy as np
38
import numexpr
39
40 75:ef24318af323 Maarten
from .Transformer import AbstractTransformer
41 318:50d30ced7cb7 Maarten
from ..utilities import CAMAException
42 5:3bfdb29d0ed0 maarten
43 318:50d30ced7cb7 Maarten
## Unity transformation function for PyCAMA
44
#
45
# No-op filter for PyCAMA.
46 120:ab92d19cd5b6 Maarten
class Unity(AbstractTransformer):
47 318:50d30ced7cb7 Maarten
    ## The actual callable interface
48
    #
49
    # This function simply returns the input array. The arguments are as defined in the superclass.
50 366:87a959e11e51 Maarten
    def __call__(self, primary,
51
                 secondary,
52
                 latitude,
53
                 longitude,
54
                 solar_zenith_angle,
55
                 viewing_zenith_angle,
56
                 solar_azimuth_angle,
57
                 viewing_azimuth_angle,
58
                 processing_quality_flags,
59
                 geolocation_flags,
60 318:50d30ced7cb7 Maarten
                 surface_classification):
61 120:ab92d19cd5b6 Maarten
        return primary
62
63 836:9e8a3379a845 maarten
64
## Select only the central part of the swath in PyCAMA
65
#
66
# Cut off the edges of the swath.
67
class SwathPartSelection(AbstractTransformer):
68
    ## The actual callable interface
69
    #
70
    # This function simply returns the input array. The arguments are as defined in the superclass.
71 838:750fb4349dae maarten
    def __init__(self, **kwargs):
72 836:9e8a3379a845 maarten
        super().__init__()
73
74
        if 'count' in kwargs:
75
            self.count = kwargs['count']
76
        else:
77
            self.count = 200
78
79
    def __call__(self, primary,
80
                 secondary,
81
                 latitude,
82
                 longitude,
83
                 solar_zenith_angle,
84
                 viewing_zenith_angle,
85
                 solar_azimuth_angle,
86
                 viewing_azimuth_angle,
87
                 processing_quality_flags,
88
                 geolocation_flags,
89
                 surface_classification):
90
91 840:a913bf0bc65f Maarten
        n_pixels = primary.shape[1]
92 836:9e8a3379a845 maarten
        offset = (n_pixels - self.count)//2
93
        try:
94
            mask = primary.mask
95
            dummy = mask[0]
96
        except (AttributeError, IndexError):
97
            primary = np.ma.asarray(primary)
98
            mask = np.zeros(primary.shape, dtype=np.bool)
99
100 839:8f07baeb93fe maarten
        mask[:, 0:offset, ...] = True
101
        mask[:, -offset:, ...] = True
102 836:9e8a3379a845 maarten
        primary.mask = mask
103
104
        return primary
105
106 318:50d30ced7cb7 Maarten
## %Unity transformation function for PyCAMA with extra debugging output.
107
#
108
# This no-op filter that prints out extra statistics on the data passing through it.
109 307:3fd65393ffee maarten
class DebugFilter(AbstractTransformer):
110 318:50d30ced7cb7 Maarten
    ## The constructor
111
    #
112
    # @param[in] name A string label for printing, to identify the output if multiple variables use this filter.
113
    # @param[in] min  Count the number of pixels with a value smaller than this minimum.
114
    # @param[in] max  Count the number of pixels with a value larger than this maximum.
115
    def __init__(self, **kwargs):
116 366:87a959e11e51 Maarten
        super().__init__()
117
118 318:50d30ced7cb7 Maarten
        if 'name' in kwargs:
119
            ## The label for printing output.
120
            self.name = kwargs['name']
121
        else:
122
            self.name = ''
123 366:87a959e11e51 Maarten
124 318:50d30ced7cb7 Maarten
        if 'min' in kwargs:
125
            ## The minimum threshold. Report the number of observations that are below this threshold.
126
            self.minimum = kwargs['min']
127
        else:
128
            self.minimum = None
129 366:87a959e11e51 Maarten
130 318:50d30ced7cb7 Maarten
        if 'max' in kwargs:
131
            ## The maximum threshold. Report the number of observations that are above this threshold.
132
            self.maximum = kwargs['max']
133
        else:
134
            self.maximum = None
135 366:87a959e11e51 Maarten
136 318:50d30ced7cb7 Maarten
    ## The actual callable interface
137
    #
138 366:87a959e11e51 Maarten
    # This function simply returns the input array after printing extra statistics
139
    # on the data passing through it. Note that standard out is used directly as the
140
    # logging system isn't available from within the filter sub-system.
141
    #
142 318:50d30ced7cb7 Maarten
    # The arguments are as defined in the superclass.
143
    #
144 366:87a959e11e51 Maarten
    def __call__(self, primary,
145
                 secondary,
146
                 latitude,
147
                 longitude,
148
                 solar_zenith_angle,
149
                 viewing_zenith_angle,
150
                 solar_azimuth_angle,
151
                 viewing_azimuth_angle,
152
                 processing_quality_flags,
153
                 geolocation_flags,
154 318:50d30ced7cb7 Maarten
                 surface_classification):
155 366:87a959e11e51 Maarten
156 307:3fd65393ffee maarten
        flattened = primary.flatten()
157
        if not hasattr(flattened, 'mask'):
158
            flattened = np.ma.asarray(flattened)
159
            flattened.mask = np.zeros(flattened.shape, dtype=np.bool)
160 366:87a959e11e51 Maarten
161 307:3fd65393ffee maarten
        if 'min' in kwargs:
162 318:50d30ced7cb7 Maarten
            select = flattened[np.logical_and(np.logical_not(flattened.mask), flattened < self.minimum)]
163 307:3fd65393ffee maarten
            if len(select) > 0:
164 366:87a959e11e51 Maarten
                print("DebugFilter({3}): {0}/{2} values smaller than minimum {1}".format(len(select),
165 318:50d30ced7cb7 Maarten
                    self.minimum, len(flattened), self.name))
166
        print("DebugFilter({1}: mimimum value {0}".format(np.min(flattened), self.name))
167 366:87a959e11e51 Maarten
168 307:3fd65393ffee maarten
        if 'max' in kwargs:
169 318:50d30ced7cb7 Maarten
            select = flattened[np.logical_and(np.logical_not(flattened.mask), flattened > self.maximum)]
170 307:3fd65393ffee maarten
            if len(select) > 0:
171 366:87a959e11e51 Maarten
                print("DebugFilter({3}): {0}/{2} values larger than maximum {1}".format(len(select),
172 318:50d30ced7cb7 Maarten
                    self.maximum, len(flattened), self.name))
173
        print("DebugFilter({1}): maximum value {0}".format(np.max(flattened), self.name))
174 366:87a959e11e51 Maarten
175 318:50d30ced7cb7 Maarten
        print("DebugFilter({0}): fill value count {1}/{2}".format(self.name, sum(flattened.mask), len(flattened)))
176 307:3fd65393ffee maarten
177
        return primary
178
179 470:ed0c5427a012 Maarten
## a surface filter
180
#
181
#  The filter removes pixels that do not match the set citeria.
182
class SurfaceFilter(AbstractTransformer):
183
    ## The constructor
184
    #
185
    #  @param select what to select: 'land' or 'water'
186
    #  @param source Allow other land/sea masks to be used. Currently accepts 'OMI', 'TROPOMI'/'S5P' (default)
187
    def __init__(self, **kwargs):
188
        super().__init__()
189
190
        if 'select' in kwargs:
191
            self.select = kwargs['select'].lower()
192
        else:
193
            self.select = 'land'
194
195
        if self.select not in ('land', 'water', 'both'):
196
            self.logger.warning("The select argument for SurfaceFilter the contains an illegal value")
197
            self.select = 'land'
198
199
        if 'source' in kwargs:
200
            self.source = kwargs['source'].upper()
201
            if self.source == 'TROPOMI':
202
                self.source = 'S5P'
203
        else:
204
            self.source = 'S5P'
205
206
        if self.source not in ('OMI', 'S5P'):
207
            self.logger.warning("The source argument for SurfaceFilter the contains an illegal value")
208
            self.source = 'S5P'
209
210
    ## land/sea filter for PyCAMA
211
    #
212
    # The arguments are as defined in the superclass.
213
    #
214
    def __call__(self,
215
                 primary,
216
                 secondary,
217
                 latitude,
218
                 longitude,
219
                 solar_zenith_angle,
220
                 viewing_zenith_angle,
221
                 solar_azimuth_angle,
222
                 viewing_azimuth_angle,
223
                 processing_quality_flags,
224
                 geolocation_flags,
225
                 surface_classification):
226
227
        if self.select == "both":
228
            return primary
229
        try:
230
            mask = primary.mask
231
            dummy = mask[0]
232
        except (AttributeError, IndexError):
233
            primary = np.ma.asarray(primary)
234
            mask = np.zeros(primary.shape, dtype=np.bool)
235
236
        if self.source == "S5P":
237
            if self.select == "land":
238
                mask = np.logical_or(mask, (surface_classification & 1) == 1)
239
            else:
240
                mask = np.logical_or(mask, (surface_classification & 1) == 0)
241
242
        elif self.source == "OMI":
243
            if self.select == "land":
244
                mask = np.logical_or(mask, (surface_classification & 0x000f) != 1)
245
            else:
246
                mask = np.logical_or(mask, (surface_classification & 0x000F) == 1)
247
248
        primary.mask = mask
249
250
        return primary
251
252 318:50d30ced7cb7 Maarten
## A threshold filter for PyCAMA
253
#
254
#  The filter removes values on the 'wrong' side of the threshold (single sided).
255 5:3bfdb29d0ed0 maarten
class ThresholdFilter(AbstractTransformer):
256 318:50d30ced7cb7 Maarten
    ## The constructor
257
    #
258
    #  @param threshold The cut off value for the filter
259 366:87a959e11e51 Maarten
    #  @param comparison Binary comparison, resolving to 'True' for elements
260 318:50d30ced7cb7 Maarten
    #    that must be removed, i.e. use '<' to remove values smaller than
261
    #    the threshold. This is a string which is used to build the comparison operation via `numexpr`.
262
    def __init__(self, **kwargs):
263 366:87a959e11e51 Maarten
        super().__init__()
264
265 5:3bfdb29d0ed0 maarten
        if 'threshold' in kwargs:
266 318:50d30ced7cb7 Maarten
            ## Threshold to apply in the filter.
267 5:3bfdb29d0ed0 maarten
            self.threshold = kwargs['threshold']
268
        else:
269 318:50d30ced7cb7 Maarten
            self.threshold = 0.2
270 366:87a959e11e51 Maarten
271 5:3bfdb29d0ed0 maarten
        if 'comparison' in kwargs:
272 318:50d30ced7cb7 Maarten
            ## Comparison symbol for the threshold.
273 5:3bfdb29d0ed0 maarten
            self.comparison = kwargs['comparison']
274
        else:
275 318:50d30ced7cb7 Maarten
            self.comparison = "<"
276 366:87a959e11e51 Maarten
277 318:50d30ced7cb7 Maarten
    ## Threshold filter for PyCAMA
278 366:87a959e11e51 Maarten
    #
279 318:50d30ced7cb7 Maarten
    # The arguments are as defined in the superclass.
280
    #
281 5:3bfdb29d0ed0 maarten
    def __call__(self,
282 366:87a959e11e51 Maarten
                 primary,
283
                 secondary,
284
                 latitude,
285
                 longitude,
286
                 solar_zenith_angle,
287
                 viewing_zenith_angle,
288
                 solar_azimuth_angle,
289
                 viewing_azimuth_angle,
290
                 processing_quality_flags,
291
                 geolocation_flags,
292 318:50d30ced7cb7 Maarten
                 surface_classification):
293 366:87a959e11e51 Maarten
294 318:50d30ced7cb7 Maarten
        # import into current namespace
295
        threshold = self.threshold
296 366:87a959e11e51 Maarten
297 5:3bfdb29d0ed0 maarten
        if secondary is not None:
298 318:50d30ced7cb7 Maarten
            mask = numexpr.evaluate("secondary {0} threshold".format(self.comparison))
299 5:3bfdb29d0ed0 maarten
        else:
300 318:50d30ced7cb7 Maarten
            mask = numexpr.evaluate("primary {0} threshold".format(self.comparison))
301 457:845c33092d0b Maarten
        n_mask = np.sum(mask)
302
        p_mask = np.sum(primary.mask)
303
304
        self.logger.debug("{0}: {1} values masked ({2} before)".format(self.__class__.__name__, n_mask, p_mask))
305 366:87a959e11e51 Maarten
306 110:3587a31f10bb Maarten
        primary.mask = np.logical_or(mask, primary.mask)
307 5:3bfdb29d0ed0 maarten
        return primary
308
309 318:50d30ced7cb7 Maarten
## A band-pass filter for PyCAMA
310
#
311
#  The filter removes values outside of the specified band.
312 312:8fc0f004eddc Maarten
class BandFilter(AbstractTransformer):
313 318:50d30ced7cb7 Maarten
    ## The constructor
314
    #
315
    # @param lower The lower threshold value for the filter. The default is \f$-\infty\f$.
316
    # @param upper The upper threshold value for the filter. The default is \f$+\infty\f$.
317 366:87a959e11e51 Maarten
    # @param parameter parameter on which to filter (string, one of 'primary', 'secondary',
318
    #            'latitude', 'longitude', 'solar_zenith_angle', 'viewing_zenith_angle',
319
    #             'solar_azimuth_angle', 'viewing_azimuth_angle', 'processing_quality_flags',
320
    #             'geolocation_flags', 'surface_classification'.
321
    #             Defaults to 'secondary if given, or primary otherwise).
322 318:50d30ced7cb7 Maarten
    #
323 366:87a959e11e51 Maarten
    # Keep in mind that the filter functions are never called on the base parameters
324 318:50d30ced7cb7 Maarten
    # (basically the list of the parameter keyword argument).
325
    def __init__(self, **kwargs):
326 366:87a959e11e51 Maarten
        super().__init__()
327
328 318:50d30ced7cb7 Maarten
        if 'lower' in kwargs and kwargs['lower'] is not None:
329
            ## Lower limit of the band filter. Values below or equal to this threshold are removed.
330 312:8fc0f004eddc Maarten
            self.lower = kwargs['lower']
331
        else:
332 318:50d30ced7cb7 Maarten
            self.lower = -np.inf
333 366:87a959e11e51 Maarten
334 318:50d30ced7cb7 Maarten
        if 'upper' in kwargs and kwargs['upper'] is not None:
335
            ## Upper limit of the band filter. Values above this threshold are removed.
336 312:8fc0f004eddc Maarten
            self.upper = kwargs['upper']
337
        else:
338 318:50d30ced7cb7 Maarten
            self.upper = np.inf
339 366:87a959e11e51 Maarten
340 318:50d30ced7cb7 Maarten
        if 'parameter' in kwargs and kwargs['parameter'] is not None:
341
            ## The parameter to apply the filter to.
342 312:8fc0f004eddc Maarten
            self.parameter = kwargs['parameter']
343
        else:
344
            self.parameter = None
345 318:50d30ced7cb7 Maarten
346
    ## A band-pass filter for PyCAMA
347 366:87a959e11e51 Maarten
    #
348 318:50d30ced7cb7 Maarten
    # The arguments are as defined in the superclass.
349 312:8fc0f004eddc Maarten
    def __call__(self,
350 366:87a959e11e51 Maarten
                 primary,
351
                 secondary,
352
                 latitude,
353
                 longitude,
354
                 solar_zenith_angle,
355
                 viewing_zenith_angle,
356
                 solar_azimuth_angle,
357
                 viewing_azimuth_angle,
358
                 processing_quality_flags,
359
                 geolocation_flags,
360 318:50d30ced7cb7 Maarten
                 surface_classification):
361 366:87a959e11e51 Maarten
362 318:50d30ced7cb7 Maarten
        lower = self.lower
363
        upper = self.upper
364 366:87a959e11e51 Maarten
365 318:50d30ced7cb7 Maarten
        if self.parameter is not None:
366 312:8fc0f004eddc Maarten
            parameter = self.parameter
367
        else:
368
            if secondary is not None:
369
                parameter = 'secondary'
370
            else:
371
                parameter = 'primary'
372 366:87a959e11e51 Maarten
373
        mask = np.logical_not(np.logical_and(numexpr.evaluate("lower <= {0}".format(parameter)),
374 318:50d30ced7cb7 Maarten
                                             numexpr.evaluate("{0} < upper".format(parameter))))
375 366:87a959e11e51 Maarten
376 312:8fc0f004eddc Maarten
        if hasattr(primary.mask, 'shape'):
377
            primary.mask = np.logical_or(mask, primary.mask)
378
        else:
379
            primary = np.ma.asarray(primary)
380
            primary.mask = mask
381 366:87a959e11e51 Maarten
382 312:8fc0f004eddc Maarten
        return primary
383 5:3bfdb29d0ed0 maarten
384 318:50d30ced7cb7 Maarten
## Filter on (warning) flags.
385
#
386
#  Error flags lead to a fill value anyway, so those are irrelevant.
387 5:3bfdb29d0ed0 maarten
class ProcessingQualityFlagsFilter(AbstractTransformer):
388 318:50d30ced7cb7 Maarten
    ## The constructor
389
    #
390 366:87a959e11e51 Maarten
    #  @param[in] mask  The warning(s) that should be treated as errors.
391
    #         A bitwise AND of the values you want to remove from the `processing_quality_flags`
392 318:50d30ced7cb7 Maarten
    #         variable in the S5P level 2 files. The default is to remove _all_ pixels with a warning.
393
    #
394 366:87a959e11e51 Maarten
    #  Note that errors are automagically filtered out, as the primary
395 318:50d30ced7cb7 Maarten
    #  already contains a fill value for those pixels.
396 366:87a959e11e51 Maarten
    #
397 318:50d30ced7cb7 Maarten
    def __init__(self, **kwargs):
398 366:87a959e11e51 Maarten
        super().__init__()
399
400 318:50d30ced7cb7 Maarten
        if 'mask' in kwargs and kwargs['mask'] is not None:
401
            ## mask, if the bitwise AND of this `mask` and the `processing_quality_flags` is not zero, the ground pixel is removed.
402 5:3bfdb29d0ed0 maarten
            self.mask = kwargs['mask']
403
        else:
404 318:50d30ced7cb7 Maarten
            self.mask = 0xffffff00
405 366:87a959e11e51 Maarten
406 318:50d30ced7cb7 Maarten
    ## Processing quality flags filter for PyCAMA
407
    #
408
    #  The arguments are as defined in the superclass.
409
    #
410 5:3bfdb29d0ed0 maarten
    def __call__(self,
411 366:87a959e11e51 Maarten
                 primary,
412
                 secondary,
413
                 latitude,
414
                 longitude,
415
                 solar_zenith_angle,
416
                 viewing_zenith_angle,
417
                 solar_azimuth_angle,
418
                 viewing_azimuth_angle,
419
                 processing_quality_flags,
420
                 geolocation_flags,
421 318:50d30ced7cb7 Maarten
                 surface_classification):
422 366:87a959e11e51 Maarten
423 292:5f823cbfd88f Maarten
        if processing_quality_flags is not None:
424 318:50d30ced7cb7 Maarten
            outmask = np.asarray((processing_quality_flags & self.mask) != 0, dtype=np.bool)
425 292:5f823cbfd88f Maarten
            primary.mask = np.logical_or(outmask, primary.mask)
426 366:87a959e11e51 Maarten
427 5:3bfdb29d0ed0 maarten
        return primary
428 20:3d6468175881 Maarten
429 318:50d30ced7cb7 Maarten
## Apply a scaling to the variable.
430
#
431
#  Scale a variable by multiplication or division with a constant.
432
#  The `numexpr` module is used for the operation.
433 20:3d6468175881 Maarten
class Multiplier(AbstractTransformer):
434 318:50d30ced7cb7 Maarten
    ## The constructor
435
    #
436
    #  @param[in] scalefactor The constant for the multiplication or division. Defaults to 100.0
437 366:87a959e11e51 Maarten
    #  @param[in] operator    The mathematical operator to use, must be a string. Defaults to '*'.
438 318:50d30ced7cb7 Maarten
    def __init__(self, **kwargs):
439 366:87a959e11e51 Maarten
        super().__init__()
440
441 20:3d6468175881 Maarten
        if 'scalefactor' in kwargs:
442 366:87a959e11e51 Maarten
            ## The constant for the multiplication or division.
443 20:3d6468175881 Maarten
            self.scalefactor = kwargs['scalefactor']
444
        else:
445 318:50d30ced7cb7 Maarten
            self.scalefactor = 100.0
446 366:87a959e11e51 Maarten
447 20:3d6468175881 Maarten
        if 'operator' in kwargs:
448 318:50d30ced7cb7 Maarten
            ## The mathematical operator to use, a string
449 20:3d6468175881 Maarten
            self.operator = kwargs['operator']
450
        else:
451 318:50d30ced7cb7 Maarten
            self.operator = "*"
452 366:87a959e11e51 Maarten
453 318:50d30ced7cb7 Maarten
    ## Multiplication filter.
454
    #
455
    #  The arguments are as defined in the superclass.
456 366:87a959e11e51 Maarten
    #
457 20:3d6468175881 Maarten
    def __call__(self,
458 366:87a959e11e51 Maarten
                 primary,
459
                 secondary,
460
                 latitude,
461
                 longitude,
462
                 solar_zenith_angle,
463
                 viewing_zenith_angle,
464
                 solar_azimuth_angle,
465
                 viewing_azimuth_angle,
466
                 processing_quality_flags,
467
                 geolocation_flags,
468 318:50d30ced7cb7 Maarten
                 surface_classification):
469
470
        scalefactor = self.scalefactor
471 20:3d6468175881 Maarten
        mask = primary.mask
472 110:3587a31f10bb Maarten
        primary = np.where(primary.mask, np.nan, primary)
473 318:50d30ced7cb7 Maarten
        primary = numexpr.evaluate("primary {0} scalefactor".format(self.operator))
474 20:3d6468175881 Maarten
        primary = np.ma.asarray(primary)
475
        primary.mask = mask
476
        return primary
477 87:0431f59fe02f Maarten
478 318:50d30ced7cb7 Maarten
479 366:87a959e11e51 Maarten
## Reduce the dimensionality of a data set by selecting a slice from the full array. This is a PyCAMA filter.
480 318:50d30ced7cb7 Maarten
#
481 366:87a959e11e51 Maarten
#  After the filtering stage, all data sets are supposed to have the same dimensionality,
482 318:50d30ced7cb7 Maarten
#  this is one of the methods to acchieve that.
483 87:0431f59fe02f Maarten
class Select(AbstractTransformer):
484 318:50d30ced7cb7 Maarten
    ## Reduction of dimensionality filter for PyCAMA
485
    #
486
    # @param dimension The dimension (index) to slice. Default = -1 (last dimension).
487
    # @param index     Which index in the slice dimension to keep. Default = 0 (first item).
488
    # @param collapse  Instead of selection use mean along given dimension. This is a boolean flag.
489
    #
490
    def __init__(self, **kwargs):
491 366:87a959e11e51 Maarten
        super().__init__()
492
493 87:0431f59fe02f Maarten
        if 'dimension' in kwargs:
494 318:50d30ced7cb7 Maarten
            ## The dimension in which to extract
495 87:0431f59fe02f Maarten
            self.dimension = kwargs['dimension']
496
        else:
497
            self.dimension = None
498 366:87a959e11e51 Maarten
499 87:0431f59fe02f Maarten
        if 'index' in kwargs:
500 318:50d30ced7cb7 Maarten
            ## the index in the dimension
501 87:0431f59fe02f Maarten
            self.index = kwargs['index']
502
        else:
503
            self.index = None
504 366:87a959e11e51 Maarten
505 87:0431f59fe02f Maarten
        if 'collapse' in kwargs:
506 318:50d30ced7cb7 Maarten
            ## Flag to indicate mean instead of slice
507 87:0431f59fe02f Maarten
            self.collapse = kwargs['collapse']
508
        else:
509
            self.collapse = None
510 366:87a959e11e51 Maarten
511 498:82f6d2bb96aa Maarten
        if 'sum' in kwargs:
512
            ## Flag to indicate mean instead of slice
513
            self.sum = kwargs['sum']
514
        else:
515
            self.sum = None
516
517
518 318:50d30ced7cb7 Maarten
    ## Select filter.
519
    #
520
    #  The arguments are as defined in the superclass.
521 366:87a959e11e51 Maarten
    #
522 87:0431f59fe02f Maarten
    def __call__(self,
523 366:87a959e11e51 Maarten
                 primary,
524
                 secondary,
525
                 latitude,
526
                 longitude,
527
                 solar_zenith_angle,
528
                 viewing_zenith_angle,
529
                 solar_azimuth_angle,
530
                 viewing_azimuth_angle,
531
                 processing_quality_flags,
532
                 geolocation_flags,
533 318:50d30ced7cb7 Maarten
                 surface_classification):
534 366:87a959e11e51 Maarten
535 318:50d30ced7cb7 Maarten
        if self.dimension is not None:
536 87:0431f59fe02f Maarten
            dimension = self.dimension
537
        else:
538
            dimension = -1
539 366:87a959e11e51 Maarten
540 318:50d30ced7cb7 Maarten
        if self.index is not None:
541 87:0431f59fe02f Maarten
            index = self.index
542
        else:
543
            index = 0
544 366:87a959e11e51 Maarten
545 318:50d30ced7cb7 Maarten
        if self.collapse is not None:
546 87:0431f59fe02f Maarten
            collapse = self.collapse
547
        else:
548
            collapse = False
549 366:87a959e11e51 Maarten
550 498:82f6d2bb96aa Maarten
        if self.sum is not None:
551
            do_sum = self.collapse
552
        else:
553
            do_sum = False
554
555 366:87a959e11e51 Maarten
        self.logger.debug("Select filter on dimension {0}, index {1}, {2}".format(dimension, index, "collapse" if collapse else "select"))
556
557 87:0431f59fe02f Maarten
        if collapse:
558 110:3587a31f10bb Maarten
            primary[np.logical_not(primary.mask)] = np.nan
559 87:0431f59fe02f Maarten
            primary = np.nanmean(primary, axis=dimension)
560
            return primary
561 366:87a959e11e51 Maarten
562 498:82f6d2bb96aa Maarten
        if do_sum:
563
            primary[np.logical_not(primary.mask)] = np.nan
564
            primary = np.nansum(primary, axis=dimension)
565
            return primary
566
567 87:0431f59fe02f Maarten
        if dimension == -1:
568 98:0fea8767b99f Maarten
            if len(primary.shape) > 1:
569 87:0431f59fe02f Maarten
                primary = primary[..., index]
570
            else:
571
                primary = primary[index]
572
        elif dimension == -2:
573
            if len(primary.shape) > 2:
574
                primary = primary[..., index, :]
575
            else:
576
                primary = primary[index, :]
577
        elif dimension == -3:
578
            if len(primary.shape) > 3:
579
                primary = primary[..., index, :, :]
580
            else:
581
                primary = primary[index, :, :]
582
        elif dimension == -4:
583
            if len(primary.shape) > 4:
584
                primary = primary[..., index, :, :, :]
585
            else:
586
                primary = primary[index, :, :, :]
587
        elif dimension == 1:
588 98:0fea8767b99f Maarten
            if len(primary.shape) > 1:
589 87:0431f59fe02f Maarten
                primary = primary[index, ...]
590
            else:
591
                primary = primary[index]
592
        elif dimension == 2:
593 98:0fea8767b99f Maarten
            if len(primary.shape) > 2:
594 87:0431f59fe02f Maarten
                primary = primary[:, index, ...]
595
            else:
596
                primary = primary[:, index]
597
        elif dimension == 3:
598 98:0fea8767b99f Maarten
            if len(primary.shape) > 3:
599 87:0431f59fe02f Maarten
                primary = primary[:, :, index, ...]
600
            else:
601
                primary = primary[:, :, index]
602
        elif dimension == 4:
603 98:0fea8767b99f Maarten
            if len(primary.shape) > 4:
604 87:0431f59fe02f Maarten
                primary = primary[:, :, :, index, ...]
605
            else:
606
                primary = primary[:, :, :, index]
607
        elif dimension == 5:
608 98:0fea8767b99f Maarten
            if len(primary.shape) > 5:
609 87:0431f59fe02f Maarten
                primary = primary[:, :, :, :, index, ...]
610
            else:
611
                primary = primary[:, :, :, :, index]
612
        elif dimension == 6:
613 98:0fea8767b99f Maarten
            if len(primary.shape) > 6:
614 87:0431f59fe02f Maarten
                primary = primary[:, :, :, :, :, index, ...]
615
            else:
616
                primary = primary[:, :, :, :, :, index]
617
        elif dimension == 7:
618 98:0fea8767b99f Maarten
            if len(primary.shape) > 7:
619 87:0431f59fe02f Maarten
                primary = primary[:, :, :, :, :, :, index, ...]
620
            else:
621
                primary = primary[:, :, :, :, :, :, index]
622
        elif dimension == 8:
623 98:0fea8767b99f Maarten
            if len(primary.shape) > 8:
624 87:0431f59fe02f Maarten
                primary = primary[:, :, :, :, :, :, :, index, ...]
625
            else:
626
                primary = primary[:, :, :, :, :, :, :, index]
627
        else:
628
            raise RuntimeError("Dimension index '{0}' not supported".format(dimension))
629 366:87a959e11e51 Maarten
630 87:0431f59fe02f Maarten
        return primary
631
632 366:87a959e11e51 Maarten
## Reduce dimensionality by integrating over a vertical column.
633 318:50d30ced7cb7 Maarten
#
634
#  This filter produces an integrated vertical column.
635
#
636 366:87a959e11e51 Maarten
#  \f[
637 318:50d30ced7cb7 Maarten
#  C = \sum_i \frac{10 N_{A} c_i p_i}{g*M}
638
#  \f]
639
#
640 366:87a959e11e51 Maarten
#  with \f$c_i\f$ the vertical profile, \f$p_i\f$ the pressure, both at level \f$i\f$,
641
#  \f$N_A\f$ Avogadro's constant (6.0221367e23), \f$g\f$ the gravitational accelleration
642 318:50d30ced7cb7 Maarten
#  (9.80665), and \f$M\f$ the molar mass of air (28.9644).
643
#
644 366:87a959e11e51 Maarten
#  The pressure levels \f$p_i\f$ use the surface pressure and two constant arrays
645 318:50d30ced7cb7 Maarten
#  \f$A_i\f$ and \f$B_i\f$:
646
#
647
#  \f[
648
#  p_i = \frac{A_i + p_s B_i}{100}
649
#  \f]
650 366:87a959e11e51 Maarten
#
651 318:50d30ced7cb7 Maarten
#  The result of this filter is \f$C\f$.
652
#
653 98:0fea8767b99f Maarten
class IntegratedColumn(AbstractTransformer):
654 318:50d30ced7cb7 Maarten
    ## The constructor
655 366:87a959e11e51 Maarten
    #
656 318:50d30ced7cb7 Maarten
    # @param[in] dimension the dimension (index) over which to integrate. Default = -1 (last dimension).
657 366:87a959e11e51 Maarten
    # @param[in] coefficients_a The constant array \f$A_i\f$.
658 87:0431f59fe02f Maarten
    def __init__(self, *args, **kwargs):
659 366:87a959e11e51 Maarten
        super().__init__()
660
661 87:0431f59fe02f Maarten
        if 'dimension' in kwargs:
662 318:50d30ced7cb7 Maarten
            ## The dimension over which to integrate
663 87:0431f59fe02f Maarten
            self.dimension = kwargs['dimension']
664
        else:
665
            self.dimension = None
666 366:87a959e11e51 Maarten
667 253:6d639d487d3c Maarten
        if 'coefficients_a' in kwargs:
668 318:50d30ced7cb7 Maarten
            ## The A constant (array) for the pressure levels in the profile. After initialization this is a string, a call to the read_additional_variables() method reads the actual data.
669 253:6d639d487d3c Maarten
            self.tm5_constant_a = kwargs['coefficients_a']
670 318:50d30ced7cb7 Maarten
        else:
671
            raise CAMAException("Required argument 'coefficients_a' for the IntegratedColumn filter is missing.")
672 366:87a959e11e51 Maarten
673 253:6d639d487d3c Maarten
        if 'coefficients_b' in kwargs:
674 318:50d30ced7cb7 Maarten
            ## The B constant (array) for the pressure levels in the profile. After initialization this is a string, a call to the read_additional_variables() method reads the actual data.
675 253:6d639d487d3c Maarten
            self.tm5_constant_b = kwargs['coefficients_b']
676 318:50d30ced7cb7 Maarten
        else:
677
            raise CAMAException("Required argument 'coefficients_b' for the IntegratedColumn filter is missing.")
678 366:87a959e11e51 Maarten
679 318:50d30ced7cb7 Maarten
    ## Read additional parameters from an input file
680 366:87a959e11e51 Maarten
    #
681
    #  Some parameters are only known at runtime, or must be read from file
682 318:50d30ced7cb7 Maarten
    #  (other than a normal variable).
683 366:87a959e11e51 Maarten
    #
684 318:50d30ced7cb7 Maarten
    #  Overridden here because we need access to the file to read `tm5_constant_a` and `tm5_constant_b`.
685
    #
686
    # @param f The pycama.File.File object.
687 366:87a959e11e51 Maarten
    #
688 318:50d30ced7cb7 Maarten
    def read_additional_variables(self, f):
689 613:1541b01eb4fa Maarten
        """
690
        Read TM5 a and b constants from file, and calculate the values at the interface.
691
        """
692 318:50d30ced7cb7 Maarten
        if hasattr(self.tm5_constant_a, 'startswith') and self.tm5_constant_a.startswith('read_from_file("'):
693 366:87a959e11e51 Maarten
            varname = self.tm5_constant_a.replace('read_from_file("', '').replace('")', '')
694
            self.logger.debug("Loading additional variable '{0}'".format(varname))
695
            var = f.find_variable(varname)
696 318:50d30ced7cb7 Maarten
            if var is None:
697
                raise CAMAException("Skipping '%s' because the variable could not be found in %s.", self.tm5_constant_a[16:-2], f.product)
698 613:1541b01eb4fa Maarten
            self.tm5_constant_a = var[0, :]
699
            self.tm5_constant_ai = [0.0]
700
            for fa in self.tm5_constant_a:
701
                self.tm5_constant_ai.append(2*fa - self.tm5_constant_ai[-1])
702
            self.tm5_constant_ai[-1] = 0.0
703
            self.tm5_constant_ai = np.asarray(self.tm5_constant_ai)
704
            self.logger.debug("hyai = [{0}]".format(", ".join([str(a) for a in self.tm5_constant_ai])))
705 366:87a959e11e51 Maarten
706 318:50d30ced7cb7 Maarten
        if hasattr(self.tm5_constant_b, 'startswith') and self.tm5_constant_b.startswith('read_from_file("'):
707 366:87a959e11e51 Maarten
            varname = self.tm5_constant_b.replace('read_from_file("', '').replace('")', '')
708
            self.logger.debug("Loading additional variable '{0}'".format(varname))
709
            var = f.find_variable(varname)
710 318:50d30ced7cb7 Maarten
            if var is None:
711
                raise CAMAException("Skipping '%s' because the variable could not be found in %s.", self.tm5_constant_b[16:-2], f.product)
712 613:1541b01eb4fa Maarten
            self.tm5_constant_b = var[0, :]
713
            self.tm5_constant_bi = [1.0]
714
            for fb in self.tm5_constant_b:
715
                self.tm5_constant_bi.append(2*fb - self.tm5_constant_bi[-1])
716
                if self.tm5_constant_bi[-1] < 1e-7: self.tm5_constant_bi[-1] = 0.0
717
            self.tm5_constant_bi = np.asarray(self.tm5_constant_bi)
718
            self.logger.debug("hybi = [{0}]".format(", ".join([str(b) for b in self.tm5_constant_bi])))
719 366:87a959e11e51 Maarten
720 318:50d30ced7cb7 Maarten
    ## The callable interface for the integration filter.
721 366:87a959e11e51 Maarten
    #
722 318:50d30ced7cb7 Maarten
    # Note that the surface pressure _must_ be supplied in the `secondary` parameter.
723
    # Otherwise the arguments are as defined in the superclass.
724
    #
725 398:a542314ec09c Maarten
    # @returns the integrated vertical column in mol/2.
726
    #
727 87:0431f59fe02f Maarten
    def __call__(self,
728 366:87a959e11e51 Maarten
                 primary,
729
                 secondary,
730
                 latitude,
731
                 longitude,
732
                 solar_zenith_angle,
733
                 viewing_zenith_angle,
734
                 solar_azimuth_angle,
735
                 viewing_azimuth_angle,
736
                 processing_quality_flags,
737
                 geolocation_flags,
738 318:50d30ced7cb7 Maarten
                 surface_classification):
739
740
        if self.dimension is not None:
741 87:0431f59fe02f Maarten
            dimension = self.dimension
742
        else:
743
            dimension = -1
744 366:87a959e11e51 Maarten
745 98:0fea8767b99f Maarten
        vcd = np.zeros(primary.shape, dtype=np.float32)
746
        dryair = 28.9644
747
        Ra = 6.0221367e23
748
        g = 9.80665
749 366:87a959e11e51 Maarten
750 412:47ffc0b476b9 Maarten
        primary[primary.mask] = np.nan
751 613:1541b01eb4fa Maarten
        p = lambda l : (self.tm5_constant_ai[l] + secondary*self.tm5_constant_bi[l])/100.0
752 412:47ffc0b476b9 Maarten
        for level in range(primary.shape[dimension]):
753 613:1541b01eb4fa Maarten
            dp = p(level) - p(level+1)
754
            vcd[:, :, level] = (10.0*Ra*primary[:, :, level]*dp)/(g*dryair)
755 366:87a959e11e51 Maarten
756 98:0fea8767b99f Maarten
        primary = np.sum(vcd, axis=dimension)
757 110:3587a31f10bb Maarten
        primary = np.ma.asarray(primary)
758
        primary.mask = np.logical_not(np.isfinite(primary))
759 398:a542314ec09c Maarten
        return primary*1.6605390404271641565e-20
760 87:0431f59fe02f Maarten
761 318:50d30ced7cb7 Maarten
## Translate a sun normalized radiance into a reflectance
762
#
763 366:87a959e11e51 Maarten
#  Or do the reverse. The conversion (either way) uses solar zenith angles up to 89 degrees.
764 318:50d30ced7cb7 Maarten
#  Larger angles are clipped.
765
#
766 87:0431f59fe02f Maarten
class SunNormalizedRadiance(AbstractTransformer):
767 318:50d30ced7cb7 Maarten
    ## The constructor
768
    #
769 366:87a959e11e51 Maarten
    #  @param[in] inverse Do the inverse conversion (from reflectance to sun
770 318:50d30ced7cb7 Maarten
    #         normalized radiance) when true.
771
    def __init__(self, **kwargs):
772 366:87a959e11e51 Maarten
        super().__init__()
773 318:50d30ced7cb7 Maarten
        ## Flag to indicate the direction of the translation.
774
        self.inverse = 'inverse' in kwargs and kwargs['inverse']
775 366:87a959e11e51 Maarten
776 318:50d30ced7cb7 Maarten
    ## The callable interface for the SunNormalizedRadiance translator.
777 366:87a959e11e51 Maarten
    #
778 318:50d30ced7cb7 Maarten
    # The arguments are as defined in the superclass.
779
    #
780 87:0431f59fe02f Maarten
    def __call__(self,
781 366:87a959e11e51 Maarten
                 primary,
782
                 secondary,
783
                 latitude,
784
                 longitude,
785
                 solar_zenith_angle,
786
                 viewing_zenith_angle,
787
                 solar_azimuth_angle,
788
                 viewing_azimuth_angle,
789
                 processing_quality_flags,
790
                 geolocation_flags,
791 318:50d30ced7cb7 Maarten
                 surface_classification):
792 366:87a959e11e51 Maarten
793 606:ede842311422 Maarten
        limit = 88.5
794
795 378:c1abd17e49c7 Maarten
        if self.inverse:
796 312:8fc0f004eddc Maarten
            try:
797 606:ede842311422 Maarten
                mu0 = np.cos(np.where(solar_zenith_angle < limit, solar_zenith_angle, limit)*np.pi/180.0)
798 312:8fc0f004eddc Maarten
            except TypeError:
799
                mu0 = np.pi
800 366:87a959e11e51 Maarten
801 312:8fc0f004eddc Maarten
            return primary*mu0/np.pi
802
        else:
803
            try:
804 606:ede842311422 Maarten
                mu0 = np.cos(np.where(solar_zenith_angle < limit, solar_zenith_angle, limit)*np.pi/180.0)
805 312:8fc0f004eddc Maarten
            except TypeError:
806
                mu0 = np.pi
807 366:87a959e11e51 Maarten
808 312:8fc0f004eddc Maarten
            return np.pi*primary/mu0
809 242:c8a4c8a0d29a Maarten
810 318:50d30ced7cb7 Maarten
811 470:ed0c5427a012 Maarten
## Snow/Ice filter
812
class SnowIceFilter(AbstractTransformer):
813
    ## The constructor
814
    #
815
    #  @param filter What to remove ("snow", "no_snow", default "snow")
816
    #  @param threshold Threshold for the ice fraction before removal (0.2)
817
    def __init__(self, **kwargs):
818
        super().__init__()
819
820
        if 'filter' in kwargs:
821
            self.filter = kwargs['filter']
822
        else:
823
            self.filter = 'snow'
824
825
        if self.filter not in ('snow', 'no_snow', 'none'):
826
            self.logger.warning("Filter must be one of ('snow', 'no_snow', 'none'), resetting to 'snow'")
827
            self.filter = 'snow'
828
829
        if 'threshold' in kwargs:
830
            self.threshold = int(100*kwargs['threshold'])
831
        else:
832
            self.threshold = 20
833
834
        self.OMI_shift = 'OMI_shift' in kwargs and kwargs['OMI_shift']
835
836
    ## The callable
837
    #
838
    # The arguments are as defined in the superclass.
839
    def __call__(self,
840
                 primary,
841
                 secondary,
842
                 latitude,
843
                 longitude,
844
                 solar_zenith_angle,
845
                 viewing_zenith_angle,
846
                 solar_azimuth_angle,
847
                 viewing_azimuth_angle,
848
                 processing_quality_flags,
849
                 geolocation_flags,
850
                 surface_classification):
851
852
        if self.filter == 'none':
853
            return primary
854
855
        try:
856
            mask = primary.mask
857
            dummy = mask.shape[1]
858
        except (AttributeError, IndexError):
859
            primary = np.ma.asarray(primary)
860
            mask = np.zeros(primary.shape, dtype=np.bool)
861
862
        total = np.sum(mask)
863
864
        if self.OMI_shift:
865
            primary = (primary >> 8) & 0x7F
866
        else:
867
            primary = primary & 0x7F
868
869
        if self.filter == 'snow':
870
            mask = np.logical_or(mask, np.logical_and(primary > self.threshold, primary <= 103))
871
            new_total = np.sum(mask)
872
            primary.mask = mask
873
            self.logger.debug("{0}: Removing {1} pixels for snow/ice".format(self.__class__.__name__,
874
                              new_total-total))
875
        else:
876
            mask = np.logical_or(mask, primary < self.threshold)
877
            new_total = np.sum(mask)
878
            primary.mask = mask
879
            self.logger.debug("{0}: Removing {1} pixels for no snow/ice".format(self.__class__.__name__, new_total-total))
880
881
882
        return primary
883
884 318:50d30ced7cb7 Maarten
## Resize the primary array to a new width of the swath.
885
#
886 242:c8a4c8a0d29a Maarten
class Resize(AbstractTransformer):
887 318:50d30ced7cb7 Maarten
    ## The constructor
888
    #
889 366:87a959e11e51 Maarten
    #  @param shape The new width (defaults to None, which will reuse the
890 318:50d30ced7cb7 Maarten
    #               width of the latitude variable).
891
    #  @param skip  The stride in which to sample (default: 1, sample all).
892 366:87a959e11e51 Maarten
    #  @param dimension  The dimension to reduce (default: 1, the `ground_pixel` dimension).
893
    #                    Can be up to 2.
894 318:50d30ced7cb7 Maarten
    def __init__(self, **kwargs):
895 366:87a959e11e51 Maarten
        super().__init__()
896
897 242:c8a4c8a0d29a Maarten
        if 'shape' in kwargs:
898 318:50d30ced7cb7 Maarten
            ## The new number of ground pixels (rows).
899 242:c8a4c8a0d29a Maarten
            self.shape = kwargs['shape']
900
        else:
901 318:50d30ced7cb7 Maarten
            self.shape = None
902 366:87a959e11e51 Maarten
903 242:c8a4c8a0d29a Maarten
        if 'skip' in kwargs:
904 318:50d30ced7cb7 Maarten
            ## The stride to use when filtering
905 242:c8a4c8a0d29a Maarten
            self.skip = kwargs['skip']
906
        else:
907
            self.skip = 1
908 366:87a959e11e51 Maarten
909 242:c8a4c8a0d29a Maarten
        if 'dimension' in kwargs:
910 318:50d30ced7cb7 Maarten
            ## The dimension to reduce.
911 242:c8a4c8a0d29a Maarten
            self.dimension = kwargs['dimension']
912
        else:
913
            self.dimension = 1
914 366:87a959e11e51 Maarten
915 318:50d30ced7cb7 Maarten
    ## The callable
916
    #
917
    # The arguments are as defined in the superclass.
918 242:c8a4c8a0d29a Maarten
    def __call__(self,
919 366:87a959e11e51 Maarten
                 primary,
920
                 secondary,
921
                 latitude,
922
                 longitude,
923
                 solar_zenith_angle,
924
                 viewing_zenith_angle,
925
                 solar_azimuth_angle,
926
                 viewing_azimuth_angle,
927
                 processing_quality_flags,
928
                 geolocation_flags,
929 318:50d30ced7cb7 Maarten
                 surface_classification):
930 366:87a959e11e51 Maarten
931 318:50d30ced7cb7 Maarten
        if self.shape is not None:
932 242:c8a4c8a0d29a Maarten
            shape = self.shape
933
        else:
934
            try:
935
                shape = latitude.shape
936
            except AttributeError:
937
                shape = 60
938 366:87a959e11e51 Maarten
939 318:50d30ced7cb7 Maarten
        if self.skip is not None:
940 242:c8a4c8a0d29a Maarten
            skip = self.skip
941
        else:
942
            skip = 1
943 366:87a959e11e51 Maarten
944 318:50d30ced7cb7 Maarten
        if self.dimension is not None:
945 242:c8a4c8a0d29a Maarten
            dimension = self.dimension
946
        else:
947
            dimension = 1
948 366:87a959e11e51 Maarten
949 242:c8a4c8a0d29a Maarten
        if dimension == 0:
950 514:597b19a4c435 Maarten
            if len(primary.shape) == 1:
951
                return np.repeat(primary.reshape((1, primary.shape[0])), shape, axis=0)
952
            elif primary.shape[0] == 1:
953
                return np.repeat(primary, shape, axis=0)
954
            else:
955
                return primary[0:shape:skip, ...]
956 242:c8a4c8a0d29a Maarten
        elif dimension == 1:
957 514:597b19a4c435 Maarten
            if len(primary.shape) == 1:
958
                return np.repeat(primary.reshape((primary.shape[0], 1)), shape, axis=1)
959
            elif len(primary.shape) == 2:
960
                if primary.shape[1] == 1:
961
                    return np.repeat(primary, shape, axis=1)
962
                else:
963
                    return primary[:, 0:shape:skip]
964 242:c8a4c8a0d29a Maarten
            else:
965 514:597b19a4c435 Maarten
                if primary.shape[1] == 1:
966
                    return np.repeat(primary, shape, axis=1)
967
                else:
968
                    return primary[:, 0:shape:skip, ...]
969 242:c8a4c8a0d29a Maarten
        elif dimension == 2:
970 514:597b19a4c435 Maarten
            if len(primary.shape) < 2:
971
                raise CAMAException("{0}: Source data has not enough dimensions".format(self.__class__.__name__))
972
            elif len(primary.shape) == 2:
973
                return np.repeat(primary.reshape((primary.shape[0], primary.shape[1], 1)), shape, axis=2)
974
            elif len(primary.shape) == 3:
975
                if primary.shape[2] == 1:
976
                    return np.repeat(primary, shape, axis=2)
977
                else:
978
                    return primary[:, :, 0:shape:skip]
979 242:c8a4c8a0d29a Maarten
            else:
980 514:597b19a4c435 Maarten
                if primary.shape[2] == 1:
981
                    return np.repeat(primary, shape, axis=2)
982
                else:
983
                    return primary[:, :, 0:shape:skip, ...]
984 592:601919372ff8 Maarten
985
986
## DLR pixel mapper
987
#
988
#  The cloud product uses both band 3, 4 and band 6.
989
#  This filter maps from band 3 to band 6 or vice versa, using a direct pixel
990
#  mapping scheme. A full colocation filter for PyCAMA is still on the to-do list.
991
#
992
class CloudProductPixelMapper(AbstractTransformer):
993
    ## The constructor
994
    #
995
    #  @param[in] target Target band for the mapping. Defaults to 6.
996
    def __init__(self, **kwargs):
997
        super().__init__()
998
        ## Flag to indicate the direction of the translation.
999
        if 'target' in kwargs and kwargs['target'] and kwargs['target'] in (3, 4, 5, 6):
1000
            self.target = 6 if kwargs['target'] in (5, 6) else 3
1001
        else:
1002
            self.target = 6
1003
1004
    ## The callable interface for the SunNormalizedRadiance translator.
1005
    #
1006
    # The arguments are as defined in the superclass.
1007
    #
1008
    def __call__(self,
1009
                 primary,
1010
                 secondary,
1011
                 latitude,
1012
                 longitude,
1013
                 solar_zenith_angle,
1014
                 viewing_zenith_angle,
1015
                 solar_azimuth_angle,
1016
                 viewing_azimuth_angle,
1017
                 processing_quality_flags,
1018
                 geolocation_flags,
1019
                 surface_classification):
1020
1021 598:e5476ccfbd0c Maarten
        if self.target == 6:
1022 592:601919372ff8 Maarten
            # cut of slice
1023 598:e5476ccfbd0c Maarten
            if primary is not None and primary.shape[0] == 450:
1024
                primary = primary[2:, ...]
1025 592:601919372ff8 Maarten
            if secondary is not None and secondary.shape[0] == 450:
1026
                secondary = secondary[2:, ...]
1027
            if latitude is not None and latitude.shape[0] == 450:
1028
                latitude = latitude[2:, ...]
1029
            if longitude is not None and longitude.shape[0] == 450:
1030
                longitude = longitude[2:, ...]
1031
            if solar_zenith_angle is not None and solar_zenith_angle.shape[0] == 450:
1032
                solar_zenith_angle = solar_zenith_angle[2:, ...]
1033
            if viewing_zenith_angle is not None and viewing_zenith_angle.shape[0] == 450:
1034
                viewing_zenith_angle = viewing_zenith_angle[2:, ...]
1035
            if solar_azimuth_angle is not None and solar_azimuth_angle.shape[0] == 450:
1036
                solar_azimuth_angle = solar_azimuth_angle[2:, ...]
1037
            if viewing_azimuth_angle is not None and viewing_azimuth_angle.shape[0] == 450:
1038
                viewing_azimuth_angle = viewing_azimuth_angle[2:, ...]
1039
            if processing_quality_flags is not None and processing_quality_flags.shape[0] == 450:
1040
                processing_quality_flags = processing_quality_flags[2:, ...]
1041
            if geolocation_flags is not None and geolocation_flags.shape[0] == 450:
1042
                geolocation_flags = geolocation_flags[2:, ...]
1043
            if surface_classification is not None and surface_classification.shape[0] == 450:
1044
                surface_classification = surface_classification[2:, ...]
1045 598:e5476ccfbd0c Maarten
        elif self.target == 3:
1046 592:601919372ff8 Maarten
            # add slice
1047 598:e5476ccfbd0c Maarten
            b = np.zeros(tuple([n if (i != 0) else 2 for i, n in enumerate(primary.shape)]), dtype=np.float64) + np.nan
1048
            if primary is not None and primary.shape[0] == 448:
1049
                primary = np.concatenate((b, primary))
1050 592:601919372ff8 Maarten
            if secondary is not None and secondary.shape[0] == 448:
1051
                secondary = np.concatenate((b, secondary))
1052
            if latitude is not None and latitude.shape[0] == 448:
1053
                latitude = np.concatenate((b, latitude))
1054
            if longitude is not None and longitude.shape[0] == 448:
1055
                longitude = np.concatenate((b, longitude))
1056
            if solar_zenith_angle is not None and solar_zenith_angle.shape[0] == 448:
1057
                solar_zenith_angle = np.concatenate((b, solar_zenith_angle))
1058
            if viewing_zenith_angle is not None and viewing_zenith_angle.shape[0] == 448:
1059
                viewing_zenith_angle = np.concatenate((b, viewing_zenith_angle))
1060
            if solar_azimuth_angle is not None and solar_azimuth_angle.shape[0] == 448:
1061
                solar_azimuth_angle = np.concatenate((b, solar_azimuth_angle))
1062
            if viewing_azimuth_angle is not None and viewing_azimuth_angle.shape[0] == 448:
1063
                viewing_azimuth_angle = np.concatenate((b, viewing_azimuth_angle))
1064
            if processing_quality_flags is not None and processing_quality_flags.shape[0] == 448:
1065
                processing_quality_flags = np.concatenate((b, processing_quality_flags))
1066
            if geolocation_flags is not None and geolocation_flags.shape[0] == 448:
1067
                geolocation_flags = np.concatenate((b, geolocation_flags))
1068
            if surface_classification is not None and surface_classification.shape[0] == 448:
1069
                surface_classification = np.concatenate((b, surface_classification))
1070
1071
        return primary