rtang/ANASEN_analysis
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

new file: MiscEdE.C plotting the radioactive beam IC PID as well as looking for timing coincideences between that and the MCP-Rf time

Browse Source 
new file:   run_misc_ede.sh same as above
	modified:   TrackRecon.C added beamE gated plots for the step aladder for 27Al
	modified:   run_17F.sh doing more files now
	modified:   run_27Al.sh scaling change to 0.9 to try it out
	modified:   run_tr.sh
	modified:   scratch/make_prettyplots.C add log scale option to the 2D plots
This commit is contained in:
Vignesh Sitaraman 2026-09-01 06:24:36 -04:00
parent 3dd678221a
commit dd6b938eb2
7 changed files with 491 additions and 25 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

359
MiscEdE.C Normal file
View File

@ -0,0 +1,359 @@
#include "TFile.h"
#include "TTree.h"
#include "TChain.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TCutG.h"
#include "TKey.h"
#include "TString.h"
#include "TSystem.h"
#include <map>
#include <vector>
#include <utility>
#include <cstdlib>
namespace
{
constexpr double kAdcMax = 16384; // 14-bit-ish ADC range for axes
struct HitData
{
double e;
unsigned long long t;
};
}
void MiscEdE(const char *f0,
const char *f1 = "", const char *f2 = "", const char *f3 = "",
const char *f4 = "", const char *f5 = "", const char *f6 = "",
const char *f7 = "", const char *f8 = "", const char *f9 = "")
{
std::vector<TString> inputs;
for (const char *f : {f0, f1, f2, f3, f4, f5, f6, f7, f8, f9})
if (f && f[0])
inputs.emplace_back(f);
if (inputs.empty())
{
printf("MiscEdE: no input files given.\n");
return;
}
TChain *tree = new TChain("tree");
for (auto &in : inputs)
{
if (gSystem->AccessPathName(in))
{
printf("MiscEdE: WARNING input not found, skipping: %s\n", in.Data());
continue;
}
tree->Add(in);
}
if (tree->GetNtrees() == 0)
{
printf("MiscEdE: no readable inputs.\n");
return;
}
// ---------------------------------------------------------
// Bash Environment Variable Extraction
// ---------------------------------------------------------
int runNumber = -1;
if (gSystem->Getenv("CURRENT_RUN")) {
runNumber = std::atoi(gSystem->Getenv("CURRENT_RUN"));
}
printf("MiscEdE: Run Number loaded from Bash: %d\n", runNumber);
TString cutDir = "~/ANASEN_analysis"; // Fallback if not run via bash
if (gSystem->Getenv("CUT_DIR")) {
cutDir = gSystem->Getenv("CUT_DIR");
}
// ---------------------------------------------------------
// Run-Dependent Channel Routing
// ---------------------------------------------------------
// Defaults (Run <= 180)
int snLolli = 405, chLolli = 9;
int snSi = 405, chSi = 11;
int snRF = 405, chRF = 15;
int snMCP = 405, chMCP = 14;
if (runNumber > 180 && runNumber <= 282) {
chRF = 10;
}
else if (runNumber > 282 && runNumber <= 322) {
snRF = 89; chRF = 0;
snMCP = 89; chMCP = 1;
}
else if (runNumber > 322) {
snRF = 89; chRF = 0;
snMCP = 89; chMCP = 2;
}
// Composite Unique IDs (UID = Board * 100 + Channel) to prevent collisions
int uidLollipop = snLolli * 100 + chLolli;
int uidSiMon = snSi * 100 + chSi;
int uidRF = snRF * 100 + chRF;
int uidMCP = snMCP * 100 + chMCP;
printf(" -> Routed Lollipop : Board %d, Ch %d (UID: %d)\n", snLolli, chLolli, uidLollipop);
printf(" -> Routed Si Monitor: Board %d, Ch %d (UID: %d)\n", snSi, chSi, uidSiMon);
printf(" -> Routed RF : Board %d, Ch %d (UID: %d)\n", snRF, chRF, uidRF);
printf(" -> Routed MCP : Board %d, Ch %d (UID: %d)\n", snMCP, chMCP, uidMCP);
// ---------------------------------------------------------
// Exact schema matching EventBuilder.cpp
// ---------------------------------------------------------
const int MAX_MULTI = 2000;
unsigned int multi = 0;
unsigned short sn[MAX_MULTI];
unsigned short ch[MAX_MULTI];
unsigned short e[MAX_MULTI];
unsigned long long e_t[MAX_MULTI]; // Timestamp branch is "e_t"
tree->SetBranchAddress("multi", &multi);
tree->SetBranchAddress("sn", sn);
tree->SetBranchAddress("ch", ch);
tree->SetBranchAddress("e", e);
tree->SetBranchAddress("e_t", e_t);
TString stem = gSystem->BaseName(inputs[0].Data());
if (stem.EndsWith(".root"))
stem.Remove(stem.Length() - 5);
TString outName = "MiscEdE_" + stem + ".root";
TFile *out = new TFile(outName, "recreate");
// ---------------------------------------------------------
// Load the 17F and 16O TCutGs from file
// ---------------------------------------------------------
TCutG *cut17F = nullptr;
TString path17F = gSystem->ExpandPathName(Form("%s/17FCut.root", cutDir.Data()));
if (!gSystem->AccessPathName(path17F))
{
TFile *fCut = TFile::Open(path17F);
if (fCut && !fCut->IsZombie())
{
for (auto keyObj : *fCut->GetListOfKeys())
{
TKey *key = (TKey *)keyObj;
if (TString(key->GetClassName()) == "TCutG")
{
cut17F = (TCutG *)key->ReadObj()->Clone("cut17F");
break;
}
}
fCut->Close();
}
}
TCutG *cut16O = nullptr;
TString path16O = gSystem->ExpandPathName(Form("%s/16OCut.root", cutDir.Data()));
if (!gSystem->AccessPathName(path16O))
{
TFile *fCut = TFile::Open(path16O);
if (fCut && !fCut->IsZombie())
{
for (auto keyObj : *fCut->GetListOfKeys())
{
TKey *key = (TKey *)keyObj;
if (TString(key->GetClassName()) == "TCutG")
{
cut16O = (TCutG *)key->ReadObj()->Clone("cut16O");
break;
}
}
fCut->Close();
}
}
out->cd();
// ---------------------------------------------------------
// Prepare Histograms
// ---------------------------------------------------------
std::map<int, TH2D*> hEvsChMap;
auto getEvsCh = [&](int boardSN) -> TH2D* {
if (hEvsChMap.count(boardSN)) return hEvsChMap[boardSN];
TH2D *h = new TH2D(Form("h2_E_vs_ch_bd%d", boardSN), Form("Board %d: E vs channel;board channel;E [ADC]", boardSN), 16, -0.5, 15.5, 800, 0, kAdcMax);
hEvsChMap[boardSN] = h;
return h;
};
TH2D *h2_EvsdT_LolliSi = new TH2D("EvsdT_Lollipop_Si", "Energy Si vs dT (Lollipop - Si);dT [ticks];Si E [ADC]", 500, 0, 2000, 400, 0, kAdcMax);
TH2D *h2_TOF_SiRF = new TH2D("TOF_Si_RF", "Time of Flight: Energy Si vs (T_Si - T_RF);T_Si - T_RF [ticks];Si E [ADC]", 500, 0, 2000, 400, 0, kAdcMax);
// Dynamic Pairwise E-dE mapped by UID
std::map<std::pair<int, int>, TH2D *> hEdE;
auto getEdE = [&](int uidA, int uidB) -> TH2D * {
if (uidA > uidB) std::swap(uidA, uidB);
auto key = std::make_pair(uidA, uidB);
auto it = hEdE.find(key);
if (it != hEdE.end()) return it->second;
TString name = Form("h2_ede_bd%d_ch%d_vs_bd%d_ch%d", uidA/100, uidA%100, uidB/100, uidB%100);
TString title = Form("E-dE: Bd %d Ch %d (x) vs Bd %d Ch %d (y);X E [ADC];Y E [ADC]", uidA/100, uidA%100, uidB/100, uidB%100);
TH2D *h = new TH2D(name, title, 400, 0, kAdcMax, 400, 0, kAdcMax);
hEdE[key] = h;
return h;
};
// 17F Gated Plots
TH2D *hEdE_17F = nullptr;
TH2D *hEdE_inv = nullptr;
TH2D *h2_EvsdT_LolliSi_17F = nullptr;
TH2D *h2_TOF_SiRF_17F = nullptr;
TH1D *h1_dT_17F = nullptr;
if (cut17F)
{
hEdE_17F = new TH2D("EdE_lollipopIC_vs_SiMonitor_17FCut", "lollipop IC vs Si monitor (17F Gated);Si monitor E [ADC];lollipop IC E [ADC]", 400, 0, kAdcMax, 400, 0, kAdcMax);
hEdE_inv = new TH2D("EdE_lollipopIC_vs_SiMonitor_17FCut_Inverse", "lollipop IC vs Si monitor (Inverse 17F Gated);Si monitor E [ADC];lollipop IC E [ADC]", 400, 0, kAdcMax, 400, 0, kAdcMax);
h2_EvsdT_LolliSi_17F = new TH2D("EvsdT_Lollipop_Si_17FCut", "Energy Si vs dT (17F Gated);dT (Lollipop - Si) [ticks];Si E [ADC]", 500, 0, 2000, 400, 0, kAdcMax);
h2_TOF_SiRF_17F = new TH2D("TOF_Si_RF_17FCut", "Si TOF vs Energy (17F Gated);T_Si - T_RF [ticks];Si E [ADC]", 500, 0, 2000, 400, 0, kAdcMax);
h1_dT_17F = new TH1D("dT_Lollipop_Si_17FCut", "dT Lollipop-Si (17F Gated);T_Lollipop - T_Si [ticks]", 500, 0, 2000);
}
// 16O Gated Plots
TH2D *hEdE_16O = nullptr;
TH2D *h2_EvsdT_LolliSi_16O = nullptr;
TH2D *h2_TOF_SiRF_16O = nullptr;
TH1D *h1_dT_16O = nullptr;
if (cut16O)
{
hEdE_16O = new TH2D("EdE_lollipopIC_vs_SiMonitor_16OCut", "lollipop IC vs Si monitor (16O Gated);Si monitor E [ADC];lollipop IC E [ADC]", 400, 0, kAdcMax, 400, 0, kAdcMax);
h2_EvsdT_LolliSi_16O = new TH2D("EvsdT_Lollipop_Si_16OCut", "Energy Si vs dT (16O Gated);dT (Lollipop - Si) [ticks];Si E [ADC]", 500, 0, 2000, 400, 0, kAdcMax);
h2_TOF_SiRF_16O = new TH2D("TOF_Si_RF_16OCut", "Si TOF vs Energy (16O Gated);T_Si - T_RF [ticks];Si E [ADC]", 500, 0, 2000, 800, 0, kAdcMax);
h1_dT_16O = new TH1D("dT_Lollipop_Si_16OCut", "dT Lollipop-Si (16O Gated);T_Lollipop - T_Si [ticks]", 500, 0, 2000);
}
Long64_t nEnt = tree->GetEntries();
printf("MiscEdE: %lld entries across %d file(s) -> %s\n", nEnt, tree->GetNtrees(), outName.Data());
Long64_t used = 0;
for (Long64_t i = 0; i < nEnt; ++i)
{
tree->GetEntry(i);
// We map hits by UID instead of standard channel to prevent overlaps between Bd405 and Bd89
std::map<int, HitData> hit;
for (unsigned int j = 0; j < multi && j < (unsigned)MAX_MULTI; ++j)
{
// Only care about MISC boards 405 and 89
if (sn[j] != 405 && sn[j] != 89)
continue;
int c = ch[j];
int uid = (sn[j] * 100) + ch[j];
double en = e[j];
unsigned long long ts = e_t[j];
auto it = hit.find(uid);
if (it == hit.end() || en > it->second.e)
hit[uid] = {en, ts};
getEvsCh(sn[j])->Fill(c, en);
}
if (hit.empty())
continue;
++used;
// ---------------------------------------------------------
// Generic Dynamic Pairwise E-dE Logic across all UIDs
// ---------------------------------------------------------
for (auto a = hit.begin(); a != hit.end(); ++a)
{
for (auto b = std::next(a); b != hit.end(); ++b)
{
if (a->second.e > 500 && b->second.e > 500)
{
getEdE(a->first, b->first)->Fill(a->second.e, b->second.e);
}
}
}
// ---------------------------------------------------------
// Core Timing and PID Gating Logic (Lollipop-Si specific)
// ---------------------------------------------------------
if (hit.count(uidLollipop) && hit.count(uidSiMon))
{
double e_lol = hit[uidLollipop].e;
double e_si = hit[uidSiMon].e;
// Calculate delta T safely across potential CAEN rollovers
double dt_lol_si = static_cast<double>(hit[uidLollipop].t) - static_cast<double>(hit[uidSiMon].t);
if (e_lol > 500 && e_si > 500)
{
h2_EvsdT_LolliSi->Fill(dt_lol_si, e_si);
bool is17F = cut17F && cut17F->IsInside(e_si, e_lol);
bool is16O = cut16O && cut16O->IsInside(e_si, e_lol);
if (is17F)
{
hEdE_17F->Fill(e_si, e_lol);
h2_EvsdT_LolliSi_17F->Fill(dt_lol_si, e_si);
h1_dT_17F->Fill(dt_lol_si);
}
else if (cut17F)
{
hEdE_inv->Fill(e_si, e_lol);
}
if (is16O)
{
hEdE_16O->Fill(e_si, e_lol);
h2_EvsdT_LolliSi_16O->Fill(dt_lol_si, e_si);
h1_dT_16O->Fill(dt_lol_si);
}
// Generate absolute TOF if the RF channel fired
if (hit.count(uidRF))
{
double tof = static_cast<double>(hit[uidSiMon].t) - static_cast<double>(hit[uidRF].t);
h2_TOF_SiRF->Fill(tof, e_si);
if (is17F)
h2_TOF_SiRF_17F->Fill(tof, e_si);
if (is16O)
h2_TOF_SiRF_16O->Fill(tof, e_si);
}
}
}
}
// ---------------------------------------------------------
// Purge Sparse E-dE Histograms (< 500 counts)
// ---------------------------------------------------------
for (auto it = hEdE.begin(); it != hEdE.end();)
{
if (it->second->GetEntries() < 500)
{
delete it->second; // Free memory and remove from ROOT file
it = hEdE.erase(it);
}
else
{
++it;
}
}
// Convenience alias: Transpose generic plot to correctly place stopping E on X
if (hEdE.count({uidLollipop, uidSiMon}))
{
TH2D *src = hEdE[{uidLollipop, uidSiMon}];
TH2D *alias = new TH2D("EdE_lollipopIC_vs_SiMonitor_Ungated",
Form("lollipop IC vs Si monitor (Ungated);Bd %d Ch %d (Si) E [ADC];Bd %d Ch %d (Lolli) E [ADC]", snSi, chSi, snLolli, chLolli),
src->GetNbinsY(), src->GetYaxis()->GetXmin(), src->GetYaxis()->GetXmax(),
src->GetNbinsX(), src->GetXaxis()->GetXmin(), src->GetXaxis()->GetXmax());
for (int ix = 1; ix <= src->GetNbinsX(); ++ix)
for (int iy = 1; iy <= src->GetNbinsY(); ++iy)
alias->SetBinContent(iy, ix, src->GetBinContent(ix, iy)); // (x,y) -> (y,x)
}
out->Write();
out->Close();
printf("MiscEdE: Processed %lld valid MISC events. Wrote %s\n", used, outName.Data());
}

29
TrackRecon.C
View File

@ -3951,6 +3951,7 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
double ebeam_kin_3774keV = invertBeamEnergyMeV(m_beam, mass_4He, m3, m4, Efix, theta * 180 / M_PI, 3.774);
double ebeam_kin_4809keV = invertBeamEnergyMeV(m_beam, mass_4He, m3, m4, Efix, theta * 180 / M_PI, 4.809);
double ebeam_kin_5614keV = invertBeamEnergyMeV(m_beam, mass_4He, m3, m4, Efix, theta * 180 / M_PI, 5.164);
double ebeam_kin_6550keV = invertBeamEnergyMeV(m_beam, mass_4He, m3, m4, Efix, theta * 180 / M_PI, 6.550);
// Gated output: only fill when this hypothesis (proton "_p" / alpha "_a") agrees
// with which side of the proton_locus gate the event fell on, so each event
@ -3962,28 +3963,38 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
auto plot_with_tag = [&](const std::string &topo)
{
std::string t = topo.empty() ? "" : ("_" + topo);
plotter->Fill1D(rx + "_Ex_from" + ejtag + t + sfx, 600, -15, 15, Ex, pmlabel);
plotter->Fill1D(rx + "_Ex_from" + ejtag + t + sfx, 600, -10, 20, Ex, pmlabel);
plotter->Fill2D(rx + "_VertexReconZ_vs_Ef" + ejtag + t + sfx, 800, -400, 400, 800, 0, ef_max, vertex_z, Efix, pmlabel);
plotter->Fill2D(rx + "_VertexReconZ_vs_Ex" + ejtag + t + sfx, 800, -400, 400, 800, -20, 20, vertex_z, Ex, pmlabel);
plotter->Fill2D(rx + "_Ex_vs_theta" + ejtag + t + sfx, 360, 0, 180, 800, -20, 20, theta * 180 / M_PI, Ex, pmlabel);
plotter->Fill2D(rx + "_VertexReconZ_vs_Ex" + ejtag + t + sfx, 800, -400, 400, 600, -10, 20, vertex_z, Ex, pmlabel);
plotter->Fill2D(rx + "_Ex_vs_theta" + ejtag + t + sfx, 360, 0, 180, 600, -10, 20, theta * 180 / M_PI, Ex, pmlabel);
if (ebeam_kin_MeV > 0.0)
plotter->Fill2D(rx + "_BeamEnergy_ETrack_vs_EKin" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_MeV, pmlabel);
if (ejtag == "_p")
{
if (beam_energy_at_vertex < 4.0 && beam_energy_at_vertex > 0.1)
plotter->Fill2D(rx + "_ETrack_vs_EKinGS" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_GS, "ETrackvsKin_assumed");
else if (beam_energy_at_vertex <= 12.0)
plotter->Fill2D(rx + "_ETrack_vs_EKin2235keV" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_2235keV, "ETrackvsKin_assumed");
else if (beam_energy_at_vertex < 24.0)
{
plotter->Fill2D(rx + "_ETrack_vs_EKin3498kev" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_3498keV, "ETrackvsKin_assumed");
plotter->Fill2D(rx + "_ETrack_vs_EKin3774keV" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_3774keV, "ETrackvsKin_assumed");
}
else if (beam_energy_at_vertex >= 36.00)
plotter->Fill2D(rx + "_ETrack_vs_EKin4809keV" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_4809keV, "ETrackvsKin_assumed");
else if (beam_energy_at_vertex > 42.0)
plotter->Fill2D(rx + "_ETrack_vs_EKin5614keV" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_5614keV, "ETrackvsKin_assumed");
else
plotter->Fill2D(rx + "_ETrack_vs_EKin6550keV" + ejtag + t + sfx, 400, 0, beamE0 * 1.5, 400, 0, beamE0 * 1.5,
beam_energy_at_vertex, ebeam_kin_6550keV, "ETrackvsKin_assumed");
}
};
@ -3999,7 +4010,7 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
}
plotter->Fill1D(rx + "_pczfix" + sfx, 600, -300, 300, pcz_fix, pmlabel);
plotter->Fill2D(rx + "_Ef_vs_theta" + ejtag + sfx, 100, 0, 180, 800, 0, ef_max, theta * 180 / M_PI, Efix, pmlabel);
plotter->Fill2D(rx + "_Ex_vs_phi" + ejtag + sfx, 180, -180, 180, 800, -20, 20, phi * 180 / M_PI, Ex, pmlabel);
plotter->Fill2D(rx + "_Ex_vs_phi" + ejtag + sfx, 180, -180, 180, 600, -10, 20, phi * 180 / M_PI, Ex, pmlabel);
if (dt_rf_mcp > -900000000)
plotter->Fill2D(rx + "_Ex_vs_TOF_rf_mcp" + ejtag + sfx, 500, -1000, 1000, 800, -20, 20, dt_rf_mcp, Ex, pmlabel);
plotter->Fill1D(rx + "_VertexReconZ" + sfx, 800, -400, 400, vertex_z, pmlabel);
@ -4029,12 +4040,10 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
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, -200, 200, 800, 0, 0.6, phi * 180 / M_PI, anodeE_MeV, pmlabel);
plotter->Fill2D(rx + "_dEgasCalib_vs_E" + ejtag + sfx + "_E<10MeV" + std::to_string(beam_energy_at_vertex < 10), 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, 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 + "_E<10MeV" + std::to_string(beam_energy_at_vertex < 10), 800, -10, 10, 800, 0, 0.6, Ex, anodeE_MeV, pmlabel);
plotter->Fill2D(rx + "_dEgasCalib_vs_Z" + ejtag + sfx + "_E<10MeV" + std::to_string(beam_energy_at_vertex < 10), 800, -400, 400, 800, 0, 0.6, vertex_z, 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);
}
}

3
run_17F.sh
View File

@ -48,7 +48,8 @@ root -q -l -b -e '.L TrackRecon.C++O'
# 3% CO2
# parallel --bar -j 6 run_once ::: {325..400}
parallel --bar -j 7 run_once ::: 351 353 355 358 359 360 362 367
parallel --bar -j 11 run_once ::: {351..400}
# parallel --bar -j 7 run_once ::: 351 353 355 358 359 360 362 367
hadd -j 4 -k ${OUT_DIR}/Output_17F.root ${OUT_DIR}/results_run*.root
unset source_vertex

2
run_27Al.sh
View File

@ -8,7 +8,7 @@ export CO2percent=3
export pressure_in_torr=250
export CATHODE_GAIN=3.0
export source_vertex=-200.0
export DEDX_SCALE=0.89
export DEDX_SCALE=0.9
export CUTLIST=cuts_list.txt
export BEAM_AXIS_X=0
export BEAM_AXIS_Y=0

93
run_misc_ede.sh Normal file
View File

@ -0,0 +1,93 @@
#!/bin/bash
# run_misc_ede.sh -- run the standalone MISC E-dE telescope plots (MiscEdE.C)
# over raw, pre-mapping event-built files. Independent of the main TrackRecon
# pipeline: no mapping/calibration/reconstruction, works on beam-monitor runs.
#
# The macro reads the flat sn/ch/e event-built tree and produces per-channel
# spectra + every pairwise E-dE 2D, so you can see all correlations and confirm
# telescopes directly from the data.
#
# Point EVT_DIR / PREFIX / SUFFIX at wherever your event-built files live and
# how they're named, then list the run numbers. Each run is processed into its
# own MiscEdE_<stem>.root; set MERGE=1 to also hadd them together.
#
# Usage:
# ./run_misc_ede.sh # uses the RUNS list below
# ./run_misc_ede.sh 351 353 355 # override runs on the command line
set -u
# ---- configure these to your layout ----
EVT_DIR="${EVT_DIR:-../ANASEN_analysis/data/17F_Data}" # where the event-built files are
PREFIX="${PREFIX:-Run_}" # filename prefix before the run number
SUFFIX="${SUFFIX:-_2000.root}" # filename suffix after the run number
PAD="${PAD:-3}" # zero-pad width for the run number (Run_053 -> 3)
OUT_DIR="${OUT_DIR:-Output_misc_ede}"
MERGE="${MERGE:-1}" # 1 = also hadd all per-run outputs
MACRO="${MACRO:-MiscEdE.C}"
CUT_DIR="${CUT_DIR:-/home/vsitaraman/ANASEN_analysis}" # Where 17FCut.root and 16OCut.root live
# ----------------------------------------
# Runs: command-line args win, else the lollipop/Si-monitor test runs below.
if [ "$#" -gt 0 ]; then
RUNS=("$@")
else
RUNS=( 52 53 54 55 63 65 68 69 72 74 75 76 78 92 168 170 226 262 324 327 337 334 369 370 371 372 373 374)
fi
mkdir -p "$OUT_DIR"
OUT_ABS="$(cd "$OUT_DIR" && pwd)"
MACRO_ABS="$(pwd)/$MACRO"
export OUT_ABS MACRO_ABS EVT_DIR PREFIX SUFFIX PAD CUT_DIR
built_path() {
local run; run=$(printf "%0${PAD}d" "$1")
echo "${EVT_DIR}/${PREFIX}${run}${SUFFIX}"
}
export -f built_path
# One run -> one MiscEdE_<stem>.root, runnable under GNU parallel.
misc_ede_one() {
local r="$1"
local infile; infile=$(built_path "$r")
if [ ! -f "$infile" ]; then
echo "SKIP: $infile not found"
return
fi
local infile_abs; infile_abs="$(cd "$(dirname "$infile")" && pwd)/$(basename "$infile")"
echo "=== MISC E-dE: run $r ($infile_abs) ==="
# EXPORT METADATA FOR C++ SCRIPT
export CURRENT_RUN="$r"
# Macro names its own output MiscEdE_<stem>.root in the CWD, so run inside
# OUT_DIR to keep outputs together. Absolute paths avoid relative-path
# surprises from the cd.
( cd "$OUT_ABS" && root -l -b -q "${MACRO_ABS}(\"${infile_abs}\")" )
}
export -f misc_ede_one
rm -f "$OUT_DIR"/*.root
JOBS="${JOBS:-7}" # match the -j 7 you use for run_once
if command -v parallel >/dev/null 2>&1; then
parallel --bar -j "$JOBS" misc_ede_one ::: "${RUNS[@]}"
else
echo "GNU parallel not found; running serially."
for r in "${RUNS[@]}"; do misc_ede_one "$r"; done
fi
if [ "$MERGE" -eq 1 ]; then
shopt -s nullglob
parts=()
for f in "$OUT_DIR"/MiscEdE_*.root; do
[ "$(basename "$f")" = "MiscEdE_ALL.root" ] && continue
parts+=("$f")
done
if [ "${#parts[@]}" -gt 1 ]; then
echo "=== merging ${#parts[@]} outputs -> ${OUT_DIR}/MiscEdE_ALL.root ==="
hadd -f "${OUT_DIR}/MiscEdE_ALL.root" "${parts[@]}"
fi
fi
echo "Done. Outputs in ${OUT_DIR}/"

3
run_tr.sh
View File

@ -55,7 +55,7 @@ if [[ 1 -eq 0 ]]; then
fi
# --- Block 3: 27Al Alpha+Gas Runs (9, 12) ---
if [[ 1 -eq 0 ]]; then
if [[ 1 -eq 1 ]]; then
export DATASET="27Al"
export PREFIX="Run_"
export OUT_DIR="Output_a"
@ -70,6 +70,7 @@ if [[ 1 -eq 0 ]]; then
unset CATHODE_GAIN
unset timecut_low
unset pressure_in_torr
exit
fi
# --- Block 4: 17F Alpha+Gas Runs (18-21) ---

3
scratch/make_prettyplots.C
View File

@ -153,6 +153,9 @@ void make_prettyplots(const char *rootFile,
h->GetYaxis()->SetTitleOffset(1.4);
h->GetXaxis()->CenterTitle(true);
h->GetYaxis()->CenterTitle(true);
// gPad->SetLogz();
h->Draw("colz");
}
else