-
Notifications
You must be signed in to change notification settings - Fork 3
Expand file tree
/
Copy pathcresis_dataaggregator.py
More file actions
142 lines (118 loc) · 6.49 KB
/
cresis_dataaggregator.py
File metadata and controls
142 lines (118 loc) · 6.49 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
import os
from collections.abc import Iterable
import numpy as np
################## NDH Tools self imports
###########################################################
from .loadmat import loadmat
from .polarstereo_fwd import polarstereo_fwd
from .savemat import savemat
from .within import within
###########################################################
def cresis_dataaggregator(filelist,remove_totaldata=0,savename='',depthcap=0,at_samples=[],at_samples_type=0):
"""
% (C) Nick Holschuh - Amherst College -- 2022 ([email protected])
%
% This function can take a list of CReSIS Radar files and extract the relevent bits
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The inputs are:
% filelist - list of filenames to read
% remove_totaldata - flag [0 or 1] which dictates whether or not to preserve the radargrams
% savename - a string for the name of the .mat file to write
% depthcap - an index for the maximum depth sample to include
% at_samples - list of lists, including the start and end along-track sample to use, or a an
% outline to act as bouunds for the aggregated data
% at_samples_type - 0 is defined indecies, 1 is an outline
%
%%%%%%%%%%%%%%%
% The outputs are:
% save_dict - dictionary containing:
% Data_Vals: The array with coordinate and profile information
% DV_Info: Metadata describing the columns in Data_Vals
% start_indecies: The index within the larger data array when each new file starts
% filenames: The original filenames included in the aggregation
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
"""
print_flag = 0
start_flag = 1
filenames = []
Aggregated_Data = []
start_indecies = [0]
Data_Vals = []
###################### We loop through the list of input files
for f_ind,fn in enumerate(filelist):
radar_data = loadmat(fn)
final_varlist = []
################## We create a Bottom object if needed
if 'Surface' not in radar_data.keys():
radar_data['Surface'] = radar_data['Latitude'].copy()*np.NaN
if 'Bottom' not in radar_data.keys():
radar_data['Bottom'] = radar_data['Surface'].copy()*np.NaN
################## Calculate Polarstereo coordinates
xy = polarstereo_fwd(radar_data['Latitude'], radar_data['Longitude'])
radar_data['x'] = xy['x']
radar_data['y'] = xy['y']
radar_data['distance'] = polarstereo_fwd(radar_data['x'],radar_data['y'])
final_varlist.append('x')
final_varlist.append('y')
################## We apply the depthcap, if prescribed
if depthcap > 0:
radar_data['Data'] = radar_data['Data'][:depthcap,:]
radar_data['Time'] = radar_data['Time'][:depthcap]
################## This object is used to subset along tack, if requested
if len(at_samples) == 0:
trace_index = np.arange(0,len(radar_data['Latitude']))
else:
if at_samples_type == 0:
trace_index = np.arange(at_samples[f_ind][0],at_samples[f_ind][1]+1)
else:
trace_index = np.where(within(np.stack([radar_data['x'],radar_data['y']]).T,at_samples))[0]
radar_data['Data'] = radar_data['Data'][:,trace_index]
################## All objects that are the same shape as latitude get subset
orig_len = len(radar_data['Latitude'])
for key_opt in radar_data.keys():
if isinstance(radar_data[key_opt],type(radar_data['Latitude'])):
shape_array = np.array(radar_data[key_opt].shape)
if len(shape_array) > 0:
if np.max(shape_array == orig_len) == 1:
final_varlist.append(key_opt)
for kk in final_varlist:
if len(radar_data[kk]) == orig_len:
radar_data[kk] = radar_data[kk][trace_index]
################## Here we extract date info from the filename
fn_parts = fn.split('/')[-1].split('.')[0].split('_')
year = int(fn_parts[1][0:4])
month = int(fn_parts[1][4:6])
day = int(fn_parts[1][6:8])
seg = int(fn_parts[2])
frm = int(fn_parts[3])
file_ind = f_ind
################## Then, we try and construct the object and concatenate everything. Some files had trouble with this
temp_Data_Vals = np.vstack((radar_data['x'],radar_data['y'],radar_data['Latitude'],radar_data['Longitude'],radar_data['Elevation'],
radar_data['Surface'],radar_data['Bottom'],radar_data['GPS_time'],trace_index,
np.ones(radar_data['Latitude'].shape)*year, np.ones(radar_data['Latitude'].shape)*month, np.ones(radar_data['Latitude'].shape)*day,
np.ones(radar_data['Latitude'].shape)*seg, np.ones(radar_data['Latitude'].shape)*frm, np.ones(radar_data['Latitude'].shape)*file_ind)).T
################## The initial file starts the objects, subsequent files add to them
if start_flag > 1:
Data_Vals = np.concatenate((Data_Vals, temp_Data_Vals), axis=0)
if print_flag == 1:
print('Completed File '+str(f_ind)+' - '+fn)
else:
start_indecies = [0]
Data_Vals = temp_Data_Vals
start_flag = start_flag+1
if print_flag == 1:
print('Started with file '+str(f_ind)+' - '+fn)
################## The start index for each new file is logged, along with the file name
start_indecies.append(start_indecies[-1] + len(trace_index))
filenames.append(fn)
################## Finally, the radargrams are appended if desired
if remove_totaldata == 0:
Aggregated_Data.append([radar_data['Data'],radar_data['distance'],radar_data['Time']])
start_indecies = start_indecies[:-1]
DV_Info = ['X Coordinate (ps)','Y Coordinate (ps)','Latitude','Longitude','Flight Elevation','Surface Pick','Bottom Pick','GPS Time','Trace Index','Year','Month','Day','Segment','Frame','File Index']
save_dict = {'Data_Vals':Data_Vals,'DV_Info':DV_Info,'start_indecies':start_indecies,'filenames':filenames,'Aggregated_Data':Aggregated_Data}
if len(savename) > 0:
savemat(save_dict,savename)
return save_dict