diff --git a/mvcm/preprocess_thresholds.py b/mvcm/preprocess_thresholds.py index dc0f1d777607451c214ca061ce67c021b96a4178..624e239d9d24b42807492257bef164ad720804e7 100644 --- a/mvcm/preprocess_thresholds.py +++ b/mvcm/preprocess_thresholds.py @@ -751,6 +751,7 @@ def vis_refl_thresholds(data, thresholds, scene, scene_idx): m128[idx] = m02[idx] cosvza = np.cos(data.sensor_zenith.values[scene_idx] * _DTR).ravel() + cosvza[cosvza == 0] = 1.0e-6 vzcpow = thresholds["VZA_correction"]["vzcpow"][0] b1_locut = b1_locut * (1.0 / np.power(cosvza, vzcpow)) diff --git a/mvcm/spectral_tests.py b/mvcm/spectral_tests.py index 86f2f3190b79bec1420d365795ca8818c4e8cb11..cbde8bb7c84d1767f6774a8c9803a2c416629a9b 100644 --- a/mvcm/spectral_tests.py +++ b/mvcm/spectral_tests.py @@ -397,33 +397,53 @@ class CloudTests: return cmin, bits - # THIS NEEDS TO BE TESTED @run_if_test_exists_for_scene("11-4um_BT_Difference_Test_Ocean") def bt_difference_11_4um_test_ocean( self, band: str, cmin: np.ndarray, bits: dict, **kwargs ) -> tuple[np.ndarray, dict]: - """Perfrom 11-4um difference spectral test over ocean.""" + """Perfrom 11-4um difference spectral test over ocean for low emissivity water clouds.""" threshold = kwargs["thresholds"] if threshold["perform"] is True and self.pixels_in_scene is True: bits["qa"][self.scene_idx] = 1 - # print(f'Testing "{self.scene_name}"\n') logger.info(f'Running test 11-4um_BT_Difference_Test_Ocean on "{self.scene_name}"\n') - thr = preproc.bt_diff_11_4um_thresholds(self.data, threshold, self.scene_idx) - # rad = self.data[band].values[self.scene_idx] + if self.scene_name in ["Ocean_Day"]: + thr = threshold["thr"] + elif self.scene_name in ["Ocean_Night"]: + c = threshold["coeffs"] + corr = threshold["corr"] + tpw = np.ravel(self.data.geos_tpw.values[self.scene_idx]) + tmp_thr = c[0] + corr + c[1] * tpw + c[2] * tpw**2 + + hi_thr = tmp_thr + threshold["hicut_coeff"][0] + midpt_thr = hi_thr + threshold["midpt_coeff"][0] + lo_thr = hi_thr + threshold["locut_coeff"][0] + thr = np.vstack([lo_thr, midpt_thr, hi_thr]) + hi_thr = tmp_thr + threshold["hicut_coeff"][1] + midpt_thr = midpt_thr + threshold["midpt_coeff"][1] + lo_thr = lo_thr + threshold["locut_coeff"][1] + thr = np.vstack( + [thr, hi_thr, midpt_thr, lo_thr, np.ones(thr.shape[1]), np.ones(thr.shape[1])] + ) + + else: + logger.error(f"Scene name not accounted for in the code: {self.scene_name}") + raise ValueError + rad = self.data.M15.values[self.scene_idx] - self.data.M12.values[self.scene_idx] - # idx = np.nonzero((rad >= thr[1, :]) & (self.data[self.scene_name] == 1)) - # tmp_bit = kwargs['test_bit'][self.scene_idx] - # tmp_bit[idx] = 1 - # kwargs['test_bit'][self.scene_idx] = tmp_bit - if np.ndim(thr) == 1: - thr = thr.reshape(thr.shape[0], 1) - bits["test"][self.scene_idx] = self.compute_test_bits(rad, thr[1, :], "<") - kwargs["confidence"][self.scene_idx] = conf.conf_test_dble(rad, thr) + if self.scene_name == "Ocean_Day": + bits["test"][self.scene_idx] = self.compute_test_bits(rad, thr[1], "<=") + kwargs["confidence"][self.scene_idx] = conf.conf_test_new(rad, thr) + elif self.scene_name == "Ocean_Night": + bits["test"][self.scene_idx] = np.where((rad >= thr[1]) & (rad <= thr[4]), 0, 1) + kwargs["confidence"][self.scene_idx] = conf.conf_test_dble(rad.ravel(), thr) + else: + logger.error(f"Scene name not accounted for in the code: {self.scene_name}") + raise ValueError cmin = np.fmin(cmin, kwargs["confidence"]) - return cmin, bits # restoral_flag(kwargs) + return cmin, bits # THIS NEEDS TO BE TESTED @run_if_test_exists_for_scene("11-4um_BT_Difference_Test_Land") @@ -461,10 +481,8 @@ class CloudTests: ) & ((sh_lake == 1) | (rad11_12 < 0)) ) - print(f"idx shape: {np.shape(idx)}") thr = preproc.land_night_thresholds(rad11_12[idx], tpw[idx], threshold, coast=False) thr[:3, idx_lake] = thr[:3, idx_lake] + 2 - print(f"thr shape: {thr.shape}") if np.shape(idx[0])[0] > 0: kwargs["confidence"][idx] = conf.conf_test_new(rad11_4[idx], thr) @@ -486,15 +504,12 @@ class CloudTests: & ((sh_lake == 1) | (rad11_12 < 0)) ) ) - print(f"idx shape: {np.shape(idx)}") - print(f"rad shape: {np.shape(rad11_12)}") if np.shape(idx[0])[0] > 0: thr = preproc.land_night_thresholds( rad11_12[idx], tpw[idx], threshold, coast=True ) thr[0:3, idx_lake] = thr[0:3, idx_lake] - 1 thr[3:6, idx_lake] = thr[3:6, idx_lake] + 1 - print(f"thr shape: {thr.shape}") kwargs["confidence"][idx] = conf.conf_test_dble(rad11_4[idx], thr[:, idx]) diff --git a/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml b/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml index 844c7f9555cf5589b9502968fbcfe5320edc6366..887ff1fdebd0ae0d98ff502553244907e41d1252 100644 --- a/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml +++ b/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml @@ -247,6 +247,11 @@ Ocean_Night: perform: True 11-4um_Oceanic_Stratus_Test: thr: [1.25, 1.0, 0.25, 1.0, 1.0] + coeffs: [-0.7093, 0.1128, -0.1567] + corr: 0.0 # no_intadj + locut_coeff: [-4.0, 4.0] # no11_4load + midpt_coeff: [-2.5, 2.5] # no11_4mdad + hicut_coeff: [-1.5, 1.5] # no11_4hiad perform: True Water_Vapor_Cloud_Test: thr: [14.0, 15.0, 16.0, 1.0, 1.0] @@ -256,6 +261,7 @@ Ocean_Night: variability: 0.40 perform: True 11-4um_BT_Difference_Test_Ocean: + thr: [1.25, 1.0, 0.25, 1.0, 1.0] coeffs: [-0.7093, 0.1128, -0.1567] corr: 0.0 # no_intadj locut_coeff: [-4.0, 4.0] # no11_4load @@ -555,6 +561,7 @@ Polar_Night_Ocean: variability: 040 perform: True 11-4um_BT_Difference_Test_Ocean: + thr: [1.25, 1.0, 0.25, 1.0, 1.0] coeffs: [-0.7093, 0.1128, -0.1567] corr: 0.0 # pno_intadj locut_coeff: [-4.0, 4.0] # pno11_4load