search

spradlin / WCSim

The WCSim GEANT4 application
0 stars 0 forks source link

Interpolation of QE functions required in WCSimDetectorConstruction::CreateCombinedPMTQE() #8

Closed spradlin closed 2 years ago

spradlin commented 2 years ago

Placeholder to collect notes and comments in advance of submitting a pull request to the main repository.

spradlin commented 2 years ago

Have completed the coding. Want to run the CI validation tests described in the Validation README.md. See https://github.com/WCSim/Validation.

spradlin commented 2 years ago

Hmm, the Validation README.md appears to be mostly about updating the references for tests rather than running the validation tests. It may be performed automatically when a PR is made.

Let me move on to collating discussion and material for an issue.

spradlin commented 2 years ago

Description of Issue

Intended behavior

When constructing a detector with multiple PMT types, WCSim invokes WCSimDetectorConstruction::CreateCombinedPMTQE() to calculate a hybrid quantum efficiency function for the combination of types. The purpose of the hybrid QE function is to allow a QE check on photons at emission (in the Stacking Action) so that they may be efficiently discarded before propagation through the bulk of the detector. Because it is not know at emission what type of PMT, if any, the photon will hit, a QE function representative of the whole detector is required in order to make such a check.

The apparent intent of the implementation of WCSimDetectorConstruction::CreateCombinedPMTQE() is to create a hybrid QE function that, for each photon wavelength, is equal to the largest QE of any used PMT type at that wavelength.

Implementation

The QE(lambda) functions for each PMT type are defined as linear interpolations of parallel arrays of lambda_i and QE_i values. The domain array of the hybrid QE function, $lambda_k$, is a union of the domain points for each of the individual PMT types. The range array, QE_k, has the maximum QE_i value for any PMT with lambda_k in its domain array.

Issue

When the lambda_i arrays of the PMT types are not identical, the hybrid QE_k array values can jump wildly among the PMT types in a way that does not represent the maximum QE efficiency at a given wavelength.

spradlin commented 2 years ago

Demonstration of issue

I ran some tests using a .mac file that uses the HyperKWithOD geometry. In the current WCSim/develop branch, this will construct an inner detector with BoxandLine20inchHQE PMTs and an outer detector with PMT8inch PMTs. It calls WCSimDetectorConstruction::CreateCombinedPMTQE() to compute a hybrid QE function for these two PMT types.

WCSimDetectorConstruction::CreateCombinedPMTQE() has debugging output that reports the input QE function arrays and the hybrid arrays that it computes.

Control

I ran my WCSim-HKOD.mac file with the current head of WCSim/develop and collected the debugging output for WCSimDetectorConstruction::CreateCombinedPMTQE() in the table below.

WL ID QE OD QE Combined
280nm 0 0 0
300nm 0.00109589 0.0139 0.0139
320nm 0.171918 0.0854 0.171918
340nm 0.349263 0.169 0.349263
360nm 0.40137 0.203 0.40137
380nm 0.428356 0.206 0.428356
400nm 0.428767 0.211 0.428767
420nm 0.410137 0.202 0.410137
440nm 0.382329 0.188 0.382329
460nm 0.341233 0.167 0.341233
480nm 0.283562 0.14 0.283562
500nm 0.240822 0.116 0.240822
520nm 0.189589 0.0806 0.189589
540nm 0.106712 0.0432 0.106712
560nm 0.0647945 0.0265 0.0647945
580nm 0.0394521 0.0146 0.0394521
600nm 0.020411 0.00756 0.020411
620nm 0.00849315 0.00508 0.00849315
640nm 0.000273973 0.00158 0.00158
660nm 0.000136986 0 0.000136986

It works as intended because the wavelength arrays of the two PMT types are initialized with identical floating point literals.

Problematic case

In order to demonstrate the potential issue, i incremented every entry in the BoxandLine20inchHQE PMT wavelength array by +5nm and reran my WCSim-HKOD.mac. The debugging output from WCSimDetectorConstruction::CreateCombinedPMTQE() for this case is summarized below.

WL ID QE OD QE Combined
280nm - 0 0
285nm 0 - 0
300nm - 0.0139 0.0139
305nm 0.00109589 - 0.00109589
320nm - 0.0854 0.0854
325nm 0.171918 - 0.171918
340nm - 0.169 0.169
345nm 0.349263 - 0.349263
360nm - 0.203 0.203
365nm 0.40137 - 0.40137
380nm - 0.206 0.206
385nm 0.428356 - 0.428356
400nm - 0.211 0.211
405nm 0.428767 - 0.428767
420nm - 0.202 0.202
425nm 0.410137 - 0.410137
440nm - 0.188 0.188
445nm 0.382329 - 0.382329
460nm - 0.167 0.167
465nm 0.341233 - 0.341233
480nm - 0.14 0.14
485nm 0.283562 - 0.283562
500nm - 0.116 0.116
505nm 0.240822 - 0.240822
520nm - 0.0806 0.0806
525nm 0.189589 - 0.189589
540nm - 0.0432 0.0432
545nm 0.106712 - 0.106712
560nm - 0.0265 0.0265
565nm 0.0647945 - 0.0647945
580nm - 0.0146 0.0146
585nm 0.0394521 - 0.0394521
600nm - 0.00756 0.00756
605nm 0.020411 - 0.020411
620nm - 0.00508 0.00508
625nm 0.00849315 - 0.00849315
640nm - 0.00158 0.00158
645nm 0.000273973 - 0.000273973
660nm - 0 0
665nm 0.000136986 - 0.000136986

The result in an interlacing of the two PMT efficiency functions. The combined QE no longer describes the maximum efficiency of any PMT. It vascillates wildly.

spradlin commented 2 years ago

Proposed fix

In MR (cite MR), i updated WCSimDetectorConstruction::CreateCombinedPMTQE() to interpolate the QE function for each PMT type in order to find the maximum QE for the given wavelength. After the update, the Problematic case above has the intended behavior. The debugging output is summarized in the table below.

WL ID QE OD QE Combined
280nm -- 0 0
285nm 0 -- 0.003475
300nm -- 0.0139 0.0139
305nm 0.00109589 -- 0.031775
320nm -- 0.0854 0.129212
325nm 0.171918 -- 0.171918
340nm -- 0.169 0.304927
345nm 0.349263 -- 0.349263
360nm -- 0.203 0.388343
365nm 0.40137 -- 0.40137
380nm -- 0.206 0.42161
385nm 0.428356 -- 0.428356
400nm -- 0.211 0.428664
405nm 0.428767 -- 0.428767
420nm -- 0.202 0.414795
425nm 0.410137 -- 0.410137
440nm -- 0.188 0.389281
445nm 0.382329 -- 0.382329
460nm -- 0.167 0.351507
465nm 0.341233 -- 0.341233
480nm -- 0.14 0.297979
485nm 0.283562 -- 0.283562
500nm -- 0.116 0.251507
505nm 0.240822 -- 0.240822
520nm -- 0.0806 0.202397
525nm 0.189589 -- 0.189589
540nm -- 0.0432 0.127432
545nm 0.106712 -- 0.106712
560nm -- 0.0265 0.075274
565nm 0.0647945 -- 0.0647945
580nm -- 0.0146 0.0457877
585nm 0.0394521 -- 0.0394521
600nm -- 0.00756 0.0251712
605nm 0.020411 -- 0.020411
620nm -- 0.00508 0.0114726
625nm 0.00849315 -- 0.00849315
640nm -- 0.00158 0.00232877
645nm 0.000273973 -- 0.001185
660nm -- 0 0.000171233
665nm 0.000136986 -- 0.000136986
spradlin commented 2 years ago

Recheck control

Verify that PR (cite PR) has the intended behavior for the control case. Running my WCSim-HKOD.mac with a clean checkout of the PR branch (no modifications to BoxandLine20inchHQE wavelengths) produces the combined QE function summarized below.

WL ID QE OD QE Combined
280nm 0 0 0
300nm 0.00109589 0.0139 0.0139
320nm 0.171918 0.0854 0.171918
340nm 0.349263 0.169 0.349263
360nm 0.40137 0.203 0.40137
380nm 0.428356 0.206 0.428356
400nm 0.428767 0.211 0.428767
420nm 0.410137 0.202 0.410137
440nm 0.382329 0.188 0.382329
460nm 0.341233 0.167 0.341233
480nm 0.283562 0.14 0.283562
500nm 0.240822 0.116 0.240822
520nm 0.189589 0.0806 0.189589
540nm 0.106712 0.0432 0.106712
560nm 0.0647945 0.0265 0.0647945
580nm 0.0394521 0.0146 0.0394521
600nm 0.020411 0.00756 0.020411
620nm 0.00849315 0.00508 0.00849315
640nm 0.000273973 0.00158 0.00158
660nm 0.000136986 0 0.000136986
spradlin commented 2 years ago

Created Issue WCSim#323 on the main repository with this information. No longer need this placeholder.