Source code for act.corrections.mpl

"""
This module contains corrections for micropulse lidars

"""
import warnings

import numpy as np
import xarray as xr


[docs]def correct_mpl( ds, co_pol_var_name='signal_return_co_pol', cross_pol_var_name='signal_return_cross_pol', co_pol_afterpuls_var_name='afterpulse_correction_co_pol', cross_pol_afterpulse_var_name='afterpulse_correction_cross_pol', overlap_corr_var_name='overlap_correction', overlap_corr_heights_var_name='overlap_correction_heights', range_bins_var_name='range_bins', height_var_name='height', ratio_var_name='cross_co_ratio', ): """ This procedure corrects MPL data: 1.) Throw out data before laser firing (heights < 0). 2.) Remove background signal. 3.) Afterpulse Correction - Subtraction of (afterpulse-darkcount). NOTE: Currently the Darkcount in VAPS is being calculated as the afterpulse at ~30km. But that might not be absolutely correct and we will likely start providing darkcount profiles ourselves along with other corrections. 4.) Range Correction. 5.) Overlap Correction (Multiply). If the height variable changes between netCDF files, Xarray will turn the height dimention variables from a 1D to a 2D array. This will cause issues with processing and other data manipulation. To fix this the 2D height will be converted to a 1D array by using the median value for each height value. Note: Deadtime and darkcount corrections are not being applied yet. Parameters ---------- ds : xarray.Dataset The Xarray Dataset containing data co_pol_var_name : str The Co-polar variable name used in Dataset cross_pol_var_name : str The Cross-polar variable name used in Dataset co_pol_afterpuls_var_name : str The Co-polar after pulse variable name used in Dataset cross_pol_afterpulse_var_name : str The Cross-polar afterpulse variable name used in Dataset overlap_corr_var_name : str The overlap correction variable name used in Dataset overlap_corr_heights_var_name : str The overlap correction height variable name used in Dataset range_bins_var_name : str The range bins variable name used in Dataset height_var_name : str The height variable name used in Dataset ratio_var_name : str The variable name to use for newly created variable of co/cross ratio. Set to None if do not want this new variables created in Dataset. Returns ------- ds : xarray.Dataset Xarray dataset containing the corrected values. The original Xarray Dataset passed in is modified. """ data_dims = ds[co_pol_var_name].dims # Overlap Correction Variable op = ds[overlap_corr_var_name].values if len(op.shape) > 1: op = op[0, :] op_height = ds[overlap_corr_heights_var_name].values if len(op_height.shape) > 1: op_height = op_height[0, :] # Check if height has dimentionality of time and height. If so reduce # height to only dimentionality of height in the dataset before removing # values less than 0. if len(ds[height_var_name].shape) > 1: reduce_dim_name = {'time'} & set(ds[height_var_name].dims) ds[height_var_name] = ds[height_var_name].reduce( func=np.median, dim=reduce_dim_name, keep_attrs=True ) # 1 - Remove negative height data ds = ds.where(ds[height_var_name].load() > 0.0, drop=True) height = ds[height_var_name].values # Get indices for calculating background if len(ds.height.shape) > 1: ind = [ds.height.shape[1] - 50, ds.height.shape[1] - 2] else: ind = [ds.height.shape[0] - 50, ds.height.shape[0] - 2] # Subset last gates into new dataset dummy = ds.isel(range_bins=xr.DataArray(np.arange(ind[0], ind[1]))) # Turn off warnings warnings.filterwarnings('ignore') # Run through co and cross pol data for corrections co_bg = dummy[co_pol_var_name] co_bg = co_bg.where(co_bg.load() > -9998.0) co_bg = co_bg.mean(dim='dim_0').values x_bg = dummy[cross_pol_var_name] x_bg = x_bg.where(x_bg.load() > -9998.0) x_bg = x_bg.mean(dim='dim_0').values # Seems to be the fastest way of removing background signal at the moment co_data = ds[co_pol_var_name].where(ds[co_pol_var_name].load() > 0).values x_data = ds[cross_pol_var_name].where(ds[cross_pol_var_name].load() > 0).values for i in range(len(ds['time'].values)): co_data[i, :] = co_data[i, :] - co_bg[i] x_data[i, :] = x_data[i, :] - x_bg[i] # After Pulse Correction Variable co_ap = ds[co_pol_afterpuls_var_name].values # Fix dimentionality if backwards co_ap_dims = ds[co_pol_afterpuls_var_name].dims if len(co_ap_dims) > 1 and co_ap_dims[::-1] == data_dims: co_ap = np.transpose(co_ap) x_ap = ds[cross_pol_afterpulse_var_name].values # Fix dimentionality if backwards x_ap_dims = ds[cross_pol_afterpulse_var_name].dims if len(x_ap_dims) > 1 and x_ap_dims[::-1] == data_dims: x_ap = np.transpose(x_ap) # Afterpulse Correction co_data = co_data - co_ap x_data = x_data - x_ap # R-Squared Correction co_data = co_data * height**2 x_data = x_data * height**2 # Overlap Correction for j in range(ds[range_bins_var_name].size): if len(height.shape) > 1: idx = (np.abs(op_height - height[0, j])).argmin() else: # Overlap Correction idx = (np.abs(op_height - height[j])).argmin() co_data[:, j] = co_data[:, j] * op[idx] x_data[:, j] = x_data[:, j] * op[idx] # Create the co/cross ratio variable if ratio_var_name is not None: ratio = (x_data / (x_data + co_data)) * 100.0 ds[ratio_var_name] = ds[co_pol_var_name].copy(data=ratio) ds[ratio_var_name].attrs['long_name'] = 'Cross-pol / Co-pol ratio * 100' ds[ratio_var_name].attrs['units'] = '1' try: del ds[ratio_var_name].attrs['ancillary_variables'] del ds[ratio_var_name].attrs['description'] except KeyError: pass # Convert data to decibels co_data = 10.0 * np.log10(co_data) x_data = 10.0 * np.log10(x_data) # Write data to Xarray dataset ds[co_pol_var_name].values = co_data ds[cross_pol_var_name].values = x_data # Update units ds[co_pol_var_name].attrs['units'] = f"10 * log10({ds[co_pol_var_name].attrs['units']})" ds[cross_pol_var_name].attrs['units'] = f"10 * log10({ds[cross_pol_var_name].attrs['units']})" return ds