From 955ae585e79702ae3189a22f046e77493c485dad Mon Sep 17 00:00:00 2001
From: Paolo Veglio <paolo.veglio@ssec.wisc.edu>
Date: Mon, 20 Jun 2022 22:13:28 +0000
Subject: [PATCH] initial draft of the function that defines the scene based on
ancillary inputs and geolocation
---
scene.py | 100 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 100 insertions(+)
create mode 100644 scene.py
diff --git a/scene.py b/scene.py
new file mode 100644
index 0000000..8e1ac1a
--- /dev/null
+++ b/scene.py
@@ -0,0 +1,100 @@
+import numpy as np
+
+# lsf: land sea flag
+_scene_list = ['ocean_day', 'ocean_night', 'land_day', 'land_night', 'snow_day', 'coast_day',
+ 'desert_day', 'antarctic_day', 'polar_day_snow', 'polar_day_desert',
+ 'polar_day_desert_coast', 'polar_day_coast', 'polar_day_land', 'polar_night_snow',
+ 'polar_night_land', 'polar_ocean_night']
+_flags = ['day', 'night', 'land', 'coast', 'sh_lake', 'sh_ocean', 'water', 'polar', 'sunglint',
+ 'greenland', 'high_elevation', 'antarctica', 'desert', 'vrused', 'map_snow', 'map_ice',
+ 'ndsi_snow', 'ice', 'new_zealand']
+
+# I'm defining here the flags for difference scenes. Eventually I want to find a better way of doing this
+land = 1
+coast = .2
+sh_lake = .3
+sh_ocean = .4
+water = 5
+polar = 60
+sunglint = 70
+day = 100
+night = 200
+
+
+def find_scene(data):
+ eco = data['eco'].values
+ lsf = data['land_water_mask'].values
+ lat = data['latitude'].values
+ lon = data['longitude'].values
+ sza = data['solar_zenith'].values
+ day = np.zeros((3232, 3200))
+ day[sza <= 85] = 1
+
+ tmp = np.ones((3232, 3200))
+
+ tmp[day == 1] = day
+ tmp[day == 0] = night
+
+ scene = {scn: np.zeros((3232, 3200)) for scn in _scene_list}
+ scene_flag = {flg: np.zeros((3232, 3200)) for flg in _flags}
+
+ scene_flag['day'][sza <= 85] = 1
+ scene_flag['visusd'][sza <= 85] = 1
+ scene_flag['night'][sza > 85] = 1
+
+ scene_flag[np.abs(lat) > 60]['polar'] = 1
+ ################## need to pass refang (once I figure out what it is) and sunglint_angle. The latter
+ # comes from the thresholds file. In the C code is snglnt_bounds[3]
+ idx = np.nonzero((scene_flag['day'] == 1) & (refang <= sunglint_angle))
+ scene_flag['sunglint'][idx] = 1
+
+ ##################################
+ # CONTINUE FROM HERE #
+ ##################################
+
+ # Force consistency between lsf and ecosystem type for water
+ idx = np.nonzero((lsf == 0) | ((lsf >= 5) & (lsf <= 7)))
+ eco[idx] = 14
+
+ idx = np.nonzero((lsf == 1) | (lsf == 4))
+ tmp[idx] = tmp[idx] + land
+
+ idx = np.nonzero(eco == 14)
+
+ # fix-up for missing ecosystem data in eastern Greenland and north-eastern Siberia.
+ # Without this, these regions become completely "coast".
+ idx1 = np.nonzero((eco == 14) & (lat < 64))
+ idx2 = np.nonzero((eco == 14) & (lat >= 67.5) & (lon > -168.5) & (lon < -40))
+ idx3 = np.nonzero((eco == 14) & (lat >= 67.5) & (lon > -30))
+ idx4 = np.nonzero((eco == 14) & (lat >= 64) & (lat < 67.5) & (lon > -168.5) & (lon < -40))
+ idx5 = np.nonzero((eco == 14) & (lat >= 64) & (lat < 67.5) & (lon > -30))
+ tmp[idx1] = tmp[idx1] + coast
+ tmp[idx2] = tmp[idx2] + coast
+ tmp[idx3] = tmp[idx3] + coast
+ tmp[idx4] = tmp[idx4] + coast
+ tmp[idx5] = tmp[idx5] + coast
+
+ tmp[lsf == 2] = tmp[lsf == 2] + land + coast
+ # need shallow lakes to be processed as "coast" for day, but not night
+ tmp[(lsf == 3) & (day == 1)] = tmp[(lsf == 3) & (day == 1)] + land + sh_lake + coast
+ tmp[(lsf > 4) & (lsf != 255)] = tmp[(lsf > 4) & (lsf != 255)] + water + sh_ocean
+
+ # need to cover the case where eco is 255 and I use te vis ratio between gbnd065 and gbnd086
+ # to determine if land or water
+
+ greenland_flag = 0
+ hi_elevation_flag = 0
+ antarctic_flag = 0
+ desert_flag = 0
+ vrused = 0 # visible ratio test used based on eco_type. Check C code
+
+ map_snow_flag = 0
+ map_ice_flag = 0
+ ndsi_snow_flag = 0
+ ice_flag = 0
+ new_zealand_flag = 0
+
+ return scene
+
+
+
--
GitLab