From 05cf73d4917a3a064120cb0dd2dfb02d4942442a Mon Sep 17 00:00:00 2001
From: tomrink <rink@ssec.wisc.edu>
Date: Thu, 14 Oct 2021 08:57:24 -0500
Subject: [PATCH] add icing altitude to tiles training/test files

---
 modules/icing/pirep_goes.py | 32 +++++++++++++++++++++++++++++---
 1 file changed, 29 insertions(+), 3 deletions(-)

diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py
index 6676e8e6..f781fec0 100644
--- a/modules/icing/pirep_goes.py
+++ b/modules/icing/pirep_goes.py
@@ -1142,6 +1142,8 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     no_ice_lon_s = []
     ice_lat_s = []
     no_ice_lat_s = []
+    ice_flt_alt_s = []
+    no_ice_flt_alt_s = []
 
     h5_s_icing = []
     h5_s_no_icing = []
@@ -1167,6 +1169,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         lats = f['latitude']
 
         icing_int = f['icing_intensity'][:]
+        flt_altitude = f['icing_altitude'][:]
 
         for i in range(num_obs):
             if check_keep_out and check_times(times[i]):
@@ -1183,6 +1186,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
             ice_time_s.append(times[i])
             ice_lon_s.append(lons[i])
             ice_lat_s.append(lats[i])
+            ice_flt_alt_s.append(flt_altitude[i])
 
         print(fname)
 
@@ -1193,6 +1197,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     ice_time_s = np.array(ice_time_s)
     ice_lon_s = np.array(ice_lon_s)
     ice_lat_s = np.array(ice_lat_s)
+    ice_flt_alt_s = np.array(ice_flt_alt_s)
     num_ice = icing_int_s.shape[0]
 
     # No icing  ------------------------------------------------------------
@@ -1206,6 +1211,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         num_obs = len(times)
         lons = f['longitude']
         lats = f['latitude']
+        flt_altitude = f['icing_altitude'][:]
 
         for i in range(num_obs):
             if check_keep_out and not check_times(times[i]):
@@ -1221,6 +1227,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
             no_ice_time_s.append(times[i])
             no_ice_lon_s.append(lons[i])
             no_ice_lat_s.append(lats[i])
+            no_ice_flt_alt_s.append(flt_altitude[i])
 
         print(fname)
 
@@ -1231,11 +1238,13 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     no_ice_time_s = np.array(no_ice_time_s)
     no_ice_lon_s = np.array(no_ice_lon_s)
     no_ice_lat_s = np.array(no_ice_lat_s)
+    no_ice_flt_alt_s = np.array(no_ice_flt_alt_s)
 
     icing_intensity = np.concatenate([icing_int_s, no_icing_int_s])
     icing_times = np.concatenate([ice_time_s, no_ice_time_s])
     icing_lons = np.concatenate([ice_lon_s, no_ice_lon_s])
     icing_lats = np.concatenate([ice_lat_s, no_ice_lat_s])
+    icing_alt = np.concatenate([ice_flt_alt_s, no_ice_flt_alt_s])
 
     data_dct = {}
     for ds_name in params:
@@ -1249,6 +1258,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     icing_times = icing_times[ds_indexes]
     icing_lons = icing_lons[ds_indexes]
     icing_lats = icing_lats[ds_indexes]
+    icing_alt = icing_alt[ds_indexes]
 
     if do_split:
         trn_idxs, tst_idxs = split_data(icing_times)
@@ -1264,6 +1274,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     trn_icing_times = icing_times[trn_idxs,]
     trn_icing_lons = icing_lons[trn_idxs,]
     trn_icing_lats = icing_lats[trn_idxs,]
+    trn_icing_alt = icing_alt[trn_idxs,]
 
     #  Data augmentation -------------------------------------------------------------
     if augment:
@@ -1272,12 +1283,14 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         trn_icing_times_aug = []
         trn_icing_lons_aug = []
         trn_icing_lats_aug = []
+        trn_icing_alt_aug = []
 
         for k in range(trn_icing_intensity.shape[0]):
             iceint = trn_icing_intensity[k]
             icetime = trn_icing_times[k]
             icelon = trn_icing_lons[k]
             icelat = trn_icing_lats[k]
+            icealt = trn_icing_alt[k]
             if iceint == 3 or iceint == 4 or iceint == 5 or iceint == 6:
                 for ds_name in params:
                     dat = trn_data_dct[ds_name]
@@ -1301,12 +1314,17 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
                 trn_icing_lats_aug.append(icelat)
                 trn_icing_lats_aug.append(icelat)
 
+                trn_icing_alt_aug.append(icealt)
+                trn_icing_alt_aug.append(icealt)
+                trn_icing_alt_aug.append(icealt)
+
         for ds_name in params:
             trn_data_dct_aug[ds_name] = np.stack(trn_data_dct_aug[ds_name])
         trn_icing_intensity_aug = np.stack(trn_icing_intensity_aug)
         trn_icing_times_aug = np.stack(trn_icing_times_aug)
         trn_icing_lons_aug = np.stack(trn_icing_lons_aug)
         trn_icing_lats_aug = np.stack(trn_icing_lats_aug)
+        trn_icing_alt_aug = np.stack(trn_icing_alt_aug)
 
         for ds_name in params:
             trn_data_dct[ds_name] = np.concatenate([trn_data_dct[ds_name], trn_data_dct_aug[ds_name]])
@@ -1314,6 +1332,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         trn_icing_times = np.concatenate([trn_icing_times, trn_icing_times_aug])
         trn_icing_lons = np.concatenate([trn_icing_lons, trn_icing_lons_aug])
         trn_icing_lats = np.concatenate([trn_icing_lats, trn_icing_lats_aug])
+        trn_icing_alt = np.concatenate([trn_icing_alt, trn_icing_alt_aug])
 
     # do sort
     ds_indexes = np.argsort(trn_icing_times)
@@ -1323,8 +1342,9 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
     trn_icing_times = trn_icing_times[ds_indexes]
     trn_icing_lons = trn_icing_lons[ds_indexes]
     trn_icing_lats = trn_icing_lats[ds_indexes]
+    trn_icing_alt = trn_icing_alt[ds_indexes]
 
-    write_file(trnfile, params, param_types, trn_data_dct, trn_icing_intensity, trn_icing_times, trn_icing_lons, trn_icing_lats)
+    write_file(trnfile, params, param_types, trn_data_dct, trn_icing_intensity, trn_icing_times, trn_icing_lons, trn_icing_lats, trn_icing_alt)
 
     if do_split:
         tst_data_dct = {}
@@ -1334,6 +1354,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         tst_icing_times = icing_times[tst_idxs,]
         tst_icing_lons = icing_lons[tst_idxs,]
         tst_icing_lats = icing_lats[tst_idxs,]
+        tst_icing_alt = icing_alt[tst_idxs,]
 
         # do sort
         ds_indexes = np.argsort(tst_icing_times)
@@ -1343,8 +1364,9 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         tst_icing_times = tst_icing_times[ds_indexes]
         tst_icing_lons = tst_icing_lons[ds_indexes]
         tst_icing_lats = tst_icing_lats[ds_indexes]
+        tst_icing_alt = tst_icing_alt[ds_indexes]
 
-        write_file(tstfile, params, param_types, tst_data_dct, tst_icing_intensity, tst_icing_times, tst_icing_lons, tst_icing_lats)
+        write_file(tstfile, params, param_types, tst_data_dct, tst_icing_intensity, tst_icing_times, tst_icing_lons, tst_icing_lats, tst_icing_alt)
 
     # --- close files
     for h5f in h5_s_icing:
@@ -1354,7 +1376,7 @@ def tile_extract(trnfile='/home/rink/tiles_l1b_train.h5', tstfile='/home/rink/ti
         h5f.close()
 
 
-def write_file(outfile, params, param_types, data_dct, icing_intensity, icing_times, icing_lons, icing_lats):
+def write_file(outfile, params, param_types, data_dct, icing_intensity, icing_times, icing_lons, icing_lats, icing_alt):
     h5f_expl = h5py.File(a_clvr_file, 'r')
     h5f_out = h5py.File(outfile, 'w')
 
@@ -1378,6 +1400,10 @@ def write_file(outfile, params, param_types, data_dct, icing_intensity, icing_ti
     lat_ds.attrs.create('units', data='degrees_north')
     lat_ds.attrs.create('long_name', data='PIREP latitude')
 
+    alt_ds = h5f_out.create_dataset('flight_altitude', data=icing_alt, dtype='f4')
+    alt_ds.attrs.create('units', data='meter')
+    alt_ds.attrs.create('long_name', data='PIREP altitude')
+
     # copy relevant attributes
     for ds_name in params:
         h5f_ds = h5f_out[ds_name]
-- 
GitLab