Skip to content
Snippets Groups Projects
Commit f81f8758 authored by Greg Quinn's avatar Greg Quinn
Browse files

Pull out use of weights

parent b7f6606c
No related branches found
No related tags found
No related merge requests found
......@@ -23,8 +23,6 @@ sd = SD(first_file)
n = sd.select('VIIRS_N').get()
rads_sum = sd.select('VIIRS_Radiance_Sum').get()
rads_squared_sum = sd.select('VIIRS_Radiance_Squared_Sum').get()
total_weights = sd.select('VIIRS_Total_Weights').get()
weighted_sum = sd.select('VIIRS_Weighted_Sum').get()
flags = sd.select('VIIRS_Flags').get()
# now fold in data from the other intermediate files
......@@ -33,8 +31,6 @@ for filename in other_files:
n += sd.select('VIIRS_N').get()
rads_sum += sd.select('VIIRS_Radiance_Sum').get()
rads_squared_sum += sd.select('VIIRS_Radiance_Squared_Sum').get()
total_weights += sd.select('VIIRS_Total_Weights').get()
weighted_sum += sd.select('VIIRS_Weighted_Sum').get()
flags |= sd.select('VIIRS_Flags').get()
# initialize output arrays
......@@ -42,19 +38,13 @@ rad_mean = np.empty(n.shape, np.float64)
rad_std = np.empty(n.shape, np.float64)
bt_mean = np.empty(n.shape, np.float64)
std_bt = np.empty(n.shape, np.float64)
weighted_mean = np.empty(n.shape, np.float64)
bt_weighted_mean = np.empty(n.shape, np.float64)
# calculate radiance mean and standard deviation
rad_mean[:] = rads_sum / n
rad_std[:] = np.sqrt((rads_squared_sum / n) - (rad_mean * rad_mean))
# calculate weighted mean
weighted_mean[:] = weighted_sum / total_weights
# now convert the radiances to brightness temperature
bt_mean[:] = viirs_bright(rad_mean)
bt_weighted_mean[:] = viirs_bright(weighted_mean)
# to get a feel for the standard deviation in terms of
# brightness temperature, see what BT difference results from
......@@ -68,18 +58,14 @@ rad_mean[missing_idxs] = -9999.0
rad_std[missing_idxs] = -9999.0
bt_mean[missing_idxs] = -9999.0
std_bt[missing_idxs] = -9999.0
weighted_mean[missing_idxs] = -9999.0
bt_weighted_mean[missing_idxs] = -9999.0
flags[missing_idxs] |= MISSING
# output to HDF
writer = HdfWriter(output_file)
writer.write('VIIRS_N', n)
writer.write('VIIRS_Radiance_Unweighted', rad_mean, fill=-9999.0)
writer.write('VIIRS_Radiance', weighted_mean, fill=-9999.0)
writer.write('VIIRS_Radiance', rad_mean, fill=-9999.0)
writer.write('VIIRS_Radiance_Std', rad_std, fill=-9999.0)
writer.write('VIIRS_Brightness_Temp_Unweighted', bt_mean, fill=-9999.0)
writer.write('VIIRS_Brightness_Temp', bt_weighted_mean, fill=-9999.0)
writer.write('VIIRS_Brightness_Temp', bt_mean, fill=-9999.0)
writer.write('VIIRS_Std_Brightness', std_bt, fill=-9999.0)
writer.write('VIIRS_Flags', flags)
......@@ -18,7 +18,7 @@ def read_viirs_radiances(sdr_filenames):
for i, sdr_filename in enumerate(sdr_filenames):
f = h5py.File(sdr_filename)
f = h5py.File(sdr_filename, 'r')
g = f['/All_Data']
g, = g.values()
......@@ -49,7 +49,6 @@ coll_sd = SD.SD(coll_filename)
num_coll = coll_sd.select('Number_Of_VIIRS_FOV')[:]
viirs_along_idx = coll_sd.select('VIIRS_Along_Track_IDX')[:]
viirs_across_idx = coll_sd.select('VIIRS_Across_Track_IDX')[:]
viirs_weight = coll_sd.select('Weights')[:]
num_cris_scans, num_cris_angs, num_cris_dets = num_coll.shape
......@@ -59,10 +58,6 @@ rad_sum_arr = np.empty([3, num_cris_scans, num_cris_angs, num_cris_dets],
np.float64)
rad_squared_sum_arr = np.empty([3, num_cris_scans, num_cris_angs, num_cris_dets],
np.float64)
total_weights_arr = np.empty([3, num_cris_scans, num_cris_angs, num_cris_dets],
np.float64)
weighted_sum_arr = np.empty([3, num_cris_scans, num_cris_angs, num_cris_dets],
np.float64)
flags_arr = np.zeros([3, num_cris_scans, num_cris_angs, num_cris_dets],
np.int8)
......@@ -73,8 +68,6 @@ for cris_row in range(num_cris_scans):
n = n_arr[:,cris_row,cris_col,cris_det]
rads_sum = rad_sum_arr[:,cris_row,cris_col,cris_det]
rads_squared_sum = rad_squared_sum_arr[:,cris_row,cris_col,cris_det]
total_weights = total_weights_arr[:,cris_row,cris_col,cris_det]
weighted_sum = weighted_sum_arr[:,cris_row,cris_col,cris_det]
flags = flags_arr[:,cris_row,cris_col,cris_det]
num = num_coll[cris_row,cris_col,cris_det]
......@@ -82,12 +75,9 @@ for cris_row in range(num_cris_scans):
n[:] = 0
rads_sum[:] = 0.0
rads_squared_sum[:] = 0.0
total_weights[:] = 0.0
weighted_sum[:] = 0.0
continue
along_idx = viirs_along_idx[cris_row,cris_col,cris_det,:num] - 1
across_idx = viirs_across_idx[cris_row,cris_col,cris_det,:num] - 1
weight = viirs_weight[cris_row,cris_col,cris_det,:num]
rads = in_rads[:,along_idx,across_idx]
# if there's any invalid data amongst the collocated VIIRS
......@@ -102,16 +92,11 @@ for cris_row in range(num_cris_scans):
rads[bad_idxs] = 0.0
rads_sum[:] = rads.sum(axis=1)
rads_squared_sum[:] = (rads * rads).sum(axis=1)
total_weights[:] = \
(np.logical_not(bad_idxs) * weight).sum(axis=1)
weighted_sum[:] = (rads * weight).sum(axis=1)
# write results to HDF
writer = HdfWriter(out_filename)
writer.write('VIIRS_N', n_arr)
writer.write('VIIRS_Radiance_Sum', rad_sum_arr)
writer.write('VIIRS_Radiance_Squared_Sum', rad_squared_sum_arr)
writer.write('VIIRS_Total_Weights', total_weights_arr)
writer.write('VIIRS_Weighted_Sum', weighted_sum_arr)
writer.write('VIIRS_Flags', flags_arr)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment