1
|
|
|
#!/usr/bin/env python |
2
|
|
|
# -*- coding: utf-8 -*- |
3
|
|
|
|
4
|
|
|
# Copyright (c) 2014-2018 Adam.Dybbroe |
5
|
|
|
|
6
|
|
|
# Author(s): |
7
|
|
|
|
8
|
|
|
# Adam.Dybbroe <[email protected]> |
9
|
|
|
|
10
|
|
|
# This program is free software: you can redistribute it and/or modify |
11
|
|
|
# it under the terms of the GNU General Public License as published by |
12
|
|
|
# the Free Software Foundation, either version 3 of the License, or |
13
|
|
|
# (at your option) any later version. |
14
|
|
|
|
15
|
|
|
# This program is distributed in the hope that it will be useful, |
16
|
|
|
# but WITHOUT ANY WARRANTY; without even the implied warranty of |
17
|
|
|
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
18
|
|
|
# GNU General Public License for more details. |
19
|
|
|
|
20
|
|
|
# You should have received a copy of the GNU General Public License |
21
|
|
|
# along with this program. If not, see <http://www.gnu.org/licenses/>. |
22
|
|
|
|
23
|
|
|
"""Testing the radiance to brightness temperature conversion""" |
24
|
|
|
|
25
|
|
|
from pyspectral.radiance_tb_conversion import RadTbConverter |
26
|
|
|
from pyspectral.radiance_tb_conversion import SeviriRadTbConverter |
27
|
|
|
from pyspectral.utils import get_central_wave |
28
|
|
|
import unittest |
29
|
|
|
import numpy as np |
30
|
|
|
from mock import patch |
31
|
|
|
|
32
|
|
|
TEST_TBS = np.array([200., 270., 300., 302., 350.], dtype='float32') |
33
|
|
|
|
34
|
|
|
TRUE_RADS = np.array([856.937353205, 117420.385297, |
35
|
|
|
479464.582505, 521412.9511, 2928735.18944], |
36
|
|
|
dtype='float64') |
37
|
|
|
TRUE_RADS_SEVIRI = np.array([2.391091e-08, |
38
|
|
|
2.559173e-06, |
39
|
|
|
9.797091e-06, |
40
|
|
|
1.061431e-05, |
41
|
|
|
5.531423e-05], |
42
|
|
|
dtype='float64') |
43
|
|
|
|
44
|
|
|
TEST_RSR = {'20': {}} |
45
|
|
|
TEST_RSR['20']['det-1'] = {} |
46
|
|
|
TEST_RSR['20']['det-1']['wavelength'] = np.array([ |
47
|
|
|
3.6123999, 3.6163599, 3.6264927, 3.6363862, 3.646468, |
48
|
|
|
3.6564937, 3.6664478, 3.6765388, 3.6865413, 3.6964585, |
49
|
|
|
3.7065142, 3.716509, 3.7264658, 3.7364102, 3.7463682, |
50
|
|
|
3.7563652, 3.7664226, 3.7763396, 3.7863384, 3.7964207, |
51
|
|
|
3.8063589, 3.8163606, 3.8264089, 3.8364836, 3.8463381, |
52
|
|
|
3.8563975, 3.8664163, 3.8763755, 3.8864797, 3.8964978, |
53
|
|
|
3.9064275, 3.9164873, 3.9264729, 3.9364026, 3.9465107, |
54
|
|
|
3.9535347], dtype='double') |
55
|
|
|
|
56
|
|
|
TEST_RSR['20']['det-1']['response'] = np.array([ |
57
|
|
|
0.01, 0.0118, 0.01987, 0.03226, 0.05028, 0.0849, |
58
|
|
|
0.16645, 0.33792, 0.59106, 0.81815, 0.96077, 0.92855, |
59
|
|
|
0.86008, 0.8661, 0.87697, 0.85412, 0.88922, 0.9541, |
60
|
|
|
0.95687, 0.91037, 0.91058, 0.94256, 0.94719, 0.94808, |
61
|
|
|
1., 0.92676, 0.67429, 0.44715, 0.27762, 0.14852, |
62
|
|
|
0.07141, 0.04151, 0.02925, 0.02085, 0.01414, 0.01], dtype='double') |
63
|
|
|
|
64
|
|
|
TEST_RSR['20']['det-1']['central_wavelength'] = get_central_wave(TEST_RSR['20']['det-1']['wavelength'], |
65
|
|
|
TEST_RSR['20']['det-1']['response']) |
66
|
|
|
|
67
|
|
|
|
68
|
|
|
SEV_RSR = {'IR3.9': {}} |
69
|
|
|
SEV_RSR['IR3.9']['det-1'] = {} |
70
|
|
|
WAVN = np.array([2083.33325195, 2091.00048828, 2098.72387695, 2106.50488281, |
71
|
|
|
2114.34375, 2122.24121094, 2130.19775391, 2138.21435547, |
72
|
|
|
2146.29101562, 2154.42944336, 2162.62963867, 2170.89257812, |
73
|
|
|
2179.21899414, 2187.609375, 2196.06494141, 2204.58569336, |
74
|
|
|
2213.17285156, 2221.82714844, 2230.54956055, 2239.34082031, |
75
|
|
|
2248.20141602, 2257.13256836, 2266.13500977, 2275.20947266, |
76
|
|
|
2284.35668945, 2293.57788086, 2302.87402344, 2312.24584961, |
77
|
|
|
2321.6940918, 2331.21972656, 2340.82421875, 2350.5078125, |
78
|
|
|
2360.27197266, 2370.11743164, 2380.0456543, 2390.05737305, |
79
|
|
|
2400.15356445, 2410.33544922, 2420.60449219, 2430.9609375, |
80
|
|
|
2441.40625, 2451.94189453, 2462.5690918, 2473.28857422, |
81
|
|
|
2484.10180664, 2495.01025391, 2506.01416016, 2517.11645508, |
82
|
|
|
2528.31713867, 2539.61816406, 2551.02050781, 2562.52587891, |
83
|
|
|
2574.13525391, 2585.8503418, 2597.67236328, 2609.60351562, |
84
|
|
|
2621.64453125, 2633.796875, 2646.06274414, 2658.44335938, |
85
|
|
|
2670.94018555, 2683.55541992, 2696.28979492, 2709.14624023, |
86
|
|
|
2722.12548828, 2735.22973633, 2748.4609375, 2761.82055664, |
87
|
|
|
2775.31103516, 2788.93359375, 2802.69042969, 2816.58422852, |
88
|
|
|
2830.61621094, 2844.78833008, 2859.10351562, 2873.56323242, |
89
|
|
|
2888.17016602, 2902.92626953, 2917.83374023, 2932.89550781, |
90
|
|
|
2948.11328125, 2963.48999023, 2979.02758789, 2994.72924805, |
91
|
|
|
3010.59741211, 3026.63452148, 3042.84326172, 3059.22705078, |
92
|
|
|
3075.78735352, 3092.52880859, 3109.45263672, 3126.56323242, |
93
|
|
|
3143.86328125, 3161.35571289, 3179.04394531, 3196.9309082, |
94
|
|
|
3215.02050781, 3233.31640625, 3251.82104492, 3270.5390625], |
95
|
|
|
dtype='float32') |
96
|
|
|
|
97
|
|
|
RESP = np.array([5.85991074e-07, 5.05963471e-05, 1.54738867e-04, |
98
|
|
|
8.75972546e-07, 2.23005936e-05, 6.17855985e-05, |
99
|
|
|
1.41724333e-04, 1.87453145e-06, 3.19355922e-06, |
100
|
|
|
1.08511595e-04, 2.12896630e-04, 5.65914146e-04, |
101
|
|
|
5.93333738e-04, 2.45316158e-04, 1.77410198e-04, |
102
|
|
|
3.18188017e-04, 5.27926895e-05, 1.41405777e-04, |
103
|
|
|
1.64295849e-03, 2.69834511e-03, 4.89762053e-03, |
104
|
|
|
2.71760323e-03, 2.49398337e-03, 4.83754929e-03, |
105
|
|
|
1.08462553e-02, 5.53890038e-03, 8.30772892e-03, |
106
|
|
|
1.33131407e-02, 2.89320182e-02, 4.69624363e-02, |
107
|
|
|
6.85162693e-02, 1.17517754e-01, 2.26854816e-01, |
108
|
|
|
3.69935125e-01, 5.16705751e-01, 6.70479536e-01, |
109
|
|
|
8.18419516e-01, 9.00036395e-01, 9.59491372e-01, |
110
|
|
|
9.60837066e-01, 9.63596582e-01, 9.77563441e-01, |
111
|
|
|
9.98380423e-01, 9.98030603e-01, 9.93735969e-01, |
112
|
|
|
9.84225452e-01, 9.98880267e-01, 1.00000000e+00, |
113
|
|
|
9.90870714e-01, 9.75207090e-01, 9.68391836e-01, |
114
|
|
|
9.73213553e-01, 9.75407243e-01, 9.57278728e-01, |
115
|
|
|
9.68693912e-01, 9.78199899e-01, 9.73649919e-01, |
116
|
|
|
9.81804073e-01, 9.71176386e-01, 9.72167253e-01, |
117
|
|
|
9.60459769e-01, 9.40638900e-01, 9.24033165e-01, |
118
|
|
|
9.16043043e-01, 8.79902899e-01, 8.11953366e-01, |
119
|
|
|
6.69838488e-01, 4.60774124e-01, 2.68200457e-01, |
120
|
|
|
1.34857073e-01, 6.40064552e-02, 3.31763141e-02, |
121
|
|
|
1.24335978e-02, 6.22070907e-03, 2.53354642e-03, |
122
|
|
|
1.81269188e-05, 4.63075470e-03, 1.78873568e-04, |
123
|
|
|
1.01367442e-03, 1.28920563e-03, 4.91134451e-05, |
124
|
|
|
6.77187869e-04, 2.44393433e-03, 2.62995227e-03, |
125
|
|
|
6.38825062e-04, 1.70478446e-03, 1.03909883e-03, |
126
|
|
|
1.27910142e-04, 2.95412028e-04, 8.80619162e-04, |
127
|
|
|
2.42782771e-04, 7.55985593e-06, 3.55220342e-04, |
128
|
|
|
8.71264958e-04, 2.01994626e-04, 8.14358555e-06, |
129
|
|
|
2.14082262e-04, 1.07610082e-04, 5.82974189e-06, |
130
|
|
|
4.16795141e-04], dtype='float32') |
131
|
|
|
|
132
|
|
|
SEV_RSR['IR3.9']['det-1']['wavenumber'] = WAVN |
133
|
|
|
SEV_RSR['IR3.9']['det-1']['response'] = RESP |
134
|
|
|
|
135
|
|
|
|
136
|
|
|
class RSRTestDataModis(object): |
137
|
|
|
|
138
|
|
|
"""RSR test data for Aqua Modis""" |
139
|
|
|
|
140
|
|
|
def __init__(self): |
141
|
|
|
"""Making a testdata set of relative spectral responses""" |
142
|
|
|
|
143
|
|
|
self.rsr = TEST_RSR |
144
|
|
|
|
145
|
|
|
|
146
|
|
|
class TestSeviriConversions(unittest.TestCase): |
147
|
|
|
|
148
|
|
|
"""Testing the conversions between radiances and brightness temperatures""" |
149
|
|
|
|
150
|
|
|
def setUp(self): |
151
|
|
|
"""Set up""" |
152
|
|
|
|
153
|
|
|
with patch('pyspectral.radiance_tb_conversion.RelativeSpectralResponse') as mymock: |
154
|
|
|
instance = mymock.return_value |
155
|
|
|
instance.rsr = SEV_RSR |
156
|
|
|
instance.unit = 'cm-1' |
157
|
|
|
instance.si_scale = 100. |
158
|
|
|
|
159
|
|
|
self.sev1 = RadTbConverter('Meteosat-9', 'seviri', 'IR3.9', |
160
|
|
|
wavespace='wavenumber') |
161
|
|
|
|
162
|
|
|
self.sev2 = SeviriRadTbConverter('Meteosat-9', 'IR3.9') |
163
|
|
|
|
164
|
|
|
def test_rad2tb(self): |
165
|
|
|
"""Unit testing the radiance to brightness temperature conversion""" |
166
|
|
|
res = self.sev1.tb2radiance(TEST_TBS, lut=False) |
167
|
|
|
self.assertTrue(np.allclose(TRUE_RADS_SEVIRI, res['radiance'])) |
168
|
|
|
|
169
|
|
|
def test_conversion_simple(self): |
170
|
|
|
"""Test the tb2radiance function to convert radiances to Tb's |
171
|
|
|
using tabulated coefficients based on a non-linear approximation |
172
|
|
|
|
173
|
|
|
""" |
174
|
|
|
retv = self.sev2.tb2radiance(TEST_TBS) |
175
|
|
|
rads = retv['radiance'] |
176
|
|
|
# Units space = wavenumber (cm-1): |
177
|
|
|
tbs = self.sev2.radiance2tb(rads) |
178
|
|
|
self.assertTrue(np.allclose(TEST_TBS, tbs)) |
179
|
|
|
|
180
|
|
|
np.random.seed() |
181
|
|
|
tbs1 = 200.0 + np.random.random(50) * 150.0 |
182
|
|
|
retv = self.sev2.tb2radiance(tbs1) |
183
|
|
|
rads = retv['radiance'] |
184
|
|
|
tbs = self.sev2.radiance2tb(rads) |
185
|
|
|
self.assertTrue(np.allclose(tbs1, tbs)) |
186
|
|
|
|
187
|
|
|
def test_conversions_methods(self): |
188
|
|
|
"""Using the two diferent conversion methods to verify that they give |
189
|
|
|
approximately the same results. Conversion from Tb's to Radiances |
190
|
|
|
only |
191
|
|
|
""" |
192
|
|
|
# Units space = wavenumber (cm-1): |
193
|
|
|
retv2 = self.sev2.tb2radiance(TEST_TBS) |
194
|
|
|
retv1 = self.sev1.tb2radiance(TEST_TBS) |
195
|
|
|
|
196
|
|
|
rads1 = retv1['radiance'] |
197
|
|
|
rads2 = retv2['radiance'] |
198
|
|
|
self.assertTrue(np.allclose(rads1, rads2)) |
199
|
|
|
|
200
|
|
|
def tearDown(self): |
201
|
|
|
"""Clean up""" |
202
|
|
|
pass |
203
|
|
|
|
204
|
|
|
|
205
|
|
|
class TestRadTbConversions(unittest.TestCase): |
206
|
|
|
|
207
|
|
|
"""Testing the conversions between radiances and brightness temperatures""" |
208
|
|
|
|
209
|
|
|
def setUp(self): |
210
|
|
|
"""Set up""" |
211
|
|
|
# mymock: |
212
|
|
|
with patch('pyspectral.radiance_tb_conversion.RelativeSpectralResponse') as mymock: |
213
|
|
|
instance = mymock.return_value |
214
|
|
|
instance.rsr = TEST_RSR |
215
|
|
|
instance.unit = '1e-6 m' |
216
|
|
|
instance.si_scale = 1e-6 |
217
|
|
|
|
218
|
|
|
self.modis = RadTbConverter('EOS-Aqua', 'modis', '20') |
219
|
|
|
self.modis2 = RadTbConverter('EOS-Aqua', 'modis', 3.75) |
220
|
|
|
|
221
|
|
|
def test_rad2tb(self): |
222
|
|
|
"""Unit testing the radiance to brightness temperature conversion""" |
223
|
|
|
res = self.modis.tb2radiance(TEST_TBS, lut=False) |
224
|
|
|
self.assertTrue(np.allclose(TRUE_RADS, res['radiance'])) |
225
|
|
|
|
226
|
|
|
res = self.modis2.tb2radiance(TEST_TBS, lut=False) |
227
|
|
|
self.assertTrue(np.allclose(TRUE_RADS, res['radiance'])) |
228
|
|
|
|
229
|
|
|
rad = res['radiance'] |
230
|
|
|
tbs = self.modis.radiance2tb(rad) |
231
|
|
|
self.assertTrue(np.allclose(TEST_TBS, tbs, atol=0.25)) |
232
|
|
|
|
233
|
|
|
res = self.modis.tb2radiance(TEST_TBS, lut=False, normalized=False) |
234
|
|
|
integral = self.modis.rsr_integral |
235
|
|
|
self.assertTrue(np.allclose(TRUE_RADS * integral, res['radiance'])) |
236
|
|
|
|
237
|
|
|
res = self.modis.tb2radiance(237., lut=False) |
238
|
|
|
self.assertAlmostEqual(16570.592171157, res['radiance']) |
239
|
|
|
|
240
|
|
|
res = self.modis.tb2radiance(277., lut=False) |
241
|
|
|
self.assertAlmostEqual(167544.3823631, res['radiance']) |
242
|
|
|
|
243
|
|
|
res = self.modis.tb2radiance(1.1, lut=False) |
244
|
|
|
self.assertAlmostEqual(0.0, res['radiance']) |
245
|
|
|
|
246
|
|
|
res = self.modis.tb2radiance(11.1, lut=False) |
247
|
|
|
self.assertAlmostEqual(0.0, res['radiance']) |
248
|
|
|
|
249
|
|
|
res = self.modis.tb2radiance(100.1, lut=False) |
250
|
|
|
self.assertAlmostEqual(5.3940515573e-06, res['radiance']) |
251
|
|
|
|
252
|
|
|
res = self.modis.tb2radiance(200.1, lut=False) |
253
|
|
|
self.assertAlmostEqual(865.09776189, res['radiance']) |
254
|
|
|
|
255
|
|
|
def tearDown(self): |
256
|
|
|
"""Clean up""" |
257
|
|
|
pass |
258
|
|
|
|
259
|
|
|
|
260
|
|
|
def suite(): |
261
|
|
|
"""The suite for test_reflectance.""" |
262
|
|
|
loader = unittest.TestLoader() |
263
|
|
|
mysuite = unittest.TestSuite() |
264
|
|
|
mysuite.addTest(loader.loadTestsFromTestCase(TestRadTbConversions)) |
265
|
|
|
mysuite.addTest(loader.loadTestsFromTestCase(TestSeviriConversions)) |
266
|
|
|
|
267
|
|
|
return mysuite |
268
|
|
|
|