snapshot...

modules/amv/caliop_clavrx_amv.py

@@ -8,6 +8,7 @@ from aeolus.datasource import get_datasource, CLAVRx_CALIPSO, get_parameters_cal
 from util.util import haversine_np
 from amv.intercompare import *
 from util.util import pressure_to_altitude
+from util.gfs_reader import *
 
 amv_file_duration = 60  # minutes
 half_width = 30  # search box centered on CALIOP profile (FGF coordinates)
@@ -356,6 +357,35 @@ new_header = ' TARG     LAT      LON     BOX SRCH    SPD      DIR      PW     LL
 num_fmts = ['%6i', '%7.2f', '%8.2f', '%4i', '%4i', '%7.2f', '%7.2f', '%8.2f', '%4i', '%7.2f', '%7.2f', '%8.2f', '%7.2f', '%4i', '%7.2f', '%4i', '%4i', '%6i', '%6i', '%9.2f']
 
 
+def get_temp_prof_s(gfs_ds, nom_time, amv_lons, amv_lats):
+    gfs_fname, _, _ = gfs_ds.get_file(nom_time, window=180.0)
+    if gfs_fname is None:
+        return None
+    print(gfs_fname)
+    gfs_xr = xr.open_dataset(gfs_fname)
+    gfs_press = gfs_xr['pressure levels']
+    gfs_press = gfs_press.values
+
+    temp_prof = get_vert_profile_s(gfs_xr, ['temperature'], amv_lons, amv_lats, method='nearest')
+    temp_prof = temp_prof.values
+    temp_prof_s = temp_prof[0, :, :]
+
+    return temp_prof_s, gfs_press
+
+
+def get_temp_prof_s_intrp(gfs_ds, nom_time, amv_lons, amv_lats):
+    gfs_fname_left, nom_time_left, gfs_fname_rght, nom_time_rght = get_bounding_gfs_files(nom_time)
+    gfs_xr_left = xr.open_dataset(gfs_fname_left)
+    gfs_xr_rght = xr.open_dataset(gfs_fname_rght)
+    gfs_press = gfs_xr_left['pressure_values'].values
+
+    temp_prof = get_vert_profile_s([gfs_xr_left, gfs_xr_rght], [nom_time_left, nom_time_rght], ['temperature'], amv_lons, amv_lats, method='nearest')
+    temp_prof = temp_prof.values
+    temp_prof_s = temp_prof[0, :, :]
+
+    return temp_prof_s, gfs_press
+
+
 def compute_and_add_geo_hgt(path_to_amvs, path_to_gfs, amv_source, band='14', out_file=None):
     gfs_ds = get_datasource(path_to_gfs, 'GFS')
     amv_ds = get_datasource(path_to_amvs, amv_source, file_time_span=20, band=band)
@@ -374,17 +404,19 @@ def compute_and_add_geo_hgt(path_to_amvs, path_to_gfs, amv_source, band='14', ou
         amv_lons = amvs_nd[:, amv_lon_idx]
         amv_lats = amvs_nd[:, amv_lat_idx]
 
-        gfs_fname, _, _ = gfs_ds.get_file(nom_time, window=180.0)
-        if gfs_fname is None:
-            return None
-        print(gfs_fname)
-        gfs_xr = xr.open_dataset(gfs_fname)
-        gfs_press = gfs_xr['pressure levels']
-        gfs_press = gfs_press.values
-
-        temp_prof = get_vert_profile_s(gfs_xr, ['temperature'], amv_lons, amv_lats, method='nearest')
-        temp_prof = temp_prof.values
-        temp_prof_s = temp_prof[0, :, :]
+        # gfs_fname, _, _ = gfs_ds.get_file(nom_time, window=180.0)
+        # if gfs_fname is None:
+        #     return None
+        # print(gfs_fname)
+        # gfs_xr = xr.open_dataset(gfs_fname)
+        # gfs_press = gfs_xr['pressure levels']
+        # gfs_press = gfs_press.values
+        #
+        # temp_prof = get_vert_profile_s(gfs_xr, ['temperature'], amv_lons, amv_lats, method='nearest')
+        # temp_prof = temp_prof.values
+        # temp_prof_s = temp_prof[0, :, :]
+
+        temp_prof_s, gfs_press = get_temp_prof_s(gfs_ds, nom_time, amv_lons, amv_lats)
 
         alt_f = []
         for k in range(num_amvs):