diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py index 7cc7fed24c6d1ab2ceb5b0892f59d7548c1f3dc5..c1e72f76fd226e84c7c041d48db8b9ff406aa391 100644 --- a/modules/icing/pirep_goes.py +++ b/modules/icing/pirep_goes.py @@ -686,8 +686,10 @@ def run_qc(filename, filename_l1b, day_night='ANY', pass_thresh_frac=0.20, icing print('num pireps all: ', len(icing_alt)) - mask, idxs, num_tested = apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, bt_11um, sol_zen, day_night=day_night) - + if icing: + mask, idxs, num_tested = apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, bt_11um, sol_zen, day_night=day_night) + else: + mask, idxs, num_tested = apply_qc_no_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) keep_idxs = [] @@ -788,6 +790,59 @@ def apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, return mask, idxs, num_tested +def apply_qc_no_icing_pireps(icing_alt, cld_top_hgt, cld_phase, cld_opd, cld_mask, bt_11um, solzen, day_night='ANY'): + opd_thick_threshold = 2 + if day_night == 'DAY': + opd_thick_threshold = 20 + closeness = 100.0 # meters + num_obs = len(icing_alt) + cld_mask = cld_mask.reshape((num_obs, -1)) + cld_top_hgt = cld_top_hgt.reshape((num_obs, -1)) + cld_phase = cld_phase.reshape((num_obs, -1)) + cld_opd = cld_opd.reshape((num_obs, -1)) + bt_11um = bt_11um.reshape((num_obs, -1)) + + mask = [] + idxs = [] + num_tested = [] + for i in range(num_obs): + if day_night == 'NIGHT': + if is_day(solzen[i,]): + continue + elif day_night == 'DAY': + if is_night(solzen[i,]): + continue + + keep_0 = np.logical_or(cld_mask[i,] == 2, cld_mask[i,] == 3) # cloudy + keep_1 = np.invert(np.isnan(cld_top_hgt[i,])) + keep_2 = np.invert(np.isnan(bt_11um[i,])) + keep_3 = np.invert(np.isnan(cld_opd[i,])) + keep = keep_0 & keep_1 & keep_2 & keep_3 + num_keep = np.sum(keep) + if num_keep == 0: + continue + + keep = np.where(keep, cld_top_hgt[i,] > icing_alt[i], False) + + # keep = np.where(keep, + # np.invert((cld_phase[i,] == 4) & + # np.logical_and(cld_top_hgt[i,]+closeness > icing_alt[i], cld_top_hgt[i,]-closeness < icing_alt[i])), + # False) + # + # keep = np.where(keep, (cld_opd[i,] >= opd_thick_threshold) & (cld_phase[i,] == 4) & (cld_top_hgt[i,] > icing_alt[i]), False) + # + # keep = np.where(keep, np.invert((cld_phase[i,] == 4) & (cld_opd[i,] < 0.1) & (cld_top_hgt[i,] > icing_alt[i])), False) + # + # keep = np.where(keep, np.invert(bt_11um[i,] > 270.0), False) + # + # keep = np.where(keep, np.invert(bt_11um[i,] < 228.0), False) + + mask.append(keep) + idxs.append(i) + num_tested.append(num_keep) + + return mask, idxs, num_tested + def fov_extract(outfile='/home/rink/fovs_l1b_out.h5', train_params=l1b_ds_list, ds_types=l1b_ds_types): ice_times = [] icing_int_s = []