From 952f68731fabfed14f38d23c6c36176dcb68effe Mon Sep 17 00:00:00 2001 From: tomrink <rink@ssec.wisc.edu> Date: Thu, 12 Nov 2020 01:16:21 -0600 Subject: [PATCH] snapshot... --- modules/aeolus/aeolus_amv.py | 64 +++++++++++++++++++++++++++++------- 1 file changed, 53 insertions(+), 11 deletions(-) diff --git a/modules/aeolus/aeolus_amv.py b/modules/aeolus/aeolus_amv.py index 327d5eeb..1bd33a5e 100644 --- a/modules/aeolus/aeolus_amv.py +++ b/modules/aeolus/aeolus_amv.py @@ -671,7 +671,7 @@ def subset_by_time(match_file, dt_str_0, dt_str_1, amv_var_names): return prof_da, amvs_da, prof_locs_da -def analyze(prof_da, amvs_da, prof_locs_da): +def analyze(prof_da, amvs_da, prof_locs_da, dist_threshold=5.0): # sort amvs by distance to profile dst = amvs_da.sel(num_params='dist_to_prof') s_i = np.argsort(dst, axis=1) @@ -686,14 +686,11 @@ def analyze(prof_da, amvs_da, prof_locs_da): for k in range(prof_da.shape[0]): prof_da[k] = prof_da[k, s_i[k, ::-1]] + # analyze single cloud layer profiles one_lyr = prof_locs_da.sel(space='num_levels') == 1 - dst = amvs_da.sel(num_params='dist_to_prof') - close = np.logical_and(dst > 0.0, dst < 5.0) - - # compare for one cloud layer profile - one_lyr_profs = prof_da[one_lyr] one_lyr_amvs = amvs_da[one_lyr] + print('number of one layer profs: ', one_lyr_profs.shape[0]) hgt_vld = one_lyr_amvs.sel(num_params='H_3D') > 0 hgt_vld = hgt_vld.values @@ -701,19 +698,64 @@ def analyze(prof_da, amvs_da, prof_locs_da): one_lyr_amvs = one_lyr_amvs[has_vld] one_lyr_profs = one_lyr_profs[has_vld] + print('number of one layer profs with at least one valid AMV height', one_lyr_profs.shape[0]) # compare profile highest cloud layer top to closest within 5km dst = one_lyr_amvs.sel(num_params='dist_to_prof') - close = (np.logical_and(dst > 0.0, dst < 2.0)).values + close = (np.logical_and(dst > 0.0, dst < dist_threshold)).values close = close.sum(1) > 0 one_lyr_amvs = one_lyr_amvs[close] one_lyr_profs = one_lyr_profs[close] + print('number of one layer profs with at least one AMV within 5km: ', one_lyr_profs.shape[0]) + cnt = 0 for k in range(one_lyr_profs.shape[0]): - dst = one_lyr_amvs[k, :, 2] - b = np.logical_and(dst > 0.0, dst < 2.0) - diffs = one_lyr_amvs[k, b, 0] - one_lyr_profs[k, 0, 1] - print(diffs.values, one_lyr_amvs[k, b, 2].values) + dst = one_lyr_amvs[k, :, ].sel(num_params='dist_to_prof') + b = np.logical_and(dst > 0.0, dst < dist_threshold) + h_3d = one_lyr_amvs[k, :].sel(num_params='H_3D') + h_3d = h_3d[b] + vld = h_3d > 0 + h_3d = h_3d[vld] + if len(h_3d) > 0: + in_lyr = np.logical_and(h_3d > one_lyr_profs[k, 0, 0], h_3d < one_lyr_profs[k, 0, 1]) + cnt += np.sum(in_lyr) + + print(cnt) + + # Do calculations for single and multi-layer combined + hgt_vld = amvs_da.sel(num_params='H_3D') > 0 + hgt_vld = hgt_vld.values + has_vld = hgt_vld.sum(1) > 0 + + amvs_da = amvs_da[has_vld] + prof_da = prof_da[has_vld] + prof_locs_da = prof_locs_da[has_vld] + print('number of profs with at least one valid AMV height', prof_da.shape[0]) + # compare profile highest cloud layer top to closest within 5km + dst = amvs_da.sel(num_params='dist_to_prof') + close = (np.logical_and(dst > 0.0, dst < dist_threshold)).values + close = close.sum(1) > 0 + + amvs_da = amvs_da[close] + prof_da = prof_da[close] + prof_locs_da = prof_locs_da[close] + print('number of profs with at least one AMV within 5km: ', prof_da.shape[0]) + + cnt = 0 + for k in range(prof_da.shape[0]): + dst = amvs_da[k, :, ].sel(num_params='dist_to_prof') + b = np.logical_and(dst > 0.0, dst < dist_threshold) + h_3d = amvs_da[k, :].sel(num_params='H_3D') + h_3d = h_3d[b] + vld = h_3d > 0 + h_3d = h_3d[vld] + if len(h_3d) > 0: + nlevs = prof_locs_da[k].values.astype(int)[3] + for j in range(nlevs): + in_lyr = np.logical_and(h_3d > prof_da[k, j, 0], h_3d < prof_da[k, j, 1]) + cnt += np.sum(in_lyr) + + print(cnt) return one_lyr_profs, one_lyr_amvs -- GitLab