Mercurial > public > atmospheric_lidar / changeset
changeset
Updated documentation and cleanup of the LidarChannel class.
Wed, 06 Dec 2017 17:14:02 +0200
- author
- Iannis <ulalume3@yahoo.com>
- date
- Wed, 06 Dec 2017 17:14:02 +0200
- changeset 93
- c352254b4650
- parent 92
- 6d26002aaeed
- child 94
- 7e2ab78dffd9
Updated documentation and cleanup of the LidarChannel class.
| atmospheric_lidar/generic.py | file | annotate | diff | comparison | revisions |
--- a/atmospheric_lidar/generic.py Wed Dec 06 13:36:42 2017 +0200 +++ b/atmospheric_lidar/generic.py Wed Dec 06 17:14:02 2017 +0200 @@ -21,7 +21,7 @@ You can override the set_PT method to define a custom procedure to get ground temperature and pressure. - The class assumes that the input files are consequtive, i.e. there are no measurements gaps. + The class assumes that the input files are consecutive, i.e. there are no measurements gaps. """ def __init__(self, file_list=None): """ @@ -237,7 +237,7 @@ def subset_by_bins(self, b_min=0, b_max=None): """ - Remove some height bins from the file. + Remove some height bins from the measurement data. This could be needed to remove acquisition artifacts at the first or last bins of the profiles. @@ -650,8 +650,26 @@ class LidarChannel: + """ + This class represents a general measurement channel, independent of the input files. + + The class assumes that the input files are consecutive, i.e. there are no measurements gaps. + """ + def __init__(self, channel_parameters): - c = 299792458 # Speed of light + """ + This is run when creating a new object. + + The input dictionary should contain 'name', 'binwidth', and 'data' keys. + + Parameters + ---------- + channel_parameters : dict + A dict containing channel parameters. + """ + # TODO: Change channel_parameters to explicit keyword arguments? + + c = 299792458. # Speed of light self.wavelength = channel_parameters['name'] self.name = str(self.wavelength) self.binwidth = float(channel_parameters['binwidth']) # in microseconds @@ -664,35 +682,80 @@ self.duration = 60 def calculate_rc(self, idx_min=4000, idx_max=5000): + """ Calculate range corrected signal. + + The background signal is estimated as the mean signal between idx_min and idx_max. + + The calculations is stored in the self.rc parameters. + + Parameters + ---------- + idx_min : int + Minimum index to calculate background signal. + idx_max : int + Maximum index to calculate background signal. + """ background = np.mean(self.matrix[:, idx_min:idx_max], axis=1) # Calculate the background from 30000m and above self.rc = (self.matrix.transpose() - background).transpose() * (self.z ** 2) def update(self): + """ + Update the time parameters and data according to the raw input data. + """ self.start_time = min(self.data.keys()) self.stop_time = max(self.data.keys()) + datetime.timedelta(seconds=self.duration) self.time = tuple(sorted(self.data.keys())) sorted_data = sorted(self.data.iteritems(), key=itemgetter(0)) self.matrix = np.array(map(itemgetter(1), sorted_data)) - def _nearest_dt(self, dtime): - margin = datetime.timedelta(seconds=300) - if ((dtime + margin) < self.start_time) | ((dtime - margin) > self.stop_time): + def _nearest_datetime(self, input_time): + """ + Find the profile datetime and index that is nearest to the given datetime. + + Parameters + ---------- + input_time : datetime.datetime + Input datetime object. + + Returns + ------- + profile_datetime : datetime + The datetime of the selected profile. + profile_idx : int + The index of the selected profile. + """ + margin = datetime.timedelta(seconds=self.duration * 5) + + if ((input_time + margin) < self.start_time) | ((input_time - margin) > self.stop_time): logging.error("Requested date not covered in this file") raise ValueError("Requested date not covered in this file") - dt = abs(self.time - np.array(dtime)) - dtmin = min(dt) + dt = np.abs(np.array(self.time) - input_time) + + profile_idx = np.argmin(dt) + profile_datetime = self.time[profile_idx] - if dtmin > datetime.timedelta(seconds=60): - logging.warning("Nearest profile more than 60 seconds away. dt = %s." % dtmin) - ind_t = np.where(dt == dtmin) - ind_a = ind_t[0] - if len(ind_a) > 1: - ind_a = ind_a[0] - chosen_time = self.time[ind_a] - return chosen_time, ind_a + dt_min = dt[profile_idx] + if dt_min > datetime.timedelta(seconds=self.duration): + logging.warning("Nearest profile more than 60 seconds away. dt = %s." % dt_min) + + return profile_datetime, profile_idx def subset_by_time(self, start_time, stop_time): + """ + Subset the channel for a specific time interval. + Parameters + ---------- + start_time : datetime + The starting datetime to subset. + stop_time : datetime + The stopping datetime to subset. + + Returns + ------- + c : LidarChannel object + A new channel object + """ time_array = np.array(self.time) condition = (time_array >= start_time) & (time_array <= stop_time) @@ -706,16 +769,32 @@ # We should use __class__ instead of class name, so that this works properly # with subclasses - # Ex: c = self.__class__(parameters_values) + # Eg: c = self.__class__(parameters_values) # This does not currently work with Licel files though + # TODO: Fix this! c = LidarChannel(parameter_values) c.data = subset_data c.update() return c def subset_by_bins(self, b_min=0, b_max=None): - """Remove some height bins from the file. This could be needed to - remove acquisition artifacts at the start or the end of the files. + """ + Remove some height bins from the channel data. + + This could be needed to remove acquisition artifacts at + the first or last bins of the profiles. + + Parameters + ---------- + b_min : int + The fist valid data bin + b_max : int or None + The last valid data bin. If equal to None, all bins are used. + + Returns + ------- + m : LidarChannel object + A new channel object """ subset_data = {} @@ -734,50 +813,106 @@ c.update() return c - def profile(self, dtime, signal_type='rc'): - t, idx = self._nearest_dt(dtime) - if signal_type == 'rc': + def get_profile(self, input_datetime, range_corrected=True): + """ Returns a single profile, that is nearest to the input datetime. + + Parameters + ---------- + input_datetime : datetime + The required datetime of the profile + range_corrected : bool + If True, the range corrected profile is returned. Else, normal signal is returned. + + Returns + ------- + profile : ndarray + The profile nearest to input_datetime. + t : datetime + The actual datetime corresponding to the profile. + """ + t, idx = self._nearest_datetime(input_datetime) + if range_corrected: data = self.rc else: data = self.matrix - prof = data[idx, :][0] - return prof, t + profile = data[idx, :][0] + + return profile, t + + def get_slice(self, start_time, stop_time, range_corrected=True): + """ Returns a slice of data, between the provided start and stop time. - def get_slice(self, starttime, endtime, signal_type='rc'): - if signal_type == 'rc': + Parameters + ---------- + start_time : datetime + The starting time of the slice + stop_time : datetime + The stoping time of the slice + range_corrected : bool + If True, the range corrected profile is returned. Else, normal signal is returned. + + Returns + ------- + slice : ndarray + The slice of profiles. + slice_datetimes : ndarray of datetime + The actual datetimes corresponding to the slice. + """ + if range_corrected: data = self.rc else: data = self.matrix - tim = np.array(self.time) - starttime = self._nearest_dt(starttime)[0] - endtime = self._nearest_dt(endtime)[0] - condition = (tim >= starttime) & (tim <= endtime) - sl = data[condition, :] - t = tim[condition] - return sl, t + + time_array = np.array(self.time) + start_time = self._nearest_datetime(start_time)[0] + stop_time = self._nearest_datetime(stop_time)[0] - def profile_for_duration(self, tim, duration=datetime.timedelta(seconds=0), signal_type='rc'): - """ Calculates the profile around a specific time (tim). """ - starttime = tim - duration / 2 - endtime = tim + duration / 2 - d, t = self.get_slice(starttime, endtime, signal_type=signal_type) - prof = np.mean(d, axis=0) - tmin = min(t) - tmax = max(t) - tav = tmin + (tmax - tmin) / 2 - return prof, (tav, tmin, tmax) + condition = (time_array >= start_time) & (time_array <= stop_time) + + slice = data[condition, :] + slice_datetimes = time_array[condition] + + return slice, slice_datetimes def average_profile(self): - """ Return the average profile (NOT range corrected) for all the duration of the measurement. """ - prof = np.mean(self.matrix, axis=0) - return prof + """ Return the average profile (NOT range corrected) for all the duration of + the measurement. + + Returns + ------- + profile : ndarray + The average profile for all data. + """ + profile = np.mean(self.matrix, axis=0) + return profile - def plot(self, figsize=(8, 4), signal_type='rc', zoom=[0, 12000, 0, -1], show_plot=True, cmap=plt.cm.jet, z0=None, - title=None, vmin=0, vmax=1.3 * 10 ** 7): - # if filename is not None: - # matplotlib.use('Agg') - + def plot(self, figsize=(8, 4), signal_type='rc', zoom=[0, 12000, 0, -1], show_plot=True, + cmap=plt.cm.jet, z0=None, title=None, vmin=0, vmax=1.3 * 10 ** 7): + """ + Plot of the channel data. + + Parameters + ---------- + figsize : tuple + (width, height) of the output figure (inches) + signal_type : str + If 'rc', the range corrected signal is ploted. Else, the raw signals are used. + zoom : list + A four elemet list defined as [xmin, xmax, ymin, ymax]. Use ymin=0, ymax=-1 to plot full range. + show_plot : bool + If True, the show_plot command is run. + cmap : cmap + An instance of a matplotlib colormap to use. + z0 : float + The ground-level altitude. If provided the plot shows altitude above sea level. + title : str + Optional title for the plot. + vmin : float + Minimum value for the color scale. + vmax : float + Maximum value for the color scale. + """ fig = plt.figure(figsize=figsize) ax1 = fig.add_subplot(111) @@ -795,7 +930,32 @@ zoom=[0, 12000, 0, -1], z0=None, add_colorbar=True, cmap_label='a.u.', cb_format=None, vmin=0, vmax=1.3 * 10 ** 7): - + """ + Draw channel data on the given axis. + + Parameters + ---------- + ax1 : axis object + The axis object to draw. + cmap : cmap + An instance of a matplotlib colormap to use. + signal_type : str + If 'rc', the range corrected signal is ploted. Else, the raw signals are used. + zoom : list + A four elemet list defined as [xmin, xmax, ymin, ymax]. Use ymin=0, ymax=-1 to plot full range. + z0 : float + The ground-level altitude. If provided the plot shows altitude above sea level. + add_colorbar : bool + If True, a colorbar will be added to the plot. + cmap_label : str + Label for the colorbar. Ignored if add_colorbar is False. + cb_format : str + + vmin : float + Minimum value for the color scale. + vmax : float + Maximum value for the color scale. + """ if signal_type == 'rc': if len(self.rc) == 0: self.calculate_rc() @@ -862,6 +1022,14 @@ date_labels=False, vmin=0, vmax=1.3 * 10 ** 7): + """ + Updated drawing routine, using pcolormesh. It is slower but provides more flexibility / precision + in time plotting. The 2 drawing routines should be merged. + + Check draw_plot method for the meaning of the input arguments. + """ + # TODO: Merge the two drawing routines. + if signal_type == 'rc': if len(self.rc) == 0: self.calculate_rc() @@ -934,6 +1102,19 @@ cb1.ax.yaxis.get_offset_text().set_fontsize(cb_font_size) def _index_at_height(self, height): + """ + Get the altitude index nearest to the specified height. + + Parameters + ---------- + height : float + Height (m) + + Returns + ------- + idx : int + Index corresponding to the provided height. + """ idx = np.array(np.abs(self.z - height).argmin()) if idx.size > 1: idx = idx[0]
