From 71679618409a5c808c66ac16dae6b3d875321c1a Mon Sep 17 00:00:00 2001
From: tomrink <rink@ssec.wisc.edu>
Date: Wed, 31 Mar 2021 14:11:41 -0500
Subject: [PATCH] snapshot...
---
modules/icing/pirep_goes.py | 34 +++++++++++++++++++++++++++++++++-
1 file changed, 33 insertions(+), 1 deletion(-)
diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py
index 122fc9a1..71d685b8 100644
--- a/modules/icing/pirep_goes.py
+++ b/modules/icing/pirep_goes.py
@@ -1,6 +1,7 @@
from icing.pireps import pirep_icing
import numpy as np
import pickle
+import matplotlib.pyplot as plt
import os
from util.util import get_time_tuple_utc, GenericException, add_time_range_to_filename, is_night
from aeolus.datasource import CLAVRx, GOESL1B
@@ -594,6 +595,35 @@ def process_1(ice_dct, no_ice_dct, neg_ice_dct):
return new_ice_dct, new_no_ice_dct, new_neg_ice_dct
+def analyze2(filename, filename_l1b):
+ f = h5py.File(filename, 'r')
+ icing_alt = f['icing_altitude'][:]
+ cld_top_hgt = f['cld_height_acha'][:, 10:30, 10:30]
+ cld_phase = f['cloud_phase'][:, 10:30, 10:30]
+ cld_opd_dc = f['cld_opd_dcomp'][:, 10:30, 10:30]
+ cld_opd = f['cld_opd_acha'][:, 10:30, 10:30]
+ solzen = f['solar_zenith_angle'][:, 10:30, 10:30]
+
+ f_l1b = h5py.File(filename_l1b, 'r')
+ bt_11um = f_l1b['temp_11_0um_nom'][:, 10:30, 10:30]
+
+ cld_opd = cld_opd.flatten()
+ cld_opd_dc = cld_opd_dc.flatten()
+ solzen = solzen.flatten()
+
+ keep1 = np.invert(np.isnan(cld_opd))
+ keep2 = np.invert(np.isnan(solzen))
+ keep = keep1 & keep2
+ cld_opd = cld_opd[np.invert(np.isnan(cld_opd))]
+ cld_opd_dc = cld_opd_dc[keep]
+ solzen = solzen[keep]
+
+ plt.hist(cld_opd, bins=20)
+ plt.show()
+ plt.hist(cld_opd_dc, bins=20)
+ plt.show()
+
+
def run_qc(filename, filename_l1b, day_night='ANY'):
f = h5py.File(filename, 'r')
icing_alt = f['icing_altitude'][:]
@@ -610,7 +640,9 @@ def run_qc(filename, filename_l1b, day_night='ANY'):
bt_11um = f_l1b['temp_11_0um_nom'][:, 10:30, 10:30]
print('num pireps all: ', len(icing_alt))
+
mask, idxs = apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, bt_11um, sol_zen, day_night=day_night)
+
print('num pireps day_night: ', len(mask), day_night)
bts_pass = []
@@ -687,7 +719,7 @@ def run_qc(filename, filename_l1b, day_night='ANY'):
def apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, bt_11um, solzen, day_night='ANY'):
opd_threshold = 2
if day_night == 'DAY':
- opd_threshold = 5
+ opd_threshold = 20
closeness = 100.0 # meters
num_obs = len(icing_alt)
cld_mask = cld_mask.reshape((num_obs, -1))
--
GitLab