PollyNET / Pollynet_Processing_Chain

NRT lidar data processing program for multiwavelength polarization Raman lidar network (PollyNET)
https://polly.tropos.de/
GNU General Public License v3.0
20 stars 8 forks source link

Improved filtering of particle depol #80

Closed martin-rdz closed 2 years ago

martin-rdz commented 3 years ago

Problem

The particle depolarization ratio estimate is rather noisy and plotted incompletely. Example from the Punta Arenas measurements is shown in the left two plots below, two features are marked at 2.0 and 5.5 km. The retrieved data as saved in the netcdf is shown on the right (blue curve in the 3rd plot and orange curve in the 4th plot).

grafik

The exploding values in the particle depol are caused by a singularity in the formula for low ratios of molecular_bsc/aersol_backscatter. When the aerosol backscatter is smaller than than the molecular backscatter, this behaviour becomes very pronounced (everything left of green curve in the plot below)

grafik

voldepol_to_pardepol

Proposed solution

Include a filter based on the volume depol, the molecular depol and the ratio molecular_bsc/aerosol bsc. An example is shown in the 4th plot (first figure, green line).

Several steps seem to work rather well, without requiring additional smoothing (only tested for Klett retrieved 532 backscatter so far)

  1. set all values of volume depol < molecular depol to molecular_depol
  2. mask values with very low molecular/aerosol bsc ratio, e.g. everything >40 [optional, but already reduces noise quite a bit]
  3. mask volume depol, where the particle depol formula causes physically unreasonable high values, e.g. particle depol should be below 0.7

Formula for step 3 is grafik with the coefficients a,b depending on the desired particle depol threshold

particle depol thres a1 a2 a3 b1 b2 b3
0.6 6 1 1 5 6 5
0.7 7 2 2 5 7 5
0.8 8 3 3 5 8 5
1.0 2 1 1 1 2 1

Example implementation

In python, as I am not fluent in matlab and not sure how masking is done in the processing chain. Sorry for not including it on my own.

# filter step 1
voldepol[voldepol < moldepol] = moldepol
# filter step 2 (optional)
aerbsc[mol_beta_att/aerbsc > 40] = np.nan

# calculate the particle depol
enum = (voldepol + 1)
denom1 = (mol_beta_att/aerbsc)
denom2 = ((moldepol - voldepol)/(1 + moldepol))
denom = (denom1*denom2 + 1)
pardepol = enum/denom - 1

# filter step 3  for different thresholds
# no pardepol larger 1
# threshold = (2*moldepol*denom1+moldepol+1)/(moldepol+2*denom1+1)
# no pardepol larger 0.8
# threshold = (8*moldepol*denom1+3*moldepol+3)/(5*moldepol+8*denom1+5)
# no pardepol larger 0.7
threshold = (7*moldepol*denom1+2*moldepol+2)/(5*moldepol+7*denom1+5)
# no pardepol larger 0.6
# threshold = (6*moldepol*denom1+1*moldepol+1)/(5*moldepol+6*denom1+5)
pardepol[voldepol>threshold] = np.nan

Likely the two thresholds (40 from step 2 and 0.7 from step 3) should be subject of the config file.

ZPYin commented 3 years ago

@martin-rdz Thanks for your informative comments.

Under current version, the mask for particle depolarization ratio is rather simple with only filtering data points with relative uncertainty larger than 30%. See below

https://github.com/PollyNET/Pollynet_Processing_Chain/blob/facccb3e8f3f0a76d6cce227ddb3576cbe706fed/lib/polly_general_func_lib/pollyxt_display_retrieving.m#L167-L168

But apparently, your method is more suitable with more features remained.

I will have a test with your method first. And hopefully implement it soon...

martin-rdz commented 3 years ago

@ZPYin Thanks for taking care 👍 Looking forward for version 2.0

I was looking into the python version of the plotting routine, but didn't grasp where the filtering was done. /lib/polly_general_func_lib/pollyxt_display_retrieving.py

Would it be possible to also filter the data that goes into the netcdf file? Or at least introduce a flag? Otherwise the end-users would have to calculate the molecular profiles thereself.

ZPYin commented 3 years ago

Some improvements have been implemented in commit: c57dbdc

Improvement of the filter of particle depolarization ratio

I increased the upper limit of the relative uncertainty for filtering noise par-depol points (from 30% to 60%), instead of implementing your criteria. Because, the uncertainty caused by particle backscatter, molecular depol and signal shot noise have been taken into account, based on the error propagation equation. Therefore, your criteria 2 can be realized by introducing a threshold for uncertainty of par-depol.

Your criteria 1 and 3 were adapted to better contrain the results, but with a different threshold. (see below)

https://github.com/PollyNET/Pollynet_Processing_Chain/blob/c57dbdc3563be43c7c06407cc9c325de6b228bc2/lib/polly_general_func_lib/pollyxt_display_retrieving.m#L172-L175

image

Uncertainty of particle depolarization in the netCDF file

To help end-users to use the par-depol in a better way, the uncertainty of par-depol was stored in the **profile** files. Besides, molecule depolarization ratio was saved in the attributes of parDepol***. see below,

image

No keywords for the threshold

Data visualization requires a lot of controls, in order to increase the readability of each figure. However, the framework of Picasso v2.0 is dedicated for data processing instead of data visualization, which makes it hard to be maintained and improved for features of data visualization. That's one of reason that I want to introduce Picasso v3.0.

I will summarize issues associated with data visualization and fix them together in the new version.

HolgerPollyNet commented 3 years ago

@Moritz-TROPOS , wanted to work on a substantial way to calculate the depol for the Polly systems. Thus, the issue might be reopened using another name.....

ZPYin commented 3 years ago

@Moritz-TROPOS , wanted to work on a substantial way to calculate the depol for the Polly systems. Thus, the issue might be reopened using another name.....

Welcome~~~ It's reopened now, feel free to add comments.

With regard to code development, please refer to issue #12 and #15. 🍻

HolgerPollyNet commented 3 years ago

any further comments or can it be closed or moved to discussions? @Moritz-TROPOS @HolgerPollyNet