From 12709d7e4cf41ea40c4233e75995965f91606da1 Mon Sep 17 00:00:00 2001
From: Paolo Veglio <paolo.veglio@ssec.wisc.edu>
Date: Thu, 19 Oct 2023 14:47:40 +0000
Subject: [PATCH] solved several issues that caused the code to crash
---
mvcm/main.py | 38 +++++-
mvcm/preprocess_thresholds.py | 132 +++++++++++---------
mvcm/read_data.py | 9 +-
mvcm/spectral_tests.py | 36 ++++--
mvcm/utility_functions.py | 2 +-
thresholds/thresholds.mvcm.snpp.v0.0.1.yaml | 12 +-
6 files changed, 148 insertions(+), 81 deletions(-)
diff --git a/mvcm/main.py b/mvcm/main.py
index 08f5f42..999f97d 100644
--- a/mvcm/main.py
+++ b/mvcm/main.py
@@ -73,6 +73,25 @@ _scene_list = [
"Night_Snow",
]
+_mod_bands = [
+ "M01",
+ "M02",
+ "M03",
+ "M04",
+ "M05",
+ "M06",
+ "M07",
+ "M08",
+ "M09",
+ "M10",
+ "M11",
+ "M12",
+ "M13",
+ "M14",
+ "M15",
+ "M16",
+]
+
def timer(func):
"""Compute elapsed time."""
@@ -184,6 +203,16 @@ def main(
use_hires = True
#################
#################
+
+ # for the time being we're not processing night granules
+ # vnp02 = xr.open_dataset(file_names["MOD02"], group="observation_data")
+ # for b in _mod_bands:
+ # try:
+ # vnp02[b]
+ # except KeyError:
+ # logging.info(f"Band {b} not found in file. No output will be written.")
+ # return
+
viirs_data = rd.get_data(satellite, sensor, file_names, hires=use_hires)
sunglint_angle = thresholds["Sun_Glint"]["bounds"][3]
@@ -253,9 +282,10 @@ def main(
# Group 1
cmin_G1, bits["01"] = MyScene.test_11um(m15_name, cmin_G1, bits["01"])
- cmin_G1, bits["02"] = MyScene.surface_temperature_test(
- "", viirs_data, cmin_G1, bits["02"]
- )
+ # this test needs to be implemented properly
+ # cmin_G1, bits["02"] = MyScene.surface_temperature_test(
+ # "", viirs_data, cmin_G1, bits["02"]
+ # )
cmin_G1, bits["03"] = MyScene.sst_test("M15", "M16", cmin_G1, bits["03"])
# Group 2
@@ -268,7 +298,7 @@ def main(
"M15-M12", cmin_G2, bits["07"]
)
cmin_G2, bits["08"] = MyScene.bt_difference_11_4um_test_land(
- "M15-M12", cmin_G2, bits["08"]
+ "M15-M16", "M15-M12", cmin_G2, bits["08"]
)
cmin_G2, bits["09"] = MyScene.oceanic_stratus_11_4um_test(
"M15-M12", cmin_G2, bits["09"]
diff --git a/mvcm/preprocess_thresholds.py b/mvcm/preprocess_thresholds.py
index bdd57b5..99b6bc0 100644
--- a/mvcm/preprocess_thresholds.py
+++ b/mvcm/preprocess_thresholds.py
@@ -1,3 +1,4 @@
+"""Preprcessing thresholds module."""
from typing import Dict
import numpy as np
@@ -12,6 +13,7 @@ _DTR = np.pi / 180
def prepare_11_12um_thresholds(thresholds: np.ndarray, dim1: int) -> Dict:
+ """Prepare 11-12um test thresholds."""
coeff_values = np.empty((dim1, 2))
coeff_values[:, 0] = np.full(dim1, thresholds["coeffs"][0])
coeff_values[:, 1] = np.full(dim1, thresholds["coeffs"][1])
@@ -39,6 +41,7 @@ def prepare_11_12um_thresholds(thresholds: np.ndarray, dim1: int) -> Dict:
def thresholds_11_12um(
data: xr.Dataset, thresholds: np.ndarray, scene: str, scene_idx: np.ndarray
) -> np.ndarray:
+ """Compute 11-12um Test thresholds."""
cosvza = np.cos(data.sensor_zenith.values[scene_idx].ravel() * _DTR)
schi = np.full(cosvza.shape, 99.0)
schi[np.abs(cosvza) > 0] = 1.0 / cosvza
@@ -114,6 +117,7 @@ def thresholds_11_12um(
def thresholds_NIR(data, thresholds, scene, test_name, scene_idx):
+ """Compute NIR test thresholds."""
sza = data.solar_zenith.values[scene_idx].ravel()
# NOTE: the visud condition in the C code is equivalent to having sza <= 85
# For the time being the visud filtering is not implemented
@@ -172,12 +176,14 @@ def thresholds_NIR(data, thresholds, scene, test_name, scene_idx):
def thresholds_surface_temperature(data, thresholds, scene_idx):
+ """Compute thresholds for surface temperature test."""
# def preproc_surf_temp(data, thresholds):
thr_sfc1 = thresholds["desert_thr"]
thr_sfc2 = thresholds["regular_thr"]
thr_df1 = thresholds["channel_diff_11-12um_thr"]
thr_df2 = thresholds["channel_diff_11-4um_thr"]
- max_vza = 70.13 # This values is set based on sensor. Check mask_processing_constants.h for MODIS value
+ max_vza = 70.13 # This values is set based on sensor.
+ # Check mask_processing_constants.h for MODIS value
df1 = (data.M15 - data.M16).values[scene_idx].ravel()
df2 = (data.M15 - data.M13).values[scene_idx].ravel()
@@ -190,7 +196,7 @@ def thresholds_surface_temperature(data, thresholds, scene_idx):
)
thresh[idx] = thr_sfc2
idx = np.where(desert_flag == 1)
- thresh[idx] == thr_sfc1
+ thresh[idx] = thr_sfc1
midpt = thresh
idx = np.where(df1 >= thr_df1[1])
@@ -210,6 +216,7 @@ def thresholds_surface_temperature(data, thresholds, scene_idx):
# This function is currently not used
def preproc_sst(data, thresholds):
+ """Compue SST test thresholds."""
m31c = data.M15 - 273.16
m32c = data.M16 - 273.16
m31c_m32c = m31c - m32c
@@ -231,6 +238,7 @@ def preproc_sst(data, thresholds):
def var_11um(data, thresholds):
+ """Compute thresholds for 11um spatial variability test."""
rad = data.M15.values
var = np.zeros((rad.shape[0], rad.shape[1], 9))
var_thr = thresholds["Daytime_Ocean_Spatial_Variability"]["dovar11"]
@@ -246,6 +254,7 @@ def var_11um(data, thresholds):
def get_b1_thresholds(data, thresholds, scene_idx):
+ """Compute b1 thresholds."""
ndvi = data.ndvi.values[scene_idx].ravel()
sctang = data.scattering_angle.values[scene_idx].ravel()
@@ -328,7 +337,8 @@ def get_b1_thresholds(data, thresholds, scene_idx):
midpt[idx] = -999
locut[idx] = -999
- # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut, np.ones(data.ndvi.shape),
+ # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut,
+ # np.ones(data.ndvi.shape),
# np.full(data.ndvi.shape, 2))),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
#
@@ -337,34 +347,27 @@ def get_b1_thresholds(data, thresholds, scene_idx):
def polar_night_thresholds(data, thresholds, scene, test_name, scene_idx):
- thresholds = thresholds[scene]
+ """Compute thresholds for polar night test."""
+ # thresholds = thresholds[scene]
if (
(test_name == "4-12um_BTD_Thin_Cirrus_Test")
and (scene in ["Land_Night", "Night_Snow"])
or (test_name == "7.3-11um_BTD_Mid_Level_Cloud_Test")
and (scene == "Land_Night")
):
- locut = thresholds[test_name]["thr"][0] * np.ones(
- data.M15.values[scene_idx].ravel()
- )
- midpt = thresholds[test_name]["thr"][1] * np.ones(
- data.M15.values[scene_idx].ravel()
- )
- hicut = thresholds[test_name]["thr"][2] * np.ones(
- data.M15.values[scene_idx].ravel()
- )
- power = thresholds[test_name]["thr"][3] * np.ones(
- data.M15.values[scene_idx].ravel()
- )
-
- # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut, np.ones(data.ndvi.shape), power)),
+ locut = thresholds["thr"][0] * np.ones((len(scene_idx[0]),))
+ midpt = thresholds["thr"][1] * np.ones((len(scene_idx[0]),))
+ hicut = thresholds["thr"][2] * np.ones((len(scene_idx[0]),))
+ power = thresholds["thr"][3] * np.ones((len(scene_idx[0]),))
+ # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut,
+ # np.ones(data.ndvi.shape), power)),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
out_thr = np.dstack((locut, midpt, hicut, np.ones(locut.shape), power))
- return out_thr.T
+ return np.squeeze(out_thr.T)
rad = data.M15.values[scene_idx].ravel()
- bt_bounds = thresholds[test_name]["bt11_bounds"]
+ bt_bounds = thresholds["bt11_bounds"]
locut, midpt = np.empty(rad.shape), np.empty(rad.shape)
hicut, power = np.empty(rad.shape), np.empty(rad.shape)
@@ -373,16 +376,16 @@ def polar_night_thresholds(data, thresholds, scene, test_name, scene_idx):
conf_range = np.empty(rad.shape)
idx = np.nonzero(rad < bt_bounds[0])
- locut[idx] = thresholds[test_name]["low"][0]
- midpt[idx] = thresholds[test_name]["low"][1]
- hicut[idx] = thresholds[test_name]["low"][2]
- power[idx] = thresholds[test_name]["low"][3]
+ locut[idx] = thresholds["low"][0]
+ midpt[idx] = thresholds["low"][1]
+ hicut[idx] = thresholds["low"][2]
+ power[idx] = thresholds["low"][3]
idx = np.nonzero(rad > bt_bounds[3])
- locut[idx] = thresholds[test_name]["high"][0]
- midpt[idx] = thresholds[test_name]["high"][1]
- hicut[idx] = thresholds[test_name]["high"][2]
- power[idx] = thresholds[test_name]["high"][3]
+ locut[idx] = thresholds["high"][0]
+ midpt[idx] = thresholds["high"][1]
+ hicut[idx] = thresholds["high"][2]
+ power[idx] = thresholds["high"][3]
# # # # #
idx = np.nonzero(
@@ -391,36 +394,36 @@ def polar_night_thresholds(data, thresholds, scene, test_name, scene_idx):
& (bt_bounds[1] == 0)
& (bt_bounds[2] == 0)
)
- lo[idx] = thresholds[test_name]["bt11_bounds"][0]
- hi[idx] = thresholds[test_name]["bt11_bounds"][3]
- lo_thr[idx] = thresholds[test_name]["mid1"][0]
- hi_thr[idx] = thresholds[test_name]["mid1"][1]
- power[idx] = thresholds[test_name]["mid1"][3]
- conf_range[idx] = thresholds[test_name]["mid1"][2]
+ lo[idx] = thresholds["bt11_bounds"][0]
+ hi[idx] = thresholds["bt11_bounds"][3]
+ lo_thr[idx] = thresholds["mid1"][0]
+ hi_thr[idx] = thresholds["mid1"][1]
+ power[idx] = thresholds["mid1"][3]
+ conf_range[idx] = thresholds["mid1"][2]
idx = np.nonzero((rad >= bt_bounds[0]) & (rad < bt_bounds[1]))
- lo[idx] = thresholds[test_name]["bt11_bounds"][0]
- hi[idx] = thresholds[test_name]["bt11_bounds"][1]
- lo_thr[idx] = thresholds[test_name]["mid1"][0]
- hi_thr[idx] = thresholds[test_name]["mid1"][1]
- power[idx] = thresholds[test_name]["mid1"][3]
- conf_range[idx] = thresholds[test_name]["mid1"][2]
+ lo[idx] = thresholds["bt11_bounds"][0]
+ hi[idx] = thresholds["bt11_bounds"][1]
+ lo_thr[idx] = thresholds["mid1"][0]
+ hi_thr[idx] = thresholds["mid1"][1]
+ power[idx] = thresholds["mid1"][3]
+ conf_range[idx] = thresholds["mid1"][2]
idx = np.nonzero((rad >= bt_bounds[1]) & (rad < bt_bounds[2]))
- lo[idx] = thresholds[test_name]["bt11_bounds"][1]
- hi[idx] = thresholds[test_name]["bt11_bounds"][2]
- lo_thr[idx] = thresholds[test_name]["mid2"][0]
- hi_thr[idx] = thresholds[test_name]["mid2"][1]
- power[idx] = thresholds[test_name]["mid2"][3]
- conf_range[idx] = thresholds[test_name]["mid2"][2]
+ lo[idx] = thresholds["bt11_bounds"][1]
+ hi[idx] = thresholds["bt11_bounds"][2]
+ lo_thr[idx] = thresholds["mid2"][0]
+ hi_thr[idx] = thresholds["mid2"][1]
+ power[idx] = thresholds["mid2"][3]
+ conf_range[idx] = thresholds["mid2"][2]
idx = np.nonzero((rad >= bt_bounds[2]) & (rad < bt_bounds[3]))
- lo[idx] = thresholds[test_name]["bt11_bounds"][2]
- hi[idx] = thresholds[test_name]["bt11_bounds"][3]
- lo_thr[idx] = thresholds[test_name]["mid3"][0]
- hi_thr[idx] = thresholds[test_name]["mid3"][1]
- power[idx] = thresholds[test_name]["mid3"][3]
- conf_range[idx] = thresholds[test_name]["mid3"][2]
+ lo[idx] = thresholds["bt11_bounds"][2]
+ hi[idx] = thresholds["bt11_bounds"][3]
+ lo_thr[idx] = thresholds["mid3"][0]
+ hi_thr[idx] = thresholds["mid3"][1]
+ power[idx] = thresholds["mid3"][3]
+ conf_range[idx] = thresholds["mid3"][2]
idx = np.nonzero(
((rad >= bt_bounds[0]) & (rad < bt_bounds[3]))
@@ -438,15 +441,17 @@ def polar_night_thresholds(data, thresholds, scene, test_name, scene_idx):
# hicut = hicut.reshape(data.M15.shape)
# power = power.reshape(data.M15.shape)
- # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut, np.ones(data.ndvi.shape), power)),
+ # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut,
+ # np.ones(data.ndvi.shape), power)),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
out_thr = np.dstack((locut, midpt, hicut, np.ones(locut.shape), power))
- return out_thr.T
+ return np.squeeze(out_thr.T)
# get_nl_thresholds
def land_night_thresholds(data, threshold, coast=True):
+ """Compyte land_night_thresholds."""
if coast is False:
lo_val = threshold["bt_diff_bounds"][0]
hi_val = threshold["bt_diff_bounds"][1]
@@ -520,6 +525,7 @@ def land_night_thresholds(data, threshold, coast=True):
def vis_refl_thresholds(data, thresholds, scene, scene_idx):
+ """Compute thresholds for visible reflectance test."""
locut, midpt, hicut = get_b1_thresholds(data, thresholds, scene_idx)
bias_adj = thresholds[scene]["Visible_Reflectance_Test"]["adj"]
ndvi = data.ndvi.values[scene_idx].ravel()
@@ -559,16 +565,19 @@ def vis_refl_thresholds(data, thresholds, scene, scene_idx):
b1_midpt = b1_midpt * (1.0 / np.power(cosvza, vzcpow))
b1_hicut = b1_hicut * (1.0 / np.power(cosvza, vzcpow))
- # out_thr = xr.DataArray(data=np.dstack((b1_locut, b1_midpt, b1_hicut, np.ones(data.ndvi.shape),
+ # out_thr = xr.DataArray(data=np.dstack((b1_locut, b1_midpt,
+ # b1_hicut, np.ones(data.ndvi.shape),
# b1_power.reshape(data.ndvi.shape))),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
- # out_rad = xr.DataArray(data=m128.reshape(data.M01.shape), dims=('number_of_lines', 'number_of_pixels'))
+ # out_rad = xr.DataArray(data=m128.reshape(data.M01.shape),
+ # dims=('number_of_lines', 'number_of_pixels'))
out_thr = np.dstack((b1_locut, b1_midpt, b1_hicut, b1_power))
out_rad = m128
return np.squeeze(out_thr.T), out_rad
def gemi_thresholds(data, thresholds, scene_name, scene_idx):
+ """Compute thresholds for GEMI test."""
ndvi = data.ndvi.values[scene_idx].ravel()
gemi_thr = np.ones((ndvi.shape[0], 5))
@@ -579,7 +588,8 @@ def gemi_thresholds(data, thresholds, scene_name, scene_idx):
idx = np.nonzero((ndvi >= 0.2) & (ndvi < 0.3))
gemi_thr[idx, :3] = thresholds["gemi2"][:3]
- # thr_out = xr.DataArray(data=np.dstack((gemi_thr[:, :, 0], gemi_thr[:, :, 1], gemi_thr[:, :, 2],
+ # thr_out = xr.DataArray(data=np.dstack((gemi_thr[:, :, 0],
+ # gemi_thr[:, :, 1], gemi_thr[:, :, 2],
# np.ones(gemi_thr[:, :, 0].shape),
# np.ones(gemi_thr[:, :, 0].shape))),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
@@ -588,6 +598,7 @@ def gemi_thresholds(data, thresholds, scene_name, scene_idx):
def bt_diff_11_4um_thresholds(data, threshold, scene_idx):
+ """Compute thresholds for 11-4um test."""
c = threshold["coeffs"]
tpw = data.geos_tpw.values[scene_idx].ravel()
@@ -616,6 +627,7 @@ def bt_diff_11_4um_thresholds(data, threshold, scene_idx):
def thresholds_1_38um_test(data, thresholds, scene_name, scene_idx):
+ """Compute thresholds for 1.38um test."""
sza = data.solar_zenith.values[scene_idx]
vza = data.sensor_zenith.values[scene_idx]
thresh = thresholds[scene_name]["1.38um_High_Cloud_Test"]
@@ -642,7 +654,8 @@ def thresholds_1_38um_test(data, thresholds, scene_name, scene_idx):
midpt = midpt * (1 / np.power(cosvza, vzcpow))
hicut = hicut * (1 / np.power(cosvza, vzcpow))
- # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut, np.ones(data.ndvi.shape),
+ # out_thr = xr.DataArray(data=np.dstack((locut, midpt, hicut,
+ # np.ones(data.ndvi.shape),
# np.ones(data.ndvi.shape))),
# dims=('number_of_lines', 'number_of_pixels', 'z'))
out_thr = np.dstack(
@@ -679,7 +692,8 @@ def thresholds_1_38um_test(data, thresholds, scene_name, scene_idx):
# - Night_Snow
# NOTE: 1.38um High Cloud Test
-# thresholds are not always computed the same way. In group 1 there's no preprocessing required,
+# thresholds are not always computed the same way. In group 1 there's no a
+# preprocessing required,
# in group 2 some calcuations are needed
# 1.
# - Land_Day
diff --git a/mvcm/read_data.py b/mvcm/read_data.py
index ec404c9..c070961 100644
--- a/mvcm/read_data.py
+++ b/mvcm/read_data.py
@@ -276,7 +276,6 @@ class ReadData(CollectInputs):
rad_data = xr.Dataset()
for band in list(l1b_data.variables):
- # if "reflectance" in l1b_data[band].long_name:
if band in _reflectance_bands:
if hasattr(l1b_data[band], "VCST_scale_factor"):
scale_factor = (
@@ -291,7 +290,6 @@ class ReadData(CollectInputs):
)
elif band in _emissive_bands:
- # elif "radiance" in l1b_data[band].long_name:
bt_lut = f"{band}_brightness_temperature_lut"
rad_data[band] = (
self.dims,
@@ -323,8 +321,13 @@ class ReadData(CollectInputs):
viirs_out: xr.Dataset
VIIRS L1b data plus band combinations for the spectral tests
"""
+ viirs_out = xr.Dataset()
if hires_data is None:
- viirs_out = viirs
+ for band in _mod_bands:
+ if band not in viirs:
+ viirs_out[band] = viirs.M15 * 0 + _bad_data
+ else:
+ viirs_out[band] = viirs[band]
if hires_data is not None:
viirs_out = xr.Dataset()
diff --git a/mvcm/spectral_tests.py b/mvcm/spectral_tests.py
index 9049e33..2f7760e 100644
--- a/mvcm/spectral_tests.py
+++ b/mvcm/spectral_tests.py
@@ -18,7 +18,7 @@ import mvcm.scene as scn
# import mvcm.utility_functions as utils
-
+_bad_data = -999.0
_RTD = 180.0 / np.pi
_DTR = np.pi / 180
_scene_list = [
@@ -224,6 +224,8 @@ class CloudTests(object):
# idx = np.nonzero((sfcdif < thr[1, :]) & (self.data[self.scene_name] == 1))
# kwargs['test_bit'][idx] = 1
+ if np.ndim(thr) == 1:
+ thr = thr.reshape(thr.shape[0], 1)
bits["test"][self.scene_idx] = self.compute_test_bits(
sfcdif, thr[:, 1], "<"
)
@@ -325,9 +327,12 @@ class CloudTests(object):
comparison = "<"
else:
comparison = "<="
+ if np.ndim(thr) == 1:
+ thr = thr.reshape(thr.shape[0], 1)
bits["test"][self.scene_idx] = self.compute_test_bits(
rad, thr[1, :], comparison
)
+
kwargs["confidence"][self.scene_idx] = conf.conf_test_new(rad, thr)
cmin = np.fmin(cmin, kwargs["confidence"])
@@ -357,6 +362,8 @@ class CloudTests(object):
# 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)
@@ -405,15 +412,17 @@ class CloudTests(object):
rad11_12 = self.data[band1].values
rad11_4 = self.data[band2].values
- thr = preproc.get_nl_thresholds(self.data, threshold, coast=False)
+ thr = preproc.land_night_thresholds(self.data, threshold, coast=False)
idx_lake = np.where((sh_lake == 1) | (rad11_12 < 0))
# test_bit[self.scene_idx] = 1 # THIS NEEDS TO BE CALCULATED PROPERLY
thr[:3, idx_lake] = thr[:3, idx_lake] + 2
temp_thr = np.squeeze(thr[:, idx[0]])
temp_rad = rad11_4[idx]
- bits["test"][idx] = self.compute_test_bits(
- temp_rad[idx], thr[1, idx], "<="
- )
+ # if np.ndim(thr) == 1:
+ # thr = thr.reshape(thr.shape[0], 1)
+ # bits["test"][idx] = self.compute_test_bits(
+ # temp_rad[idx], thr[1, idx], "<="
+ # )
kwargs["confidence"][idx] = conf.conf_test_new(temp_rad, temp_thr)
thr = preproc.land_night_thresholds(self.data, threshold, coast=True)
@@ -636,6 +645,8 @@ class CloudTests(object):
# 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_new(rad, thr)
@@ -769,6 +780,8 @@ class CloudTests(object):
# tmp_bit = kwargs['test_bit'][self.scene_name]
# 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_new(rad, thr)
@@ -799,6 +812,8 @@ class CloudTests(object):
# tmp_bit = kwargs['test_bit'][self.scene_name]
# 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"][scene_idx] = self.compute_test_bits(
rad, thr[1, :], ">", scene_idx=scene_idx
)
@@ -844,6 +859,8 @@ class CloudTests(object):
# 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"][scene_idx] = self.compute_test_bits(
rad, thr[1, :], "<=", scene_idx=scene_idx
)
@@ -864,7 +881,7 @@ class CloudTests(object):
bits["qa"][self.scene_idx] = 1
thr = preproc.polar_night_thresholds(
self.data,
- self.thresholds,
+ threshold,
self.scene_name,
"4-12um_BTD_Thin_Cirrus_Test",
self.scene_idx,
@@ -873,13 +890,16 @@ class CloudTests(object):
logger.info(
f'Running test 4-12um_BTD_Thin_Cirrus_Test on "{self.scene_name}"\n'
)
-
# idx = np.nonzero((rad <= thr[1, :]) &
# (self.data[self.scene_name].values[self.scene_idx] == 1))
# tmp_bit = kwargs['test_bit'][self.scene_idx]
# tmp_bit[idx] = 1
# kwargs['test_bit'][self.scene_idx] = tmp_bit
- bits["test"][self.scene_idx] = self.compute_test_bits(rad, thr[1, :], "<=")
+ # 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_new(rad, thr)
cmin = np.fmin(cmin, kwargs["confidence"])
diff --git a/mvcm/utility_functions.py b/mvcm/utility_functions.py
index 63033e7..ea817ed 100644
--- a/mvcm/utility_functions.py
+++ b/mvcm/utility_functions.py
@@ -98,4 +98,4 @@ def combine_indices(index1: tuple, index2: tuple, shape: tuple) -> tuple:
intersect_idx = list(set(idx1).intersection(set(idx2)))
- return tuple(np.unravel_index(intersect_idx, shape))
+ return tuple(np.unravel_index(np.array(intersect_idx, dtype=np.intp), shape))
diff --git a/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml b/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml
index e9e2195..c641564 100644
--- a/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml
+++ b/thresholds/thresholds.mvcm.snpp.v0.0.1.yaml
@@ -296,7 +296,7 @@ Polar_Night_Land:
mid2: [0.00, 0.00, 0.00, 0.0]
mid3: [0.00, 0.00, 0.00, 0.0]
high: [1.00, 0.70, 0.40, 1.0]
- bt11_bounds: [235, 0, 0, 265]
+ bt11_bounds: [235, 0.0, 0.0, 265]
perform: True
7.3-11um_BTD_Mid_Level_Cloud_Test:
low: [-1.00, 0.00, 1.00, 1.0, 1.0] # pn_7_11l
@@ -311,7 +311,7 @@ Polar_Night_Land:
mid2: [0.00, 0.00, 0.00, 0.0]
mid3: [0.00, 0.00, 0.00, 0.0]
high: [1.00, 0.70, 0.40, 1.0]
- bt11_bounds: [235.0, 000.0, 000.0, 265.0]
+ bt11_bounds: [235.0, 0.0, 0.0, 265.0]
perform: True
# pnl_11_4_pfm : 1.0
# pnl_7_11_pfm : 1.0
@@ -345,7 +345,8 @@ Polar_Day_Desert:
11-4um_Oceanic_Stratus_Test:
# thr: [-22.0, -20.0, -18.0, -2.0, 0.0, 2.0, 1.0, 1.0]
thr: [2.0, 0.0, -2.0, -18.0, -20.0, -22.0, 1.0, 1.0]
- perform: True
+ perform: True # I disable this test while I figure out if I need
+ bt_cutoff: 320.0 # the bt_cutoff parameter
Water_Vapor_High_Clouds_Test: [215.0, 220.0, 225.0, 1.0, 1.0]
Visible_Reflectance_Test:
thr: [0.326, 0.288, 0.250, 1.00, 1.0]
@@ -360,7 +361,6 @@ Polar_Day_Desert:
gemi2: [0.310, 0.335, 0.360, 1.00]
perform: True
pds_ref3_tpw : -10.00
- pdsbt1 : 320.0
Polar_Day_Desert_Coast:
11-12um_Cirrus_Test:
@@ -371,7 +371,8 @@ Polar_Day_Desert_Coast:
11-4um_Oceanic_Stratus_Test:
# thr: [-23.0, -21.0, -19.0, -2.0, 0.0, 2.0, 1.0, 1.0]
thr: [2.0, 0.0, -2.0, -19.0, -21.0, -23.0, 1.0, 1.0]
- perform: True
+ perform: True # see Polar_Day_Desert
+ bt_cutoff: 320.0
Water_Vapor_High_Clouds_Test: [215.0, 220.0, 225.0, 1.0, 1.0]
Visible_Reflectance_Test:
thr: [0.326, 0.288, 0.250, 1.00, 1.0]
@@ -381,7 +382,6 @@ Polar_Day_Desert_Coast:
thr: [0.0375, 0.0250, 0.0125, 1.00, 1.0]
perform: True
pds_ref3_tpw_c : -10.00
- pdsbt1_c : 320.0
Polar_Day_Snow:
11-12um_Cirrus_Test:
--
GitLab