modified: GainMatchSX3.C

modified:   GainMatchSX3Front.C
This commit is contained in:
Vignesh Sitaraman 2025-08-20 15:25:40 -07:00
parent b44ffd7fdf
commit 06fbc1afd9
2 changed files with 133 additions and 36 deletions

View File

@ -40,8 +40,8 @@ bool frontGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}};
// ==== Configuration Flags ==== // ==== Configuration Flags ====
const bool interactiveMode = false; // If true: show canvas + wait for user const bool interactiveMode = false; // If true: show canvas + wait for user
const bool verboseFit = true; // If true: print fit summary and chi² const bool verboseFit = true; // If true: print fit summary and chi²
const bool drawCanvases = true; // If false: canvases won't be drawn at all const bool drawCanvases = true; // If false: canvases won't be drawn at all
void GainMatchSX3::Begin(TTree * /*tree*/) void GainMatchSX3::Begin(TTree * /*tree*/)
{ {
@ -207,14 +207,14 @@ Bool_t GainMatchSX3::Process(Long64_t entry)
{ {
auto key = std::make_tuple(sx3.id[i], sx3ChBk, sx3ChUp, sx3ChDn); auto key = std::make_tuple(sx3.id[i], sx3ChBk, sx3ChUp, sx3ChDn);
comboCounts[key]++; comboCounts[key]++;
// If we have a valid front and back channel, fill the histograms
hSX3->Fill(sx3ChDn+4, sx3ChBk);
hSX3->Fill(sx3ChUp, sx3ChBk);
// Fill the histogram for the front vs back
hSX3FvsB->Fill(sx3EUp + sx3EDn, sx3EBk);
} }
// If we have a valid front and back channel, fill the histograms
hSX3->Fill(sx3ChDn, sx3ChBk);
hSX3->Fill(sx3ChUp, sx3ChBk);
// Fill the histogram for the front vs back
hSX3FvsB->Fill(sx3EUp + sx3EDn, sx3EBk);
for (int i = 0; i < sx3.multi; i++) for (int i = 0; i < sx3.multi; i++)
{ {
@ -322,7 +322,7 @@ void GainMatchSX3::Terminate()
// if (TMath::Abs(f.GetParameter(0) - 1) > 3.0) // if (TMath::Abs(f.GetParameter(0) - 1) > 3.0)
// continue; // continue;
const double fixedError = 20.0; // in ADC channels const double fixedError = 10.0; // in ADC channels
std::vector<double> xVals, yVals, exVals, eyVals; std::vector<double> xVals, yVals, exVals, eyVals;
@ -334,7 +334,7 @@ void GainMatchSX3::Terminate()
xVals.push_back(x); xVals.push_back(x);
yVals.push_back(y); yVals.push_back(y);
// exVals.push_back(fixedError); // error in front energy exVals.push_back(fixedError); // error in front energy
eyVals.push_back(fixedError); // error in back energy eyVals.push_back(fixedError); // error in back energy
} }
@ -382,22 +382,36 @@ void GainMatchSX3::Terminate()
gainArray[id][bk][u][d] = f.GetParameter(0); gainArray[id][bk][u][d] = f.GetParameter(0);
gainValid[id][bk][u][d] = true; gainValid[id][bk][u][d] = true;
// }
// Check if the gain is valid for this detector, back, up, and down // // Output results
if (gainValid[id][bk][u][d]) // for (int id = 0; id < MAX_DET; ++id)
// {
// for (int bk = 0; bk < MAX_BK; ++bk)
// {
// for (int u = 0; u < MAX_UP; ++u)
// {
// for (int d = 0; d < MAX_DOWN; ++d)
// {
// // Check if the gain is valid for this detector, back, up, and down
// if (gainValid[id][bk][u][d])
// {
if (TMath::Abs(gainArray[id][u][d][bk] - 1) < 0.3)
{ {
if (TMath::Abs(gainArray[id][u][d][bk] - 1) < 0.3) printf("Gain match Det%d Up%dDn%d Backs%d → %.4f \n", id, u, d, bk, gainArray[id][u][d][bk]);
{ outFile << id << " " << bk << " " << u << " " << d << " " << gainArray[id][u][d][bk] << std::endl;
printf("Gain match Det%d Up%dDn%d Back%d → %.4f \n", id, u, d, bk, gainArray[id][u][d][bk]); }
outFile << id << " " << bk << " " << u << " " << d << " " << gainArray[id][u][d][bk] << std::endl; else if (gainArray[id][u][d][bk] != 0)
} {
else if (gainArray[id][u][d][bk] != 0) std::cerr << "Warning: Gain value out of range for Det " << id << " Up " << u << " Dn " << d << " Back " << bk << ": "
{ << gainArray[id][u][d][bk] << std::endl;
std::cerr << "Warning: Gain value out of range for Det " << id << " Up " << u << " Dn " << d << " Back " << bk << ": "
<< gainArray[id][u][d][bk] << std::endl;
}
} }
} }
// }
// }
// }
// }
// }
// for (int bk = 0; bk < MAX_BK; ++bk) // for (int bk = 0; bk < MAX_BK; ++bk)
// { // {
@ -411,7 +425,7 @@ void GainMatchSX3::Terminate()
// === Create histograms === // === Create histograms ===
TH2F *hFVB = new TH2F("hFVB", "Corrected Up+Dn vs Corrected Back;Corrected Back E;Up+Dn E", TH2F *hFVB = new TH2F("hFVB", "Corrected Up+Dn vs Corrected Back;Corrected Back E;Up+Dn E",
400, 0, 16000, 400, 0, 16000); 600, 0, 16000, 600, 0, 16000);
TH2F *hAsym = new TH2F("hAsym", "Up vs Dn dvide corrected back;Up/Back E;Dn/Back E", TH2F *hAsym = new TH2F("hAsym", "Up vs Dn dvide corrected back;Up/Back E;Dn/Back E",
400, 0.0, 1.0, 400, 0.0, 1.0); 400, 0.0, 1.0, 400, 0.0, 1.0);

View File

@ -41,11 +41,13 @@ double backGain[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{0}}}};
bool backGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}}; bool backGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}};
double frontGain[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{0}}}}; double frontGain[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{0}}}};
bool frontGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}}; bool frontGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}};
double uvdslope[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{0}}}};
// ==== Configuration Flags ==== // ==== Configuration Flags ====
const bool interactiveMode = false; // If true: show canvas + wait for user const bool interactiveMode = false; // If true: show canvas + wait for user
const bool interactiveMode1 = true; // If true: show canvas + wait for user
const bool verboseFit = true; // If true: print fit summary and chi² const bool verboseFit = true; // If true: print fit summary and chi²
const bool drawCanvases = true; // If false: canvases won't be drawn at all const bool drawCanvases = false; // If false: canvases won't be drawn at all
const bool drawCanvases1 = true; // If false: canvases won't be drawn at all
void GainMatchSX3Front::Begin(TTree * /*tree*/) void GainMatchSX3Front::Begin(TTree * /*tree*/)
{ {
@ -195,15 +197,16 @@ Bool_t GainMatchSX3Front::Process(Long64_t entry)
sx3EBk = sx3.e[index]; sx3EBk = sx3.e[index];
} }
} }
// If we have a valid front and back channel, fill the histograms
hSX3->Fill(sx3ChDn + 4, sx3ChBk);
hSX3->Fill(sx3ChUp, sx3ChBk);
// Fill the histogram for the front vs back
hSX3FvsB->Fill(sx3EUp + sx3EDn, sx3EBk);
for (int i = 0; i < sx3.multi; i++) for (int i = 0; i < sx3.multi; i++)
{ {
// If we have a valid front and back channel, fill the histograms
hSX3->Fill(sx3ChDn + 4, sx3ChBk);
hSX3->Fill(sx3ChUp, sx3ChBk);
// Fill the histogram for the front vs back
hSX3FvsB->Fill(sx3EUp + sx3EDn, sx3EBk);
if (sx3.e[i] > 100) // && sx3.id[i] == 4) if (sx3.e[i] > 100) // && sx3.id[i] == 4)
{ {
// back gain correction // back gain correction
@ -330,9 +333,9 @@ void GainMatchSX3Front::Terminate()
int ndf = f.GetNDF(); int ndf = f.GetNDF();
double reducedChi2 = (ndf != 0) ? chi2 / ndf : -1; double reducedChi2 = (ndf != 0) ? chi2 / ndf : -1;
std::cout << Form("Det%d U%d D%d B%d → Gain: %.4f | χ²/ndf = %.2f/%d = %.2f", // std::cout << Form("Det%d U%d D%d B%d → Gain: %.4f | χ²/ndf = %.2f/%d = %.2f",
id, u, d, bk, f.GetParameter(0), chi2, ndf, reducedChi2) // id, u, d, bk, f.GetParameter(0), chi2, ndf, reducedChi2)
<< std::endl; // << std::endl;
} }
if (interactiveMode) if (interactiveMode)
@ -354,7 +357,87 @@ void GainMatchSX3Front::Terminate()
frontGainValid[id][bk][u][d] = true; frontGainValid[id][bk][u][d] = true;
outFile << id << " " << bk << " " << u << " " << d << " " << frontGain[id][bk][u][d] << std::endl; outFile << id << " " << bk << " " << u << " " << d << " " << frontGain[id][bk][u][d] << std::endl;
printf("Front gain Det%d Back%d Up%dDn%d → %.4f\n", id, bk, u, d, frontGain[id][bk][u][d]); // printf("Front gain Det%d Back%d Up%dDn%d → %.4f\n", id, bk, u, d, frontGain[id][bk][u][d]);
}
for (const auto &kv : dataPoints)
{
auto [id, bk, u, d] = kv.first;
const auto &pts = kv.second;
std::vector<double> uvals, dvals;
if (pts.size() < 5)
continue;
std::vector<double> uE, dE, udE, corrBkE;
for (const auto &pr : pts)
{
double eBkCorr, eUp, eDn;
std::tie(eBkCorr, eUp, eDn) = pr;
if ((eBkCorr < 100) || (eUp < 100) || (eDn < 100))
continue; // Skip if any energy is zero
eUp *= frontGain[id][bk][u][d];
eDn *= frontGain[id][bk][u][d];
uE.push_back(eUp / eBkCorr);
dE.push_back(eDn / eBkCorr);
corrBkE.push_back(eBkCorr);
}
for (size_t i = 0; i < uE.size(); ++i)
{
uvals.push_back(uE[i]);
dvals.push_back(dE[i]);
}
TGraph g1(uvals.size(), uvals.data(), dvals.data());
g1.SetMarkerStyle(20);
g1.SetMarkerColor(kBlue);
g1.SetTitle(Form("Det %d: U%d D%d B%d", id, u, d, bk));
// Fit function (straight line through origin)
TF1 f1("f1", "[0]*x", 0, 16000);
f1.SetParameter(0, -1.0);
if (drawCanvases1)
{
TCanvas *c1 = new TCanvas(Form("c_%d_%d_%d_%d", id, bk, u, d), "Fit", 800, 600);
g1.SetTitle(Form("Detector %d: U%d D%d B%d", id, u, d, bk));
g1.SetMarkerStyle(20);
g1.SetMarkerColor(kBlue);
g1.Draw("AP");
g1.Fit(&f1, interactiveMode ? "Q" : "QNR"); // 'R' avoids refit, 'N' skips drawing
if (verboseFit)
{
double chi2 = f1.GetChisquare();
int ndf1 = f1.GetNDF();
double reducedChi2 = (ndf1 != 0) ? chi2 / ndf1 : -1;
std::cout << Form("Det%d U%d D%d B%d → Gain: %.4f | χ²/ndf = %.2f/%d = %.2f",
id, u, d, bk, f1.GetParameter(0), chi2, ndf1, reducedChi2)
<< std::endl;
}
if (interactiveMode1)
{
c1->Update();
gPad->WaitPrimitive();
}
else
{
c1->Close(); // Optionally avoid clutter in batch
}
}
else
{
g1.Fit(&f1, "QNR");
}
// Save slope
uvdslope[id][bk][u][d] = f1.GetParameter(0);
printf("UvD slope Det%d Back%d Up%dDn%d → %.4f\n", id, bk, u, d, uvdslope[id][bk][u][d]);
} }
outFile.close(); outFile.close();