modified: Armory/ClassPW.h changed because lowest cluster was getting deleted leading to loss of data
modified: TrackRecon.C loss/ halucination for pc calib due to ADC ceiling being placed on gainmatched files instead of the raw modified: pc_energy_calibration.dat calibration fie,updated to get the "real" calibrated signal as proposed by the Alpha calibration run
This commit is contained in:
parent
ea79180595
commit
c9339a1165
|
|
@ -322,14 +322,14 @@ PW::Make_Clusters(const std::unordered_map<int, std::tuple<int, double, double>>
|
|||
}
|
||||
|
||||
if (wireClusters.size() > 1)
|
||||
{ // Deal with wraparound if required
|
||||
auto first_cluster = wireClusters.front(); // front and back provide references to the elements themselves. less copy, can modify etc
|
||||
auto last_cluster = wireClusters.back();
|
||||
{
|
||||
auto &first_cluster = wireClusters.front();
|
||||
auto &last_cluster = wireClusters.back();
|
||||
if (std::get<0>(last_cluster.back()) == 23 && std::get<0>(first_cluster.front()) == 0)
|
||||
{
|
||||
last_cluster.insert(last_cluster.end(), first_cluster.begin(), first_cluster.end());
|
||||
wireClusters.erase(wireClusters.begin()); // invalidates the references above; not used after this
|
||||
}
|
||||
wireClusters.erase(wireClusters.begin()); // canonically, erase() needs an iterator, hence begin() not front()
|
||||
// TODO: Can also deal with 'gaps' of missing wires similarly. end of one segment and beginning of another segment will be separated by missing wire --> combine the two
|
||||
// TODO: Also needs some development regarding the time-correlation. Don't put wires in the same cluster if they aren't time coincident
|
||||
}
|
||||
|
|
@ -621,12 +621,12 @@ inline double PW::GetZ0()
|
|||
// the ring radius. The anode waist was fit to the measured anode crossover points; the
|
||||
// other families share the flare 'c' and scale their waist with the ring radius.
|
||||
//
|
||||
// Radial ordering, inner to outer: guard (ring 32) < anode (ring 37) < cathode (ring 43).
|
||||
// Radial ordering, inner to outer: guard (ring 32) < anode (ring 37) < cathode (ring 42).
|
||||
// Charge collection spans guard -> cathode. The anodes in the middle are the readout, the
|
||||
// cathodes see the induced/mirror charge, and the guard wires are field-shaping, not read out.
|
||||
const double kPCHyperbC = 301.895;
|
||||
const double kPCAnodeWaist = 32.0429; // ring radius 37 mm -- readout
|
||||
const double kPCCathodeWaist = 32.0429 * 42.0 / 37.0; // ring radius 43 mm -- induced charge
|
||||
const double kPCCathodeWaist = 32.0429 * 42.0 / 37.0; // ring radius 42 mm -- induced charge
|
||||
const double kPCGuardWaist = 32.0429 * 32.0 / 37.0; // ring radius 32 mm -- not read out
|
||||
|
||||
// Intersection of the segment x1 -> x1+dx with the one-sheet hyperboloid of waist a, flare c.
|
||||
|
|
|
|||
443
TrackRecon.C
443
TrackRecon.C
|
|
@ -38,20 +38,22 @@ Int_t colors[40] = {
|
|||
#include <utility>
|
||||
#include <stdexcept>
|
||||
#include <algorithm>
|
||||
#include <limits>
|
||||
|
||||
// --- Analysis Control Flags ---
|
||||
bool process_alpha_proton_scattering = false,
|
||||
doMiscHistograms = true,
|
||||
doRawHistos = false,
|
||||
doPCSX3ClusterAnalysis = false,
|
||||
doPCQQQClusterAnalysis = false,
|
||||
doOldAnalysis = false,
|
||||
BenchMark = false,
|
||||
onewire_analysis = true,
|
||||
doSingles = true,
|
||||
diagnostic_eplots = true,
|
||||
diagnostic_tplots = true,
|
||||
doPCSX3ClusterAnalysis = true,
|
||||
doPCQQQClusterAnalysis = true,
|
||||
BenchMark = true,
|
||||
onewire_analysis = true,
|
||||
doOldAnalysis = false,
|
||||
reactiondata = false,
|
||||
doPCEnergyCalibration = true,
|
||||
doPCEnergyCalibration = false,
|
||||
ta_foil_run = false,
|
||||
source_run = false;
|
||||
|
||||
|
|
@ -448,11 +450,15 @@ public:
|
|||
double EnergySum = -1;
|
||||
double rawEnergy1 = -1; // pre-calibration Energy1 (apSumE), for cfrac -- MeV-scale Energy1 is wrong for this
|
||||
double rawEnergy2 = -1; // pre-calibration Energy2 (cpMaxE), for cfrac
|
||||
double Energy1MeV = -1;
|
||||
double Energy2MeV = -1;
|
||||
double Time1 = -1;
|
||||
double Time2 = -1;
|
||||
int Anodech = -1;
|
||||
int Cathodech = -1;
|
||||
|
||||
bool siMatched = false;
|
||||
|
||||
// misc elements;
|
||||
int multi1 = -1, multi2 = -1;
|
||||
};
|
||||
|
|
@ -487,6 +493,8 @@ bool pcEnergyCalibLoaded = false;
|
|||
// impact on downstream histograms can be compared against the default (off).
|
||||
static const std::set<int> badAnodeWires = {6, 12, 19, 21, 22, 23};
|
||||
bool excludeBadAnodeWires = false; // set in Begin() from DISABLE_BAD_ANODE_WIRES
|
||||
bool pcUniqueCathode = true; // set in Begin() from PC_UNIQUE_CATHODE
|
||||
bool pcA1C0PerCluster = true; // set in Begin() from PC_A1C0_PER_CLUSTER
|
||||
inline bool isAnodeWireExcluded(int wire)
|
||||
{
|
||||
return excludeBadAnodeWires && badAnodeWires.count(wire) > 0;
|
||||
|
|
@ -1042,6 +1050,7 @@ inline void pcEnergyCalibrationAccumulate(const std::vector<Event> &PC_Events,
|
|||
considerSi(sx3event, TMath::Pi() / 3.0);
|
||||
}
|
||||
}
|
||||
constexpr double kPcCalApSiWindowMeV = 0.5;
|
||||
|
||||
inline void pcEnergyCalibrationAccumulateProton(const std::vector<Event> &PC_Events, const std::vector<Event> &QQQ_Events, const std::vector<Event> &SX3_Events)
|
||||
{
|
||||
|
|
@ -1092,13 +1101,21 @@ inline void pcEnergyCalibrationAccumulateProton(const std::vector<Event> &PC_Eve
|
|||
if (predicted_alpha_E <= 0.0)
|
||||
return;
|
||||
|
||||
double total_cm = pathLengthCm(vertex, sievent.pos);
|
||||
double Efix_alpha = evalEloss(MeV_to_cm_spl, cm_to_MeV_spl, sievent.Energy1, total_cm);
|
||||
const double dEfix = Efix_alpha - predicted_alpha_E;
|
||||
plotter->Fill1D("PCCal_ap_EfixAlpha_minus_pred", 400, -8, 8, Efix_alpha > 0.0 ? dEfix : -7.99, "PCCalib");
|
||||
plotter->Fill2D("PCCal_ap_EfixAlpha_minus_pred_vs_ESi", 200, 0, 10, 400, -8, 8, sievent.Energy1,
|
||||
Efix_alpha > 0.0 ? dEfix : -7.99, "PCCalib");
|
||||
if (Efix_alpha <= 0.0 || TMath::Abs(dEfix) > kPcCalApSiWindowMeV)
|
||||
return;
|
||||
|
||||
// pcCollectionPath: guard_cm = si->guard, cathode_cm = si->cathode (both from the si end).
|
||||
// Crossing order from the beam axis: vertex -> guard -> cathode -> si, so measured from
|
||||
// the vertex, dist_to_entry = total - guard_cm < dist_to_exit = total - cathode_cm.
|
||||
PCCollect pc = pcCollectionPath(vertex, sievent.pos);
|
||||
if (!pc.ok)
|
||||
return;
|
||||
double total_cm = pathLengthCm(vertex, sievent.pos);
|
||||
double dist_to_entry = total_cm - pc.guard_cm; // vertex -> guard wires, cm
|
||||
double dist_to_exit = total_cm - pc.cathode_cm; // vertex -> cathode, cm
|
||||
if (!std::isfinite(dist_to_entry) || dist_to_entry <= 0.0 ||
|
||||
|
|
@ -1302,7 +1319,7 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
|||
continue;
|
||||
}
|
||||
auto det = Fsx3.at(id);
|
||||
if (det.valid && diagnostic_eplots)
|
||||
if (det.valid && diagnostic_eplots && doSingles)
|
||||
{
|
||||
// std::cout << det.frontEL << " " << det.frontEL*sx3RightGain[id][det.stripF] << std::endl;
|
||||
plotter->Fill2D("be_vs_x_sx3_id_" + std::to_string(id) + "_f" + std::to_string(det.stripF) + "_b" + std::to_string(det.stripB), 200, -1, 1, 800, 0, 8192, det.frontX, det.backE, "evsx");
|
||||
|
|
@ -1675,6 +1692,95 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
|||
std::vector<std::vector<std::tuple<int, double, double>>> aClusters = pwinstance.Make_Clusters(aWireEvents);
|
||||
std::vector<std::vector<std::tuple<int, double, double>>> cClusters = pwinstance.Make_Clusters(cWireEvents);
|
||||
|
||||
auto pushAnodeOnly = [&](const std::vector<std::tuple<int, double, double>> &aCl) -> bool
|
||||
{
|
||||
if (aCl.size() < 1 || aCl.size() > 2)
|
||||
return false;
|
||||
|
||||
auto aPw = pwinstance.GetPseudoWire(aCl, "ANODE");
|
||||
auto apwire = std::get<0>(aPw);
|
||||
double apSumE = std::get<1>(aPw);
|
||||
double apTSMaxE = std::get<3>(aPw);
|
||||
int anodeIdx = std::get<0>(aCl[0]); // representative wire index (tag/sanity-check only,
|
||||
if (anodeIdx < 0 || anodeIdx >= 24) // not assumed to be "the" wire for A2C0's 2-wire cluster)
|
||||
return false;
|
||||
|
||||
const Event *bestSi = nullptr;
|
||||
bool bestIsQQQ = true;
|
||||
double bestDphi = 1e9;
|
||||
auto consider = [&](const std::vector<Event> &sis, bool isQQQ)
|
||||
{
|
||||
for (const auto &si : sis)
|
||||
{
|
||||
if (!siPcCoincident(si.Time1, apTSMaxE))
|
||||
continue;
|
||||
TVector3 pc = pwinstance.getClosestWirePosAtWirePhi(apwire, si.pos.Phi());
|
||||
double dphi = TMath::Abs(si.pos.DeltaPhi(pc));
|
||||
double phi_win = isQQQ ? TMath::Pi() / 4.0 : TMath::Pi() / 3.0; // per-detector, as elsewhere
|
||||
if (dphi <= phi_win && dphi < bestDphi)
|
||||
{
|
||||
bestDphi = dphi;
|
||||
bestSi = &si;
|
||||
bestIsQQQ = isQQQ;
|
||||
}
|
||||
}
|
||||
};
|
||||
consider(QQQ_Events, true);
|
||||
consider(SX3_Events, false);
|
||||
if (!bestSi)
|
||||
return false; // no phi reference -> no anode-only position can be built
|
||||
|
||||
bool isA2C0 = (aCl.size() == 2);
|
||||
TVector3 pc = isA2C0 ? a2c0_wirePos(apwire, bestSi->pos.Phi(), bestIsQQQ)
|
||||
: a1c0_wirePos(apwire, bestSi->pos.Phi(), bestIsQQQ); // same z reference as the benchmark
|
||||
|
||||
Event PCEventRaw(pc, apSumE, -1.0, apTSMaxE, -1.0);
|
||||
PCEventRaw.multi1 = static_cast<int>(aCl.size());
|
||||
PCEventRaw.multi2 = 0;
|
||||
PCEventRaw.Anodech = anodeIdx;
|
||||
PCEventRaw.Cathodech = -1;
|
||||
PCEventRaw.siMatched = true; // true by construction -- bestSi is required above
|
||||
|
||||
double anodeCalibSum = 0.0; // per-wire, see Energy1MeV on Event
|
||||
for (const auto &w : aCl)
|
||||
{
|
||||
int wi = std::get<0>(w);
|
||||
if (wi >= 0 && wi < 24)
|
||||
anodeCalibSum += pcEnergySlope[wi] * std::get<1>(w);
|
||||
}
|
||||
PCEventRaw.Energy1MeV = anodeCalibSum; // Energy2MeV stays -1: no cathode
|
||||
PC_Events.push_back(PCEventRaw);
|
||||
|
||||
if (pcEnergyCalibLoaded)
|
||||
{
|
||||
Event ev(pc, anodeCalibSum, -1.0, apTSMaxE, -1.0);
|
||||
ev.multi1 = static_cast<int>(aCl.size());
|
||||
ev.multi2 = 0; // no cathode -> a1c0/a2c0 topology in pcCalibratedHistograms
|
||||
ev.Anodech = anodeIdx;
|
||||
ev.Cathodech = -1;
|
||||
ev.siMatched = true;
|
||||
PC_Events_calibrated.push_back(ev);
|
||||
}
|
||||
return true;
|
||||
};
|
||||
|
||||
auto anySiCoincident = [&](const TVector3 &pcPos, double pcTime) -> bool
|
||||
{
|
||||
auto scan = [&](const std::vector<Event> &sis, bool isQQQ)
|
||||
{
|
||||
double phi_win = isQQQ ? TMath::Pi() / 4.0 : TMath::Pi() / 3.0;
|
||||
for (const auto &si : sis)
|
||||
{
|
||||
if (!siPcCoincident(si.Time1, pcTime))
|
||||
continue;
|
||||
if (TMath::Abs(si.pos.DeltaPhi(pcPos)) <= phi_win)
|
||||
return true;
|
||||
}
|
||||
return false;
|
||||
};
|
||||
return scan(QQQ_Events, true) || scan(SX3_Events, false);
|
||||
};
|
||||
|
||||
for (const auto &aCluster : aClusters)
|
||||
{
|
||||
if (clusterHasExcludedAnode(aCluster))
|
||||
|
|
@ -1711,7 +1817,7 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
|||
double ahi_MeV = std::max(ae0_MeV, ae1_MeV);
|
||||
double aratio_MeV = alo_MeV / ahi_MeV;
|
||||
plotter->Fill1D("A2_anode_ratio_calib", 120, 0, 1.2, aratio_MeV, "hGMPC");
|
||||
plotter->Fill2D("A1_vs_A2_calib", 800, 0, 0.6, 800, 0, 0.6, ae0_MeV, ae1_MeV, "hGMPC");
|
||||
plotter->Fill2D("A1_vs_A2_calib", 800, 0, 0.4, 800, 0, 0.4, ae0_MeV, ae1_MeV, "hGMPC");
|
||||
plotter->Fill2D("A2_anode_ratio_calib_vs_lowerIndex", 24, 0, 24, 120, 0, 1.2,
|
||||
std::min(wi0, wi1), aratio_MeV, "hGMPC");
|
||||
}
|
||||
|
|
@ -1722,6 +1828,71 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
|||
{
|
||||
plotter->Fill1D("Raw_A1_AnodeSum", 800, 0, 40000, std::get<1>(aCluster[0]), "hGMPC");
|
||||
}
|
||||
auto pushPaired = [&](const std::vector<std::tuple<int, double, double>> &cCluster)
|
||||
{
|
||||
auto [crossover, alpha, apSumE, cpSumE, apMaxE, cpMaxE, apTSMaxE, cpTSMaxE] = pwinstance.FindCrossoverProperties(aCluster, cCluster);
|
||||
Event PCEvent(crossover, apSumE, cpMaxE, cpSumE, apTSMaxE, cpTSMaxE); // run12 shows cathode-max and anode-sum provide best dE signals.
|
||||
// std::cout << apSumE << " " << crossover.Perp() << " " << apMaxE << " " << apTSMaxE << std::endl;
|
||||
PCEvent.multi1 = aCluster.size();
|
||||
PCEvent.multi2 = cCluster.size();
|
||||
PCEvent.Anodech = std::get<0>(aCluster[0]);
|
||||
PCEvent.Cathodech = std::get<0>(cCluster[0]);
|
||||
PCEvent.siMatched = anySiCoincident(crossover, apTSMaxE);
|
||||
|
||||
// Per-wire calibrated energies, computed once here and carried on PCEvent
|
||||
// itself so PC_Events consumers don't re-derive them from Anodech/Cathodech
|
||||
// (see Energy1MeV/Energy2MeV on Event). With no table loaded pcEnergySlope
|
||||
// is all 1.0, so these reduce to apSumE / cpMaxE.
|
||||
double anodeCalibSum = 0.0;
|
||||
for (const auto &w : aCluster)
|
||||
{
|
||||
int wi = std::get<0>(w);
|
||||
if (wi >= 0 && wi < 24)
|
||||
anodeCalibSum += pcEnergySlope[wi] * std::get<1>(w);
|
||||
}
|
||||
// Cathode uses the single max-energy wire (cpMaxE). That wire is NOT
|
||||
// necessarily cCluster[0], which is all PCEvent.Cathodech records, so
|
||||
// pcEnergySlope[24 + Cathodech] was applying the wrong wire's constant to
|
||||
// cpMaxE for every multi-wire cathode cluster -- i.e. for A1C2, the primary
|
||||
// topology. GetPseudoWire tracks the max energy but not its index, so find
|
||||
// it here rather than change that signature for its five call sites.
|
||||
int cMaxWire = PCEvent.Cathodech;
|
||||
double cMaxE = -1.0;
|
||||
for (const auto &w : cCluster)
|
||||
{
|
||||
if (std::get<1>(w) > cMaxE)
|
||||
{
|
||||
cMaxE = std::get<1>(w);
|
||||
cMaxWire = std::get<0>(w);
|
||||
}
|
||||
}
|
||||
const bool cMaxWireOk = (cMaxWire >= 0 && cMaxWire < 24);
|
||||
PCEvent.Energy1MeV = anodeCalibSum;
|
||||
PCEvent.Energy2MeV = cMaxWireOk ? pcEnergySlope[24 + cMaxWire] * cpMaxE : -1.0;
|
||||
PC_Events.push_back(PCEvent);
|
||||
|
||||
if (pcEnergyCalibLoaded)
|
||||
{
|
||||
Event PCEventCalibrated = PCEvent;
|
||||
PCEventCalibrated.rawEnergy1 = PCEvent.Energy1; // stash BEFORE overwriting -- see rawEnergy1/2 comment on Event
|
||||
PCEventCalibrated.rawEnergy2 = PCEvent.Energy2;
|
||||
PCEventCalibrated.Energy1 = anodeCalibSum;
|
||||
PCEventCalibrated.Energy2 = cMaxWireOk ? PCEvent.Energy2MeV : cpMaxE;
|
||||
PC_Events_calibrated.push_back(PCEventCalibrated);
|
||||
}
|
||||
};
|
||||
|
||||
// --- Cathode selection. Previously this loop pushed inside the body, once
|
||||
// per crossing cathode cluster, so an anode cluster crossing N cathode
|
||||
// clusters landed N times in PC_Events/PC_Events_calibrated while the
|
||||
// anode-only topologies could only ever land once -- the source of the
|
||||
// a1c1/a1c2 vs a1c0 count asymmetry in the Calib_AnodeE* family. Now the
|
||||
// candidates are enumerated first and at most one is pushed. ---
|
||||
const std::vector<std::tuple<int, double, double>> *bestC = nullptr;
|
||||
double bestCathodeE = -1.0;
|
||||
int bestCathodeWire = std::numeric_limits<int>::max();
|
||||
int nCathodeMatches = 0;
|
||||
|
||||
for (const auto &cCluster : cClusters)
|
||||
{
|
||||
if (aCluster.size() == 0)
|
||||
|
|
@ -1730,129 +1901,41 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
|||
continue;
|
||||
// both have at least 1, here. Keep the a1, c1 events
|
||||
auto [crossover, alpha, apSumE, cpSumE, apMaxE, cpMaxE, apTSMaxE, cpTSMaxE] = pwinstance.FindCrossoverProperties(aCluster, cCluster);
|
||||
if (alpha != 9999999 && apSumE != -1)
|
||||
{
|
||||
// Event PCEvent(crossover,apMaxE,cpMaxE,apTSMaxE,cpTSMaxE);
|
||||
// Event PCEvent(crossover,apSumE,cpSumE,apTSMaxE,cpTSMaxE);
|
||||
Event PCEvent(crossover, apSumE, cpMaxE, cpSumE, apTSMaxE, cpTSMaxE); // run12 shows cathode-max and anode-sum provide best dE signals.
|
||||
// std::cout << apSumE << " " << crossover.Perp() << " " << apMaxE << " " << apTSMaxE << std::endl;
|
||||
PCEvent.multi1 = aCluster.size();
|
||||
PCEvent.multi2 = cCluster.size();
|
||||
PCEvent.Anodech = std::get<0>(aCluster[0]);
|
||||
PCEvent.Cathodech = std::get<0>(cCluster[0]);
|
||||
PC_Events.push_back(PCEvent);
|
||||
if (!(alpha != 9999999 && apSumE != -1))
|
||||
continue; // no valid crossover -- was the silent `else ;` branch
|
||||
++nCathodeMatches;
|
||||
|
||||
if (pcEnergyCalibLoaded)
|
||||
{
|
||||
Event PCEventCalibrated = PCEvent;
|
||||
PCEventCalibrated.rawEnergy1 = PCEvent.Energy1; // stash BEFORE overwriting -- see rawEnergy1/2 comment on Event
|
||||
PCEventCalibrated.rawEnergy2 = PCEvent.Energy2;
|
||||
double anodeCalibSum = 0.0;
|
||||
for (const auto &w : aCluster)
|
||||
{
|
||||
int wi = std::get<0>(w);
|
||||
if (wi >= 0 && wi < 24)
|
||||
anodeCalibSum += pcEnergySlope[wi] * std::get<1>(w);
|
||||
}
|
||||
PCEventCalibrated.Energy1 = anodeCalibSum;
|
||||
// Cathode uses the single max-energy wire (cpMaxE). That wire is NOT
|
||||
// necessarily cCluster[0], which is all PCEvent.Cathodech records, so
|
||||
// pcEnergySlope[24 + Cathodech] was applying the wrong wire's constant to
|
||||
// cpMaxE for every multi-wire cathode cluster -- i.e. for A1C2, the primary
|
||||
// topology. GetPseudoWire tracks the max energy but not its index, so find
|
||||
// it here rather than change that signature for its five call sites.
|
||||
int cMaxWire = PCEvent.Cathodech;
|
||||
double cMaxE = -1.0;
|
||||
for (const auto &w : cCluster)
|
||||
{
|
||||
if (std::get<1>(w) > cMaxE)
|
||||
{
|
||||
cMaxE = std::get<1>(w);
|
||||
cMaxWire = std::get<0>(w);
|
||||
}
|
||||
}
|
||||
PCEventCalibrated.Energy2 = (cMaxWire >= 0 && cMaxWire < 24)
|
||||
? pcEnergySlope[24 + cMaxWire] * cpMaxE
|
||||
: cpMaxE;
|
||||
PC_Events_calibrated.push_back(PCEventCalibrated);
|
||||
}
|
||||
}
|
||||
else
|
||||
if (!pcUniqueCathode)
|
||||
{
|
||||
; // std::cout << "AAAA " << std::endl;
|
||||
pushPaired(cCluster); // legacy: one event per crossing pair
|
||||
continue;
|
||||
}
|
||||
|
||||
int firstWire = std::get<0>(cCluster[0]);
|
||||
if (cpMaxE > bestCathodeE || (cpMaxE == bestCathodeE && firstWire < bestCathodeWire))
|
||||
{
|
||||
bestC = &cCluster;
|
||||
bestCathodeE = cpMaxE;
|
||||
bestCathodeWire = firstWire;
|
||||
}
|
||||
}
|
||||
|
||||
plotter->Fill1D("Calib_nCathodeMatch_per_anode", 10, 0, 10, nCathodeMatches, "hCalibPC");
|
||||
|
||||
if (pcUniqueCathode && bestC)
|
||||
pushPaired(*bestC);
|
||||
|
||||
if (pcA1C0PerCluster && nCathodeMatches == 0)
|
||||
pushAnodeOnly(aCluster);
|
||||
}
|
||||
|
||||
if (cClusters.empty())
|
||||
if (!pcA1C0PerCluster && cClusters.empty())
|
||||
{
|
||||
for (const auto &aCl : aClusters)
|
||||
{
|
||||
if (aCl.size() < 1 || aCl.size() > 2) // A1C0 (1 wire) or A2C0 (2 wires) --
|
||||
continue; // reaction_ax_core / protonAlphaElastic_core's
|
||||
// a1c0 convention, one wire wider for A2C0.
|
||||
if (clusterHasExcludedAnode(aCl))
|
||||
continue;
|
||||
auto aPw = pwinstance.GetPseudoWire(aCl, "ANODE");
|
||||
auto apwire = std::get<0>(aPw);
|
||||
double apSumE = std::get<1>(aPw);
|
||||
double apTSMaxE = std::get<3>(aPw);
|
||||
int anodeIdx = std::get<0>(aCl[0]); // representative wire index (tag/sanity-check only,
|
||||
if (anodeIdx < 0 || anodeIdx >= 24) // not assumed to be "the" wire for A2C0's 2-wire cluster)
|
||||
continue;
|
||||
|
||||
const Event *bestSi = nullptr;
|
||||
bool bestIsQQQ = true;
|
||||
double bestDphi = 1e9;
|
||||
auto consider = [&](const std::vector<Event> &sis, bool isQQQ)
|
||||
{
|
||||
for (const auto &si : sis)
|
||||
{
|
||||
if (!siPcCoincident(si.Time1, apTSMaxE))
|
||||
continue;
|
||||
TVector3 pc = pwinstance.getClosestWirePosAtWirePhi(apwire, si.pos.Phi());
|
||||
double dphi = TMath::Abs(si.pos.DeltaPhi(pc));
|
||||
double phi_win = isQQQ ? TMath::Pi() / 4.0 : TMath::Pi() / 3.0; // per-detector, as elsewhere
|
||||
if (dphi <= phi_win && dphi < bestDphi)
|
||||
{
|
||||
bestDphi = dphi;
|
||||
bestSi = &si;
|
||||
bestIsQQQ = isQQQ;
|
||||
}
|
||||
}
|
||||
};
|
||||
consider(QQQ_Events, true);
|
||||
consider(SX3_Events, false);
|
||||
if (!bestSi)
|
||||
continue;
|
||||
|
||||
bool isA2C0 = (aCl.size() == 2);
|
||||
TVector3 pc = isA2C0 ? a2c0_wirePos(apwire, bestSi->pos.Phi(), bestIsQQQ)
|
||||
: a1c0_wirePos(apwire, bestSi->pos.Phi(), bestIsQQQ); // same z reference as the benchmark
|
||||
|
||||
Event PCEventRaw(pc, apSumE, -1.0, apTSMaxE, -1.0);
|
||||
PCEventRaw.multi1 = static_cast<int>(aCl.size());
|
||||
PCEventRaw.multi2 = 0;
|
||||
PCEventRaw.Anodech = anodeIdx;
|
||||
PCEventRaw.Cathodech = -1;
|
||||
PC_Events.push_back(PCEventRaw);
|
||||
|
||||
if (pcEnergyCalibLoaded)
|
||||
{
|
||||
double anodeCalibSum = 0.0;
|
||||
for (const auto &w : aCl)
|
||||
{
|
||||
int wi = std::get<0>(w);
|
||||
if (wi >= 0 && wi < 24)
|
||||
anodeCalibSum += pcEnergySlope[wi] * std::get<1>(w);
|
||||
}
|
||||
Event ev(pc, anodeCalibSum, -1.0, apTSMaxE, -1.0);
|
||||
ev.multi1 = static_cast<int>(aCl.size());
|
||||
ev.multi2 = 0; // no cathode -> a1c0/a2c0 topology in pcCalibratedHistograms
|
||||
ev.Anodech = anodeIdx;
|
||||
ev.Cathodech = -1;
|
||||
PC_Events_calibrated.push_back(ev);
|
||||
}
|
||||
pushAnodeOnly(aCl);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -2141,6 +2224,15 @@ void TrackRecon::Terminate()
|
|||
std::string outname = "pc_calib_raw/points_" + tag + ".dat";
|
||||
std::ofstream outfile(outname);
|
||||
outfile << std::scientific << std::setprecision(6);
|
||||
|
||||
const double kPcRawOverflowADC = 64000.0;
|
||||
for (int wire = 0; wire < 48; ++wire)
|
||||
{
|
||||
double thr = pcSlope[wire] * kPcRawOverflowADC + pcIntercept[wire];
|
||||
if (wire >= 24)
|
||||
thr *= cathode_gain;
|
||||
outfile << "# overflow " << wire << " " << thr << "\n";
|
||||
}
|
||||
long long nPoints = 0;
|
||||
for (int wire = 0; wire < 48; ++wire)
|
||||
{
|
||||
|
|
@ -2459,35 +2551,40 @@ void pcCalibratedHistograms(HistPlotter *plotter, const std::vector<Event> &QQQ_
|
|||
const std::string topo = "_a" + std::to_string(pcevent.multi1) + "c" + std::to_string(pcevent.multi2);
|
||||
const bool hasCathode = (pcevent.Cathodech >= 0);
|
||||
if (hasCathode)
|
||||
plotter->Fill2D("Calib_AnodeE_vs_CathodeE_a1c1andup", 800, 0, 0.6, 800, 0, 0.6, pcevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_AnodeE_vs_CathodeE_a1c1andup", 800, 0, 0.4, 800, 0, 0.4, pcevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
for (const std::string &t : {std::string(""), topo})
|
||||
{
|
||||
plotter->Fill2D("Calib_AnodeE_vs_AnodeIndex" + t, 24, 0, 24, 800, 0, 0.6, pcevent.Anodech, pcevent.Energy1, "hCalibPC");
|
||||
plotter->Fill1D("Calib_AnodeE" + t, 800, 0, 0.6, pcevent.Energy1, "hCalibPC");
|
||||
plotter->Fill2D("Calib_AnodeE_vs_AnodeIndex" + t, 24, 0, 24, 800, 0, 0.4, pcevent.Anodech, pcevent.Energy1, "hCalibPC");
|
||||
plotter->Fill1D("Calib_AnodeE" + t, 800, 0, 0.4, pcevent.Energy1, "hCalibPC");
|
||||
if (pcevent.siMatched)
|
||||
{
|
||||
plotter->Fill1D("Calib_AnodeE" + t + "_siGated", 800, 0, 0.4, pcevent.Energy1, "hCalibPC");
|
||||
plotter->Fill2D("Calib_AnodeE_vs_AnodeIndex" + t + "_siGated", 24, 0, 24, 800, 0, 0.4, pcevent.Anodech, pcevent.Energy1, "hCalibPC");
|
||||
}
|
||||
if (hasCathode)
|
||||
{
|
||||
plotter->Fill2D("Calib_CathodeE_vs_CathodeIndex" + t, 24, 0, 24, 800, 0, 0.6, pcevent.Cathodech, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill1D("Calib_CathodeE" + t, 800, 0, 0.6, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_AnodeE_vs_CathodeE" + t, 800, 0, 0.6, 800, 0, 0.6, pcevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_CathodeE_vs_CathodeIndex" + t, 24, 0, 24, 800, 0, 0.4, pcevent.Cathodech, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill1D("Calib_CathodeE" + t, 800, 0, 0.4, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_AnodeE_vs_CathodeE" + t, 800, 0, 0.4, 800, 0, 0.4, pcevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
}
|
||||
|
||||
for (const auto &qqqevent : QQQ_Events)
|
||||
{
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_QQQE" + t, 400, 0, 10, 800, 0, 0.6, qqqevent.Energy1, pcevent.Energy1, "hCalibPC");
|
||||
// if (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
||||
// plotter->Fill2D("Calib_dE_AnodeE_vs_QQQE" + t + "_anode" + pad2(pcevent.Anodech),
|
||||
// 400, 0, 10, 800, 0, 0.6, qqqevent.Energy1, pcevent.Energy1, "EdE_wire");
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_QQQE" + t, 400, 0, 10, 800, 0, 0.4, qqqevent.Energy1, pcevent.Energy1, "hCalibPC"); // 0.4 clipped the alpha blob
|
||||
if (pcevent.Anodech >= 0 && pcevent.Anodech < 24 && doSingles)
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_QQQE" + t + "_anode" + pad2(pcevent.Anodech),
|
||||
400, 0, 10, 800, 0, 0.4, qqqevent.Energy1, pcevent.Energy1, "EdE_wire");
|
||||
if (hasCathode)
|
||||
plotter->Fill2D("Calib_dE_CathodeE_vs_QQQE" + t, 400, 0, 10, 800, 0, 0.6, qqqevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dE_CathodeE_vs_QQQE" + t, 400, 0, 10, 800, 0, 0.4, qqqevent.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
}
|
||||
for (const auto &sx3event : SX3_Events)
|
||||
{
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_SX3E" + t, 400, 0, 10, 800, 0, 0.6, sx3event.Energy1, pcevent.Energy1, "hCalibPC");
|
||||
// if (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
||||
// plotter->Fill2D("Calib_dE_AnodeE_vs_SX3E" + t + "_anode" + pad2(pcevent.Anodech),
|
||||
// 400, 0, 10, 800, 0, 0.6, sx3event.Energy1, pcevent.Energy1, "EdE_wire");
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_SX3E" + t, 400, 0, 10, 800, 0, 0.4, sx3event.Energy1, pcevent.Energy1, "hCalibPC"); // 0.4 clipped the alpha blob
|
||||
if (pcevent.Anodech >= 0 && pcevent.Anodech < 24 && doSingles)
|
||||
plotter->Fill2D("Calib_dE_AnodeE_vs_SX3E" + t + "_anode" + pad2(pcevent.Anodech),
|
||||
400, 0, 10, 800, 0, 0.4, sx3event.Energy1, pcevent.Energy1, "EdE_wire");
|
||||
if (hasCathode)
|
||||
plotter->Fill2D("Calib_dE_CathodeE_vs_SX3E" + t, 400, 0, 10, 800, 0, 0.6, sx3event.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dE_CathodeE_vs_SX3E" + t, 400, 0, 10, 800, 0, 0.4, sx3event.Energy1, pcevent.Energy2, "hCalibPC");
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -2541,11 +2638,11 @@ void pcCalibratedHistograms(HistPlotter *plotter, const std::vector<Event> &QQQ_
|
|||
|
||||
double Egu_p = evalEloss(MeV_to_cm_p_spl, cm_to_MeVp_spl, qqqevent.Energy1, pcc.guard_cm);
|
||||
double Eca_p = evalEloss(MeV_to_cm_p_spl, cm_to_MeVp_spl, qqqevent.Energy1, pcc.cathode_cm);
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asProton" + topo, 400, 0, 0.6, 400, 0, 0.6, pcevent.Energy1, Egu_p - Eca_p, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asProton" + topo, 400, 0, 0.4, 400, 0, 0.4, pcevent.Energy1, Egu_p - Eca_p, "hCalibPC");
|
||||
|
||||
double Egu_a = evalEloss(MeV_to_cm_spl, cm_to_MeV_spl, qqqevent.Energy1, pcc.guard_cm);
|
||||
double Eca_a = evalEloss(MeV_to_cm_spl, cm_to_MeV_spl, qqqevent.Energy1, pcc.cathode_cm);
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asAlpha" + topo, 400, 0, 0.6, 400, 0, 0.6, pcevent.Energy1, Egu_a - Eca_a, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asAlpha" + topo, 400, 0, 0.4, 400, 0, 0.4, pcevent.Energy1, Egu_a - Eca_a, "hCalibPC");
|
||||
}
|
||||
for (const auto &sx3event : SX3_Events)
|
||||
{
|
||||
|
|
@ -2566,11 +2663,11 @@ void pcCalibratedHistograms(HistPlotter *plotter, const std::vector<Event> &QQQ_
|
|||
|
||||
double Egu_p = evalEloss(MeV_to_cm_p_spl, cm_to_MeVp_spl, sx3event.Energy1, pcc.guard_cm);
|
||||
double Eca_p = evalEloss(MeV_to_cm_p_spl, cm_to_MeVp_spl, sx3event.Energy1, pcc.cathode_cm);
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asProton" + topo, 400, 0, 0.6, 400, 0, 0.6, pcevent.Energy1, Egu_p - Eca_p, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asProton" + topo, 400, 0, 0.4, 400, 0, 0.4, pcevent.Energy1, Egu_p - Eca_p, "hCalibPC");
|
||||
|
||||
double Egu_a = evalEloss(MeV_to_cm_spl, cm_to_MeV_spl, sx3event.Energy1, pcc.guard_cm);
|
||||
double Eca_a = evalEloss(MeV_to_cm_spl, cm_to_MeV_spl, sx3event.Energy1, pcc.cathode_cm);
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asAlpha" + topo, 400, 0, 0.6, 400, 0, 0.6, pcevent.Energy1, Egu_a - Eca_a, "hCalibPC");
|
||||
plotter->Fill2D("Calib_dEgasPred_vs_dEgasCalib_asAlpha" + topo, 400, 0, 0.4, 400, 0, 0.4, pcevent.Energy1, Egu_a - Eca_a, "hCalibPC");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
@ -3948,7 +4045,7 @@ void protonAlphaElastic_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
if (TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) > phi_win)
|
||||
continue;
|
||||
plotter->Fill2D(rx + "_Ex_vs_dT" + ejtag + sfx, 500, -2000, 2000, 600, -10, 20, (sievent.Time1 - pcevent.Time1), Ex, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.6, (sievent.Time1 - pcevent.Time1), anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.4, (sievent.Time1 - pcevent.Time1), anodeE_MeV, pmlabel);
|
||||
}
|
||||
|
||||
// Ground-state beam-energy consistency check -- elastic scattering has
|
||||
|
|
@ -3984,7 +4081,7 @@ void protonAlphaElastic_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
PCCollect pcc = pcCollectionPath(r_rhoMin_fix, sievent.pos);
|
||||
if (pcc.ok && anodeE_MeV > 0.0)
|
||||
{
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ex" + ejtag + t + sfx, 800, -10, 10, 800, 0, 0.6, Ex, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ex" + ejtag + t + sfx, 800, -10, 10, 800, 0, 0.4, Ex, anodeE_MeV, pmlabel);
|
||||
}
|
||||
};
|
||||
std::string topo1;
|
||||
|
|
@ -4005,19 +4102,19 @@ void protonAlphaElastic_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
double E_gu = evalEloss(ej_fwd, ej_inv, sievent.Energy1, pcc.guard_cm);
|
||||
double E_ca = evalEloss(ej_fwd, ej_inv, sievent.Energy1, pcc.cathode_cm);
|
||||
double dE_pred = E_gu - E_ca;
|
||||
plotter->Fill2D(rx + "_dEgas_vs_Ef" + ejtag + sfx, 400, 0, 10, 400, 0, 0.6, Efix, dE_pred, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgas_vs_Ef" + ejtag + sfx, 400, 0, 10, 400, 0, 0.4, Efix, dE_pred, pmlabel);
|
||||
if (anodeE_MeV > 0.0)
|
||||
{
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ef" + ejtag + sfx, 400, 0, 10, 800, 0, 0.6, Efix, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx, 400, 0, 10, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 0.6, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_theta" + ejtag + sfx, 100, 0, 180, 800, 0, 0.6, theta * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_phi" + ejtag + sfx, 100, -200, 200, 800, 0, 0.6, sievent.pos.Phi() * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ef" + ejtag + sfx, 400, 0, 10, 800, 0, 0.4, Efix, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx, 400, 0, 10, 800, 0, 0.4, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 0.4, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_theta" + ejtag + sfx, 100, 0, 180, 800, 0, 0.4, theta * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_phi" + ejtag + sfx, 100, -200, 200, 800, 0, 0.4, sievent.pos.Phi() * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
if (anodeCh >= 0 && anodeCh < 24)
|
||||
// plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx + "_anode" + pad2(anodeCh),
|
||||
// 400, 0, 10, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Z" + ejtag + sfx, 800, -400, 400, 800, 0, 0.6, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasPred_vs_dEgasCalib" + ejtag + sfx, 800, 0, 0.6, 800, 0, 0.6, anodeE_MeV, dE_pred, pmlabel);
|
||||
// 400, 0, 10, 800, 0, 0.4, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Z" + ejtag + sfx, 800, -400, 400, 800, 0, 0.4, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasPred_vs_dEgasCalib" + ejtag + sfx, 800, 0, 0.4, 800, 0, 0.4, anodeE_MeV, dE_pred, pmlabel);
|
||||
}
|
||||
}
|
||||
};
|
||||
|
|
@ -4050,9 +4147,7 @@ void protonAlphaElastic_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
if (hasCathode)
|
||||
plotter->Fill1D(rx + "_dt_Cathode" + sfx, 600, -2000, 2000, pcevent.Time2 - sievent.Time1, misclabel);
|
||||
|
||||
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
||||
? pcEnergySlope[pcevent.Anodech] * pcevent.Energy1
|
||||
: -1.0;
|
||||
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24) ? pcevent.Energy1MeV : -1.0;
|
||||
SiPcPid pid = classifyByAnodeDe(anodeE_MeV);
|
||||
if (pid == SiPcPid::kUnknown)
|
||||
plotter->Fill1D(rx + "_pidUnknown" + sfx, 2, 0, 2, 1.0, misclabel);
|
||||
|
|
@ -4424,8 +4519,8 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
if (TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) > phi_win)
|
||||
continue;
|
||||
plotter->Fill2D(rx + "_Ex_vs_dT" + ejtag + sfx, 500, -2000, 2000, 600, -10, 20, (sievent.Time1 - pcevent.Time1), Ex, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.6, (sievent.Time1 - pcevent.Time1), anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalibCathode_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.6, (sievent.Time1 - pcevent.Time1), cathodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.4, (sievent.Time1 - pcevent.Time1), anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalibCathode_vs_dT" + ejtag + sfx, 500, -2000, 2000, 800, 0, 0.4, (sievent.Time1 - pcevent.Time1), cathodeE_MeV, pmlabel);
|
||||
}
|
||||
|
||||
if (dt_rf_mcp > -900000000)
|
||||
|
|
@ -4441,25 +4536,25 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
double E_gu = evalEloss(ej_fwd, ej_inv, sievent.Energy1, pcc.guard_cm);
|
||||
double E_ca = evalEloss(ej_fwd, ej_inv, sievent.Energy1, pcc.cathode_cm);
|
||||
double dE_pred = E_gu - E_ca;
|
||||
plotter->Fill2D(rx + "_dEgas_vs_Ef" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.6, Efix, dE_pred, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgas_vs_Ef" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.4, Efix, dE_pred, pmlabel);
|
||||
if (anodeE_MeV >= 0.0)
|
||||
{
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + sfx, 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, folderPrefix + "EdEComparison");
|
||||
plotter->Fill2D(rx + "_dEgasCalib*sintheta_vs_E" + sfx, 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, anodeE_MeV * sin(theta), folderPrefix + "EdEComparison");
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ef" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.6, Efix, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_EBeam" + ejtag + sfx, 400, 0, beamE0 * 1.5, 800, 0, 0.6, beam_energy_at_vertex, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + sfx, 400, 0, ef_max, 800, 0, 0.4, sievent.Energy1, anodeE_MeV, folderPrefix + "EdEComparison");
|
||||
plotter->Fill2D(rx + "_dEgasCalib*sintheta_vs_E" + sfx, 400, 0, ef_max, 800, 0, 0.4, sievent.Energy1, anodeE_MeV * sin(theta), folderPrefix + "EdEComparison");
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ef" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.4, Efix, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_EBeam" + ejtag + sfx, 400, 0, beamE0 * 1.5, 800, 0, 0.4, beam_energy_at_vertex, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasRaw_vs_EBeam" + ejtag + sfx, 400, 0, beamE0 * 1.5, 800, 0, 20000, beam_energy_at_vertex, anodeE, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 0.6, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx, 400, 0, ef_max, 800, 0, 0.4, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 0.4, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasRaw_vs_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 20000, vertex_z, anodeE, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasRaw_vs_theta" + ejtag + sfx, 180, 0, 180, 800, 0, 20000, theta * 180 / M_PI, anodeE, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_theta" + ejtag + sfx, 360, 0, 180, 800, 0, 0.6, theta * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_phi" + ejtag + sfx, 90, -180, 180, 800, 0, 0.6, phi * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_theta" + ejtag + sfx, 360, 0, 180, 800, 0, 0.4, theta * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_phi" + ejtag + sfx, 90, -180, 180, 800, 0, 0.4, phi * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||
// if (anodeCh >= 0)
|
||||
// plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx + "_anode" + pad2(anodeCh), 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ex" + ejtag + sfx, 600, -10, 20, 800, 0, 0.6, Ex, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Z" + ejtag + sfx, 800, -400, 400, 800, 0, 0.6, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasPred_vs_dEgasCalib" + ejtag + sfx, 800, 0, 0.6, 800, 0, 0.6, anodeE_MeV, dE_pred, pmlabel);
|
||||
// plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx + "_anode" + pad2(anodeCh), 400, 0, ef_max, 800, 0, 0.4, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Ex" + ejtag + sfx, 600, -10, 20, 800, 0, 0.4, Ex, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasCalib_vs_Z" + ejtag + sfx, 800, -400, 400, 800, 0, 0.4, vertex_z, anodeE_MeV, pmlabel);
|
||||
plotter->Fill2D(rx + "_dEgasPred_vs_dEgasCalib" + ejtag + sfx, 800, 0, 0.4, 800, 0, 0.4, anodeE_MeV, dE_pred, pmlabel);
|
||||
}
|
||||
}
|
||||
};
|
||||
|
|
@ -4491,12 +4586,8 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
|||
bool timecut = siPcCoincident(sievent.Time1, pcevent.Time1);
|
||||
if (!(phicut && timecut))
|
||||
continue;
|
||||
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
||||
? pcEnergySlope[pcevent.Anodech] * pcevent.Energy1
|
||||
: -1.0;
|
||||
double cathodeE_MeV = (pcevent.Cathodech >= 0 && pcevent.Cathodech < 24)
|
||||
? pcEnergySlope[24 + pcevent.Cathodech] * pcevent.Energy2
|
||||
: -1.0;
|
||||
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24) ? pcevent.Energy1MeV : -1.0;
|
||||
double cathodeE_MeV = (pcevent.Cathodech >= 0 && pcevent.Cathodech < 24) ? pcevent.Energy2MeV : -1.0;
|
||||
|
||||
const bool isA2 = (pcevent.multi1 == 2);
|
||||
const std::string mg = isA2 ? "a2c1c2" : "a1c1c2";
|
||||
|
|
|
|||
|
|
@ -1,48 +1,48 @@
|
|||
0 7.569000E-06 0.000000E+00 1
|
||||
1 1.283082E-05 0.000000E+00 1
|
||||
2 1.283082E-05 0.000000E+00 1
|
||||
3 1.283082E-05 0.000000E+00 1
|
||||
4 1.283082E-05 0.000000E+00 1
|
||||
5 1.496409E-05 0.000000E+00 2
|
||||
6 2.036882E-05 0.000000E+00 1
|
||||
7 1.124799E-05 0.000000E+00 1
|
||||
8 1.123432E-05 0.000000E+00 1
|
||||
9 3.333333E-01 0.000000E+00 3
|
||||
10 1.292001E-05 0.000000E+00 1
|
||||
11 8.787520E-06 0.000000E+00 1
|
||||
12 3.333333E-01 0.000000E+00 3
|
||||
13 7.867300E-06 0.000000E+00 1
|
||||
14 6.707640E-06 0.000000E+00 1
|
||||
15 8.515706E-06 0.000000E+00 1
|
||||
16 1.041633E-05 0.000000E+00 1
|
||||
17 1.041633E-05 0.000000E+00 1
|
||||
18 1.283082E-05 0.000000E+00 1
|
||||
19 9.043133E-06 0.000000E+00 2
|
||||
20 1.046743E-05 0.000000E+00 1
|
||||
21 1.046743E-05 0.000000E+00 2
|
||||
22 1.047560E-05 0.000000E+00 2
|
||||
23 1.283082E-05 0.000000E+00 2
|
||||
24 1.008216E-05 0.000000E+00 1
|
||||
25 1.244045E-05 0.000000E+00 1
|
||||
26 1.162146E-05 0.000000E+00 1
|
||||
27 1.837130E-05 0.000000E+00 1
|
||||
28 1.490705E-05 0.000000E+00 1
|
||||
29 1.218593E-05 0.000000E+00 1
|
||||
30 1.193124E-05 0.000000E+00 1
|
||||
31 1.169320E-05 0.000000E+00 1
|
||||
32 1.204609E-05 0.000000E+00 1
|
||||
33 9.622476E-06 0.000000E+00 1
|
||||
34 1.010416E-05 0.000000E+00 1
|
||||
35 9.946100E-06 0.000000E+00 1
|
||||
36 1.245206E-05 0.000000E+00 1
|
||||
37 2.004189E-05 0.000000E+00 2
|
||||
38 1.101200E-05 0.000000E+00 1
|
||||
39 2.013155E-05 0.000000E+00 1
|
||||
40 1.503905E-05 0.000000E+00 1
|
||||
41 1.695527E-05 0.000000E+00 1
|
||||
42 1.255106E-05 0.000000E+00 1
|
||||
43 1.885628E-05 0.000000E+00 1
|
||||
44 3.987430E-05 0.000000E+00 2
|
||||
45 9.909176E-06 0.000000E+00 1
|
||||
46 9.432880E-06 0.000000E+00 1
|
||||
47 1.103287E-05 0.000000E+00 1
|
||||
0 1.000000e+00 0.000000e+00 0
|
||||
1 6.627072e-06 0.000000e+00 7
|
||||
2 7.082738e-06 0.000000e+00 7
|
||||
3 7.410826e-06 0.000000e+00 7
|
||||
4 7.098589e-06 0.000000e+00 7
|
||||
5 9.176313e-06 0.000000e+00 8
|
||||
6 8.113908e-06 0.000000e+00 7
|
||||
7 8.211761e-06 0.000000e+00 7
|
||||
8 8.498843e-06 0.000000e+00 7
|
||||
9 6.849563e-06 0.000000e+00 6
|
||||
10 6.936034e-06 0.000000e+00 7
|
||||
11 5.811927e-06 0.000000e+00 7
|
||||
12 1.000000e+00 0.000000e+00 0
|
||||
13 5.764281e-06 0.000000e+00 7
|
||||
14 6.766942e-06 0.000000e+00 7
|
||||
15 7.071716e-06 0.000000e+00 8
|
||||
16 6.849563e-06 0.000000e+00 6
|
||||
17 6.731738e-06 0.000000e+00 7
|
||||
18 6.930737e-06 0.000000e+00 7
|
||||
19 1.000000e+00 0.000000e+00 0
|
||||
20 5.465672e-06 0.000000e+00 7
|
||||
21 5.837881e-06 0.000000e+00 7
|
||||
22 9.659138e-06 0.000000e+00 9
|
||||
23 1.072884e-05 0.000000e+00 9
|
||||
24 5.817766e-06 0.000000e+00 7
|
||||
25 6.345783e-06 0.000000e+00 7
|
||||
26 5.519944e-06 0.000000e+00 7
|
||||
27 6.197149e-06 0.000000e+00 7
|
||||
28 5.479187e-06 0.000000e+00 7
|
||||
29 6.309151e-06 0.000000e+00 7
|
||||
30 6.548096e-06 0.000000e+00 7
|
||||
31 6.844764e-06 0.000000e+00 7
|
||||
32 6.543819e-06 0.000000e+00 7
|
||||
33 6.372436e-06 0.000000e+00 7
|
||||
34 6.446273e-06 0.000000e+00 7
|
||||
35 6.440614e-06 0.000000e+00 7
|
||||
36 6.339770e-06 0.000000e+00 7
|
||||
37 1.000000e+00 0.000000e+00 0
|
||||
38 5.996007e-06 0.000000e+00 7
|
||||
39 6.386000e-06 0.000000e+00 7
|
||||
40 1.159312e-05 0.000000e+00 6
|
||||
41 1.165020e-05 0.000000e+00 6
|
||||
42 1.135841e-05 0.000000e+00 6
|
||||
43 1.045216e-05 0.000000e+00 6
|
||||
44 1.017683e-05 0.000000e+00 7
|
||||
45 8.145607e-06 0.000000e+00 7
|
||||
46 7.644835e-06 0.000000e+00 7
|
||||
47 9.592766e-06 0.000000e+00 7
|
||||
|
|
|
|||
Loading…
Reference in New Issue
Block a user