--- a/generic.py Mon Nov 24 11:35:17 2014 +0200
+++ b/generic.py Mon Nov 24 11:35:40 2014 +0200
@@ -5,6 +5,7 @@
# Science imports
import numpy as np
import matplotlib as mpl
+from matplotlib.ticker import ScalarFormatter
from matplotlib import pyplot as plt
import netCDF4 as netcdf
@@ -151,6 +152,20 @@
m.update()
return m
+ def subset_by_bins(self, b_min = 0, b_max = None):
+ """Remove some height bins from the file. This could be needed to
+ remove aquisition artifacts at the start or the end of the files.
+ """
+
+ m = self.__class__() # Create an object of the same type as this one.
+
+ for (channel_name, channel) in self.channels.items():
+ m.channels[channel_name] = channel.subset_by_bins(b_min, b_max)
+
+ m.update()
+
+ return m
+
def r_plot(self):
#Make a basic plot of the data.
#Should include some dictionary with params to make plot stable.
@@ -390,13 +405,13 @@
return t_mean
-class Lidar_channel:
+class LidarChannel:
- def __init__(self,channel_parameters):
- c = 299792458 #Speed of light
+ def __init__(self, channel_parameters):
+ c = 299792458 # Speed of light
self.wavelength = channel_parameters['name']
self.name = str(self.wavelength)
- self.binwidth = float(channel_parameters['binwidth']) # in microseconds
+ self.binwidth = float(channel_parameters['binwidth']) # in microseconds
self.data = {}
self.resolution = self.binwidth * c / 2
self.z = np.arange(len(channel_parameters['data'])) * self.resolution + self.resolution/2.0 # Change: add half bin in the z
@@ -405,7 +420,7 @@
self.duration = 60
def calculate_rc(self):
- background = np.mean(self.matrix[:,4000:], axis = 1) #Calculate the background from 30000m and above
+ background = np.mean(self.matrix[:,4000:], axis = 1) # Calculate the background from 30000m and above
self.rc = (self.matrix.transpose()- background).transpose() * (self.z **2)
@@ -442,16 +457,39 @@
subset_data = dict([(c_time, self.data[c_time]) for c_time in subset_time])
#Create a list with the values needed by channel's __init__()
- parameters_values = {'name': self.wavelength,
+ parameter_values = {'name': self.wavelength,
'binwidth': self.binwidth,
'data': subset_data[subset_time[0]],}
+
+ # We should use __class__ instead of class name, so that this works properly
+ # with subclasses
+ # Ex: c = self.__class__(parameters_values)
+ # This does not currently work with Licel files though
+ 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 aquisition artifacts at the start or the end of the files.
+ """
+
+ subset_data = {}
+
+ for ctime, cdata in self.data.items():
+ subset_data[ctime] = cdata[b_min:b_max]
+
+ #Create a list with the values needed by channel's __init__()
+ parameters_values = {'name': self.wavelength,
+ 'binwidth': self.binwidth,
+ 'data': subset_data[subset_data.keys()[0]],} # We just need any array with the correct length
- c = Lidar_channel(parameters_values)
+ c = LidarChannel(parameters_values)
c.data = subset_data
c.update()
return c
-
def profile(self,dtime, signal_type = 'rc'):
t, idx = self._nearest_dt(dtime)
if signal_type == 'rc':
@@ -512,7 +550,10 @@
plt.show()
#plt.close() ???
- def draw_plot(self,ax1, cmap = plt.cm.jet, signal_type = 'rc', zoom = [0,12000,0,-1], z0 = None, cmap_label = 'a.u.', cb_format = None, vmin = 0, vmax = 1.3 * 10 ** 7):
+ def draw_plot(self,ax1, cmap = plt.cm.jet, signal_type = 'rc',
+ zoom = [0,12000,0,-1], z0 = None,
+ add_colorbar = True, cmap_label = 'a.u.', cb_format = None,
+ vmin = 0, vmax = 1.3 * 10 ** 7):
if signal_type == 'rc':
if len(self.rc) == 0:
@@ -561,19 +602,101 @@
#vmax = data[:,10:400].max() * 0.9,
extent = [ts1,ts2,self.z[zoom[0]]/1000.0 + z0/1000., self.z[hmax_idx]/1000.0 + z0/1000.],
)
+
+ if add_colorbar:
+ if cb_format:
+ cb1 = plt.colorbar(im1, format = cb_format)
+ else:
+ cb1 = plt.colorbar(im1)
+ cb1.ax.set_ylabel(cmap_label)
+
+ # Make the ticks of the colorbar smaller, two points smaller than the default font size
+ cb_font_size = mpl.rcParams['font.size'] - 2
+ for ticklabels in cb1.ax.get_yticklabels():
+ ticklabels.set_fontsize(cb_font_size)
+ cb1.ax.yaxis.get_offset_text().set_fontsize(cb_font_size)
+
+
+ def draw_plot_new(self, ax1, cmap = plt.cm.jet, signal_type = 'rc',
+ zoom = [0, 12000, 0, None], z0 = None,
+ add_colorbar = True, cmap_label = 'a.u.',
+ cb_format = None, power_limits = (-2, 2),
+ date_labels = False,
+ vmin = 0, vmax = 1.3 * 10 ** 7):
- if cb_format:
- cb1 = plt.colorbar(im1, format = cb_format)
+ if signal_type == 'rc':
+ if len(self.rc) == 0:
+ self.calculate_rc()
+ data = self.rc
else:
- cb1 = plt.colorbar(im1)
- cb1.ax.set_ylabel(cmap_label)
+ data = self.matrix
+
+ hmax_idx = self.index_at_height(zoom[1])
+ hmin_idx = self.index_at_height(zoom[0])
+
+ # If z0 is given, then the plot is a.s.l.
+ if z0:
+ ax1.set_ylabel('Altitude a.s.l. [km]')
+ else:
+ ax1.set_ylabel('Altitude a.g.l. [km]')
+ z0 = 0
+
+ ax1.set_xlabel('Time UTC [hh:mm]')
+ #y axis in km, xaxis /2 to make 30s measurements in minutes. Only for 1064
+ #dateFormatter = mpl.dates.DateFormatter('%H.%M')
+ #hourlocator = mpl.dates.HourLocator()
+
- # Make the ticks of the colorbar smaller, two points smaller than the default font size
- cb_font_size = mpl.rcParams['font.size'] - 2
- for ticklabels in cb1.ax.get_yticklabels():
- ticklabels.set_fontsize(cb_font_size)
- cb1.ax.yaxis.get_offset_text().set_fontsize(cb_font_size)
+ if date_labels:
+ dayFormatter = mpl.dates.DateFormatter('%H:%M\n%d/%m/%Y')
+ daylocator = mpl.dates.AutoDateLocator(minticks = 3, maxticks = 8, interval_multiples=True)
+ ax1.axes.xaxis.set_major_formatter(dayFormatter)
+ ax1.axes.xaxis.set_major_locator(daylocator)
+ else:
+ hourFormatter = mpl.dates.DateFormatter('%H:%M')
+ hourlocator = mpl.dates.AutoDateLocator(minticks = 3, maxticks = 8, interval_multiples=True)
+ ax1.axes.xaxis.set_major_formatter(hourFormatter)
+ ax1.axes.xaxis.set_major_locator(hourlocator)
+
+ # Get the values of the time axis
+ dt = datetime.timedelta(seconds = self.duration)
+
+ time_cut = self.time[zoom[2]:zoom[3]]
+ time_last = time_cut[-1] + dt # The last element needed for pcolormesh
+
+ time_all = time_cut + (time_last,)
+
+ t_axis = mpl.dates.date2num(time_all)
+
+ # Get the values of the z axis
+ z_cut = self.z[hmin_idx:hmax_idx] - self.resolution / 2.
+ z_last = z_cut[-1] + self.resolution
+
+ z_axis = np.append(z_cut, z_last) / 1000. + z0 / 1000. # Convert to km
+
+ # Plot
+ im1 = ax1.pcolormesh(t_axis, z_axis, data.T[hmin_idx:hmax_idx, zoom[2]:zoom[3]],
+ cmap = cmap,
+ vmin = vmin,
+ vmax = vmax,
+ )
+
+ if add_colorbar:
+ if cb_format:
+ cb1 = plt.colorbar(im1, format = cb_format)
+ else:
+ cb_formatter = ScalarFormatter()
+ cb_formatter.set_powerlimits(power_limits)
+ cb1 = plt.colorbar(im1, format = cb_formatter)
+ cb1.ax.set_ylabel(cmap_label)
+
+ # Make the ticks of the colorbar smaller, two points smaller than the default font size
+ cb_font_size = mpl.rcParams['font.size'] - 2
+ for ticklabels in cb1.ax.get_yticklabels():
+ ticklabels.set_fontsize(cb_font_size)
+ cb1.ax.yaxis.get_offset_text().set_fontsize(cb_font_size)
+
def index_at_height(self, height):
idx = np.array(np.abs(self.z - height).argmin())
if idx.size > 1:
