|
1 # General imports |
|
2 import datetime |
|
3 from operator import itemgetter |
|
4 |
|
5 # Science imports |
|
6 import numpy as np |
|
7 import matplotlib as mpl |
|
8 from matplotlib import pyplot as plt |
|
9 import netCDF4 as netcdf |
|
10 |
|
11 # CNR-IMAA specific imports |
|
12 import milos |
|
13 |
|
14 |
|
15 netcdf_format = 'NETCDF3_CLASSIC' # choose one of 'NETCDF3_CLASSIC', 'NETCDF3_64BIT', 'NETCDF4_CLASSIC' and 'NETCDF4' |
|
16 |
|
17 class BaseLidarMeasurement(): |
|
18 """ This is the general measurement object. |
|
19 It is meant to become a general measurement object |
|
20 independent of the input files. |
|
21 |
|
22 Each subclass should implement the following: |
|
23 * the import_file method. |
|
24 * set the "extra_netcdf_parameters" variable to a dictionary that includes the appropriate parameters. |
|
25 |
|
26 You can override the get_PT method to define a custom procedure to get ground temperature and pressure. |
|
27 The one implemented by default is by using the MILOS meteorological station data. |
|
28 |
|
29 """ |
|
30 |
|
31 def __init__(self, filelist= None): |
|
32 self.info = {} |
|
33 self.dimensions = {} |
|
34 self.variables = {} |
|
35 self.channels = {} |
|
36 self.attributes = {} |
|
37 self.files = [] |
|
38 self.dark_measurement = None |
|
39 if filelist: |
|
40 self.import_files(filelist) |
|
41 |
|
42 def import_files(self,filelist): |
|
43 for f in filelist: |
|
44 self.import_file(f) |
|
45 self.update() |
|
46 |
|
47 def import_file(self,filename): |
|
48 raise NotImplementedError('Importing files should be defined in the instrument-specific subclass.') |
|
49 |
|
50 def update(self): |
|
51 ''' |
|
52 Update the the info, variables and dimensions of the lidar measurement based |
|
53 on the information found in the channels. |
|
54 |
|
55 Reading of the scan_angles parameter is not implemented. |
|
56 ''' |
|
57 |
|
58 # Initialize |
|
59 start_time =[] |
|
60 stop_time = [] |
|
61 points = [] |
|
62 all_time_scales = [] |
|
63 channel_name_list = [] |
|
64 |
|
65 # Get the information from all the channels |
|
66 for channel_name, channel in self.channels.items(): |
|
67 channel.update() |
|
68 start_time.append(channel.start_time) |
|
69 stop_time.append(channel.stop_time) |
|
70 points.append(channel.points) |
|
71 all_time_scales.append(channel.time) |
|
72 channel_name_list.append(channel_name) |
|
73 |
|
74 # Find the unique time scales used in the channels |
|
75 time_scales = set(all_time_scales) |
|
76 |
|
77 # Update the info dictionary |
|
78 self.info['start_time'] = min(start_time) |
|
79 self.info['stop_time'] = max(stop_time) |
|
80 self.info['duration'] = self.info['stop_time'] - self.info['start_time'] |
|
81 |
|
82 # Update the dimensions dictionary |
|
83 self.dimensions['points'] = max(points) |
|
84 self.dimensions['channels'] = len(self.channels) |
|
85 # self.dimensions['scan angles'] = 1 |
|
86 self.dimensions['nb_of_time_scales'] = len(time_scales) |
|
87 |
|
88 # Update the variables dictionary |
|
89 # Write time scales in seconds |
|
90 raw_Data_Start_Time = [] |
|
91 raw_Data_Stop_Time = [] |
|
92 |
|
93 for current_time_scale in list(time_scales): |
|
94 raw_start_time = np.array(current_time_scale) - min(start_time) # Time since start_time |
|
95 raw_start_in_seconds = np.array([t.seconds for t in raw_start_time]) # Convert in seconds |
|
96 raw_Data_Start_Time.append(raw_start_in_seconds) # And add to the list |
|
97 # Check if this time scale has measurements every 30 or 60 seconds. |
|
98 |
|
99 duration = self._get_duration(raw_start_in_seconds) |
|
100 |
|
101 raw_stop_in_seconds = raw_start_in_seconds + duration |
|
102 raw_Data_Stop_Time.append(raw_stop_in_seconds) |
|
103 |
|
104 self.variables['Raw_Data_Start_Time']= raw_Data_Start_Time |
|
105 self.variables['Raw_Data_Stop_Time']= raw_Data_Stop_Time |
|
106 |
|
107 # Make a dictionary to match time scales and channels |
|
108 channel_timescales = [] |
|
109 for (channel_name, current_time_scale) in zip(channel_name_list, all_time_scales): |
|
110 # The following lines are PEARL specific. The reason they are here is not clear. |
|
111 # if channel_name =='1064BLR': |
|
112 # channel_name = '1064' |
|
113 for (ts,n) in zip(time_scales, range(len(time_scales))): |
|
114 if current_time_scale == ts: |
|
115 channel_timescales.append([channel_name,n]) |
|
116 self.variables['id_timescale'] = dict(channel_timescales) |
|
117 |
|
118 def _get_duration(self, raw_start_in_seconds): |
|
119 ''' Return the duration for a given time scale. In some files (ex. Licel) this |
|
120 can be specified from the files themselves. In others this must be guessed. |
|
121 |
|
122 ''' |
|
123 # The old method, kept here for reference |
|
124 #dt = np.mean(np.diff(raw_start_in_seconds)) |
|
125 #for d in duration_list: |
|
126 # if abs(dt - d) <15: #If the difference of measuremetns is 10s near the(30 or 60) seconds |
|
127 # duration = d |
|
128 |
|
129 duration = np.diff(raw_start_in_seconds)[0] |
|
130 |
|
131 return duration |
|
132 |
|
133 def subset_by_channels(self, channel_subset): |
|
134 ''' Get a measurement object containing only the channels with names |
|
135 contained in the channel_sublet list ''' |
|
136 |
|
137 m = self.__class__() # Create an object of the same type as this one. |
|
138 m.channels = dict([(channel, self.channels[channel]) for channel |
|
139 in channel_subset]) |
|
140 m.update() |
|
141 return m |
|
142 |
|
143 def subset_by_time(self, start_time, stop_time): |
|
144 |
|
145 if start_time > stop_time: |
|
146 raise ValueError('Stop time should be after start time') |
|
147 |
|
148 if (start_time < self.info['start_time']) or (stop_time > self.info['stop_time']): |
|
149 raise ValueError('The time interval specified is not part of the measurement') |
|
150 |
|
151 m = self.__class__() # Create an object of the same type as this one. |
|
152 for (channel_name, channel) in self.channels.items(): |
|
153 m.channels[channel_name] = channel.subset_by_time(start_time, stop_time) |
|
154 m.update() |
|
155 return m |
|
156 |
|
157 def r_plot(self): |
|
158 #Make a basic plot of the data. |
|
159 #Should include some dictionary with params to make plot stable. |
|
160 pass |
|
161 |
|
162 def r_pdf(self): |
|
163 # Create a pdf report using a basic plot and meta-data. |
|
164 pass |
|
165 |
|
166 def save(self): |
|
167 #Save the current state of the object to continue the analysis later. |
|
168 pass |
|
169 |
|
170 def get_PT(self): |
|
171 ''' Sets the pressure and temperature at station level . |
|
172 The results are stored in the info dictionary. |
|
173 ''' |
|
174 |
|
175 self.info['Temperature'] = 10.0 |
|
176 self.info['Pressure'] = 930.0 |
|
177 |
|
178 |
|
179 def save_as_netcdf(self, filename): |
|
180 """Saves the measurement in the netcdf format as required by the SCC. |
|
181 Input: filename |
|
182 """ |
|
183 params = self.extra_netcdf_parameters |
|
184 needed_parameters = ['Measurement_ID', 'Temperature', 'Pressure'] |
|
185 |
|
186 for parameter in needed_parameters: |
|
187 stored_value = self.info.get(parameter, None) |
|
188 if stored_value is None: |
|
189 raise ValueError('A value needs to be specified for %s' % parameter) |
|
190 |
|
191 |
|
192 dimensions = {'points': 1, |
|
193 'channels': 1, |
|
194 'time': None, |
|
195 'nb_of_time_scales': 1, |
|
196 'scan_angles': 1,} # Mandatory dimensions. Time bck not implemented |
|
197 |
|
198 global_att = {'Measurement_ID': None, |
|
199 'RawData_Start_Date': None, |
|
200 'RawData_Start_Time_UT': None, |
|
201 'RawData_Stop_Time_UT': None, |
|
202 'RawBck_Start_Date': None, |
|
203 'RawBck_Start_Time_UT': None, |
|
204 'RawBck_Stop_Time_UT': None, |
|
205 'Sounding_File_Name': None, |
|
206 'LR_File_Name': None, |
|
207 'Overlap_File_Name': None, |
|
208 'Location': None, |
|
209 'System': None, |
|
210 'Latitude_degrees_north': None, |
|
211 'Longitude_degrees_east': None, |
|
212 'Altitude_meter_asl': None} |
|
213 |
|
214 channel_variables = \ |
|
215 {'channel_ID': (('channels', ), 'i'), |
|
216 'Background_Low': (('channels', ), 'd'), |
|
217 'Background_High': (('channels', ), 'd'), |
|
218 'LR_Input': (('channels', ), 'i'), |
|
219 'DAQ_Range': (('channels', ), 'd'), |
|
220 'Depolarization_Factor': (('channels', ), 'd'), } |
|
221 |
|
222 |
|
223 channels = self.channels.keys() |
|
224 |
|
225 input_values = dict(self.dimensions, **self.variables) |
|
226 |
|
227 # Add some mandatory global attributes |
|
228 input_values['Measurement_ID'] = self.info['Measurement_ID'] |
|
229 input_values['RawData_Start_Date'] = '\'%s\'' % self.info['start_time'].strftime('%Y%m%d') |
|
230 input_values['RawData_Start_Time_UT'] = '\'%s\'' % self.info['start_time'].strftime('%H%M%S') |
|
231 input_values['RawData_Stop_Time_UT'] = '\'%s\'' % self.info['stop_time'].strftime('%H%M%S') |
|
232 |
|
233 # Add some optional global attributes |
|
234 input_values['System'] = params.general_parameters['System'] |
|
235 input_values['Latitude_degrees_north'] = params.general_parameters['Latitude_degrees_north'] |
|
236 input_values['Longitude_degrees_east'] = params.general_parameters['Longitude_degrees_east'] |
|
237 input_values['Altitude_meter_asl'] = params.general_parameters['Altitude_meter_asl'] |
|
238 |
|
239 # Open a netCDF4 file |
|
240 f = netcdf.Dataset(filename,'w', format = netcdf_format) # the format is specified in the begining of the file |
|
241 |
|
242 # Create the dimensions in the file |
|
243 for (d,v) in dimensions.iteritems(): |
|
244 v = input_values.pop(d, v) |
|
245 f.createDimension(d,v) |
|
246 |
|
247 # Create global attributes |
|
248 for (attrib,value) in global_att.iteritems(): |
|
249 val = input_values.pop(attrib,value) |
|
250 if val: |
|
251 exec('f.%s = %s' % (attrib,val)) |
|
252 |
|
253 """ Variables """ |
|
254 # Write the values of fixes channel parameters |
|
255 for (var,t) in channel_variables.iteritems(): |
|
256 temp_v = f.createVariable(var,t[1],t[0]) |
|
257 for (channel, n) in zip(channels, range(len(channels))): |
|
258 temp_v[n] = params.channel_parameters[channel][var] |
|
259 |
|
260 # Write the id_timescale values |
|
261 temp_id_timescale = f.createVariable('id_timescale','i',('channels',)) |
|
262 for (channel, n) in zip(channels, range(len(channels))): |
|
263 temp_id_timescale[n] = self.variables['id_timescale'][channel] |
|
264 |
|
265 # Laser pointing angle |
|
266 temp_v = f.createVariable('Laser_Pointing_Angle','d',('scan_angles',)) |
|
267 temp_v[:] = params.general_parameters['Laser_Pointing_Angle'] |
|
268 |
|
269 # Molecular calculation |
|
270 temp_v = f.createVariable('Molecular_Calc','i') |
|
271 temp_v[:] = params.general_parameters['Molecular_Calc'] |
|
272 |
|
273 # Laser pointing angles of profiles |
|
274 temp_v = f.createVariable('Laser_Pointing_Angle_of_Profiles','i',('time','nb_of_time_scales')) |
|
275 for (time_scale,n) in zip(self.variables['Raw_Data_Start_Time'], |
|
276 range(len(self.variables['Raw_Data_Start_Time']))): |
|
277 temp_v[:len(time_scale), n] = 0 # The lidar has only one laser pointing angle |
|
278 |
|
279 # Raw data start/stop time |
|
280 temp_raw_start = f.createVariable('Raw_Data_Start_Time','i',('time','nb_of_time_scales')) |
|
281 temp_raw_stop = f.createVariable('Raw_Data_Stop_Time','i',('time','nb_of_time_scales')) |
|
282 for (start_time, stop_time,n) in zip(self.variables['Raw_Data_Start_Time'], |
|
283 self.variables['Raw_Data_Stop_Time'], |
|
284 range(len(self.variables['Raw_Data_Start_Time']))): |
|
285 temp_raw_start[:len(start_time),n] = start_time |
|
286 temp_raw_stop[:len(stop_time),n] = stop_time |
|
287 |
|
288 #Laser shots |
|
289 temp_v = f.createVariable('Laser_Shots','i',('time','channels')) |
|
290 for (channel,n) in zip(channels, range(len(channels))): |
|
291 time_length = len(self.variables['Raw_Data_Start_Time'][self.variables['id_timescale'][channel]]) |
|
292 temp_v[:time_length, n] = params.channel_parameters[channel]['Laser_Shots'] |
|
293 |
|
294 #Raw lidar data |
|
295 temp_v = f.createVariable('Raw_Lidar_Data','d',('time', 'channels','points')) |
|
296 for (channel,n) in zip(channels, range(len(channels))): |
|
297 c = self.channels[channel] |
|
298 temp_v[:len(c.time),n, :c.points] = c.matrix |
|
299 |
|
300 self.add_dark_measurements_to_netcdf(f, channels) |
|
301 |
|
302 #Pressure at lidar station |
|
303 temp_v = f.createVariable('Pressure_at_Lidar_Station','d') |
|
304 temp_v[:] = self.info['Pressure'] |
|
305 |
|
306 #Temperature at lidar station |
|
307 temp_v = f.createVariable('Temperature_at_Lidar_Station','d') |
|
308 temp_v[:] = self.info['Temperature'] |
|
309 |
|
310 self.save_netcdf_extra(f) |
|
311 f.close() |
|
312 |
|
313 def add_dark_measurements_to_netcdf(self, f, channels): |
|
314 |
|
315 # Get dark measurements. If it is not given in self.dark_measurement |
|
316 # try to get it using the get_dark_measurements method. If none is found |
|
317 # return without adding something. |
|
318 if self.dark_measurement is None: |
|
319 self.dark_measurement = self.get_dark_measurements() |
|
320 |
|
321 if self.dark_measurement is None: |
|
322 return |
|
323 |
|
324 dark_measurement = self.dark_measurement |
|
325 |
|
326 # Calculate the length of the time_bck dimensions |
|
327 number_of_profiles = [len(c.time) for c in dark_measurement.channels.values()] |
|
328 max_number_of_profiles = np.max(number_of_profiles) |
|
329 |
|
330 # Create the dimension |
|
331 f.createDimension('time_bck', max_number_of_profiles) |
|
332 |
|
333 # Save the dark measurement data |
|
334 temp_v = f.createVariable('Background_Profile','d',('time_bck', 'channels', 'points')) |
|
335 for (channel,n) in zip(channels, range(len(channels))): |
|
336 c = dark_measurement.channels[channel] |
|
337 temp_v[:len(c.time),n, :c.points] = c.matrix |
|
338 |
|
339 # Dark profile start/stop time |
|
340 temp_raw_start = f.createVariable('Raw_Bck_Start_Time','i',('time','nb_of_time_scales')) |
|
341 temp_raw_stop = f.createVariable('Raw_Bck_Stop_Time','i',('time','nb_of_time_scales')) |
|
342 for (start_time, stop_time,n) in zip(dark_measurement.variables['Raw_Data_Start_Time'], |
|
343 dark_measurement.variables['Raw_Data_Stop_Time'], |
|
344 range(len(dark_measurement.variables['Raw_Data_Start_Time']))): |
|
345 temp_raw_start[:len(start_time),n] = start_time |
|
346 temp_raw_stop[:len(stop_time),n] = stop_time |
|
347 |
|
348 # Dark measurement start/stop time |
|
349 f.RawBck_Start_Date = dark_measurement.info['start_time'].strftime('%Y%m%d') |
|
350 f.RawBck_Start_Time_UT = dark_measurement.info['start_time'].strftime('%H%M%S') |
|
351 f.RawBck_Stop_Time_UT = dark_measurement.info['stop_time'].strftime('%H%M%S') |
|
352 |
|
353 |
|
354 |
|
355 def save_netcdf_extra(self, f): |
|
356 pass |
|
357 |
|
358 def _gettime(self, date_str, time_str): |
|
359 t = datetime.datetime.strptime(date_str+time_str,'%d/%m/%Y%H.%M.%S') |
|
360 return t |
|
361 |
|
362 def plot(self): |
|
363 for channel in self.channels: |
|
364 self.channels[channel].plot(show_plot = False) |
|
365 plt.show() |
|
366 |
|
367 def get_dark_measurements(self): |
|
368 return None |
|
369 |
|
370 |
|
371 class Lidar_channel: |
|
372 |
|
373 def __init__(self,channel_parameters): |
|
374 c = 299792458 #Speed of light |
|
375 self.wavelength = channel_parameters['name'] |
|
376 self.name = str(self.wavelength) |
|
377 self.binwidth = float(channel_parameters['binwidth']) # in microseconds |
|
378 self.data = {} |
|
379 self.resolution = self.binwidth * c / 2 |
|
380 self.z = np.arange(len(channel_parameters['data'])) * self.resolution + self.resolution/2.0 # Change: add half bin in the z |
|
381 self.points = len(channel_parameters['data']) |
|
382 self.rc = [] |
|
383 self.duration = 60 |
|
384 |
|
385 def calculate_rc(self): |
|
386 background = np.mean(self.matrix[:,4000:], axis = 1) #Calculate the background from 30000m and above |
|
387 self.rc = (self.matrix.transpose()- background).transpose() * (self.z **2) |
|
388 |
|
389 |
|
390 def update(self): |
|
391 self.start_time = min(self.data.keys()) |
|
392 self.stop_time = max(self.data.keys()) + datetime.timedelta(seconds = self.duration) |
|
393 self.time = tuple(sorted(self.data.keys())) |
|
394 sorted_data = sorted(self.data.iteritems(), key=itemgetter(0)) |
|
395 self.matrix = np.array(map(itemgetter(1),sorted_data)) |
|
396 |
|
397 def _nearest_dt(self,dtime): |
|
398 margin = datetime.timedelta(seconds = 300) |
|
399 if ((dtime + margin) < self.start_time)| ((dtime - margin) > self.stop_time): |
|
400 print "Requested date not covered in this file" |
|
401 raise |
|
402 dt = abs(self.time - np.array(dtime)) |
|
403 dtmin = min(dt) |
|
404 |
|
405 if dtmin > datetime.timedelta(seconds = 60): |
|
406 print "Nearest profile more than 60 seconds away. dt = %s." % dtmin |
|
407 ind_t = np.where(dt == dtmin) |
|
408 ind_a= ind_t[0] |
|
409 if len(ind_a) > 1: |
|
410 ind_a = ind_a[0] |
|
411 chosen_time = self.time[ind_a] |
|
412 return chosen_time, ind_a |
|
413 |
|
414 def subset_by_time(self, start_time, stop_time): |
|
415 |
|
416 time_array = np.array(self.time) |
|
417 condition = (time_array >= start_time) & (time_array <= stop_time) |
|
418 |
|
419 subset_time = time_array[condition] |
|
420 subset_data = dict([(c_time, self.data[c_time]) for c_time in subset_time]) |
|
421 |
|
422 #Create a list with the values needed by channel's __init__() |
|
423 parameters_values = {'name': self.wavelength, |
|
424 'binwidth': self.binwidth, |
|
425 'data': subset_data[subset_time[0]],} |
|
426 |
|
427 c = Lidar_channel(parameters_values) |
|
428 c.data = subset_data |
|
429 c.update() |
|
430 return c |
|
431 |
|
432 |
|
433 def profile(self,dtime, signal_type = 'rc'): |
|
434 t, idx = self._nearest_dt(dtime) |
|
435 if signal_type == 'rc': |
|
436 data = self.rc |
|
437 else: |
|
438 data = self.matrix |
|
439 |
|
440 prof = data[idx,:][0] |
|
441 return prof, t |
|
442 |
|
443 def get_slice(self, starttime, endtime, signal_type = 'rc'): |
|
444 if signal_type == 'rc': |
|
445 data = self.rc |
|
446 else: |
|
447 data = self.matrix |
|
448 tim = np.array(self.time) |
|
449 starttime = self._nearest_dt(starttime)[0] |
|
450 endtime = self._nearest_dt(endtime)[0] |
|
451 condition = (tim >= starttime) & (tim <= endtime) |
|
452 sl = data[condition, :] |
|
453 t = tim[condition] |
|
454 return sl,t |
|
455 |
|
456 def av_profile(self, tim, duration = datetime.timedelta(seconds = 0), signal_type = 'rc'): |
|
457 starttime = tim - duration/2 |
|
458 endtime = tim + duration/2 |
|
459 d,t = self.get_slice(starttime, endtime, signal_type = signal_type) |
|
460 prof = np.mean(d, axis = 0) |
|
461 tmin = min(t) |
|
462 tmax = max(t) |
|
463 tav = tmin + (tmax-tmin)/2 |
|
464 return prof,(tav, tmin,tmax) |
|
465 |
|
466 def plot(self, signal_type = 'rc', filename = None, zoom = [0,12000,0,-1], show_plot = True, cmap = plt.cm.jet): |
|
467 #if filename is not None: |
|
468 # matplotlib.use('Agg') |
|
469 |
|
470 fig = plt.figure() |
|
471 ax1 = fig.add_subplot(111) |
|
472 self.draw_plot(ax1, cmap = cmap, signal_type = signal_type, zoom = zoom) |
|
473 ax1.set_title("%s signal - %s" % (signal_type.upper(), self.name)) |
|
474 |
|
475 if filename is not None: |
|
476 pass |
|
477 #plt.savefig(filename) |
|
478 else: |
|
479 if show_plot: |
|
480 plt.show() |
|
481 #plt.close() ??? |
|
482 |
|
483 def draw_plot(self,ax1, cmap = plt.cm.jet, signal_type = 'rc', zoom = [0,12000,0,-1]): |
|
484 |
|
485 if signal_type == 'rc': |
|
486 if len(self.rc) == 0: |
|
487 self.calculate_rc() |
|
488 data = self.rc |
|
489 else: |
|
490 data = self.matrix |
|
491 |
|
492 hmax_idx = self.index_at_height(zoom[1]) |
|
493 |
|
494 ax1.set_ylabel('Altitude (km)') |
|
495 ax1.set_xlabel('Time UTC') |
|
496 #y axis in km, xaxis /2 to make 30s measurements in minutes. Only for 1064 |
|
497 #dateFormatter = mpl.dates.DateFormatter('%H.%M') |
|
498 #hourlocator = mpl.dates.HourLocator() |
|
499 |
|
500 #dayFormatter = mpl.dates.DateFormatter('\n\n%d/%m') |
|
501 #daylocator = mpl.dates.DayLocator() |
|
502 hourFormatter = mpl.dates.DateFormatter('%H.%M') |
|
503 hourlocator = mpl.dates.AutoDateLocator(interval_multiples=True) |
|
504 |
|
505 |
|
506 #ax1.axes.xaxis.set_major_formatter(dayFormatter) |
|
507 #ax1.axes.xaxis.set_major_locator(daylocator) |
|
508 ax1.axes.xaxis.set_major_formatter(hourFormatter) |
|
509 ax1.axes.xaxis.set_major_locator(hourlocator) |
|
510 |
|
511 |
|
512 ts1 = mpl.dates.date2num(self.start_time) |
|
513 ts2 = mpl.dates.date2num(self.stop_time) |
|
514 |
|
515 |
|
516 im1 = ax1.imshow(data.transpose()[zoom[0]:hmax_idx,zoom[2]:zoom[3]], |
|
517 aspect = 'auto', |
|
518 origin = 'lower', |
|
519 cmap = cmap, |
|
520 #vmin = 0, |
|
521 vmin = data[:,10:400].max() * 0.1, |
|
522 #vmax = 1.4*10**7, |
|
523 vmax = data[:,10:400].max() * 0.9, |
|
524 extent = [ts1,ts2,self.z[zoom[0]]/1000.0, self.z[hmax_idx]/1000.0], |
|
525 ) |
|
526 |
|
527 cb1 = plt.colorbar(im1) |
|
528 cb1.ax.set_ylabel('a.u.') |
|
529 |
|
530 def index_at_height(self, height): |
|
531 idx = np.array(np.abs(self.z - height).argmin()) |
|
532 if idx.size >1: |
|
533 idx =idx[0] |
|
534 return idx |
|
535 |
|
536 def netcdf_from_files(LidarClass, filename, files, channels, measurement_ID): |
|
537 #Read the lidar files and select channels |
|
538 temp_m = LidarClass(files) |
|
539 m = temp_m.subset_by_channels(channels) |
|
540 m.get_PT() |
|
541 m.info['Measurement_ID'] = measurement_ID |
|
542 m.save_as_netcdf(filename) |
|
543 |