From c018caf9db13aaee9e74ea87264c945450a2ec25 Mon Sep 17 00:00:00 2001
From: Paolo Veglio <paolo.veglio@ssec.wisc.edu>
Date: Mon, 14 Nov 2022 22:35:40 +0000
Subject: [PATCH] Surface_Temperature_Test roughly implemented. code runs but
 function and implementation need to be polished... a lot

---
 main.py                          |  4 ++--
 preprocess_thresholds.py         | 39 ++++++++++++++++++++++++++++++++
 tests.py                         |  9 +++++++-
 thresholds.mvcm.snpp.v0.0.1.yaml |  8 ++++++-
 4 files changed, 56 insertions(+), 4 deletions(-)

diff --git a/main.py b/main.py
index e4eb5aa..32e01b3 100644
--- a/main.py
+++ b/main.py
@@ -98,7 +98,7 @@ def main(*, data_path=_datapath, mod02=_fname_mod02, mod03=_fname_mod03,
     perform = {'11um BT Test': False,
                'CO2 High Clouds Test': False,
                'Water Vapor High Clouds Test': False,
-               'Surface Temperature Test': False,
+               'Surface Temperature Test': True,
                'SST Test': False,
                '8.6-11um BT Difference Test': False,
                '11-12um BTD Thin Cirrus Test': False,
@@ -107,7 +107,7 @@ def main(*, data_path=_datapath, mod02=_fname_mod02, mod03=_fname_mod03,
                'Water Vapor Cloud Test': False,
                '11um BT Variability Test': False,
                '11-4um BTD Oceanic Stratus': False,
-               'NIR Reflectance Test': True,
+               'NIR Reflectance Test': False,
                'Vis/NIR Ratio Test': False,
                '1.6um or 2.1um NIR Reflectance Test': False,
                'Visible Reflectance Test': False,
diff --git a/preprocess_thresholds.py b/preprocess_thresholds.py
index 6c95d1d..f4b7e59 100644
--- a/preprocess_thresholds.py
+++ b/preprocess_thresholds.py
@@ -105,6 +105,45 @@ def preproc_nir(data, thresholds, scene):
     return corr_thr
 
 
+def preproc_surf_temp(data, thresholds):
+    thr_sfc1 = thresholds['Surface_Temperature_Test_1']
+    thr_sfc2 = thresholds['Surface_Temperature_Test_2']
+    thr_df1 = thresholds['Surface_Temperature_Test_df1']
+    thr_df2 = thresholds['Surface_Temperature_Test_df2']
+    max_vza = 70.13  # This values is set based on sensor. Check mask_processing_constants.h for MODIS value
+
+    rs = np.prod(data.M15.shape)
+    df1 = (data.M15 - data.M16).values.reshape(rs)
+    df2 = (data.M15 - data.M13).values.reshape(rs)
+    desert_flag = data.Desert.values.reshape(rs)
+    thresh = np.ones((rs, )) * thr_sfc1
+
+    idx = np.where((df1 >= thr_df1[0]) | ((df1 < thr_df1[0]) & ((df2 <= thr_df2[0]) | (df2 >= thr_df2[1]))))
+    thresh[idx] = thr_sfc2
+    idx = np.where(desert_flag == 1)
+    thresh[idx] == thr_sfc1
+
+    midpt = thresh
+    idx = np.where(df1 >= thr_df1[1])
+    midpt[idx] = thresh[idx] + 2.0*df1[idx]
+
+    corr = np.power(data.sensor_zenith.values/max_vza, 4) * 3.0
+    midpt = midpt.reshape(corr.shape) + corr
+    locut = midpt + 2.0
+    hicut = midpt - 2.0
+
+    thr_out = xr.DataArray(data=np.dstack((locut, midpt, hicut, np.ones(locut.shape), np.ones(locut.shape))),
+                           dims=('number_of_lines', 'number_of_pixels', 'z'))
+
+    return thr_out
+
+
+def get_b1_thresholds():
+    # fill_ndvi[0] is for fill_ndvi1
+    # fill_ndvi[1] is for fill_ndvi2
+    pass
+
+
 # NOTE: 11-12um Cirrus Test
 # hicut is computed in different ways depending on the scene
 # 1. midpt - adj
diff --git a/tests.py b/tests.py
index 0508626..0e72173 100644
--- a/tests.py
+++ b/tests.py
@@ -8,6 +8,11 @@ import conf
 import conf_xr
 import preprocess_thresholds as pt
 
+import importlib
+
+# this is used for testing, eventually we want to remove it
+importlib.reload(pt)
+
 
 # ############## GROUP 1 TESTS ############## #
 
@@ -266,7 +271,9 @@ class CloudTests:
             thr_xr['threshold'] = pt.preproc(self.data, self.thresholds[self.scene_name])
             thr = np.ones((5,))  # This is only temporary to force the logic of the code
                                  # I need to find a better solution at some point
-
+        elif test_name == 'Surface_Temperature_Test':
+            thr_xr['threshold'] = pt.preproc_surf_temp(self.data, self.thresholds[self.scene_name])
+            thr = np.ones((5,))
         elif test_name == 'SST_Test':
             thr_xr['threshold'] = (('number_of_lines', 'number_of_pixels', 'z'),
                                    np.ones((self.data[band].shape[0], self.data[band].shape[1], 5))*thr)
diff --git a/thresholds.mvcm.snpp.v0.0.1.yaml b/thresholds.mvcm.snpp.v0.0.1.yaml
index 654002e..cf79141 100644
--- a/thresholds.mvcm.snpp.v0.0.1.yaml
+++ b/thresholds.mvcm.snpp.v0.0.1.yaml
@@ -179,7 +179,13 @@ Polar_Night_Land:
     adj: []    # I NEED TO WORK ON THIS
   pnlbt1         : 270.0
   pnlbt2         : 270.0
-  Surface_Temperature_Test_pfm: 1.0
+  Surface_Temperature_Test_df1: [0.0, 1.0]
+  Surface_Temperature_Test_df2: [-0.5, 1.0]
+  Surface_Temperature_Test_difference: [-0.2, 1.0, -0.5, 1.0]  # <- this merges the previous two arrays
+  Surface_Temperature_Test_1: 20.0  #   | might be worth figuring out if we can
+  Surface_Temperature_Test_2: 12.0  #   | merge these three coefficients
+  Surface_Temperature_Test_pfm: 1.0 # __|
+  Surface_Temperature_Test: [20.0, 12.0, 1.0]  # <- First attempt to merge the three values above
   pnl_11_4_pfm   : 1.0
   pnl_7_11_pfm   : 1.0
   pnl_4_12_pfm   : 1.0
-- 
GitLab