diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py index 4fc488ed9475cb9f4ade15cb994f2f2b3a697ce6..bc09c855803a04e93c60cc21c4255093bb9fc053 100644 --- a/modules/icing/pirep_goes.py +++ b/modules/icing/pirep_goes.py @@ -567,70 +567,6 @@ def pirep_info(pirep_dct): return flt_lvl_s, ice_intensity_s, lat_s, lon_s -# Keep for example -# def analyze(ice_dct, no_ice_dct): -# -# last_file = None -# ice_files = [] -# ice_times = [] -# for ts in list(ice_dct.keys()): -# try: -# ds = get_goes_datasource(ts) -# goes_file, t_0, _ = ds.get_file(ts) -# if goes_file is not None and goes_file != last_file: -# ice_files.append(goes_file) -# ice_times.append(t_0) -# last_file = goes_file -# except Exception: -# continue -# -# last_file = None -# no_ice_files = [] -# no_ice_times = [] -# for ts in list(no_ice_dct.keys()): -# try: -# ds = get_goes_datasource(ts) -# goes_file, t_0, _ = ds.get_file(ts) -# if goes_file is not None and goes_file != last_file: -# no_ice_files.append(goes_file) -# no_ice_times.append(t_0) -# last_file = goes_file -# except Exception: -# continue -# -# ice_times = np.array(ice_times) -# no_ice_times = np.array(no_ice_times) -# -# itrsct_vals, comm1, comm2 = np.intersect1d(no_ice_times, ice_times, return_indices=True) -# -# ice_indexes = np.arange(len(ice_times)) -# -# ucomm2 = np.setxor1d(comm2, ice_indexes) -# np.random.seed(42) -# np.random.shuffle(ucomm2) -# ucomm2 = ucomm2[0:8000] -# -# files_comm = [] -# for i in comm2: -# files_comm.append(ice_files[i]) -# -# files_extra = [] -# times_extra = [] -# for i in ucomm2: -# files_extra.append(ice_files[i]) -# times_extra.append(ice_times[i]) -# -# files = files_comm + files_extra -# times = itrsct_vals.tolist() + times_extra -# times = np.array(times) -# -# sidxs = np.argsort(times) -# for i in sidxs: -# filename = os.path.split(files[i])[1] -# so = re.search('_s\\d{11}', filename) -# dt_str = so.group() -# print(dt_str[2:]) - lon_space_hdeg = np.linspace(-180, 180, 721) lat_space_hdeg = np.linspace(-90, 90, 361) @@ -1139,141 +1075,6 @@ def run_qc(filename, filename_l1b, day_night='ANY', pass_thresh_frac=0.20, cloud return len(icing_alt), len(keep_idxs) -# def apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_geo_dz, cld_phase, cld_opd, cld_mask, bt_11um, solzen, satzen, day_night='ANY'): -# -# if day_night == 'DAY': -# opd_thick_threshold = 20 -# opd_thin_threshold = 1 -# elif day_night == 'NIGHT' or day_night == 'ANY': -# opd_thick_threshold = 2 -# opd_thin_threshold = 0.1 -# -# closeness_top = 100.0 # meters -# closeness_bot = 100.0 -# max_depth = 3000.0 -# -# num_obs = len(icing_alt) -# cld_mask = cld_mask.reshape((num_obs, -1)) -# cld_top_hgt = cld_top_hgt.reshape((num_obs, -1)) -# cld_geo_dz = cld_geo_dz.reshape((num_obs, -1)) -# bt_11um = bt_11um.reshape((num_obs, -1)) -# -# mask = [] -# idxs = [] -# num_tested = [] -# for i in range(num_obs): -# if not check_oblique(satzen[i,]): -# continue -# if day_night == 'NIGHT' and not is_night(solzen[i,]): -# continue -# elif day_night == 'DAY' and not is_day(solzen[i,]): -# continue -# elif day_night == 'ANY': -# pass -# # if not (is_day(solzen[i,]) or 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_geo_dz[i,])) -# keep = keep_0 & keep_1 & keep_2 & keep_3 -# num_keep = np.sum(keep) -# if num_keep == 0: -# continue -# -# # Test 1 -# keep = np.where(keep, (cld_top_hgt[i,] + closeness_top) > icing_alt[i], False) -# # keep = np.where(keep, (cld_top_hgt[i,] - cld_geo_dz[i,] - closeness_bot) < icing_alt[i], False) -# keep = np.where(keep, (cld_top_hgt[i,] - max_depth) < icing_alt[i], False) -# -# # # Test2 -# # 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) -# # -# # # Test4 -# # keep = np.where(keep, np.invert((cld_phase[i,] == 4) & (cld_opd[i,] < opd_thin_threshold)), False) -# # -# # # Test5 and Test6 -# # keep = np.where(keep, np.logical_and(bt_11um[i,] > 228.0, bt_11um[i,] < 270.0), False) -# # -# # # Test3 -# # keep = np.where(keep, (cld_opd[i,] >= opd_thick_threshold) & (cld_phase[i,] == 4), False) -# -# mask.append(keep) -# idxs.append(i) -# num_tested.append(num_keep) -# -# return mask, idxs, num_tested - - -# def apply_qc_no_icing_pireps(icing_alt, cld_top_hgt, cld_geo_dz, cld_phase, cld_opd, cld_mask, bt_11um, solzen, satzen, day_night='ANY'): -# -# if day_night == 'DAY': -# opd_thick_threshold = 20 -# opd_thin_threshold = 1 -# elif day_night == 'NIGHT' or day_night == 'ANY': -# opd_thick_threshold = 2 -# opd_thin_threshold = 0.1 -# -# closeness_top = 100.0 # meters -# closeness_bot = 200.0 -# max_depth = 3000.0 -# -# num_obs = len(icing_alt) -# cld_mask = cld_mask.reshape((num_obs, -1)) -# cld_top_hgt = cld_top_hgt.reshape((num_obs, -1)) -# cld_geo_dz = cld_geo_dz.reshape((num_obs, -1)) -# bt_11um = bt_11um.reshape((num_obs, -1)) -# -# mask = [] -# idxs = [] -# num_tested = [] -# for i in range(num_obs): -# if not check_oblique(satzen[i,]): -# continue -# if day_night == 'NIGHT' and not is_night(solzen[i,]): -# continue -# elif day_night == 'DAY' and not is_day(solzen[i,]): -# continue -# elif day_night == 'ANY': -# pass -# # if not (is_day(solzen[i,]) or 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_geo_dz[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,] + closeness_top) > icing_alt[i], False) -# # keep = np.where(keep, (cld_top_hgt[i,] - cld_geo_dz[i,] - closeness_bot) < icing_alt[i], False) -# keep = np.where(keep, (cld_top_hgt[i,] - max_depth) < icing_alt[i], False) -# -# # keep = np.where(keep, np.logical_and(bt_11um[i,] > 228.0, bt_11um[i,] < 273.0), False) -# -# # Test6 -# # keep = np.where(keep, (bt_11um[i,] < 228.0), False) -# -# # Test1 -# # keep = np.where(keep, cld_top_hgt[i,] > icing_alt[i], False) -# -# # Test3 -# # keep = np.where(keep, np.invert((cld_opd[i,] >= opd_thick_threshold) & (cld_phase[i,] == 4) & (cld_top_hgt[i,] > icing_alt[i])), False) -# -# mask.append(keep) -# idxs.append(i) -# num_tested.append(num_keep) -# -# return mask, idxs, num_tested - - def apply_qc_icing_pireps(icing_alt, cld_top_hgt, cld_geo_dz, cld_phase, cld_opd, cld_mask, bt_11um, solzen, satzen, cld_top_temp, day_night='ANY', cloudy_frac=0.5): if day_night == 'DAY':