diff --git a/modules/icing/pirep_goes.py b/modules/icing/pirep_goes.py
index 54ace0c9be3c99dd34e3dfc04743088e043ffd73..42f9f70a698f71a8e6323c0c744745b62d4a23ad 100644
--- a/modules/icing/pirep_goes.py
+++ b/modules/icing/pirep_goes.py
@@ -1083,6 +1083,7 @@ def tile_extract(outfile='/home/rink/tiles_l1b_out.h5', train_params=l1b_ds_list
     ice_lat_s = np.array(ice_lat_s)
     num_ice = icing_int_s.shape[0]
 
+    # No icing  ------------------------------------------------------------
     num_no_ice = 0
     for fidx in range(len(no_icing_files)):
         fname = no_icing_files[fidx]
@@ -1189,6 +1190,48 @@ def tile_extract(outfile='/home/rink/tiles_l1b_out.h5', train_params=l1b_ds_list
     h5f_expl.close()
 
 
+def write_file(outfile, train_params, data_dct, icing_intensity, icing_times, icing_lons, icing_lats):
+    h5f_expl = h5py.File(a_clvr_file, 'r')
+    h5f_out = h5py.File(outfile, 'w')
+
+    for idx, ds_name in enumerate(train_params):
+        dt = ds_types[idx]
+        data = data_dct[ds_name]
+        h5f_out.create_dataset(ds_name, data=data, dtype=dt)
+
+    icing_int_ds = h5f_out.create_dataset('icing_intensity', data=icing_intensity, dtype='i4')
+    icing_int_ds.attrs.create('long_name', data='From PIREP. -1:No Icing, 1:Trace, 2:Light, 3:Light Moderate, 4:Moderate, 5:Moderate Severe, 6:Severe')
+
+    time_ds = h5f_out.create_dataset('time', data=icing_times, dtype='f4')
+    time_ds.attrs.create('units', data='seconds since 1970-1-1 00:00:00')
+    time_ds.attrs.create('long_name', data='PIREP time')
+
+    lon_ds = h5f_out.create_dataset('longitude', data=icing_lons, dtype='f4')
+    lon_ds.attrs.create('units', data='degrees_east')
+    lon_ds.attrs.create('long_name', data='PIREP longitude')
+
+    lat_ds = h5f_out.create_dataset('latitude', data=icing_lats, dtype='f4')
+    lat_ds.attrs.create('units', data='degrees_north')
+    lat_ds.attrs.create('long_name', data='PIREP latitude')
+
+    # copy relevant attributes
+    for ds_name in train_params:
+        h5f_ds = h5f_out[ds_name]
+        h5f_ds.attrs.create('standard_name', data=h5f_expl[ds_name].attrs.get('standard_name'))
+        h5f_ds.attrs.create('long_name', data=h5f_expl[ds_name].attrs.get('long_name'))
+        h5f_ds.attrs.create('units', data=h5f_expl[ds_name].attrs.get('units'))
+        attr = h5f_expl[ds_name].attrs.get('actual_range')
+        if attr is not None:
+            h5f_ds.attrs.create('actual_range', data=attr)
+        attr = h5f_expl[ds_name].attrs.get('flag_values')
+        if attr is not None:
+            h5f_ds.attrs.create('flag_values', data=attr)
+
+    # --- close files
+    h5f_out.close()
+    h5f_expl.close()
+
+
 def run_mean_std(check_cloudy=False, no_icing_to_icing_ratio=5):
     ds_list = ['cld_height_acha', 'cld_geo_thick', 'cld_press_acha',
                'supercooled_cloud_fraction', 'cld_temp_acha', 'cld_opd_acha',