Skip to content
Snippets Groups Projects
Commit c46a70fa authored by Paolo Veglio's avatar Paolo Veglio
Browse files

results finally look good for Ocean Day using mostly only hires channels (except for 11-12um)

parent 0954b891
No related branches found
No related tags found
No related merge requests found
......@@ -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()
......@@ -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'),
......
......@@ -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'),
......
......@@ -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)
......
......@@ -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,
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment