XENON1T / pax

The XENON1T raw data processor [deprecated]
BSD 3-Clause "New" or "Revised" License
16 stars 15 forks source link

S1 correction map: interpolation and scale factor #229

Closed sanderbreur closed 8 years ago

sanderbreur commented 9 years ago

The s1_lce map (commit f431059a0972c6dbc658ed92f7a8b3f9502d60e3) was intoduced (Thanks @JelleAalbers !) to correct the S1 depending on position (X,R). This is needed to do S1 energy estimation in xenon100 analysis with PAX.

Two issues are found when comparing xerawdp and pax data (radon events):

The xerawdp and PAX S1 corrections in R,Z space: position with s1 correction - xerawdp position with s1 correction - without scale factor - pax

After applying a scale factor of 1.65: position with s1 correction - pax

In the uncorrected S1 and R we now see the same shapes: (Xerawdp has more data points) s1 correction vs uncorrected s1 - xerawdp s1 correction vs uncorrected s1 - pax

s1 correction vs r - xerawdp s1 correction vs r - pax

In Z we see a problem with interpolation ('binning'): s1 correction vs z - xerawdp s1 correction vs drift time - pax

andrewgrahambrown commented 9 years ago

Is it possible to draw: "S1 Corrected PAX vs. S1 corrected Xerawdp"?

If not trivial due to different trees for each variable, you could output event-number : S1 corrected as a text file, then merge the tables...

On Mon, Sep 14, 2015 at 5:41 PM, Sander Breur notifications@github.com wrote:

The s1_lce map (commit f431059 https://github.com/XENON1T/pax/commit/f431059a0972c6dbc658ed92f7a8b3f9502d60e3) was intoduced (Thanks @JelleAalbers https://github.com/JelleAalbers !) to correct the S1 depending on position (X,R). This is needed to do S1 energy estimation in xenon100 analysis with PAX.

Two issues are found when comparing xerawdp and pax data (radon events):

  • There is a general scale factor difference between xerawdp and pax S1 correction. (~1.65)
  • The method of interpolating from points outside the region of known values (convex hull of bin centers of the original map) uses a nearest neighbor approach instead of a real interpolation. This effects only events close to the edge, but the effect is shown in the plots below.

The xerawdp and PAX S1 corrections in R,Z space: [image: position with s1 correction - xerawdp] https://cloud.githubusercontent.com/assets/6515893/9854155/a2371e82-5b07-11e5-8bca-a0976f4feefa.png [image: position with s1 correction - without scale factor - pax] https://cloud.githubusercontent.com/assets/6515893/9854158/a3da6c80-5b07-11e5-94ac-c0fb381d7fd5.png

After applying a scale factor of 1.65: [image: position with s1 correction - pax] https://cloud.githubusercontent.com/assets/6515893/9854160/a7a1d394-5b07-11e5-9c31-85e57bbf0abb.png

In the uncorrected S1 and R we now see the same shapes: (Xerawdp has more data points) [image: s1 correction vs uncorrected s1 - xerawdp] https://cloud.githubusercontent.com/assets/6515893/9854174/b07ddff8-5b07-11e5-803b-3e860e915aed.png [image: s1 correction vs uncorrected s1 - pax] https://cloud.githubusercontent.com/assets/6515893/9854177/b2dcb4cc-5b07-11e5-9cd9-308843db30e0.png

[image: s1 correction vs r - xerawdp] https://cloud.githubusercontent.com/assets/6515893/9854179/b7bcede0-5b07-11e5-8181-5cd0928ba59b.png [image: s1 correction vs r - pax] https://cloud.githubusercontent.com/assets/6515893/9854181/b9a5f4c6-5b07-11e5-9540-99e7b10d8f17.png

In Z we see a problem with interpolation ('binning'): [image: s1 correction vs z - xerawdp] https://cloud.githubusercontent.com/assets/6515893/9854191/c77ce884-5b07-11e5-907c-e76d61e468d1.png [image: s1 correction vs drift time - pax] https://cloud.githubusercontent.com/assets/6515893/9854193/c9b68740-5b07-11e5-975f-4bd91f0c45dd.png

— Reply to this email directly or view it on GitHub https://github.com/XENON1T/pax/issues/229.

JelleAalbers commented 9 years ago

I hope the interpolation problem is fixed now. I guess Xerawdp normalizes the map to give 1 at a specific point in the detector?

tunnell commented 9 years ago

it's 1 at center.

tunnell commented 9 years ago

Wait, I thought I worked with @sanderbreur on this and we found out the issue? The light yield (constant) and correction map (3D) are applied as two steps in xerawdp.

JelleAalbers commented 9 years ago

Yes, you're right; around line 500 of SignalCorrections.cpp. They normalize the correction map so that the average correction is 1 (where each bin of the map counts equally). I've changed the map in the same way now.

JelleAalbers commented 9 years ago

Sander, could you repeat the nice checks you did before? Then we'll know if this is solved.

sanderbreur commented 9 years ago

Good point, will do.

S

On Thu, 22 Oct 2015 11:35 Jelle Aalbers notifications@github.com wrote:

Assigned #229 https://github.com/XENON1T/pax/issues/229 to @sanderbreur https://github.com/sanderbreur.

— Reply to this email directly or view it on GitHub https://github.com/XENON1T/pax/issues/229#event-442449882.

Shintler commented 8 years ago

Hi Paxers!

I was testing the S1 LY correction map in my analysis. For this I used Cecilias map directly and "pythonized" the correction function from xerawdp. This works, but runs very slow (~3h for correcting 20k events). I got in contact with Jelle, since I knew he was implementing the correction map in PAX. My very naive approach was to just C&P the class and use it within my code. After some issues, I got it working but the results were completely different. So I started to investigate why:

Using all my findings (and Jelles hints), I can reproduce and verify your implementation of the correctoin map:

xerawdp_map Fig. 1: Raw map as stored and used in XeRawDP

pax_map Fig. 2: Modified PAX map (modifications see points above)

ly_comparison Fig. 3: Comparison of map entries

difference_mps Fig. 4: Differences of the map entreis

JelleAalbers commented 8 years ago

Thanks for checking this Stefan! You're right, we changed the units and the coordinates. In the latest pax Ly is again stored rather than 1/Ly, and the normalization is taken into account (see #229).

This shows the map in pax uses the same data as Xerawdp does, but we still don't know if the interpolation routine is working correctly; for this we'd need something like Sander's last two plots.

sanderbreur commented 8 years ago

Ok I just reprocessed some high energy radon events and made the comparison plots for XerawDP (27k events) and PAX (23k events):

The changes made a good improvement. The following features are still different:

A difference unrelated to the S1 LCE correction:

Position dependent correction:

position with s1 correction - xerawdp

position with s1 correction - pax

Correction vs uncor S1:

s1 correction vs uncorrected s1 - xerawdp

s1 correction vs uncorrected s1 - pax

\ Correction vs Radius:** s1 correction vs r - xerawdp

s1 correction vs r - pax

Correction vs Depth: s1 correction vs z - xerawdp

s1 correction vs drift time - pax

Shintler commented 8 years ago

I also did another test to compare the interpolations: I generate ~21k events uniformly within the FV (34Kg). After that I overplot the correction of PAX and the pythonized XeRawDP against the z position: ly_comparison_zposition So there is a small difference, with the biggest discrepancies in the bulge or step pattern present in the PAX interpolation. To quantify the difference, I plot the difference of the interpolations for each position normalized to the XeRawDP result: difference_hist

tunnell commented 8 years ago

Thanks for the plot. How should we proceed with these differences

richardsaldanha commented 8 years ago

Hi,

Not sure if this is very relevant, but it could be a nice check on the new S1 correction interpolation in pax.

I found a bug in the S1 3d interpolation calculated by xerawdp for events around the theta = 0 line (see https://xecluster.lngs.infn.it/dokuwiki/doku.php?id=xenon:saldanha:xe100:kr83m#discontinuity_in_ly_correction_interpolation_bug )

The effect is small, but you should see a difference in the estimated correction for xerawdp and pax (or any interpolater without this bug) for events near this position.

Cheers, Richard

On 25 Nov 2015, at 10:12, Christopher Tunnell wrote:

Thanks for the plot. How should we proceed with these differences

— Reply to this email directly or view it on GitHubhttps://github.com/XENON1T/pax/issues/229#issuecomment-159636549.

tunnell commented 8 years ago

Is there anything action-able here @sanderbreur? It seems like there is a percent or so difference, but won't this be smaller than other smearings? It's not clear which interpolation is correct, right?

sanderbreur commented 8 years ago

The general code works, @jaalbers would you change the interpolation code after seeing the difference seen by Stefan?

If not then we can close, we documented the differences now and when the Xe1T maps are checked they can always look trough these results.

Op di 5 jan. 2016 om 12:01 schreef Christopher Tunnell < notifications@github.com>:

Is there anything action-able here @sanderbreur https://github.com/sanderbreur? It seems like there is a percent or so difference, but won't this be smaller than other smearings? It's not clear which interpolation is correct, right?

— Reply to this email directly or view it on GitHub https://github.com/XENON1T/pax/issues/229#issuecomment-168972174.

JelleAalbers commented 8 years ago

The step pattern is due to using 1/distance weighted average interpolation/extrapolation. Cecilia tried this for XENON100/Xerawdp as well, see:

https://xecluster.lngs.infn.it/dokuwiki/doku.php?id=xenon:xenon100:analysis:cecilia:lycorrectionmap#update_8_11_2011

She concluded switching to linear interpolation is better. However, scipy.interpolate.LinearNDInterpolator gives nan if you're outside the convex hull of interpolation points, which doesn't cover the entire TPC. My first solution (before this issue) was to use LinearNDInterpolator inside the hull and NearestNDINterpolator (which does extrapolate) outside. However, there's quite a lot of mass at high radii, so you get the very heavy step pattern seen in Sander's first results.

I think a 1-2% error on the LCE is well within the error on the current light yield map (looking also at the differences with Richard's new map), so for now this is fine. However, we should probably revisit the issue of interpolation for XENON1T maps: 1/d weighted average may only work if the map is sufficiently smooth & highly sampled. If we can't get such a map, we may have to first upsample it using a more intelligent interpolation, then continue to use the 'dumb' interpolation in pax.

tunnell commented 8 years ago

Revisit once we get XENON1T LCE maps.