From 19c337b65388e4ebd03034d22a232188dd840f0d Mon Sep 17 00:00:00 2001 From: Sean Pike Date: Thu, 22 May 2025 11:19:27 -0400 Subject: [PATCH 1/5] CHG: Calculate X and Y position using energy-weighted strip positions --- include/MModuleDepthCalibration2024.h | 2 + src/MModuleDepthCalibration2024.cxx | 67 +++++++++++++++++++++++---- 2 files changed, 60 insertions(+), 9 deletions(-) diff --git a/include/MModuleDepthCalibration2024.h b/include/MModuleDepthCalibration2024.h index 9923efac..a928bcf8 100644 --- a/include/MModuleDepthCalibration2024.h +++ b/include/MModuleDepthCalibration2024.h @@ -99,6 +99,8 @@ class MModuleDepthCalibration2024 : public MModule vector norm_pdf(vector x, double mu, double sigma); //! Adds a Depth-to-CTD relation bool AddDepthCTD(vector depthvec, vector> ctdarr, int DetID, unordered_map>& DepthGrid, unordered_map>>& CTDMap); + // Calculate the Energy-weighted X and Y strip position + bool CalculateEnergyWeightedPosition(vector XStrips, vector YStrips, double& WeightedXStripID, double& WeightedYStripID); //! Determine the Grade (geometry of charge sharing) of the Hit int GetHitGrade(MHit* H); //! Load in the specified coefficients file diff --git a/src/MModuleDepthCalibration2024.cxx b/src/MModuleDepthCalibration2024.cxx index 64890753..e6770dea 100644 --- a/src/MModuleDepthCalibration2024.cxx +++ b/src/MModuleDepthCalibration2024.cxx @@ -265,13 +265,25 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) int YStripID = YSH->GetStripID(); int pixel_code = 10000*DetID + 100*XStripID + YStripID; - // TODO: Calculate X and Y positions more rigorously using charge sharing. + double WeightedXStripID = 0; + double WeightedYStripID = 0; + + bool WeightedPosResult = CalculateEnergyWeightedPosition(XStrips, YStrips, WeightedXStripID, WeightedYStripID); + // Somewhat confusing notation: XStrips run parallel to X-axis, so we calculate X position with YStrips. - double Xpos = m_YPitches[DetID]*((m_NYStrips[DetID]/2.0) - ((double)YStripID)); - double Ypos = m_XPitches[DetID]*((m_NXStrips[DetID]/2.0) - ((double)XStripID)); - // cout << "X position " << Xpos << endl; - // cout << "Y position " << Ypos << endl; - double Zpos = 0.0; + + double XPos = 0.0; + double YPos = 0.0; + double ZPos = 0.0; + + // Try to calculate the energy-weighted strip positions. If that doesn't work, place the Hit in the middle of the dominant strips. + if ( WeightedPosResult == true ) { + XPos = m_YPitches[DetID]*((m_NYStrips[DetID]/2.0) - (WeightedYStripID)); + YPos = m_XPitches[DetID]*((m_NXStrips[DetID]/2.0) - (WeightedXStripID)); + } else { + XPos = m_YPitches[DetID]*((m_NYStrips[DetID]/2.0) - ((double)YStripID)); + YPos = m_XPitches[DetID]*((m_NXStrips[DetID]/2.0) - ((double)XStripID)); + } double Xsigma = m_YPitches[DetID]/sqrt(12.0); double Ysigma = m_XPitches[DetID]/sqrt(12.0); @@ -373,17 +385,17 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) } Zsigma = sqrt(depth_var/prob_sum); - Zpos = mean_depth - (m_Thicknesses[DetID]/2.0); + ZPos = mean_depth - (m_Thicknesses[DetID]/2.0); // Add the depth to the GUI histogram. if (Event->IsStripPairingIncomplete()==false) { - m_ExpoDepthCalibration->AddDepth(DetID, Zpos); + m_ExpoDepthCalibration->AddDepth(DetID, ZPos); } m_NoError+=1; } } - LocalPosition.SetXYZ(Xpos, Ypos, Zpos); + LocalPosition.SetXYZ(XPos, YPos, ZPos); LocalOrigin.SetXYZ(0.0,0.0,0.0); // cout << m_DetectorNames[DetID] << endl; GlobalPosition = m_Detectors[DetID]->GetSensitiveVolume(0)->GetPositionInWorldVolume(LocalPosition); @@ -583,6 +595,43 @@ bool MModuleDepthCalibration2024::LoadSplinesFile(MString FName) } +bool MModuleDepthCalibration2024::CalculateEnergyWeightedPosition(vector XStrips, vector YStrips, double& WeightedXStripID, double& WeightedYStripID) { + // Determine the weighted strip ID positions based on charge sharing. + WeightedXStripID = 0; + WeightedYStripID = 0; + double TotalXEnergy = 0; + double TotalYEnergy = 0; + + // If one side or another doesn't have strips, abort. + if ( (XStrips.size() == 0) || (YStrips.size()==0) ) { + return false; + } + + // Calculate the weighted values based on energy per strip. + for ( unsigned int i=0; iGetEnergy(); + int XSHID = XStrips[i]->GetStripID(); + TotalXEnergy += XSHEnergy; + WeightedXStripID += XSHEnergy*XSHID; + } + for ( unsigned int i=0; iGetEnergy(); + int YSHID = YStrips[i]->GetStripID(); + TotalYEnergy += YSHEnergy; + WeightedYStripID += YSHEnergy*YSHID; + } + + // If one side or another doesn't have energy, abort. + if ( (TotalXEnergy == 0) || (TotalYEnergy == 0) ) { + return false; + } + + WeightedXStripID /= TotalXEnergy; + WeightedYStripID /= TotalYEnergy; + + return true; +} + int MModuleDepthCalibration2024::GetHitGrade(MHit* H){ // Function for choosing which Depth-to-CTD relation to use for a given event. // At time of writing, intention is to choose a CTD based on sub-pixel region determined via charge sharing (Event "grade"). From 1899e354e3339800b8d9a55658513f3f78d1f263 Mon Sep 17 00:00:00 2001 From: Sean Pike Date: Thu, 18 Sep 2025 14:08:44 -0700 Subject: [PATCH 2/5] CHG: Update spline handling to ensure the full detector is sampled and the depth grid is evenly spaced --- include/MModuleDepthCalibration2024.h | 5 +- src/MModuleDepthCalibration2024.cxx | 211 ++++++++++++++++---------- 2 files changed, 137 insertions(+), 79 deletions(-) diff --git a/include/MModuleDepthCalibration2024.h b/include/MModuleDepthCalibration2024.h index 21f48291..17237b4f 100644 --- a/include/MModuleDepthCalibration2024.h +++ b/include/MModuleDepthCalibration2024.h @@ -95,10 +95,12 @@ class MModuleDepthCalibration2024 : public MModule vector GetDepth(int DetID); //! Retrieve the appropriate CTD values given the DetID and Grade vector GetCTD(int DetID, int Grade); + //! Retrieve the appropriate depth-to-CTD spline given the DetID and Grade + TSpline3* GetSpline(int DetID, int Grade); //! Normal distribution vector norm_pdf(vector x, double mu, double sigma); //! Adds a Depth-to-CTD relation - bool AddDepthCTD(vector depthvec, vector> ctdarr, int DetID, unordered_map>& DepthGrid, unordered_map>>& CTDMap); + bool AddDepthCTD(vector Depth, vector> CTDArr, int DetID, unordered_map>& DepthGrid, unordered_map>>& CTDMap, unordered_map>& SplineMap, unsigned int NPoints); //! Determine the Grade (geometry of charge sharing) of the Hit int GetHitGrade(MHit* H); //! Load in the specified coefficients file @@ -145,6 +147,7 @@ class MModuleDepthCalibration2024 : public MModule // The CTD Map maps each detector (int) to a 2D array of CTD values. unordered_map>> m_CTDMap; unordered_map> m_DepthGrid; + unordered_map> m_SplineMap; bool m_SplinesFileIsLoaded; bool m_CoeffsFileIsLoaded; diff --git a/src/MModuleDepthCalibration2024.cxx b/src/MModuleDepthCalibration2024.cxx index 135db6dc..be3dcc7a 100644 --- a/src/MModuleDepthCalibration2024.cxx +++ b/src/MModuleDepthCalibration2024.cxx @@ -246,8 +246,8 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) // Calculate the position. If error is thrown, record and no depth. // Take a Hit and separate its activated X- and Y-strips into separate vectors. - std::vector LVStrips; - std::vector HVStrips; + vector LVStrips; + vector HVStrips; for( unsigned int j = 0; j < H->GetNStripHits(); ++j){ MStripHit* SH = H->GetStripHit(j); if( SH->IsLowVoltageStrip() ) LVStrips.push_back(SH); else HVStrips.push_back(SH); @@ -278,6 +278,9 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) vector* Coeffs = GetPixelCoeffs(PixelCode); + vector CTDVec = GetCTD(DetID, Grade); + vector DepthVec = GetDepth(DetID); + // TODO: For Card Cage, may need to add noise double LVTiming = LVSH->GetTiming(); double HVTiming = HVSH->GetTiming(); @@ -288,34 +291,29 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) H->SetNoDepth(); Event->SetDepthCalibrationIncomplete(); ++m_Error1; - } - // If there isn't timing information, set no depth. - // Alex's old comments suggest assigning the event to the middle of the detector and the position resolution to be large. - else if( (LVTiming < 1.0E-6) || (HVTiming < 1.0E-6) ){ + } else if (CTDVec.size() == 0) { + cout << "Empty CTD vector" << endl; + H->SetNoDepth(); + Event->SetDepthCalibrationIncomplete(); + } else if (DepthVec.size() == 0) { + cout << "Empty Depth vector" << endl; + H->SetNoDepth(); + Event->SetDepthCalibrationIncomplete(); + } else if ((LVTiming < 1.0E-6) || (HVTiming < 1.0E-6)) { ++m_Error3; H->SetNoDepth(); Event->SetDepthCalibrationIncomplete(); - } - - // If there are coefficients and timing information is loaded, try calculating the CTD and depth - else { - - vector ctdvec = GetCTD(DetID, Grade); - vector depthvec = GetDepth(DetID); - - if ( ctdvec.size() == 0){ - cout << "Empty CTD vector" << endl; - H->SetNoDepth(); - Event->SetDepthCalibrationIncomplete(); - } + } else { + + // If there are coefficients and timing information is loaded, try calculating the CTD and depth double CTD = (HVTiming - LVTiming); // Confirmed that this matches SP's python code. CTD_s = (CTD - Coeffs->at(1))/(Coeffs->at(0)); //apply inverse stretch and offset - double Xmin = * std::min_element(ctdvec.begin(), ctdvec.end()); - double Xmax = * std::max_element(ctdvec.begin(), ctdvec.end()); + double Xmin = * std::min_element(CTDVec.begin(), CTDVec.end()); + double Xmax = * std::max_element(CTDVec.begin(), CTDVec.end()); double noise = GetTimingNoiseFWHM(PixelCode, H->GetEnergy()); @@ -328,9 +326,9 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) // If the CTD is in range, calculate the depth else { - // Calculate the probability given timing noise of CTD_s corresponding to the values of depth in depthvec + // Calculate the probability given timing noise of CTD_s corresponding to the values of depth in DepthVec // Utlize symmetry of the normal distribution. - vector prob_dist = norm_pdf(ctdvec, CTD_s, noise/2.355); + vector prob_dist = norm_pdf(CTDVec, CTD_s, noise/2.355); // Weight the depth by probability double prob_sum = 0.0; @@ -339,20 +337,20 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) } //double prob_sum = std::accumulate(prob_dist.begin(), prob_dist.end(), 0); double weighted_depth = 0.0; - for( unsigned int k = 0; k < depthvec.size(); ++k ){ - weighted_depth += prob_dist[k]*depthvec[k]; + for( unsigned int k = 0; k < DepthVec.size(); ++k ){ + weighted_depth += prob_dist[k]*DepthVec[k]; } // Calculate the expectation value of the depth double mean_depth = weighted_depth/prob_sum; // Calculate the standard deviation of the depth double depth_var = 0.0; - for( unsigned int k=0; kIsStripPairingIncomplete()==false) { @@ -498,24 +496,17 @@ vector MModuleDepthCalibration2024::norm_pdf(vector x, double mu bool MModuleDepthCalibration2024::LoadSplinesFile(MString FName) { - //when invert flag is set to true, the splines returned are CTD->Depth - // Previously saved cathode and anode timing in addition to CTD. This may be redundant, commenting out for now. // Input spline files should have the following format: // ### DetID, HV, Temperature, Photopeak Energy (TODO: More? Fewer?) // depth, ctd0, ctd1, ctd2.... (Basically, allow for CTDs for different subpixel regions) // '' '' '' - MFile F; + MFile F; if( F.Open(FName) == false ){ cout << "ERROR in MModuleDepthCalibration2024::LoadSplinesFile: failed to open splines file." << endl; return false; } - // vector depthvec, ctdvec, anovec, catvec; - vector depthvec; - vector> ctdarr; - for( unsigned int i=0; i < 5; ++i ){ - vector temp_vec; - ctdarr.push_back(temp_vec); - } + vector DepthVec; + vector> CTDArr; bool Result = true; MString line; int DetID = 0; @@ -524,33 +515,28 @@ bool MModuleDepthCalibration2024::LoadSplinesFile(MString FName) if (line.BeginsWith("#")) { // If we've reached a new ctd spline then record the previous one in the m_SplineMaps and start a new one. vector tokens = line.Tokenize(" "); - if (depthvec.size() > 0) { - Result &= AddDepthCTD(depthvec, ctdarr, DetID, m_DepthGrid, m_CTDMap); - } - depthvec.clear(); ctdarr.clear(); - for (unsigned int i=0; i < 5; ++i) { - vector temp_vec; - ctdarr.push_back(temp_vec); + if (DepthVec.size() > 0) { + Result &= AddDepthCTD(DepthVec, CTDArr, DetID, m_DepthGrid, m_CTDMap, m_SplineMap, 1000); } + DepthVec.clear(); CTDArr.clear(); DetID = tokens[1].ToInt(); } else { vector tokens = line.Tokenize(","); - depthvec.push_back(tokens[0].ToDouble()); - + DepthVec.push_back(tokens[0].ToDouble()); // Multiple CTDs allowed. for (unsigned int i = 0; i < (tokens.size() - 1); ++i) { - ctdarr[i].push_back(tokens[1+i].ToDouble()); - } - // Fill in the higher grades with the GRADE=0 CTD if there are none listed in the file. - for (unsigned int i=tokens.size()-1; i<5; ++i) { - ctdarr[i].push_back(tokens[1].ToDouble()); + while (i>=CTDArr.size()) { + vector TempVec; + CTDArr.push_back(TempVec); + } + CTDArr[i].push_back(tokens[1+i].ToDouble()); } } } } //make last spline - if (depthvec.size() > 0) { - Result &= AddDepthCTD(depthvec, ctdarr, DetID, m_DepthGrid, m_CTDMap); + if (DepthVec.size() > 0) { + Result &= AddDepthCTD(DepthVec, CTDArr, DetID, m_DepthGrid, m_CTDMap, m_SplineMap, 1000); } return Result; @@ -672,44 +658,92 @@ int MModuleDepthCalibration2024::GetHitGrade(MHit* H){ return return_value; } -bool MModuleDepthCalibration2024::AddDepthCTD(vector depthvec, vector> ctdarr, int DetID, unordered_map>& DepthGrid, unordered_map>>& CTDMap){ +bool MModuleDepthCalibration2024::AddDepthCTD(vector Depth, vector> CTDArr, int DetID, unordered_map>& DepthGrid, unordered_map>>& CTDMap, unordered_map>& SplineMap, unsigned int NPoints) +{ // Saves a CTD array, basically allowing for multiple CTDs as a function of depth - // depthvec: list of simulated depth values - // ctdarr: vector of vectors of simulated CTD values. Each vector of CTDs must be the same length as depthvec + // Depth: list of simulated depth values + // CTDArr: vector of vectors of simulated CTD values. Each vector of CTDs must be the same length as Depth // DetID: An integer which specifies which detector. // CTDMap: unordered map into which the array of CTDs should be placed // TODO: Possible energy dependence of CTD? - // TODO: Depth values need to be evenly spaced. Check this when reading the files in. // Check to make sure things look right. // First check that the CTDs all have the right length. - for (unsigned int i = 0; i < ctdarr.size(); ++i) { - if ((ctdarr[i].size() != depthvec.size()) && (ctdarr[i].size() > 0)) { + + for (unsigned int i = 0; i < CTDArr.size(); ++i) { + if ((CTDArr[i].size() != Depth.size()) && (CTDArr[i].size() > 0)) { cout<<"ERROR in MModuleDepthCalibration2024::AddDepthCTD: The number of values in the CTD list is not equal to the number of depth values."< 0.01) { + double MaxDepth = * std::max_element(Depth.begin(), Depth.end()); + double MinDepth = * std::min_element(Depth.begin(), Depth.end()); + if (fabs((MaxDepth-MinDepth) - m_Thicknesses[DetID]) > 0.01) { cout<<"ERROR in MModuleDepthCalibration2024::AddDepthCTD: The thickness of detector "<0) { + cout<<"MModuleDepthCalibration2024::AddSplines: Splines already added for DetID "< TempVec; + SplineMap[DetID] = TempVec; } + + double dx_low = Depth[1] - Depth[0]; + double newx_low = Depth[0] - (dx_low/2.0); + Depth.insert(Depth.begin(), newx_low); + + size_t N = Depth.size(); + double dx_high = Depth[N-1] - Depth[N-2]; + double newx_high = Depth[N-1] + (dx_high/2.0); + Depth.push_back(newx_high); + + if (NPoints < N+1) { + NPoints = N+1; + } + + vector NewDepth; + for (unsigned int i=0; i> NewCTDArr; + for (unsigned int i=0; i CTD = CTDArr[i]; + vector NewCTD; - //Now make sure the values for the depth start with 0.0. - if (mindepth != 0.0) { - cout<<"MModuleDepthCalibration2024::AddDepthCTD: The minimum depth is not zero. Editing the depth vector."<Eval(NewDepth[d])); + } + NewCTDArr.push_back(NewCTD); } - CTDMap[DetID] = ctdarr; - DepthGrid[DetID] = depthvec; + CTDMap[DetID] = NewCTDArr; return true; } @@ -729,7 +763,7 @@ vector MModuleDepthCalibration2024::GetCTD(int DetID, int Grade) return (m_CTDMap[DetID][Grade]); } else { - cout << "MModuleDepthCalibration2024::GetCTD: No CTD map is loaded for Grade " << Grade << ". Returning Grade 0 CTD." << endl; + // cout << "MModuleDepthCalibration2024::GetCTD: No CTD map is loaded for Grade " << Grade << ". Returning Grade 0 CTD." << endl; return (m_CTDMap[DetID][0]); } } else { @@ -757,6 +791,32 @@ vector MModuleDepthCalibration2024::GetDepth(int DetID) } } + +TSpline3* MModuleDepthCalibration2024::GetSpline(int DetID, int Grade) +{ + // Retrieves the appropriate depth->CTD spline given the Detector ID and Event Grade passed + + if( !m_SplinesFileIsLoaded ){ + cout << "MModuleDepthCalibration2024::GetSpline: cannot return Depth to CTD spline because the file was not loaded." << endl; + return nullptr; + } + + // If there is a spline for the given detector, return it. + // If the Grade is larger than the number of splines stored, then just return Grade 0 spline. + if( m_SplineMap.count(DetID) > 0 ){ + if ( ((int)m_SplineMap[DetID].size()) > Grade) { + return (m_SplineMap[DetID][Grade]); + } + else { + // cout << "MModuleDepthCalibration2024::GetSpline: No spline is loaded for Grade " << Grade << ". Returning Grade 0 spline." << endl; + return (m_SplineMap[DetID][0]); + } + } else { + cout << "MModuleDepthCalibration2024::GetSpline: No spline is loaded for Det " << DetID << "." << endl; + return nullptr; + } +} + void MModuleDepthCalibration2024::ShowOptionsGUI() { // Show the options GUI - or do nothing @@ -820,11 +880,6 @@ void MModuleDepthCalibration2024::Finalize() cout << "Number of hits with no strip hits on one or both sides: " << m_ErrorSH << endl; cout << "Number of hits with null strip hits: " << m_ErrorNullSH << endl; cout << "Number of hits 0 energy on a strip hit: " << m_ErrorNoE << endl; - /* - TFile* rootF = new TFile("EHist.root","recreate"); - rootF->WriteTObject( EHist ); - rootF->Close(); - */ } From 22eb08e2cc461441c7df826ddaf10cc85b8bce4f Mon Sep 17 00:00:00 2001 From: Sean Pike Date: Tue, 23 Sep 2025 15:12:06 -0700 Subject: [PATCH 3/5] Cleanup: added comments and cleaned up extraneous print statements --- src/MModuleDepthCalibration2024.cxx | 19 +++++++++++++------ 1 file changed, 13 insertions(+), 6 deletions(-) diff --git a/src/MModuleDepthCalibration2024.cxx b/src/MModuleDepthCalibration2024.cxx index be3dcc7a..6697a1b6 100644 --- a/src/MModuleDepthCalibration2024.cxx +++ b/src/MModuleDepthCalibration2024.cxx @@ -325,6 +325,9 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) } // If the CTD is in range, calculate the depth + // Rather than plugging CTD into a spline to get depth, use the depth-CTD relation to calculate a probability-weighted depth value. + // This way we can avoid problems like non-monotonicity or assigning depth to events "outside" the detector + // Note that this requires that we don't massively overestimate the timing noise else { // Calculate the probability given timing noise of CTD_s corresponding to the values of depth in DepthVec // Utlize symmetry of the normal distribution. @@ -335,7 +338,6 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) for( unsigned int k=0; k < prob_dist.size(); ++k ){ prob_sum += prob_dist[k]; } - //double prob_sum = std::accumulate(prob_dist.begin(), prob_dist.end(), 0); double weighted_depth = 0.0; for( unsigned int k = 0; k < DepthVec.size(); ++k ){ weighted_depth += prob_dist[k]*DepthVec[k]; @@ -488,7 +490,6 @@ vector MModuleDepthCalibration2024::norm_pdf(vector x, double mu vector result; for( unsigned int i=0; i Depth, vector Depth, vector 0.01) { @@ -687,14 +690,16 @@ bool MModuleDepthCalibration2024::AddDepthCTD(vector Depth, vector0) { - cout<<"MModuleDepthCalibration2024::AddSplines: Splines already added for DetID "< TempVec; SplineMap[DetID] = TempVec; } + // Add a new point to the depth vector on the bottom and top of the detector to make sure we're sampling the full depth double dx_low = Depth[1] - Depth[0]; double newx_low = Depth[0] - (dx_low/2.0); Depth.insert(Depth.begin(), newx_low); @@ -704,16 +709,20 @@ bool MModuleDepthCalibration2024::AddDepthCTD(vector Depth, vector NewDepth; for (unsigned int i=0; i> NewCTDArr; for (unsigned int i=0; i CTD = CTDArr[i]; @@ -763,7 +772,6 @@ vector MModuleDepthCalibration2024::GetCTD(int DetID, int Grade) return (m_CTDMap[DetID][Grade]); } else { - // cout << "MModuleDepthCalibration2024::GetCTD: No CTD map is loaded for Grade " << Grade << ". Returning Grade 0 CTD." << endl; return (m_CTDMap[DetID][0]); } } else { @@ -808,7 +816,6 @@ TSpline3* MModuleDepthCalibration2024::GetSpline(int DetID, int Grade) return (m_SplineMap[DetID][Grade]); } else { - // cout << "MModuleDepthCalibration2024::GetSpline: No spline is loaded for Grade " << Grade << ". Returning Grade 0 spline." << endl; return (m_SplineMap[DetID][0]); } } else { From 0096dd6ce1514279eade274c2b9b3feb24534995 Mon Sep 17 00:00:00 2001 From: Sean Pike Date: Tue, 23 Sep 2025 15:35:09 -0700 Subject: [PATCH 4/5] BUG: Clear variables on Finalize --- src/MModuleDepthCalibration2024.cxx | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/src/MModuleDepthCalibration2024.cxx b/src/MModuleDepthCalibration2024.cxx index 6697a1b6..46de19ae 100644 --- a/src/MModuleDepthCalibration2024.cxx +++ b/src/MModuleDepthCalibration2024.cxx @@ -888,6 +888,19 @@ void MModuleDepthCalibration2024::Finalize() cout << "Number of hits with null strip hits: " << m_ErrorNullSH << endl; cout << "Number of hits 0 energy on a strip hit: " << m_ErrorNoE << endl; + // Clean up maps and vectors + m_Coeffs.clear(); + m_Thicknesses.clear(); + m_NXStrips.clear(); + m_NYStrips.clear(); + m_XPitches.clear(); + m_YPitches.clear(); + m_Detectors.clear(); + m_CTDMap.clear(); + m_DepthGrid.clear(); + m_SplineMap.clear(); + m_DetectorIDs.clear(); + } From 6a81a656fdbf2649c6e69492963da55245dc0d9c Mon Sep 17 00:00:00 2001 From: Sean Pike Date: Wed, 24 Sep 2025 14:31:34 -0700 Subject: [PATCH 5/5] CHG: Provide optional checkbox for energy-weighted X and Y calculation --- include/MGUIOptionsDepthCalibration2024.h | 2 ++ include/MModuleDepthCalibration2024.h | 5 +++++ src/MGUIOptionsDepthCalibration2024.cxx | 8 +++++++- src/MModuleDepthCalibration2024.cxx | 24 ++++++++++++++++------- 4 files changed, 31 insertions(+), 8 deletions(-) diff --git a/include/MGUIOptionsDepthCalibration2024.h b/include/MGUIOptionsDepthCalibration2024.h index 15c72ed6..363b7022 100644 --- a/include/MGUIOptionsDepthCalibration2024.h +++ b/include/MGUIOptionsDepthCalibration2024.h @@ -77,6 +77,8 @@ class MGUIOptionsDepthCalibration2024 : public MGUIOptions //! Check button if working with the Card Cage at UCSD TGCheckButton* m_UCSDOverride; + //! Check button if weighting X and Y by energy + TGCheckButton* m_WeightedXY; #ifdef ___CLING___ public: diff --git a/include/MModuleDepthCalibration2024.h b/include/MModuleDepthCalibration2024.h index 1a1551ff..11582dd1 100644 --- a/include/MModuleDepthCalibration2024.h +++ b/include/MModuleDepthCalibration2024.h @@ -77,6 +77,10 @@ class MModuleDepthCalibration2024 : public MModule //! Get whether the data came from the card cage at UCSD bool GetUCSDOverride() const {return m_UCSDOverride;} + //! Set whether X and Y positions should be calculated by weighting strips by energy + void SetWeightedXY( bool WeightedXY ) {m_WeightedXY = WeightedXY;} + //! Get whether X and Y positions should be calculated by weighting strips by energy + bool GetWeightedXY() const {return m_WeightedXY;} //! Read the XML configuration bool ReadXmlConfiguration(MXmlNode* Node); @@ -155,6 +159,7 @@ class MModuleDepthCalibration2024 : public MModule // boolean for use with the card cage at UCSD since it tags all events as detector 11 bool m_UCSDOverride; + bool m_WeightedXY; diff --git a/src/MGUIOptionsDepthCalibration2024.cxx b/src/MGUIOptionsDepthCalibration2024.cxx index cd6c07ec..e7c97f26 100644 --- a/src/MGUIOptionsDepthCalibration2024.cxx +++ b/src/MGUIOptionsDepthCalibration2024.cxx @@ -80,8 +80,13 @@ void MGUIOptionsDepthCalibration2024::Create() m_UCSDOverride = new TGCheckButton(m_OptionsFrame, "Check this box if you're using the card cage at UCSD", 1); m_UCSDOverride->SetOn(dynamic_cast(m_Module)->GetUCSDOverride()); + TGLayoutHints* Label3Layout = new TGLayoutHints(kLHintsTop | kLHintsCenterX | kLHintsExpandX, 10, 10, 10, 10); + m_OptionsFrame->AddFrame(m_UCSDOverride, Label3Layout); + + m_WeightedXY = new TGCheckButton(m_OptionsFrame, "Check this box to weight the X and Y positions by energy deposited.", 1); + m_WeightedXY->SetOn(dynamic_cast(m_Module)->GetWeightedXY()); TGLayoutHints* Label4Layout = new TGLayoutHints(kLHintsTop | kLHintsCenterX | kLHintsExpandX, 10, 10, 10, 10); - m_OptionsFrame->AddFrame(m_UCSDOverride, Label4Layout); + m_OptionsFrame->AddFrame(m_WeightedXY, Label4Layout); PostCreate(); } @@ -128,6 +133,7 @@ bool MGUIOptionsDepthCalibration2024::OnApply() dynamic_cast(m_Module)->SetCoeffsFileName(m_CoeffsFileSelector->GetFileName()); dynamic_cast(m_Module)->SetSplinesFileName(m_SplinesFileSelector->GetFileName()); dynamic_cast(m_Module)->SetUCSDOverride(m_UCSDOverride->IsOn()); + dynamic_cast(m_Module)->SetWeightedXY(m_WeightedXY->IsOn()); return true; } diff --git a/src/MModuleDepthCalibration2024.cxx b/src/MModuleDepthCalibration2024.cxx index 3d7ee353..022ab1c8 100644 --- a/src/MModuleDepthCalibration2024.cxx +++ b/src/MModuleDepthCalibration2024.cxx @@ -266,13 +266,17 @@ bool MModuleDepthCalibration2024::AnalyzeEvent(MReadOutAssembly* Event) int HVStripID = HVSH->GetStripID(); int PixelCode = 10000*DetID + 100*LVStripID + HVStripID; - double WeightedLVStripID = 0; - double WeightedHVStripID = 0; - bool WeightedPosResult = CalculateEnergyWeightedPosition(LVStrips, HVStrips, WeightedLVStripID, WeightedHVStripID); - // Try to calculate the energy-weighted strip positions. If that doesn't work, place the Hit in the middle of the dominant strips. - if ( WeightedPosResult == true ) { - LVStripID = WeightedLVStripID; - HVStripID = WeightedHVStripID; + if (m_WeightedXY==true) { + + // Try to calculate the energy-weighted strip positions. If that doesn't work, place the Hit in the middle of the dominant strips. + double WeightedLVStripID = 0; + double WeightedHVStripID = 0; + bool WeightedPosResult = CalculateEnergyWeightedPosition(LVStrips, HVStrips, WeightedLVStripID, WeightedHVStripID); + if (WeightedPosResult == true) { + LVStripID = WeightedLVStripID; + HVStripID = WeightedHVStripID; + } + } // TODO: Calculate X and Y positions more rigorously using charge sharing. @@ -898,6 +902,11 @@ bool MModuleDepthCalibration2024::ReadXmlConfiguration(MXmlNode* Node) m_UCSDOverride = (bool) UCSDOverrideNode->GetValueAsInt(); } + MXmlNode* WeightedXYNode = Node->GetNode("WeightedXY"); + if( WeightedXYNode != NULL ){ + m_WeightedXY = (bool) WeightedXYNode->GetValueAsInt(); + } + return true; } @@ -912,6 +921,7 @@ MXmlNode* MModuleDepthCalibration2024::CreateXmlConfiguration() new MXmlNode(Node, "CoeffsFileName", m_CoeffsFile); new MXmlNode(Node, "SplinesFileName", m_SplinesFile); new MXmlNode(Node, "UCSDOverride",(unsigned int) m_UCSDOverride); + new MXmlNode(Node, "WeightedXY",(unsigned int) m_WeightedXY); return Node; }