Skip to content
Snippets Groups Projects
Commit 5bfd6a3b authored by tomrink's avatar tomrink
Browse files

initial commit, just reorganizing stuff

parent 97a1dd13
Branches
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 644 additions and 0 deletions
modules/icing/util.py 0 → 100644
+ 644
− 0
View file @ 5bfd6a3b
import numpy as np
import deeplearning.icing_fcn as icing_fcn
import deeplearning.icing_cnn as icing_cnn
from icing.pirep_goes import setup, time_filter_3
from util.util import get_time_tuple_utc, is_day, check_oblique, homedir, write_icing_file_nc4,\
make_for_full_domain_predict, 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
from aeolus.datasource import CLAVRx, CLAVRx_VIIRS, GOESL1B, CLAVRx_H08
import h5py
import datetime
def get_training_parameters(day_night='DAY', l1b_andor_l2='both'):
if day_night == 'DAY':
train_params_l2 = ['cld_height_acha', 'cld_geo_thick', 'cld_temp_acha', 'cld_press_acha', 'supercooled_cloud_fraction',
'cld_emiss_acha', 'conv_cloud_fraction', 'cld_reff_dcomp', 'cld_opd_dcomp', 'iwc_dcomp', 'lwc_dcomp']
train_params_l1b = ['temp_10_4um_nom', 'temp_11_0um_nom', 'temp_12_0um_nom', 'temp_13_3um_nom', 'temp_3_75um_nom',
'temp_6_2um_nom', 'temp_6_7um_nom', 'temp_7_3um_nom', 'temp_8_5um_nom', 'temp_9_7um_nom',
'refl_0_47um_nom', 'refl_0_65um_nom', 'refl_0_86um_nom', 'refl_1_38um_nom', 'refl_1_60um_nom']
else:
train_params_l2 = ['cld_height_acha', 'cld_geo_thick', 'cld_temp_acha', 'cld_press_acha', 'supercooled_cloud_fraction',
'cld_emiss_acha', 'conv_cloud_fraction', 'cld_reff_acha', 'cld_opd_acha']
train_params_l1b = ['temp_10_4um_nom', 'temp_11_0um_nom', 'temp_12_0um_nom', 'temp_13_3um_nom', 'temp_3_75um_nom',
'temp_6_2um_nom', 'temp_6_7um_nom', 'temp_7_3um_nom', 'temp_8_5um_nom', 'temp_9_7um_nom']
if l1b_andor_l2 == 'both':
train_params = train_params_l1b + train_params_l2
elif l1b_andor_l2 == 'l1b':
train_params = train_params_l1b
elif l1b_andor_l2 == 'l2':
train_params = train_params_l2
return train_params
flt_level_ranges = {k: None for k in range(5)}
flt_level_ranges[0] = [0.0, 2000.0]
flt_level_ranges[1] = [2000.0, 4000.0]
flt_level_ranges[2] = [4000.0, 6000.0]
flt_level_ranges[3] = [6000.0, 8000.0]
flt_level_ranges[4] = [8000.0, 15000.0]
def run_make_images(clvrx_dir='/Users/tomrink/data/clavrx/RadC/', ckpt_dir_s_path='/Users/tomrink/tf_model/', prob_thresh=0.5, satellite='GOES16', domain='CONUS',
extent=[-105, -70, 15, 50],
pirep_file='/Users/tomrink/data/pirep/pireps_202109200000_202109232359.csv',
obs_lons=None, obs_lats=None, obs_times=None, obs_alt=None, flight_level=None,
use_flight_altitude=False, day_night='DAY', l1b_andor_l2='l2'):
if pirep_file is not None:
ice_dict, no_ice_dict, neg_ice_dict = setup(pirep_file)
if satellite == 'H08':
clvrx_ds = CLAVRx_H08(clvrx_dir)
else:
clvrx_ds = CLAVRx(clvrx_dir)
clvrx_files = clvrx_ds.flist
alt_lo, alt_hi = 0.0, 15000.0
if flight_level is not None:
alt_lo, alt_hi = flt_level_ranges[flight_level]
train_params = get_training_parameters(day_night=day_night, l1b_andor_l2=l1b_andor_l2)
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, ll, cc = make_for_full_domain_predict(h5f, name_list=train_params, satellite=satellite, domain=domain)
num_elems, num_lines = len(cc), len(ll)
dto, _ = get_time_tuple_utc(ts)
dto_0 = dto - datetime.timedelta(minutes=30)
dto_1 = dto + datetime.timedelta(minutes=30)
ts_0 = dto_0.timestamp()
ts_1 = dto_1.timestamp()
if pirep_file is not None:
_, keep_lons, keep_lats, _ = time_filter_3(ice_dict, ts_0, ts_1, alt_lo, alt_hi)
elif obs_times is not None:
keep = np.logical_and(obs_times >= ts_0, obs_times < ts_1)
keep = np.where(keep, np.logical_and(obs_alt >= alt_lo, obs_alt < alt_hi), False)
keep_lons = obs_lons[keep]
keep_lats = obs_lats[keep]
else:
keep_lons = None
keep_lats = None
ice_lons, ice_lats, preds_2d = icing_cnn.run_evaluate_static_avg(data_dct, ll, cc, ckpt_dir_s_path=ckpt_dir_s_path,
flight_level=flight_level, prob_thresh=prob_thresh,
satellite=satellite, domain=domain,
use_flight_altitude=use_flight_altitude)
make_icing_image(h5f, None, ice_lons, ice_lats, clvrx_str_time, satellite, domain,
ice_lons_vld=keep_lons, ice_lats_vld=keep_lats, extent=extent)
# preds_2d_dct, probs_2d_dct = run_evaluate_static(data_dct, num_lines, num_elems, day_night=day_night,
# ckpt_dir_s_path=ckpt_dir_s_path, prob_thresh=prob_thresh,
# flight_levels=[0],
# use_flight_altitude=use_flight_altitude)
#
# make_icing_image(None, probs_2d_dct[0], None, None, clvrx_str_time, satellite, domain,
# ice_lons_vld=keep_lons, ice_lats_vld=keep_lats, extent=extent)
h5f.close()
print('Done: ', clvrx_str_time)
def run_icing_predict(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=True, 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, ll, cc = make_for_full_domain_predict(h5f, name_list=train_params, satellite=satellite, domain=domain, res_fac=res_fac)
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)
# 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 = icing_cnn.run_evaluate_static(day_grd_dct, num_day_tiles, 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 = icing_cnn.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_fcn(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, offset=8)
num_elems = len(cc)
num_lines = len(ll)
if fidx == 0:
nav = get_navigation(satellite, domain)
lons_2d, lats_2d, x_rad, y_rad = get_lon_lat_2d_mesh(nav, ll, cc)
day_idxs = solzen < 80.0
num_day_tiles = np.sum(day_idxs)
nght_idxs = solzen > 100.0
num_nght_tiles = np.sum(nght_idxs)
# 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:
preds_day_dct, probs_day_dct = icing_fcn.run_evaluate_static(data_dct, 1, 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)
for flvl in flight_levels:
preds = preds_day_dct[flvl]
probs = probs_day_dct[flvl]
fd_preds = preds_2d_dct[flvl]
fd_probs = probs_2d_dct[flvl]
fd_preds[day_idxs] = preds[day_idxs]
fd_probs[day_idxs] = probs[day_idxs]
if (day_night == 'AUTO' or day_night == 'NIGHT') and num_nght_tiles > 0:
preds_nght_dct, probs_nght_dct = icing_fcn.run_evaluate_static_fcn(data_dct, 1, 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)
for flvl in flight_levels:
preds = preds_nght_dct[flvl]
probs = probs_nght_dct[flvl]
fd_preds = preds_2d_dct[flvl]
fd_probs = probs_2d_dct[flvl]
fd_preds[nght_idxs] = preds[nght_idxs]
fd_probs[nght_idxs] = probs[nght_idxs]
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',
l1b_andor_l2='BOTH', use_flight_altitude=True, res_fac=1,
extent=[-105, -70, 15, 50],
pirep_file='/Users/tomrink/data/pirep/pireps_202109200000_202109232359.csv',
obs_lons=None, obs_lats=None, obs_times=None, obs_alt=None, flight_level=None):
if use_flight_altitude is True:
flight_levels = [0, 1, 2, 3, 4]
else:
flight_levels = [0]
if pirep_file is not None:
ice_dict, no_ice_dict, neg_ice_dict = setup(pirep_file)
alt_lo, alt_hi = 0.0, 15000.0
if flight_level is not None:
alt_lo, alt_hi = flt_level_ranges[flight_level]
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, ll, cc = make_for_full_domain_predict(h5f, name_list=train_params, satellite=satellite, domain=domain, res_fac=res_fac)
if fidx == 0:
num_elems = len(cc)
num_lines = len(ll)
nav = get_navigation(satellite, domain)
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)
# 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 = icing_cnn.run_evaluate_static(day_grd_dct, num_day_tiles, 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 = icing_cnn.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))
dto, _ = get_time_tuple_utc(ts)
dto_0 = dto - datetime.timedelta(minutes=30)
dto_1 = dto + datetime.timedelta(minutes=30)
ts_0 = dto_0.timestamp()
ts_1 = dto_1.timestamp()
if pirep_file is not None:
_, keep_lons, keep_lats, _ = time_filter_3(ice_dict, ts_0, ts_1, alt_lo, alt_hi)
elif obs_times is not None:
keep = np.logical_and(obs_times >= ts_0, obs_times < ts_1)
keep = np.where(keep, np.logical_and(obs_alt >= alt_lo, obs_alt < alt_hi), False)
keep_lons = obs_lons[keep]
keep_lats = obs_lats[keep]
else:
keep_lons = None
keep_lats = None
prob_s = []
for flvl in flight_levels:
probs = probs_2d_dct[flvl]
prob_s.append(probs)
prob_s = np.stack(prob_s, axis=-1)
max_prob = np.max(prob_s, axis=2)
max_prob = np.where(max_prob < 0.5, np.nan, max_prob)
make_icing_image(h5f, max_prob, None, None, clvrx_str_time, satellite, domain,
ice_lons_vld=keep_lons, ice_lats_vld=keep_lats, extent=extent)
print('Done: ', clvrx_str_time)
h5f.close()
def run_icing_predict_image_fcn(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=True, res_fac=1,
extent=[-105, -70, 15, 50],
pirep_file='/Users/tomrink/data/pirep/pireps_202109200000_202109232359.csv',
obs_lons=None, obs_lats=None, obs_times=None, obs_alt=None, flight_level=None):
if use_flight_altitude is True:
flight_levels = [0, 1, 2, 3, 4]
else:
flight_levels = [0]
if pirep_file is not None:
ice_dict, no_ice_dict, neg_ice_dict = setup(pirep_file)
alt_lo, alt_hi = 0.0, 15000.0
if flight_level is not None:
alt_lo, alt_hi = flt_level_ranges[flight_level]
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')
dto, _ = get_time_tuple_utc(ts)
dto_0 = dto - datetime.timedelta(minutes=30)
dto_1 = dto + datetime.timedelta(minutes=30)
ts_0 = dto_0.timestamp()
ts_1 = dto_1.timestamp()
if pirep_file is not None:
_, keep_lons, keep_lats, _ = time_filter_3(ice_dict, ts_0, ts_1, alt_lo, alt_hi)
elif obs_times is not None:
keep = np.logical_and(obs_times >= ts_0, obs_times < ts_1)
keep = np.where(keep, np.logical_and(obs_alt >= alt_lo, obs_alt < alt_hi), False)
keep_lons = obs_lons[keep]
keep_lats = obs_lats[keep]
else:
keep_lons = None
keep_lats = None
data_dct, solzen, satzen, ll, cc = prepare_evaluate(h5f, name_list=train_params, satellite=satellite, domain=domain, offset=8)
num_elems = len(cc)
num_lines = len(ll)
if fidx == 0:
nav = get_navigation(satellite, domain)
lons_2d, lats_2d, x_rad, y_rad = get_lon_lat_2d_mesh(nav, ll, cc)
day_idxs = solzen < 80.0
num_day_tiles = np.sum(day_idxs)
nght_idxs = solzen > 100.0
num_nght_tiles = np.sum(nght_idxs)
# 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:
preds_day_dct, probs_day_dct = icing_fcn.run_evaluate_static(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)
for flvl in flight_levels:
preds = preds_day_dct[flvl]
probs = probs_day_dct[flvl]
fd_preds = preds_2d_dct[flvl]
fd_probs = probs_2d_dct[flvl]
fd_preds[day_idxs] = preds[day_idxs]
fd_probs[day_idxs] = probs[day_idxs]
if (day_night == 'AUTO' or day_night == 'NIGHT') and num_nght_tiles > 0:
preds_nght_dct, probs_nght_dct = icing_fcn.run_evaluate_static_fcn(data_dct, 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)
for flvl in flight_levels:
preds = preds_nght_dct[flvl]
probs = probs_nght_dct[flvl]
fd_preds = preds_2d_dct[flvl]
fd_probs = probs_2d_dct[flvl]
fd_preds[nght_idxs] = preds[nght_idxs]
fd_probs[nght_idxs] = probs[nght_idxs]
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))
prob_s = []
for flvl in flight_levels:
probs = probs_2d_dct[flvl]
prob_s.append(probs)
prob_s = np.stack(prob_s, axis=-1)
max_prob = np.max(prob_s, axis=2)
max_prob = np.where(max_prob < 0.5, np.nan, max_prob)
make_icing_image(h5f, max_prob, None, None, clvrx_str_time, satellite, domain,
ice_lons_vld=keep_lons, ice_lats_vld=keep_lats, extent=extent)
print('Done: ', clvrx_str_time)
h5f.close()
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment