diff --git a/hires_module.py b/hires_module.py index 4b6f815122099f3738ac12071864d86c4ce6c2af..d7ef7e9f0461af4d814daca7e184fdd58c169085 100644 --- a/hires_module.py +++ b/hires_module.py @@ -112,15 +112,22 @@ def main(*, viirs_data['M02'] = (('number_of_lines_2', 'number_of_pixels_2'), viirs_lores.M02.values.repeat(2, 0).repeat(2, 1)) - viirs_data['M07'] = (('number_of_lines_2', 'number_of_pixels_2'), - viirs_lores.M07.values.repeat(2, 0).repeat(2, 1)) + # viirs_data['M07'] = (('number_of_lines_2', 'number_of_pixels_2'), + # viirs_lores.M07.values.repeat(2, 0).repeat(2, 1)) viirs_data['M01'] = (('number_of_lines_2', 'number_of_pixels_2'), np.zeros(viirs_data.I01.shape)) + viirs_data['M07'] = viirs_data.I02 viirs_data['M10'] = viirs_data.I03 viirs_data['M12'] = viirs_data.I04 + viirs_data['M14'] = (('number_of_lines_2', 'number_of_pixels_2'), + viirs_lores.M14.values.repeat(2, 0).repeat(2, 1)) viirs_data['M15'] = (('number_of_lines_2', 'number_of_pixels_2'), viirs_lores.M15.values.repeat(2, 0).repeat(2, 1)) viirs_data['M16'] = (('number_of_lines_2', 'number_of_pixels_2'), viirs_lores.M16.values.repeat(2, 0).repeat(2, 1)) + viirs_data['M15-M16'] = (('number_of_lines_2', 'number_of_pixels_2'), + viirs_data.M15.values - viirs_data.M16.values) + viirs_data['M14-M15'] = (('number_of_lines_2', 'number_of_pixels_2'), + viirs_data.M14.values - viirs_data.M15.values) cmin_G1 = np.ones(viirs_data.M15.shape) cmin_G2 = np.ones(viirs_data.M15.shape) @@ -131,18 +138,20 @@ def main(*, cmin_G1, bit1 = MyScene.test_11um('I05', cmin_G1, test_name='11um_Test') cmin_G1, bit2 = MyScene.sst_test('M15', 'M16', cmin_G1, test_name='SST_Test') + cmin_G2, bit4 = MyScene.test_11_12um_diff('M15-M16', cmin_G2, test_name='11-12um_Cirrus_Test') cmin_G2, bit3 = MyScene.oceanic_stratus_11_4um_test('I05-I04', cmin_G2, test_name='11-4um_Oceanic_Stratus_Test') + # cmin_G2, bit_x = MyScene.bt_diff_86_11um('M14-M15', cmin_G2, test_name='8.6-11um_Test') - cmin_G3, bit4 = MyScene.nir_reflectance_test('M07', cmin_G3, test_name='NIR_Reflectance_Test') + cmin_G3, bit4_1 = MyScene.nir_reflectance_test('M07', cmin_G3, test_name='NIR_Reflectance_Test') cmin_G3, bit5 = MyScene.vis_nir_ratio_test('I02-I01ratio', cmin_G3, test_name='Vis/NIR_Ratio_Test') - visnir, bit5 = MyScene.vis_nir_ratio_test('I02-I01ratio', np.ones(viirs_data.M15.shape), test_name='Vis/NIR_Ratio_Test') + # visnir, bit5 = MyScene.vis_nir_ratio_test('I02-I01ratio', np.ones(viirs_data.M15.shape), test_name='Vis/NIR_Ratio_Test') cmin_G3, bit6 = MyScene.nir_reflectance_test('M10', cmin_G3, test_name='1.6_2.1um_NIR_Reflectance_Test') confidence = np.power(cmin_G1 * cmin_G2 * cmin_G3, 1/3) temp_confidence = np.power(cmin_G1 * cmin_G2 * cmin_G3, 1/3) total_bit = np.zeros(confidence.shape) - total_bit = bit1 + bit2 + total_bit = bit1 + bit2 + bit4 sunglint_angle = thresholds['Sun_Glint']['bounds'][3] scene_flags = scn.find_scene(viirs_data, sunglint_angle) @@ -155,19 +164,18 @@ def main(*, idx = np.nonzero((scene_flags['water'] == 1) & (scene_flags['sunglint'] == 1) & (scene_flags['uniform'] == 1) & (confidence <= 0.95)) confidence[idx] = restoral.sunglint(viirs_data, thresholds['Sun_Glint'], total_bit, temp_confidence)[idx] - # total_bit = bit1 total_confidence = np.sqrt(confidence * mvcm_conf.repeat(2, 0).repeat(2, 1)) date = mod02.split('.')[1] time = mod02.split('.')[2] - # np.savez(f'{data_path}/outputs/pyMVCM_{date}.{time}', - np.savez(f'{data_path}/pyMVCM_{date}.{time}', + # np.savez(f'{data_path}/pyMVCM_{date}.{time}', + np.savez(f'{data_path}/outputs/pyMVCM_{date}.{time}', confidence=total_confidence, py_conf=confidence, mvcm_conf=mvcm_conf, lat=viirs_data.latitude.values, lon=viirs_data.longitude.values, - ocean_day_mask=viirs_data.Ocean_Day.values, visnir=visnir, - c1=cmin_G1, c2=cmin_G2, c3=cmin_G3) + ocean_day_mask=viirs_data.Ocean_Day.values) + # c1=cmin_G1, c2=cmin_G2, c3=cmin_G3) # return confidence # , conf, total_bit @@ -190,4 +198,4 @@ def read_bits(byte, bits): if __name__ == "__main__": - main_fcn() + main() diff --git a/read_data.py b/read_data.py index 70a0dff197704b1dfa9f9f27f7264964b346cd04..24ec026fecf7b2fc9096db09577f3351b054f2bb 100644 --- a/read_data.py +++ b/read_data.py @@ -10,6 +10,7 @@ from typing import Dict _DTR = np.pi/180. _RTD = 180./np.pi +_bad_data = -999.0 def read_data(sensor: str, @@ -163,7 +164,23 @@ def get_data(file_names: Dict[str, str], r3 = r1 / r2 gemi = r3 * (1.0 - 0.25*r3) - ((m01 - 0.125) / (1.0 - m01)) - # Compute channel differences and ratios that are used in the tests + idx = np.nonzero((viirs_data.M05.values < -99) | (viirs_data.M05.values > 2)) + viirs_data['M05'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M07.values < -99) | (viirs_data.M07.values > 2)) + viirs_data['M07'].values[idx] = _bad_data + + idx = np.nonzero((viirs_data.M12.values < 0) | (viirs_data.M12.values > 1000)) + viirs_data['M12'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M13.values < 0) | (viirs_data.M13.values > 1000)) + viirs_data['M13'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M14.values < 0) | (viirs_data.M14.values > 1000)) + viirs_data['M14'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M15.values < 0) | (viirs_data.M15.values > 1000)) + viirs_data['M15'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M16.values < 0) | (viirs_data.M16.values > 1000)) + viirs_data['M16'].values[idx] = _bad_data + + # Compute channel differences and ratios that are used in the tests viirs_data['M14-M15'] = (('number_of_lines', 'number_of_pixels'), viirs_data.M14.values - viirs_data.M15.values) viirs_data['M15-M16'] = (('number_of_lines', 'number_of_pixels'), @@ -183,6 +200,23 @@ def get_data(file_names: Dict[str, str], viirs_data = read_data('viirs', f'{img02}', f'{img03}') viirs_data['M05'] = viirs_data.I01 viirs_data['M07'] = viirs_data.I02 + + idx = np.nonzero((viirs_data.M05.values < -99) | (viirs_data.M05.values > 2)) + viirs_data['M05'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.M07.values < -99) | (viirs_data.M07.values > 2)) + viirs_data['M07'].values[idx] = _bad_data + + idx = np.nonzero((viirs_data.I01.values < 0) | (viirs_data.I01.values > 1000)) + viirs_data['I01'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.I02.values < 0) | (viirs_data.I02.values > 1000)) + viirs_data['I02'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.I03.values < 0) | (viirs_data.I03.values > 1000)) + viirs_data['I03'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.I04.values < 0) | (viirs_data.I04.values > 1000)) + viirs_data['I04'].values[idx] = _bad_data + idx = np.nonzero((viirs_data.I05.values < 0) | (viirs_data.I05.values > 1000)) + viirs_data['I05'].values[idx] = _bad_data + viirs_data = read_ancillary_data(file_names, viirs_data, resolution=2) viirs_data['I05-I04'] = (('number_of_lines_2', 'number_of_pixels_2'), diff --git a/restoral.py b/restoral.py index df3f404ca23d4f0988cbd2023441e3a2c1e2fff9..e89e3b501632a09240da8718c46e8bfae10dd96a 100644 --- a/restoral.py +++ b/restoral.py @@ -22,9 +22,15 @@ def sunglint(viirs_data, threshold, bit, conf): m20 = viirs_data.M12.values m31 = viirs_data.M15.values m02 = viirs_data.M07.values + m01 = viirs_data.M05.values irclr = np.ones(viirs_data.M02.shape) - irclr[bit < 2] = 0 + idx = np.nonzero((m02 > 1.3) & (m01 <= 2.0)) + m02[idx] = 1.3 + idx = np.nonzero((m02 < -99) | (m02 > 1.3)) + m02[idx] = _bad_data + + irclr[bit < 3] = 0 var = spatial_var(m31, 0.4) reg_var_mean = sliding_window_view(np.pad(m02, [1, 1], mode='constant'), diff --git a/sort_inputs.py b/sort_inputs.py index 282fb914e4b8f273be8bac163ad576e5863802e7..1734c2013ac08813f27869db725f30c9a89a5ef1 100644 --- a/sort_inputs.py +++ b/sort_inputs.py @@ -69,10 +69,10 @@ def create_filenames(vnp02mod_fname): def call_mvcm(): - # flist = np.sort(glob('/ships19/hercules/pveglio/mvcm_viirs_hires/VNP02MOD*.nc')) - flist = np.sort(glob('/ships19/hercules/pveglio/mvcm_viirs_hires/VNP02MOD.A2022173.1454*.nc')) + flist = np.sort(glob('/ships19/hercules/pveglio/mvcm_viirs_hires/VNP02MOD*.nc')) + # flist = np.sort(glob('/ships19/hercules/pveglio/mvcm_viirs_hires/VNP02MOD.A2022173.1454*.nc')) - for fname in flist: + for fname in flist[:210]: print(f'Processing {fname}...') fnames_for_mvcm = create_filenames(fname) diff --git a/spectral_tests.py b/spectral_tests.py index b09130a70226f6b09103a83e30322ed03b3e6048..772448beaa4e3a5b6696cff501ea236883cc0724 100644 --- a/spectral_tests.py +++ b/spectral_tests.py @@ -114,7 +114,8 @@ class CloudTests(object): cmin = np.fmin(cmin, confidence) - return cmin, np.abs(1-test_bit)*qa_bit + # return cmin, np.abs(1-test_bit)*qa_bit + return sfcdif, test_bit @run_if_test_exists_for_scene def bt_diff_86_11um(self, @@ -165,7 +166,8 @@ class CloudTests(object): cmin = np.fmin(cmin, confidence) test_bit = test_bit.reshape(self.data[band].shape) - return cmin, np.abs(1-test_bit)*qa_bit + # return cmin, np.abs(1-test_bit)*qa_bit + return cmin, test_bit @run_if_test_exists_for_scene def oceanic_stratus_11_4um_test(self,