Trapping correction module and app

[sophieehaight]

Contains a new module for applying depth-based charge trapping correction to individual hit energies along with modules with GUI options and GUI expos. The branch also includes an app for characterizing trapping with Cs-137 data. The trapping correction module requires a csv parameter file as input:
detector_0_trapping_parameters.csv
This file can also be found in resource/dee

[fhagemann]

This is what I'm getting if I'm plotting the content of the file (parsing the base CCE curve, and applying the calibration parameters from the first line):

[fhagemann]

It looks like the units of depth_val (mm) and m_Depths (cm) do not match..

[fhagemann]

I looked at applying the inverse of the charge trapping correction to the DEE, and with the latest file, I would say that it gives reasonable results.

[fhagemann]

Can we also remove HasCalibratedTiming() == true from the trapping app and replace by something like HasFastTiming() == true && HasNearestNeighbor() == false (in an effort to completely get rid of HasCalibratedTiming)?

[fhagemann]

Here is my feedback from just looking at the code.
I will also try to compile and run both the module and the app to provide more feedback on that as well! :)

The example file in resource/dee and the results after applying the inverse charge trapping in the DEE looks reasonable:

Welcome to Codecov 🎉

Once you merge this PR into your default branch, you're all set! Codecov will compare coverage reports and display results in all future pull requests.

ℹ️ You can also turn on project coverage checks and project coverage reporting on Pull Request comment

Thanks for integrating Codecov - We've got you covered ☂️

Final code-style review, will give it a run to see that it compiles.
But then this should be good to go (tackling reading trapping parameters for multiple detectors in a follow-up PR)

Let's be consistent with XML tag naming here (see #193)

Suggested change

	m_XmlTag = "TrappingCorrection";
	m_XmlTag = "XmlTagTrappingCorrection";

I'm actually not seeing this convention in the depth or energy calibration modules? It looks like the XML tag is just the name of the module

This is tackled in #197 once it is merged

Once I change the name of the xml tag, nuclearizer can no longer recognize the module. Is there somewhere else in assembly I need to update the xml tag? I was looking but I couldn't find anything

This will affect reading in a nuclearizer config file (*.cfg). In those, you might need to update all occurences of TrappingCorrection with XmlTrappingCorrection.

So either update the .nuclearizer.cfg in your home directory, or the respective .cfg that you load when running nuclearizer.

Put this in c_Info?

I use this output to evaluate how well the trapping correction is working. If I change the verbosity to info, the comand line becomes overloaded (I believe w info from the depth cal) and nuclearizer crashes. I agree this info doesn't need to be repeated so I can change the verbosity for the following comments but I do need it read out somewhere.

I think we probably need to clean up the outputs for the info verbosity right now because right now it's outputting so much info that I can't receive the info I'm trying to get from the trapping correction.

@zoglauer maybe we can introduce the chatty verbosity from megalib also to nuclearizer, and move some of the info messages to chatty?

yeah that could be a good idea!

This would be good from my side. Just gathering here the open TODOs (for follow-up PRs):

• Allow for a charge trapping parameters file for multiple detectors, at not just one
• Think about flushing output to c_Info without overloading the output
• Maybe find an alternative to Interpolate to avoid code duplication (maybe from an existing package or another module --> the depth calibration also used some sort of interpolation algorithm to interpolate the CTD-to-depth spline).

Some code-style and some more major comments.

One general comment:
you are using a lot of snake-case when assigning variable names, we might want to use camel-case consistently here (e.g. SimCCESortedHoles instead of sim_cce_sorted_h)

This seems to be tailored to Cs137. Is this what we want? What if we pass data from other isotopes through MModuleTrappingCorrection?

I set it up in commits 8f30f9b, b95dc1e, and f191838 (took a couple tries to get my axes to update correctly here). I let the user define the energy of the photopeak they are tryting to fit on the nuclearizer expo. The trapping correction module then takes that input and uses it to define the Mu and min and max energies for the photopeak fits measuring resolutions pre and post trapping correction. I kept the default energy for Cs-137 as a benchmark.

Your last pushes reverted back to the fixed window around 662keV.
Is this intended?

my mistake! I was testing the fit with and without the dynamic fitting window. This is fixed in 6321870

I've lost track. Is this ready for a final review?

Do you also plan to add the individual detector files to resource/calibration/COSISMEX/unit?

yes, I will do this today

If the last changes conclude this PR, we might want to clean up the commit history and reduce the ~100 commits to a few major commits (adding the module, adding the GUI options, adding example files, ...), as this PR changes only a handful of files and adds a handful to the resource folder.

Happy to help here!

Comments on the Cs137 trapping app only:

Is there ever a scenario where this range might no be enough?
For payload calibrations, we will have slightly higher temperatures than in unit-level calibrations, resulting in lower charge carrier mobilities, that might push the CTD out of the [-250,250] range (maybe?)

Perhaps on the AC side we would see CTD >250. On the DC side most detectors can only get the required number of counts in a Cs-137 photopeak at CTD = -200 so it is not an issue. Right now the keeping CTD < +250 allows me to place an upper limit on the CTD value when map depth values to CTD values for the CCE parameterization.

Getting data at CTD >250 is something I want to take a look at in the future but right now we get plenty of counts on that side of the detector. I'll add a TO DO comment for this but for the purposes of this pull request I think we can leave these limits as is.

What is the difference between TACCalFile and TACCalibrationFile?
Edit: looks like m_TACCalibrationFile is used for what used to be the TAC cut file, so this can go:

Suggested change

	MString m_TACCalFile;
	MString m_TACCalibrationFile;
	MString m_TACCalFile;

I don't think that we need TAC cut files anymore with the latest changes in #194

Suggested change

	if ((Option == "-i") || (Option == "-o") || (Option == "--emin") || (Option == "--emax") || (Option == "--tcal") || (Option == "--tcut") || (Option == "-m") || (Option == "--ctdmin") || (Option == "--ctdmax")) {
	if ((Option == "-i") || (Option == "-o") || (Option == "--emin") || (Option == "--emax") || (Option == "--tcal") || (Option == "-m") || (Option == "--ctdmin") || (Option == "--ctdmax")) {

Suggested change

	if (Option == "--tcut") {
	m_TACCalibrationFile = argv[++i];
	cout<<"Accepting file name: "<<m_TACCalibrationFile<<endl;
	}

Do you require any sort of energy (slow) thresholds to be set here?

Yes, I set a slow threshold at 30 keV. I had originally had this set locally in the ecal module but I will add it to the app.

This seems like it's an identical copy of MModuleDepthCalibration::GetDominantStrip, which is a public function, so you can probably use that here instead:

1. Load #include "MModuleDepthCalibration.h"
2. In Analyze, define (with all other modules) MModuleDepthCalibration* Depth = new MModuleDepthCalibration, maybe with a comment that you're not running the module, just using its functions here.
3. Call Depth->GetDominantStrip where you just have GetDominantStrip.
4. Remove this function here (plus its declaration in the class).

You could be more explicit here:

Suggested change

	MGlobal::Initialize("Standalone", "a standalone example program");
	MGlobal::Initialize("TrappingCorrectionCs137", "A standalone program to determine CTD-dependent shifts of the Cs137 peak");

or something like that.

And some general comment on the file format (header lines starting with # or not)

The example files here in resource/dee and in resource/calibration are a bit inconsistent when it comes to the line defining the headers z_depth (in terms of if they should start with a # or not).

I would prefer the line to always start with a #,

I guess this file won't work because it doesn't start with ### 0 ?

Here and elsewhere in resource/calibration:

Suggested change

	z_depth_cm,e_CCE_HV,h_CCE_HV,e_CCE_LV,h_CCE_LV
	# z_depth_cm,e_CCE_HV,h_CCE_HV,e_CCE_LV,h_CCE_LV

This could go if we force the header lines to start with #.
Or is there any reason you chose not to?

The missing MGUIOptionsTrappingCorrection files were added

This is a ton of good work, @sophieehaight !

I've done a review of just the Nuclearizer App for now, I'll try to get through the rest tomorrow.

It looks like there are a few things that you need to update to account for recent merged PRs as identified in Felix's comments and the ones below.

I've run the app, but it took a bit of editing to make things work. Do you usually run with or without the "-n" option for nearest neighbors?

At first I was getting strange results, namely a double peaked spectra.

The double peak goes away when I don't use the -n option, but the resulting photopeak is much broader and shifted higher in energy, similar to what Julian has shown when we add NN strip hits without an appropriate energy threshold in at higher energy).

I then tried with a SlowThresholdCut of 15 keV and the double peak structure goes away. How were you doing this without the SlowThresholdCut?

//! The user interface for the trapping correction

yes, I am using the -n option to exclude nearest neighbors and I am applying a slowthreshold cut at 30 keV. I had originally just updated my Ecal module locally to apply the slow threshold cut but I am going to code the cut into the trapping app itself when I define the instance of the energy calibration. Will push that version soon!

Add you name :)

Missing the DetID map in the comment:

// [CTD][DetID][HV Energy, LV Energy] for full detector values

Add comment:

"// Enable sum of weights squared for TH1 by default"

Why *2 for the MaxEnergy? Maybe add a comment in the code to explain this logic?

HasCalibratedTiming is depreciated after PR #171. See #171 (comment)

I had to change this line to get any data in the histograms.

Suggested change

	if ((LVSH->HasCalibratedTiming() == true) && (HVSH->HasCalibratedTiming() == true) && (LVSH != nullptr) && (HVSH != nullptr)) {
	if ((LVSH->HasFastTiming() == true) && (HVSH->HasFastTiming() == true) && (LVSH != nullptr) && (HVSH != nullptr)) {

Change this comment to be descriptive of the following functionality. I.e. "// Define an output file to save the photopeak fitting results"

"// Close output fit parameter file"

I needed to add the following lines to include a SlowThreshold to remove the double peak spectra I was getting

Suggested change

	EnergyCalibrator->SetFileName(m_EcalFile);
	EnergyCalibrator->SetFileName(m_EcalFile);
	EnergyCalibrator->SetSlowThresholdCutMode(MSlowThresholdCutModes::e_Fixed);
	EnergyCalibrator->SetSlowThresholdCutFixedValue(15);

Are you ever expecting to combing data from different runs?

If this .txt file option is only meant to load multiple files from the same run, then the MModuleLoaderMeasurementsHDF already has the SetLoadContinuationFiles option that should handle this.

Having the txt file as an input for loading the HDF files allows me to control how many files I'm using to characterize the data. I have more data (30+ files of Cs data) for some of the detectors so I've been using only 20 data files via the txt file option to characterize individual detectors.