diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py index 8a63918997c15af5e4e9d6d380d67a94487955d1..c7d836f8468e504e76113cae906ac15bdebcab62 100644 --- a/modules/icing/pirep_goes.py +++ b/modules/icing/pirep_goes.py @@ -1,3 +1,4 @@ +from deeplearning.icing_fcn import run_evaluate_static_fcn from icing.pireps import pirep_icing import numpy as np import pickle @@ -5,7 +6,7 @@ import matplotlib.pyplot as plt import os from util.util import get_time_tuple_utc, GenericException, add_time_range_to_filename, is_night, is_day, \ check_oblique, get_timestamp, homedir, write_icing_file_nc4, make_for_full_domain_predict, \ - get_indexes_within_threshold + get_indexes_within_threshold, prepare_evaluate from util.plot import make_icing_image from util.geos_nav import get_navigation, get_lon_lat_2d_mesh from util.setup import model_path_day, model_path_night @@ -2307,6 +2308,143 @@ def run_icing_predict(clvrx_dir='/Users/tomrink/data/clavrx/RadC/', output_dir=h h5f.close() +def run_icing_predict_new(clvrx_dir='/Users/tomrink/data/clavrx/RadC/', output_dir=homedir, + day_model_path=model_path_day, night_model_path=model_path_night, + prob_thresh=0.5, satellite='GOES16', domain='CONUS', day_night='AUTO', + l1b_andor_l2='both', use_flight_altitude=False, res_fac=1, use_nan=False): + if use_flight_altitude is True: + flight_levels = [0, 1, 2, 3, 4] + else: + flight_levels = [0] + + day_train_params = get_training_parameters(day_night='DAY', l1b_andor_l2=l1b_andor_l2) + nght_train_params = get_training_parameters(day_night='NIGHT', l1b_andor_l2=l1b_andor_l2) + + if day_night == 'AUTO': + train_params = list(set(day_train_params + nght_train_params)) + elif day_night == 'DAY': + train_params = day_train_params + elif day_night == 'NIGHT': + train_params = nght_train_params + + if satellite == 'H08': + clvrx_ds = CLAVRx_H08(clvrx_dir) + else: + clvrx_ds = CLAVRx(clvrx_dir) + clvrx_files = clvrx_ds.flist + + for fidx, fname in enumerate(clvrx_files): + h5f = h5py.File(fname, 'r') + dto = clvrx_ds.get_datetime(fname) + ts = dto.timestamp() + clvrx_str_time = dto.strftime('%Y-%m-%d_%H:%M') + + data_dct, solzen, satzen, ll, cc = prepare_evaluate(h5f, name_list=train_params, satellite=satellite, domain=domain) + num_elems = len(cc) + num_lines = len(ll) + + # Don't need these yet + # if fidx == 0: + # num_elems = len(cc) + # num_lines = len(ll) + # nav = get_navigation(satellite, domain) + # lons_2d, lats_2d, x_rad, y_rad = get_lon_lat_2d_mesh(nav, ll, cc, offset=int(8 / res_fac)) + + # ancil_data_dct, _, _ = make_for_full_domain_predict(h5f, name_list= + # ['solar_zenith_angle', 'sensor_zenith_angle', 'cld_height_acha', 'cld_geo_thick'], + # satellite=satellite, domain=domain, res_fac=res_fac) + # + # satzen = ancil_data_dct['sensor_zenith_angle'] + # solzen = ancil_data_dct['solar_zenith_angle'] + # day_idxs = [] + # nght_idxs = [] + # for j in range(num_lines): + # for i in range(num_elems): + # k = i + j*num_elems + # if not check_oblique(satzen[k]): + # continue + # if is_day(solzen[k]): + # day_idxs.append(k) + # else: + # nght_idxs.append(k) + # + # num_tiles = num_lines * num_elems + # num_day_tiles = len(day_idxs) + # num_nght_tiles = len(nght_idxs) + num_day_tiles = 1 + + # initialize output arrays + probs_2d_dct = {flvl: None for flvl in flight_levels} + preds_2d_dct = {flvl: None for flvl in flight_levels} + for flvl in flight_levels: + fd_preds = np.zeros(num_lines * num_elems, dtype=np.int8) + fd_preds[:] = -1 + fd_probs = np.zeros(num_lines * num_elems, dtype=np.float32) + fd_probs[:] = -1.0 + preds_2d_dct[flvl] = fd_preds + probs_2d_dct[flvl] = fd_probs + + if (day_night == 'AUTO' or day_night == 'DAY') and num_day_tiles > 0: + + # day_data_dct = {name: [] for name in day_train_params} + # for name in day_train_params: + # for k in day_idxs: + # day_data_dct[name].append(data_dct[name][k]) + # day_grd_dct = {name: None for name in day_train_params} + # for ds_name in day_train_params: + # day_grd_dct[ds_name] = np.stack(day_data_dct[ds_name]) + + preds_day_dct, probs_day_dct = run_evaluate_static_fcn(data_dct, day_model_path, + day_night='DAY', l1b_or_l2=l1b_andor_l2, + prob_thresh=prob_thresh, + use_flight_altitude=use_flight_altitude, + flight_levels=flight_levels) + day_idxs = np.array(day_idxs) + for flvl in flight_levels: + day_preds = preds_day_dct[flvl] + day_probs = probs_day_dct[flvl] + fd_preds = preds_2d_dct[flvl] + fd_probs = probs_2d_dct[flvl] + fd_preds[day_idxs] = day_preds[:] + fd_probs[day_idxs] = day_probs[:] + + # if (day_night == 'AUTO' or day_night == 'NIGHT') and num_nght_tiles > 0: + # + # nght_data_dct = {name: [] for name in nght_train_params} + # for name in nght_train_params: + # for k in nght_idxs: + # nght_data_dct[name].append(data_dct[name][k]) + # nght_grd_dct = {name: None for name in nght_train_params} + # for ds_name in nght_train_params: + # nght_grd_dct[ds_name] = np.stack(nght_data_dct[ds_name]) + # + # preds_nght_dct, probs_nght_dct = run_evaluate_static(nght_grd_dct, num_nght_tiles, night_model_path, + # day_night='NIGHT', l1b_or_l2=l1b_andor_l2, prob_thresh=prob_thresh, + # use_flight_altitude=use_flight_altitude, + # flight_levels=flight_levels) + # nght_idxs = np.array(nght_idxs) + # for flvl in flight_levels: + # nght_preds = preds_nght_dct[flvl] + # nght_probs = probs_nght_dct[flvl] + # fd_preds = preds_2d_dct[flvl] + # fd_probs = probs_2d_dct[flvl] + # fd_preds[nght_idxs] = nght_preds[:] + # fd_probs[nght_idxs] = nght_probs[:] + + for flvl in flight_levels: + fd_preds = preds_2d_dct[flvl] + fd_probs = probs_2d_dct[flvl] + preds_2d_dct[flvl] = fd_preds.reshape((num_lines, num_elems)) + probs_2d_dct[flvl] = fd_probs.reshape((num_lines, num_elems)) + + # write_icing_file_nc4(clvrx_str_time, output_dir, preds_2d_dct, probs_2d_dct, + # x_rad, y_rad, lons_2d, lats_2d, cc, ll, + # satellite=satellite, domain=domain, use_nan=use_nan, prob_thresh=prob_thresh) + + print('Done: ', clvrx_str_time) + h5f.close() + + def run_icing_predict_image(clvrx_dir='/Users/tomrink/data/clavrx/RadC/', output_dir=homedir, day_model_path=model_path_day, night_model_path=model_path_night, prob_thresh=0.5, satellite='GOES16', domain='CONUS', day_night='AUTO',