modified: TrackRecon.C put in timing gates and phi gates uniformly
modified: .gitignore deleted: Analysis.C deleted: Analyzer.C deleted the older analysis formats deleted: Analyzer.h deleted: Timing_Summary_Matplotlib.png deleted: WedgeChVAnode.jpeg modified: run_tr.sh modified: scratch/make_prettyplots.C added rebinning capability and log z scale option
This commit is contained in:
parent
ec43a81030
commit
38c8c7e0e0
3
.gitignore
vendored
3
.gitignore
vendored
|
|
@ -36,4 +36,5 @@ Armory/EventBuilder
|
||||||
EventBuilder
|
EventBuilder
|
||||||
anasen_fem/anode_to_cathode_1d_17.43.csv
|
anasen_fem/anode_to_cathode_1d_17.43.csv
|
||||||
pc_calib_raw/
|
pc_calib_raw/
|
||||||
AutoFit.C
|
AutoFit.C
|
||||||
|
scratch/images/WedgeChVAnode.jpeg
|
||||||
|
|
|
||||||
20
Analysis.C
20
Analysis.C
|
|
@ -1,20 +0,0 @@
|
||||||
#include "TChain.h"
|
|
||||||
#include <iostream>
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
void Analysis(int start, int end) {
|
|
||||||
// Create a TChain
|
|
||||||
TChain *chain = new TChain("tree");
|
|
||||||
for(int i = start; i < end+1; i++) {
|
|
||||||
chain->Add(Form("data/root_data/Run_%03d_mapped.root", i));
|
|
||||||
}
|
|
||||||
|
|
||||||
// Process the chain using Analyzer.C+
|
|
||||||
chain->Process("Analyzer.C+");
|
|
||||||
}
|
|
||||||
|
|
||||||
// Define a macro with the same name as the script
|
|
||||||
void Analysis() {
|
|
||||||
Analysis(72, 194); // Adjust the range if needed
|
|
||||||
}
|
|
||||||
877
Analyzer.C
877
Analyzer.C
|
|
@ -1,877 +0,0 @@
|
||||||
#define Analyzer_cxx
|
|
||||||
|
|
||||||
#include "Analyzer.h"
|
|
||||||
#include "Armory/ClassSX3.h"
|
|
||||||
#include "Armory/ClassPW.h"
|
|
||||||
|
|
||||||
#include <TH2.h>
|
|
||||||
#include <TStyle.h>
|
|
||||||
#include <TCanvas.h>
|
|
||||||
#include <TMath.h>
|
|
||||||
#include "TVector3.h"
|
|
||||||
|
|
||||||
#include <fstream>
|
|
||||||
#include <iostream>
|
|
||||||
#include <sstream>
|
|
||||||
#include <map>
|
|
||||||
#include <utility>
|
|
||||||
#include <algorithm>
|
|
||||||
|
|
||||||
TH2F *hsx3IndexVE;
|
|
||||||
TH2F *hqqqIndexVE;
|
|
||||||
TH2F *hpcIndexVE;
|
|
||||||
|
|
||||||
TH2F *hpcIndexVE_GM;
|
|
||||||
TH2F *hsx3Coin;
|
|
||||||
TH2F *hqqqCoin;
|
|
||||||
TH2F *hpcCoin;
|
|
||||||
TH2F *hAVCcoin;
|
|
||||||
|
|
||||||
TH2F *hqqqPolar;
|
|
||||||
TH2F *hsx3VpcIndex;
|
|
||||||
TH2F *hqqqVpcIndex;
|
|
||||||
TH2F *hqqqVpcE;
|
|
||||||
TH2F *hsx3VpcE;
|
|
||||||
TH2F *hanVScatsum;
|
|
||||||
TH2F *hanVScatsum_a[24];
|
|
||||||
TH1F *hPC_E[48];
|
|
||||||
TH1F *hCat4An;
|
|
||||||
TH1F *hCat0An;
|
|
||||||
TH1F *hAnodehits;
|
|
||||||
TH2F *hNosvAe;
|
|
||||||
|
|
||||||
int padID = 0;
|
|
||||||
|
|
||||||
SX3 sx3_contr;
|
|
||||||
PW pw_contr;
|
|
||||||
PW pwinstance;
|
|
||||||
TVector3 hitPos;
|
|
||||||
// TVector3 anodeIntersection;
|
|
||||||
std::map<int, std::pair<double, double>> slopeInterceptMap;
|
|
||||||
// SX3 Calibration Arrays
|
|
||||||
const int MAX_DET = 24;
|
|
||||||
const int MAX_UP = 4;
|
|
||||||
const int MAX_DOWN = 4;
|
|
||||||
const int MAX_BK = 4;
|
|
||||||
double backGain[MAX_DET][MAX_BK] = {{0}};
|
|
||||||
bool backGainValid[MAX_DET][MAX_BK] = {{false}};
|
|
||||||
double frontGain[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{0}}}};
|
|
||||||
bool frontGainValid[MAX_DET][MAX_BK][MAX_UP][MAX_DOWN] = {{{{false}}}};
|
|
||||||
|
|
||||||
// QQQ Calibration Arrays
|
|
||||||
const int MAX_QQQ = 4;
|
|
||||||
const int MAX_RING = 16;
|
|
||||||
const int MAX_WEDGE = 16;
|
|
||||||
double qqqGain[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{0}}};
|
|
||||||
bool qqqGainValid[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{false}}};
|
|
||||||
double qqqCalib[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{0}}};
|
|
||||||
bool qqqCalibValid[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{false}}};
|
|
||||||
|
|
||||||
bool HitNonZero;
|
|
||||||
bool sx3ecut;
|
|
||||||
bool qqqEcut;
|
|
||||||
|
|
||||||
TH1F *hZProj;
|
|
||||||
TH1F *hPCZProj;
|
|
||||||
|
|
||||||
void Analyzer::Begin(TTree * /*tree*/)
|
|
||||||
{
|
|
||||||
TString option = GetOption();
|
|
||||||
|
|
||||||
hsx3IndexVE = new TH2F("hsx3IndexVE", "SX3 index vs Energy; sx3 index ; Energy", 24 * 12, 0, 24 * 12, 400, 0, 5000);
|
|
||||||
hqqqIndexVE = new TH2F("hqqqIndexVE", "QQQ index vs Energy; QQQ index ; Energy", 4 * 2 * 16, 0, 4 * 2 * 16, 400, 0, 5000);
|
|
||||||
hpcIndexVE = new TH2F("hpcIndexVE", "PC index vs Energy; PC index ; Energy", 2 * 24, 0, 2 * 24, 400, 0, 16000);
|
|
||||||
hpcIndexVE_GM = new TH2F("hpcIndexVE_GM", "PC index vs Energy; PC index ; Energy", 2 * 24, 0, 2 * 24, 400, 0, 16000);
|
|
||||||
|
|
||||||
hsx3Coin = new TH2F("hsx3Coin", "SX3 Coincident", 24 * 12, 0, 24 * 12, 24 * 12, 0, 24 * 12);
|
|
||||||
hqqqCoin = new TH2F("hqqqCoin", "QQQ Coincident", 4 * 2 * 16, 0, 4 * 2 * 16, 4 * 2 * 16, 0, 4 * 2 * 16);
|
|
||||||
hpcCoin = new TH2F("hpcCoin", "PC Coincident", 2 * 24, 0, 2 * 24, 2 * 24, 0, 2 * 24);
|
|
||||||
hAVCcoin = new TH2F("hAVCcoin", "Anode vs Cathode Coincident", 24, 0, 24, 24, 0, 24);
|
|
||||||
|
|
||||||
hqqqPolar = new TH2F("hqqqPolar", "QQQ Polar ID", 16 * 4, -TMath::Pi(), TMath::Pi(), 16, 10, 50);
|
|
||||||
|
|
||||||
hsx3VpcIndex = new TH2F("hsx3Vpcindex", "sx3 vs pc; sx3 index; pc index", 24 * 12, 0, 24 * 12, 48, 0, 48);
|
|
||||||
hqqqVpcIndex = new TH2F("hqqqVpcindex", "qqq vs pc; qqq index; pc index", 4 * 2 * 16, 0, 4 * 2 * 16, 48, 0, 48);
|
|
||||||
|
|
||||||
hqqqVpcE = new TH2F("hqqqVpcEnergy", "qqq vs pc; qqq energy; pc energy", 400, 0, 5000, 800, 0, 16000);
|
|
||||||
|
|
||||||
hsx3VpcE = new TH2F("hsx3VpcEnergy", "sx3 vs pc; sx3 energy; pc energy", 400, 0, 5000, 800, 0, 16000);
|
|
||||||
|
|
||||||
hZProj = new TH1F("hZProj", "Z Projection", 1200, -600, 600);
|
|
||||||
hPCZProj = new TH1F("hPCZProj", "PC Z Projection", 600, -300, 300);
|
|
||||||
|
|
||||||
hanVScatsum = new TH2F("hanVScatsum", "Anode vs Cathode Sum; Anode E; Cathode E", 400, 0, 16000, 400, 0, 20000);
|
|
||||||
hCat4An = new TH1F("hCat4An", "Number of Cathodes/Anode", 24, 0, 24);
|
|
||||||
hCat0An = new TH1F("hCat0An", "Number of Cathodes without Anode", 24, 0, 24);
|
|
||||||
hAnodehits = new TH1F("hAnodehits", "Number of Anode hits", 24, 0, 24);
|
|
||||||
hNosvAe = new TH2F("hnosvAe", "Number of Cathodes/Anode vs Anode Energy", 20, 0, 20, 400, 0, 16000);
|
|
||||||
// for (int i = 0; i < 24; i++)
|
|
||||||
// {
|
|
||||||
// TString histName = Form("hAnodeVsCathode_%d", i);
|
|
||||||
// TString histTitle = Form("Anode %d vs Cathode Sum; Anode E; Cathode Sum E", i);
|
|
||||||
// hanVScatsum_a[i] = new TH2F(histName, histTitle, 400, 0, 16000, 400, 0, 20000);
|
|
||||||
// }
|
|
||||||
// for (int i = 0; i < 48; i++)
|
|
||||||
// {
|
|
||||||
// TString histName = Form("hCathode_%d", i);
|
|
||||||
// TString histTitle = Form("Cathode_E_%d;", i);
|
|
||||||
// hPC_E[i] = new TH1F(histName, histTitle, 3200, 0, 32000);
|
|
||||||
// }
|
|
||||||
sx3_contr.ConstructGeo();
|
|
||||||
pw_contr.ConstructGeo();
|
|
||||||
|
|
||||||
std::ifstream inputFile("slope_intercept_results.txt");
|
|
||||||
|
|
||||||
if (inputFile.is_open())
|
|
||||||
{
|
|
||||||
std::string line;
|
|
||||||
int index;
|
|
||||||
double slope, intercept;
|
|
||||||
while (std::getline(inputFile, line))
|
|
||||||
{
|
|
||||||
std::stringstream ss(line);
|
|
||||||
ss >> index >> slope >> intercept;
|
|
||||||
// wires 37, 39, 44 have fit data that is incorrect or not present, they have thus been set to 1,0 (slope, intercept) for convenience
|
|
||||||
// wire 19 the 4th point was genereated using the slope of the line produced uising the other 3 points from the wire 1 vs wire 19 plot
|
|
||||||
if (index >= 0 && index <= 47)
|
|
||||||
{
|
|
||||||
slopeInterceptMap[index] = std::make_pair(slope, intercept);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
inputFile.close();
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
std::cerr << "Error opening slope_intercept.txt" << std::endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
std::string filename = "sx3_GainMatchback.txt";
|
|
||||||
|
|
||||||
std::ifstream infile(filename);
|
|
||||||
if (!infile.is_open())
|
|
||||||
{
|
|
||||||
std::cerr << "Error opening " << filename << "!" << std::endl;
|
|
||||||
return;
|
|
||||||
}
|
|
||||||
|
|
||||||
int id, bk;
|
|
||||||
double gain;
|
|
||||||
while (infile >> id >> bk >> gain)
|
|
||||||
{
|
|
||||||
backGain[id][bk] = gain;
|
|
||||||
if (backGain[id][bk] > 0)
|
|
||||||
backGainValid[id][bk] = true;
|
|
||||||
else
|
|
||||||
backGainValid[id][bk] = false;
|
|
||||||
}
|
|
||||||
|
|
||||||
infile.close();
|
|
||||||
std::cout << "Loaded back gains from " << filename << std::endl;
|
|
||||||
|
|
||||||
std::string filename1 = "sx3_GainMatchfront.txt";
|
|
||||||
|
|
||||||
std::ifstream infile1(filename1);
|
|
||||||
if (!infile1.is_open())
|
|
||||||
{
|
|
||||||
std::cerr << "Error opening " << filename1 << "!" << std::endl;
|
|
||||||
return;
|
|
||||||
}
|
|
||||||
|
|
||||||
int idf, bkf, uf, df;
|
|
||||||
double fgain;
|
|
||||||
while (infile1 >> idf >> bkf >> uf >> df >> fgain)
|
|
||||||
{
|
|
||||||
frontGain[idf][bkf][uf][df] = fgain;
|
|
||||||
frontGainValid[idf][bkf][uf][df] = true;
|
|
||||||
}
|
|
||||||
|
|
||||||
// QQQ Gain Matching and Calibration
|
|
||||||
// ----------------------- Load QQQ Gains
|
|
||||||
{
|
|
||||||
std::string filename = "qqq_GainMatch.txt";
|
|
||||||
std::ifstream infile(filename);
|
|
||||||
if (!infile.is_open())
|
|
||||||
{
|
|
||||||
std::cerr << "Error opening " << filename << "!" << std::endl;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
int det, ring, wedge;
|
|
||||||
double gainw, gainr;
|
|
||||||
while (infile >> det >> ring >> wedge >> gainw >> gainr)
|
|
||||||
{
|
|
||||||
qqqGain[det][ring][wedge] = gainw;
|
|
||||||
// qqqrGain[det][ring][wedge] = gainr;
|
|
||||||
qqqGainValid[det][ring][wedge] = (gainw > 0);
|
|
||||||
// qqqrGainValid[det][ring][wedge] = (gainr > 0);
|
|
||||||
}
|
|
||||||
infile.close();
|
|
||||||
std::cout << "Loaded QQQ gains from " << filename << std::endl;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
// ----------------------- Load QQQ Calibrations
|
|
||||||
{
|
|
||||||
std::string filename = "qqq_Calib.txt";
|
|
||||||
std::ifstream infile(filename);
|
|
||||||
if (!infile.is_open())
|
|
||||||
{
|
|
||||||
std::cerr << "Error opening " << filename << "!" << std::endl;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
int det, ring, wedge;
|
|
||||||
double slope;
|
|
||||||
while (infile >> det >> ring >> wedge >> slope)
|
|
||||||
{
|
|
||||||
qqqCalib[det][ring][wedge] = slope;
|
|
||||||
qqqCalibValid[det][ring][wedge] = (slope > 0);
|
|
||||||
}
|
|
||||||
infile.close();
|
|
||||||
std::cout << "Loaded QQQ calibrations from " << filename << std::endl;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
Bool_t Analyzer::Process(Long64_t entry)
|
|
||||||
{
|
|
||||||
|
|
||||||
// if ( entry > 100 ) return kTRUE;
|
|
||||||
|
|
||||||
hitPos.Clear();
|
|
||||||
HitNonZero = false;
|
|
||||||
|
|
||||||
// if( entry > 1) return kTRUE;
|
|
||||||
// printf("################### ev : %llu \n", entry);
|
|
||||||
|
|
||||||
b_sx3Multi->GetEntry(entry);
|
|
||||||
b_sx3ID->GetEntry(entry);
|
|
||||||
b_sx3Ch->GetEntry(entry);
|
|
||||||
b_sx3E->GetEntry(entry);
|
|
||||||
b_sx3T->GetEntry(entry);
|
|
||||||
b_qqqMulti->GetEntry(entry);
|
|
||||||
b_qqqID->GetEntry(entry);
|
|
||||||
b_qqqCh->GetEntry(entry);
|
|
||||||
b_qqqE->GetEntry(entry);
|
|
||||||
b_qqqT->GetEntry(entry);
|
|
||||||
b_pcMulti->GetEntry(entry);
|
|
||||||
b_pcID->GetEntry(entry);
|
|
||||||
b_pcCh->GetEntry(entry);
|
|
||||||
b_pcE->GetEntry(entry);
|
|
||||||
b_pcT->GetEntry(entry);
|
|
||||||
|
|
||||||
sx3.CalIndex();
|
|
||||||
qqq.CalIndex();
|
|
||||||
pc.CalIndex();
|
|
||||||
|
|
||||||
// sx3.Print();
|
|
||||||
|
|
||||||
// ########################################################### Raw data
|
|
||||||
// //======================= SX3
|
|
||||||
sx3ecut = false;
|
|
||||||
std::vector<std::pair<int, int>> ID; // first = id, 2nd = index
|
|
||||||
for (int i = 0; i < sx3.multi; i++)
|
|
||||||
{
|
|
||||||
ID.push_back(std::pair<int, int>(sx3.id[i], i));
|
|
||||||
hsx3IndexVE->Fill(sx3.index[i], sx3.e[i]);
|
|
||||||
|
|
||||||
if (sx3.e[i] > 100)
|
|
||||||
{
|
|
||||||
sx3ecut = true;
|
|
||||||
}
|
|
||||||
|
|
||||||
for (int j = i + 1; j < sx3.multi; j++)
|
|
||||||
{
|
|
||||||
hsx3Coin->Fill(sx3.index[i], sx3.index[j]);
|
|
||||||
}
|
|
||||||
|
|
||||||
for (int j = 0; j < pc.multi; j++)
|
|
||||||
{
|
|
||||||
if (pc.index[j] < 24 && pc.e[j] > 100)
|
|
||||||
{
|
|
||||||
hsx3VpcIndex->Fill(sx3.index[i], pc.index[j]);
|
|
||||||
// if( sx3.ch[index] > 8 ){
|
|
||||||
// hsx3VpcE->Fill( sx3.e[i], pc.e[j] );
|
|
||||||
// }
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
if (ID.size() > 0)
|
|
||||||
{
|
|
||||||
std::sort(ID.begin(), ID.end(), [](const std::pair<int, int> &a, const std::pair<int, int> &b)
|
|
||||||
{ return a.first < b.first; });
|
|
||||||
// printf("##############################\n");
|
|
||||||
// for( size_t i = 0; i < ID.size(); i++) printf("%zu | %d %d \n", i, ID[i].first, ID[i].second );
|
|
||||||
|
|
||||||
std::vector<std::pair<int, int>> sx3ID;
|
|
||||||
sx3ID.push_back(ID[0]);
|
|
||||||
bool found = false;
|
|
||||||
for (size_t i = 1; i < ID.size(); i++)
|
|
||||||
{
|
|
||||||
if (ID[i].first == sx3ID.back().first)
|
|
||||||
{
|
|
||||||
sx3ID.push_back(ID[i]);
|
|
||||||
if (sx3ID.size() >= 3)
|
|
||||||
{
|
|
||||||
found = true;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
if (!found)
|
|
||||||
{
|
|
||||||
sx3ID.clear();
|
|
||||||
sx3ID.push_back(ID[i]);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
// printf("---------- sx3ID Multi : %zu \n", sx3ID.size());
|
|
||||||
|
|
||||||
if (found)
|
|
||||||
{
|
|
||||||
int sx3ChUp, sx3ChDn, sx3ChBk;
|
|
||||||
float sx3EUp, sx3EDn;
|
|
||||||
// printf("------ sx3 ID : %d, multi: %zu\n", sx3ID[0].first, sx3ID.size());
|
|
||||||
for (size_t i = 0; i < sx3ID.size(); i++)
|
|
||||||
{
|
|
||||||
int index = sx3ID[i].second;
|
|
||||||
// printf(" %zu | index %d | ch : %d, energy : %d \n", i, index, sx3.ch[index], sx3.e[index]);
|
|
||||||
|
|
||||||
if (sx3.ch[index] < 8)
|
|
||||||
{
|
|
||||||
if (sx3.ch[index] % 2 == 0)
|
|
||||||
{
|
|
||||||
sx3ChDn = sx3.ch[index];
|
|
||||||
sx3EDn = sx3.e[index];
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
sx3ChUp = sx3.ch[index];
|
|
||||||
sx3EUp = sx3.e[index];
|
|
||||||
}
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
sx3ChBk = sx3.ch[index];
|
|
||||||
}
|
|
||||||
for (int j = 0; j < pc.multi; j++)
|
|
||||||
{
|
|
||||||
// hsx3VpcIndex->Fill( sx3.index[i], pc.index[j] );
|
|
||||||
if (sx3.ch[index] > 8)
|
|
||||||
{
|
|
||||||
hsx3VpcE->Fill(sx3.e[i], pc.e[j]);
|
|
||||||
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
// sx3_contr.CalSX3Pos(sx3ID[0].first, sx3ChUp, sx3ChDn, sx3ChBk, sx3EUp, sx3EDn);
|
|
||||||
// hitPos = sx3_contr.GetHitPos();
|
|
||||||
// HitNonZero = true;
|
|
||||||
// hitPos.Print();
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
// //======================= QQQ
|
|
||||||
|
|
||||||
qqqEcut = false;
|
|
||||||
for (int i = 0; i < qqq.multi; i++)
|
|
||||||
{
|
|
||||||
// for( int j = 0; j < pc.multi; j++){
|
|
||||||
// if(pc.index[j]==4){
|
|
||||||
hqqqIndexVE->Fill(qqq.index[i], qqq.e[i]);
|
|
||||||
// }
|
|
||||||
// printf("QQQ ID : %d, ch : %d, e : %d \n", qqq.id[i], qqq.ch[i], qqq.e[i]);
|
|
||||||
if (qqq.e[i] > 100)
|
|
||||||
{
|
|
||||||
qqqEcut = true;
|
|
||||||
}
|
|
||||||
// }
|
|
||||||
for (int j = 0; j < qqq.multi; j++)
|
|
||||||
{
|
|
||||||
if (j == i)
|
|
||||||
continue;
|
|
||||||
hqqqCoin->Fill(qqq.index[i], qqq.index[j]);
|
|
||||||
}
|
|
||||||
|
|
||||||
for (int k = 0; k < pc.multi; k++)
|
|
||||||
{
|
|
||||||
if (pc.index[k] < 24 && pc.e[k] > 50)
|
|
||||||
{
|
|
||||||
hqqqVpcE->Fill(qqq.e[i], pc.e[k]);
|
|
||||||
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
|
|
||||||
hqqqVpcIndex->Fill(qqq.index[i], pc.index[k]);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
// }
|
|
||||||
|
|
||||||
for (int j = i + 1; j < qqq.multi; j++)
|
|
||||||
{
|
|
||||||
// if( qqq.used[i] == true ) continue;
|
|
||||||
|
|
||||||
// if( qqq.id[i] == qqq.id[j] && (16 - qqq.ch[i]) * (16 - qqq.ch[j]) < 0 ){ // must be same detector and wedge and ring
|
|
||||||
if (qqq.id[i] == qqq.id[j])
|
|
||||||
{ // must be same detector
|
|
||||||
|
|
||||||
if (qqq.e[i] > 100)
|
|
||||||
qqqEcut = true;
|
|
||||||
if (qqq.id[i] == qqq.id[j])
|
|
||||||
{
|
|
||||||
int chWedge = -1;
|
|
||||||
int chRing = -1;
|
|
||||||
float eWedgeRaw = 0.0;
|
|
||||||
float eWedge = 0.0;
|
|
||||||
float eWedgeMeV = 0.0;
|
|
||||||
float eRingRaw = 0.0;
|
|
||||||
float eRing = 0.0;
|
|
||||||
float eRingMeV = 0.0;
|
|
||||||
// plug in gains
|
|
||||||
if (qqq.ch[i] < 16 && qqq.ch[j] >= 16 && /*qqqrGainValid[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16] &&*/ qqqGainValid[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16])
|
|
||||||
{
|
|
||||||
chWedge = qqq.ch[i];
|
|
||||||
eWedgeRaw = qqq.e[i];
|
|
||||||
eWedge = qqq.e[i] * qqqGain[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16];
|
|
||||||
// printf("Wedge E: %.2f Gain: %.4f \n", eWedge, qqqGain[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16]);
|
|
||||||
chRing = qqq.ch[j] - 16;
|
|
||||||
eRingRaw = qqq.e[j];
|
|
||||||
eRing = qqq.e[j]; //* qqqrGain[qqq.id[j]][qqq.ch[j]][qqq.ch[i]-16];
|
|
||||||
}
|
|
||||||
else if (qqq.ch[j] < 16 && qqq.ch[i] >= 16 /* && qqqrGainValid[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16] */ && qqqGainValid[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16])
|
|
||||||
{
|
|
||||||
chWedge = qqq.ch[j];
|
|
||||||
eWedge = qqq.e[j] * qqqGain[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16];
|
|
||||||
eWedgeRaw = qqq.e[j];
|
|
||||||
chRing = qqq.ch[i] - 16;
|
|
||||||
eRing = qqq.e[i]; // * qqqrGain[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16];
|
|
||||||
eRingRaw = qqq.e[i];
|
|
||||||
}
|
|
||||||
else
|
|
||||||
continue;
|
|
||||||
// plug in calibrations
|
|
||||||
if (qqqCalibValid[qqq.id[i]][chRing][chWedge])
|
|
||||||
{
|
|
||||||
eWedgeMeV = eWedge * qqqCalib[qqq.id[i]][chRing][chWedge] / 1000;
|
|
||||||
eRingMeV = eRing * qqqCalib[qqq.id[i]][chRing][chWedge] / 1000;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
continue;
|
|
||||||
|
|
||||||
// printf(" ID : %d , chWedge : %d, chRing : %d \n", qqq.id[i], chWedge, chRing);
|
|
||||||
|
|
||||||
double theta = -TMath::Pi() / 2 + 2 * TMath::Pi() / 16 / 4. * (qqq.id[i] * 16 + chWedge + 0.5);
|
|
||||||
double rho = 50. + 40. / 16. * (chRing + 0.5);
|
|
||||||
// if(qqq.e[i]>50){
|
|
||||||
hqqqPolar->Fill(theta, rho);
|
|
||||||
// }
|
|
||||||
// qqq.used[i] = true;
|
|
||||||
// qqq.used[j] = true;
|
|
||||||
|
|
||||||
if (!HitNonZero)
|
|
||||||
{
|
|
||||||
double x = rho * TMath::Cos(theta);
|
|
||||||
double y = rho * TMath::Sin(theta);
|
|
||||||
hitPos.SetXYZ(x, y, 23 + 75 + 30);
|
|
||||||
HitNonZero = true;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
// //======================= PC
|
|
||||||
|
|
||||||
// Calculate the crossover points and put them into an array
|
|
||||||
|
|
||||||
pwinstance.ConstructGeo();
|
|
||||||
Coord Crossover[24][24][2];
|
|
||||||
TVector3 a, c, diff;
|
|
||||||
double a2, ac, c2, adiff, cdiff, denom, alpha, beta;
|
|
||||||
int index = 0;
|
|
||||||
for (int i = 0; i < pwinstance.An.size(); i++)
|
|
||||||
{
|
|
||||||
a = pwinstance.An[i].first - pwinstance.An[i].second;
|
|
||||||
|
|
||||||
for (int j = 0; j < pwinstance.Ca.size(); j++)
|
|
||||||
{
|
|
||||||
// Ok so this method uses what is essentially the solution of 2 equations to find the point of intersection between the anode and cathode wires
|
|
||||||
// here a and c are the vectors of the anode and cathode wires respectively
|
|
||||||
// diff is the perpendicular vector between the anode and cathode wires
|
|
||||||
// The idea behind this is to then find the scalars alpha and beta that give a ratio between 0 and -1,
|
|
||||||
|
|
||||||
c = pwinstance.Ca[j].first - pwinstance.Ca[j].second;
|
|
||||||
diff = pwinstance.An[i].first - pwinstance.Ca[j].first;
|
|
||||||
a2 = a.Dot(a);
|
|
||||||
c2 = c.Dot(c);
|
|
||||||
ac = a.Dot(c);
|
|
||||||
adiff = a.Dot(diff);
|
|
||||||
cdiff = c.Dot(diff);
|
|
||||||
denom = a2 * c2 - ac * ac;
|
|
||||||
alpha = (ac * cdiff - c2 * adiff) / denom;
|
|
||||||
beta = (a2 * cdiff - ac * adiff) / denom;
|
|
||||||
|
|
||||||
Crossover[i][j][0].x = pwinstance.An[i].first.X() + alpha * a.X();
|
|
||||||
Crossover[i][j][0].y = pwinstance.An[i].first.Y() + alpha * a.Y();
|
|
||||||
Crossover[i][j][0].z = pwinstance.An[i].first.Z() + alpha * a.Z();
|
|
||||||
if (Crossover[i][j][0].z < -190 || Crossover[i][j][0].z > 190)
|
|
||||||
{
|
|
||||||
Crossover[i][j][0].z = 9999999;
|
|
||||||
}
|
|
||||||
// placeholder variable Crossover[i][j][1].x has nothing to do with the geometry of the crossover and is being used to store the alpha value-
|
|
||||||
//-so that it can be used to sort "good" hits later
|
|
||||||
Crossover[i][j][1].x = alpha;
|
|
||||||
Crossover[i][j][1].y = 0;
|
|
||||||
// if(i==0){
|
|
||||||
// printf("CID, Crossover z and alpha are : %d %f %f \n", j, Crossover[i][j][0].z, Crossover[i][j][1].x /*this is alpha*/);
|
|
||||||
// }
|
|
||||||
// }
|
|
||||||
// }
|
|
||||||
}
|
|
||||||
}
|
|
||||||
// printf("Anode and cathode indices, alpha, denom, andiff, cndiff : %d %d %f %f %f %f\n", i, j, alpha, denom, adiff, cdiff);
|
|
||||||
|
|
||||||
// anodeIntersection.Clear();
|
|
||||||
for (int i = 0; i < pc.multi; i++)
|
|
||||||
{
|
|
||||||
|
|
||||||
if (pc.e[i] > 100)
|
|
||||||
{
|
|
||||||
hpcIndexVE->Fill(pc.index[i], pc.e[i]); // non gain matched energy
|
|
||||||
}
|
|
||||||
|
|
||||||
// Gain Matching of PC wires
|
|
||||||
if (pc.index[i] >= 0 && pc.index[i] < 48)
|
|
||||||
{
|
|
||||||
// printf("index: %d, Old cathode energy: %d \n", pc.index[i],pc.e[i]);
|
|
||||||
auto it = slopeInterceptMap.find(pc.index[i]);
|
|
||||||
if (it != slopeInterceptMap.end())
|
|
||||||
{
|
|
||||||
double slope = it->second.first;
|
|
||||||
double intercept = it->second.second;
|
|
||||||
// printf("slope: %f, intercept:%f\n" ,slope, intercept);
|
|
||||||
pc.e[i] = slope * pc.e[i] + intercept;
|
|
||||||
// printf("index: %d, New cathode energy: %d \n",pc.index[i], pc.e[i]);
|
|
||||||
}
|
|
||||||
hpcIndexVE_GM->Fill(pc.index[i], pc.e[i]);
|
|
||||||
// hPC_E[pc.index[i]]->Fill(pc.e[i]); // gain matched energy per channel
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
std::vector<std::pair<int, double>> anodeHits = {};
|
|
||||||
std::vector<std::pair<int, double>> cathodeHits = {};
|
|
||||||
std::vector<std::pair<int, double>> corrcatMax = {};
|
|
||||||
std::vector<std::pair<int, double>> corrcatnextMax = {};
|
|
||||||
std::vector<std::pair<int, double>> commcat = {};
|
|
||||||
int aID = 0;
|
|
||||||
int cID = 0;
|
|
||||||
float aE = 0;
|
|
||||||
float cE = 0;
|
|
||||||
float aESum = 0;
|
|
||||||
float cESum = 0;
|
|
||||||
float aEMax = 0;
|
|
||||||
float cEMax = 0;
|
|
||||||
float aEnextMax = 0;
|
|
||||||
float cEnextMax = 0;
|
|
||||||
int aIDMax = 0;
|
|
||||||
int cIDMax = 0;
|
|
||||||
int aIDnextMax = 0;
|
|
||||||
int cIDnextMax = 0;
|
|
||||||
|
|
||||||
// Define the excluded SX3 and QQQ channels
|
|
||||||
// std::unordered_set<int> excludeSX3 = {34, 35, 36, 37, 61, 62, 67, 73, 74, 75, 76, 77, 78, 79, 80, 93, 97, 100, 103, 108, 109, 110, 111, 112};
|
|
||||||
// std::unordered_set<int> excludeQQQ = {0, 17, 109, 110, 111, 112, 113, 119, 127, 128};
|
|
||||||
// inCuth=false;
|
|
||||||
// inCutl=false;
|
|
||||||
// inPCCut=false;
|
|
||||||
for (int i = 0; i < pc.multi; i++)
|
|
||||||
{
|
|
||||||
if (pc.e[i] > 100 /*&& pc.multi < 7*/)
|
|
||||||
{
|
|
||||||
// creating a vector of pairs of anode and cathode hits
|
|
||||||
if (pc.index[i] < 24)
|
|
||||||
{
|
|
||||||
anodeHits.push_back(std::pair<int, double>(pc.index[i], pc.e[i]));
|
|
||||||
}
|
|
||||||
else if (pc.index[i] >= 24)
|
|
||||||
{
|
|
||||||
cathodeHits.push_back(std::pair<int, double>(pc.index[i] - 24, pc.e[i]));
|
|
||||||
}
|
|
||||||
|
|
||||||
for (int j = i + 1; j < pc.multi; j++)
|
|
||||||
{
|
|
||||||
// if(PCCoinc_cut1->IsInside(pc.index[i], pc.index[j]) || PCCoinc_cut2->IsInside(pc.index[i], pc.index[j])){
|
|
||||||
// // hpcCoin->Fill(pc.index[i], pc.index[j]);
|
|
||||||
// inPCCut = true;
|
|
||||||
// }
|
|
||||||
hpcCoin->Fill(pc.index[i], pc.index[j]);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
// sorting the anode and cathode hits in descending order of energy
|
|
||||||
std::sort(anodeHits.begin(), anodeHits.end(), [](const std::pair<int, double> &a, const std::pair<int, double> &b)
|
|
||||||
{ return a.second > b.second; });
|
|
||||||
std::sort(cathodeHits.begin(), cathodeHits.end(), [](const std::pair<int, double> &a, const std::pair<int, double> &b)
|
|
||||||
{ return a.second > b.second; });
|
|
||||||
|
|
||||||
bool SiPCflag;
|
|
||||||
|
|
||||||
corrcatMax.clear();
|
|
||||||
if (anodeHits.size() >= 1 && cathodeHits.size() > 1)
|
|
||||||
{
|
|
||||||
if (((TMath::TanH(hitPos.Y() / hitPos.X())) > (TMath::TanH(a.Y() / a.X()) - TMath::PiOver4())) || ((TMath::TanH(hitPos.Y() / hitPos.X())) < (TMath::TanH(a.Y() / a.X()) + TMath::PiOver4())))
|
|
||||||
{
|
|
||||||
|
|
||||||
for (const auto &anode : anodeHits)
|
|
||||||
{
|
|
||||||
aID = anode.first;
|
|
||||||
aE = anode.second;
|
|
||||||
aESum += aE;
|
|
||||||
if (aE > aEMax)
|
|
||||||
{
|
|
||||||
aEMax = aE;
|
|
||||||
aIDMax = aID;
|
|
||||||
}
|
|
||||||
if (aE > aEnextMax && aE < aEMax)
|
|
||||||
{
|
|
||||||
aEnextMax = aE;
|
|
||||||
aIDnextMax = aID;
|
|
||||||
}
|
|
||||||
// for(const auto &cat : cathodeHits){
|
|
||||||
// hAVCcoin->Fill(aID, cat.first);
|
|
||||||
// }
|
|
||||||
}
|
|
||||||
|
|
||||||
// std::cout << " Anode iD : " << aIDMax << " Energy : " << aEMax << std::endl;
|
|
||||||
|
|
||||||
// printf("aID : %d, aE : %f, cE : %f\n", aID, aE, cE);
|
|
||||||
for (const auto &cathode : cathodeHits)
|
|
||||||
{
|
|
||||||
cID = cathode.first;
|
|
||||||
cE = cathode.second;
|
|
||||||
// std::cout << "Cathode ID : " << cID << " Energy : " << cE << std::endl;
|
|
||||||
|
|
||||||
hAVCcoin->Fill(aIDMax, cID);
|
|
||||||
|
|
||||||
// This section of code is used to find the cathodes are correlated with the max and next max anodes, as well as to figure out if there are any common cathodes
|
|
||||||
// the anodes are correlated with the cathodes +/-3 from the anode number in the reverse order
|
|
||||||
|
|
||||||
for (int j = -4; j < 3; j++)
|
|
||||||
{
|
|
||||||
if ((aIDMax + 24 + j) % 24 == 23 - cID)
|
|
||||||
/* the 23-cID is used to accomodate for the fact that the order of the cathodes was reversed relative top the physical geometry */
|
|
||||||
// if (Crossover[aIDMax][cID][0].z != 9999999)
|
|
||||||
{
|
|
||||||
corrcatMax.push_back(std::pair<int, double>(cID, cE));
|
|
||||||
cESum += cE;
|
|
||||||
// printf("Max Anode : %d Correlated Cathode : %d Anode Energy : %f z value : %f \n", aIDMax, cID, cESum, Crossover[aIDMax][cID][1].z /*prints alpha*/);
|
|
||||||
// std::cout << " Cathode iD : " << cID << " Energy : " << cE << std::endl;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
TVector3 anodeIntersection;
|
|
||||||
anodeIntersection.Clear();
|
|
||||||
// Implementing a method for PC reconstruction using a single Anode event
|
|
||||||
// if (anodeHits.size() == 1)
|
|
||||||
{
|
|
||||||
float x, y, z = 0;
|
|
||||||
for (const auto &corr : corrcatMax)
|
|
||||||
{
|
|
||||||
if (cESum > 0)
|
|
||||||
{
|
|
||||||
x += (corr.second) / cESum * Crossover[aIDMax][corr.first][0].x;
|
|
||||||
y += (corr.second) / cESum * Crossover[aIDMax][corr.first][0].y;
|
|
||||||
z += (corr.second) / cESum * Crossover[aIDMax][corr.first][0].z;
|
|
||||||
// printf("Max Anode : %d Correlated Cathode : %d cathode Energy : %f cESum Energy : %f z value : %f \n", aIDMax, corr.first, corr.second, cESum, Crossover[aIDMax][corr.first][1].z /*prints alpha*/);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
printf("Warning: No valid cathode hits to correlate with anode %d! \n", aIDMax);
|
|
||||||
}
|
|
||||||
// printf("EventID : %llu, Max Anode : %d Cathode: %d PC X and Y : (%f, %f) \n", entry, aIDMax, cID, Crossover[aIDMax][cID][0].x, Crossover[aIDMax][cID][0].y);
|
|
||||||
// for (int i = 0; i < sx3.multi; i++)
|
|
||||||
// {
|
|
||||||
// printf("EventID : %llu, HitPos X, Y, Z: %f %f %f SX3ID : %d %d \n", entry, hitPos.X(), hitPos.Y(), hitPos.Z(), sx3.id[i], sx3.ch[i]);
|
|
||||||
// }
|
|
||||||
|
|
||||||
// for (int i = 0; i < qqq.multi; i++)
|
|
||||||
// {
|
|
||||||
// printf("Max Anode : %d Cathode: %d PC X and Y : %f %f \n", aIDMax, cID, Crossover[aIDMax][cID][0].x, Crossover[aIDMax][cID][0].y);
|
|
||||||
// printf("HitPos X, Y, Z, QQQID : %f %f %f %d \n", hitPos.X(), hitPos.Y(), hitPos.Z(), qqq.id[i]);
|
|
||||||
// }
|
|
||||||
}
|
|
||||||
anodeIntersection = TVector3(x, y, z);
|
|
||||||
// std::cout << "Anode Intersection " << anodeIntersection.Z() << " " << x << " " << y << " " << z << std::endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (anodeIntersection.Z() != 0)
|
|
||||||
{
|
|
||||||
hPCZProj->Fill(anodeIntersection.Z());
|
|
||||||
}
|
|
||||||
// Filling the PC Z projection histogram
|
|
||||||
// std::cout << anodeIntersection.Z() << std::endl;
|
|
||||||
// hPCZProj->Fill(anodeIntersection.Z());
|
|
||||||
|
|
||||||
// }
|
|
||||||
|
|
||||||
// inCuth = false;
|
|
||||||
// inCutl = false;
|
|
||||||
// inPCCut = false;
|
|
||||||
// for(int j=i+1;j<pc.multi;j++){
|
|
||||||
// if(PCCoinc_cut1->IsInside(pc.index[i], pc.index[j]) || PCCoinc_cut2->IsInside(pc.index[i], pc.index[j])){
|
|
||||||
// // hpcCoin->Fill(pc.index[i], pc.index[j]);
|
|
||||||
// inPCCut = true;
|
|
||||||
// }
|
|
||||||
// hpcCoin->Fill(pc.index[i], pc.index[j]);
|
|
||||||
// }
|
|
||||||
|
|
||||||
// Check if the accumulated energies are within the defined ranges
|
|
||||||
// if (AnCatSum_high && AnCatSum_high->IsInside(aESum, cESum)) {
|
|
||||||
// inCuth = true;
|
|
||||||
// }
|
|
||||||
// if (AnCatSum_low && AnCatSum_low->IsInside(aESum, cESum)) {
|
|
||||||
// inCutl = true;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// Fill histograms based on the cut conditions
|
|
||||||
// if (inCuth && inPCCut) {
|
|
||||||
// hanVScatsum_hcut->Fill(aESum, cESum);
|
|
||||||
// }
|
|
||||||
// if (inCutl && inPCCut) {
|
|
||||||
// hanVScatsum_lcut->Fill(aESum, cESum);
|
|
||||||
// }
|
|
||||||
// for(auto anode : anodeHits){
|
|
||||||
|
|
||||||
// float aE = anode.second;
|
|
||||||
// aESum += aE;
|
|
||||||
// if(inPCCut){
|
|
||||||
hanVScatsum->Fill(aEMax, cESum);
|
|
||||||
// }
|
|
||||||
|
|
||||||
// if (sx3ecut)
|
|
||||||
// {
|
|
||||||
hCat4An->Fill(corrcatMax.size());
|
|
||||||
hNosvAe->Fill(corrcatMax.size(), aEMax);
|
|
||||||
hAnodehits->Fill(anodeHits.size());
|
|
||||||
// }
|
|
||||||
|
|
||||||
// }
|
|
||||||
if (anodeHits.size() < 1)
|
|
||||||
{
|
|
||||||
hCat0An->Fill(cathodeHits.size());
|
|
||||||
}
|
|
||||||
|
|
||||||
if (HitNonZero && anodeIntersection.Z() != 0)
|
|
||||||
{
|
|
||||||
pw_contr.CalTrack2(hitPos, anodeIntersection);
|
|
||||||
hZProj->Fill(pw_contr.GetZ0());
|
|
||||||
}
|
|
||||||
|
|
||||||
// ########################################################### Track constrcution
|
|
||||||
|
|
||||||
// ############################## DO THE KINEMATICS
|
|
||||||
|
|
||||||
return kTRUE;
|
|
||||||
}
|
|
||||||
|
|
||||||
void Analyzer::Terminate()
|
|
||||||
{
|
|
||||||
|
|
||||||
// gStyle->SetOptStat("neiou");
|
|
||||||
// TCanvas *canvas = new TCanvas("cANASEN", "ANASEN", 2000, 2000);
|
|
||||||
// canvas->Divide(3, 3);
|
|
||||||
|
|
||||||
// // hsx3VpcIndex->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-1
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// hsx3IndexVE->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-2
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// hqqqIndexVE->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-3
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// hpcIndexVE->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-4
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// hsx3Coin->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-5
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// canvas->cd(padID)->SetLogz(true);
|
|
||||||
|
|
||||||
// hqqqCoin->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-6
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// hpcCoin->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-7
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// // hsx3VpcIndex ->Draw("colz");
|
|
||||||
// hsx3VpcE->Draw("colz");
|
|
||||||
|
|
||||||
// //=============================================== pad-8
|
|
||||||
// padID++;
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
|
|
||||||
// // hqqqVpcIndex ->Draw("colz");
|
|
||||||
|
|
||||||
// hqqqVpcE->Draw("colz");
|
|
||||||
// //=============================================== pad-9
|
|
||||||
// padID++;
|
|
||||||
|
|
||||||
// // canvas->cd(padID)->DrawFrame(-50, -50, 50, 50);
|
|
||||||
// // hqqqPolar->Draw("same colz pol");
|
|
||||||
|
|
||||||
// canvas->cd(padID);
|
|
||||||
// canvas->cd(padID)->SetGrid(1);
|
|
||||||
// // hZProj->Draw();
|
|
||||||
// hanVScatsum->Draw("colz");
|
|
||||||
|
|
||||||
// // TFile *outRoot = new TFile("Histograms.root", "RECREATE");
|
|
||||||
|
|
||||||
// // if (!outRoot->IsOpen())
|
|
||||||
// // {
|
|
||||||
// // std::cerr << "Error opening file for writing!" << std::endl;
|
|
||||||
// // return;
|
|
||||||
// // }
|
|
||||||
|
|
||||||
// // // Loop through histograms and write them to the ROOT file
|
|
||||||
// // for (int i = 0; i < 48; i++)
|
|
||||||
// // {
|
|
||||||
// // if (hPC_E[i] != nullptr)
|
|
||||||
// // {
|
|
||||||
// // hPC_E[i]->Write(); // Write histogram to file
|
|
||||||
// // }
|
|
||||||
// // }
|
|
||||||
|
|
||||||
// // outRoot->Close();
|
|
||||||
}
|
|
||||||
129
Analyzer.h
129
Analyzer.h
|
|
@ -1,129 +0,0 @@
|
||||||
#ifndef Analyzer_h
|
|
||||||
#define Analyzer_h
|
|
||||||
|
|
||||||
#include <TROOT.h>
|
|
||||||
#include <TChain.h>
|
|
||||||
#include <TFile.h>
|
|
||||||
#include <TSelector.h>
|
|
||||||
|
|
||||||
#include "Armory/ClassDet.h"
|
|
||||||
|
|
||||||
class Analyzer : public TSelector {
|
|
||||||
public :
|
|
||||||
TTree *fChain; //!pointer to the analyzed TTree or TChain
|
|
||||||
|
|
||||||
// Fixed size dimensions of array or collections stored in the TTree if any.
|
|
||||||
|
|
||||||
// Declaration of leaf types
|
|
||||||
Det sx3;
|
|
||||||
Det qqq;
|
|
||||||
Det pc ;
|
|
||||||
Det misc;
|
|
||||||
|
|
||||||
ULong64_t evID;
|
|
||||||
UInt_t run;
|
|
||||||
|
|
||||||
// List of branches
|
|
||||||
TBranch *b_eventID; //!
|
|
||||||
TBranch *b_run; //!
|
|
||||||
TBranch *b_sx3Multi; //!
|
|
||||||
TBranch *b_sx3ID; //!
|
|
||||||
TBranch *b_sx3Ch; //!
|
|
||||||
TBranch *b_sx3E; //!
|
|
||||||
TBranch *b_sx3T; //!
|
|
||||||
TBranch *b_qqqMulti; //!
|
|
||||||
TBranch *b_qqqID; //!
|
|
||||||
TBranch *b_qqqCh; //!
|
|
||||||
TBranch *b_qqqE; //!
|
|
||||||
TBranch *b_qqqT; //!
|
|
||||||
TBranch *b_pcMulti; //!
|
|
||||||
TBranch *b_pcID; //!
|
|
||||||
TBranch *b_pcCh; //!
|
|
||||||
TBranch *b_pcE; //!
|
|
||||||
TBranch *b_pcT; //!
|
|
||||||
TBranch *b_miscMulti; //!
|
|
||||||
TBranch *b_miscID; //!
|
|
||||||
TBranch *b_miscCh; //!
|
|
||||||
TBranch *b_miscE; //!
|
|
||||||
TBranch *b_miscT; //!
|
|
||||||
TBranch *b_miscTf; //!
|
|
||||||
|
|
||||||
|
|
||||||
Analyzer(TTree * /*tree*/ =0) : fChain(0) { }
|
|
||||||
virtual ~Analyzer() { }
|
|
||||||
virtual Int_t Version() const { return 2; }
|
|
||||||
virtual void Begin(TTree *tree);
|
|
||||||
virtual void SlaveBegin(TTree *tree);
|
|
||||||
virtual void Init(TTree *tree);
|
|
||||||
virtual Bool_t Notify();
|
|
||||||
virtual Bool_t Process(Long64_t entry);
|
|
||||||
virtual Int_t GetEntry(Long64_t entry, Int_t getall = 0) { return fChain ? fChain->GetTree()->GetEntry(entry, getall) : 0; }
|
|
||||||
virtual void SetOption(const char *option) { fOption = option; }
|
|
||||||
virtual void SetObject(TObject *obj) { fObject = obj; }
|
|
||||||
virtual void SetInputList(TList *input) { fInput = input; }
|
|
||||||
virtual TList *GetOutputList() const { return fOutput; }
|
|
||||||
virtual void SlaveTerminate();
|
|
||||||
virtual void Terminate();
|
|
||||||
|
|
||||||
ClassDef(Analyzer,0);
|
|
||||||
};
|
|
||||||
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#ifdef Analyzer_cxx
|
|
||||||
void Analyzer::Init(TTree *tree){
|
|
||||||
|
|
||||||
// Set branch addresses and branch pointers
|
|
||||||
if (!tree) return;
|
|
||||||
fChain = tree;
|
|
||||||
fChain->SetMakeClass(1);
|
|
||||||
|
|
||||||
fChain->SetBranchAddress("evID", &evID, &b_eventID);
|
|
||||||
fChain->SetBranchAddress("run", &run, &b_run);
|
|
||||||
|
|
||||||
sx3.SetDetDimension(24,12);
|
|
||||||
qqq.SetDetDimension(4,32);
|
|
||||||
pc.SetDetDimension(2,24);
|
|
||||||
|
|
||||||
fChain->SetBranchAddress("sx3Multi", &sx3.multi, &b_sx3Multi);
|
|
||||||
fChain->SetBranchAddress("sx3ID", &sx3.id, &b_sx3ID);
|
|
||||||
fChain->SetBranchAddress("sx3Ch", &sx3.ch, &b_sx3Ch);
|
|
||||||
fChain->SetBranchAddress("sx3E", &sx3.e, &b_sx3E);
|
|
||||||
fChain->SetBranchAddress("sx3T", &sx3.t, &b_sx3T);
|
|
||||||
fChain->SetBranchAddress("qqqMulti", &qqq.multi, &b_qqqMulti);
|
|
||||||
fChain->SetBranchAddress("qqqID", &qqq.id, &b_qqqID);
|
|
||||||
fChain->SetBranchAddress("qqqCh", &qqq.ch, &b_qqqCh);
|
|
||||||
fChain->SetBranchAddress("qqqE", &qqq.e, &b_qqqE);
|
|
||||||
fChain->SetBranchAddress("qqqT", &qqq.t, &b_qqqT);
|
|
||||||
fChain->SetBranchAddress("pcMulti", &pc.multi, &b_pcMulti);
|
|
||||||
fChain->SetBranchAddress("pcID", &pc.id, &b_pcID);
|
|
||||||
fChain->SetBranchAddress("pcCh", &pc.ch, &b_pcCh);
|
|
||||||
fChain->SetBranchAddress("pcE", &pc.e, &b_pcE);
|
|
||||||
fChain->SetBranchAddress("pcT", &pc.t, &b_pcT);
|
|
||||||
fChain->SetBranchAddress("miscMulti", &misc.multi, &b_miscMulti);
|
|
||||||
fChain->SetBranchAddress("miscID", &misc.id, &b_miscID);
|
|
||||||
fChain->SetBranchAddress("miscCh", &misc.ch, &b_miscCh);
|
|
||||||
fChain->SetBranchAddress("miscE", &misc.e, &b_miscE);
|
|
||||||
fChain->SetBranchAddress("miscT", &misc.t, &b_miscT);
|
|
||||||
// fChain->SetBranchAddress("miscF", &misc.tf, &b_miscTf);
|
|
||||||
|
|
||||||
|
|
||||||
}
|
|
||||||
|
|
||||||
Bool_t Analyzer::Notify(){
|
|
||||||
|
|
||||||
return kTRUE;
|
|
||||||
}
|
|
||||||
|
|
||||||
void Analyzer::SlaveBegin(TTree * /*tree*/){
|
|
||||||
|
|
||||||
TString option = GetOption();
|
|
||||||
|
|
||||||
}
|
|
||||||
|
|
||||||
void Analyzer::SlaveTerminate(){
|
|
||||||
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
#endif // #ifdef Analyzer_cxx
|
|
||||||
Binary file not shown.
|
Before Width: | Height: | Size: 70 KiB |
54
TrackRecon.C
54
TrackRecon.C
|
|
@ -42,8 +42,8 @@ Int_t colors[40] = {
|
||||||
bool process_alpha_proton_scattering = false,
|
bool process_alpha_proton_scattering = false,
|
||||||
doMiscHistograms = true,
|
doMiscHistograms = true,
|
||||||
doRawHistos = false,
|
doRawHistos = false,
|
||||||
doPCSX3ClusterAnalysis = true,
|
doPCSX3ClusterAnalysis = false,
|
||||||
doPCQQQClusterAnalysis = true,
|
doPCQQQClusterAnalysis = false,
|
||||||
doOldAnalysis = false,
|
doOldAnalysis = false,
|
||||||
BenchMark = false,
|
BenchMark = false,
|
||||||
onewire_analysis = true,
|
onewire_analysis = true,
|
||||||
|
|
@ -196,7 +196,8 @@ inline void fillBeamProfile(HistPlotter *plotter, const TVector3 &vertex,
|
||||||
plotter->Fill1D(one + "pocaDist", 400, 0, 100, pocaDist, folder);
|
plotter->Fill1D(one + "pocaDist", 400, 0, 100, pocaDist, folder);
|
||||||
const double zLo = -440.0, zHi = 40.0;
|
const double zLo = -440.0, zHi = 40.0;
|
||||||
const int nSlice = 16;
|
const int nSlice = 16;
|
||||||
const double sliceW = (zHi - zLo) / nSlice; // 40 mm
|
const double sliceW = (zHi - zLo) / nSlice; // 30 mm
|
||||||
|
|
||||||
if (vz >= zLo && vz < zHi)
|
if (vz >= zLo && vz < zHi)
|
||||||
{
|
{
|
||||||
int is = static_cast<int>((vz - zLo) / sliceW);
|
int is = static_cast<int>((vz - zLo) / sliceW);
|
||||||
|
|
@ -1296,20 +1297,20 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
||||||
continue;
|
continue;
|
||||||
}
|
}
|
||||||
auto det = Fsx3.at(id);
|
auto det = Fsx3.at(id);
|
||||||
// if (det.valid)
|
if (det.valid && diagnostic_eplots)
|
||||||
// {
|
{
|
||||||
// // std::cout << det.frontEL << " " << det.frontEL*sx3RightGain[id][det.stripF] << std::endl;
|
// 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");
|
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");
|
||||||
// // plotter->Fill2D("unmatched_be_vs_x_sx3_id_" + std::to_string(id), 200, -1, 1, 800, 0, 4096, det.frontX, det.backE, "evsx");
|
plotter->Fill2D("unmatched_be_vs_x_sx3_id_" + std::to_string(id), 200, -1, 1, 800, 0, 4096, det.frontX, det.backE, "evsx");
|
||||||
// // plotter->Fill2D("unmatched_be_vs_x_sx3", 200, -1, 1, 800, 0, 4096, det.frontX, det.backE, "evsx");
|
plotter->Fill2D("unmatched_be_vs_x_sx3", 200, -1, 1, 800, 0, 4096, det.frontX, det.backE, "evsx");
|
||||||
// // plotter->Fill2D("matched_be_vs_x_sx3", 200, -60, 60, 800, 0, 8192, det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx");
|
plotter->Fill2D("matched_be_vs_x_sx3", 200, -60, 60, 800, 0, 8192, det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx");
|
||||||
// // plotter->Fill2D("matched_be_vs_x_sx3_id_" + std::to_string(id), 200, -60, 60, 800, 0, 8192, det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx");
|
plotter->Fill2D("matched_be_vs_x_sx3_id_" + std::to_string(id), 200, -60, 60, 800, 0, 8192, det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx");
|
||||||
|
|
||||||
// // plotter->Fill2D("matched_be_vs_x_sx3_id_" + std::to_string(id) + "_f" + std::to_string(det.stripF), 200, -60, 60, 800, 0, 8192,
|
plotter->Fill2D("matched_be_vs_x_sx3_id_" + std::to_string(id) + "_f" + std::to_string(det.stripF), 200, -60, 60, 800, 0, 8192,
|
||||||
// // det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx_matched");
|
det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF], det.backE * sx3BackGain[id][det.stripF][det.stripB], "evsx_matched");
|
||||||
// // plotter->Fill2D("fe_vs_x_sx3_id_"+std::to_string(id)+"_f"+std::to_string(det.stripF)+"_"+std::to_string(det.stripB),200,-1,1,800,0,4096,det.frontX,det.backE,"evsx");
|
plotter->Fill2D("fe_vs_x_sx3_id_"+std::to_string(id)+"_f"+std::to_string(det.stripF)+"_"+std::to_string(det.stripB),200,-1,1,800,0,4096,det.frontX,det.backE,"evsx");
|
||||||
// // plotter->Fill2D("l_vs_r_sx3_id_" + std::to_string(id) + "_f" + std::to_string(det.stripF), 800, 0, 4096, 800, 0, 4096, det.frontEL, det.frontER, "l_vs_r");
|
plotter->Fill2D("l_vs_r_sx3_id_" + std::to_string(id) + "_f" + std::to_string(det.stripF), 800, 0, 4096, 800, 0, 4096, det.frontEL, det.frontER, "l_vs_r");
|
||||||
// }
|
}
|
||||||
if (det.valid && (id == 9 || id == 7 || id == 1 || id == 3) && det.stripF != DEFAULT_NULL && det.stripB != DEFAULT_NULL)
|
if (det.valid && (id == 9 || id == 7 || id == 1 || id == 3) && det.stripF != DEFAULT_NULL && det.stripB != DEFAULT_NULL)
|
||||||
{
|
{
|
||||||
double z = det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF];
|
double z = det.frontX * sx3FrontGain[id][det.stripF] + sx3FrontOffset[id][det.stripF];
|
||||||
|
|
@ -1334,7 +1335,7 @@ Bool_t TrackRecon::Process(Long64_t entry)
|
||||||
plotter->Fill2D("sx3backs_gm", 100, 0, 100, 800, 0, 8192, det.stripB + 4 * id, backE, "hCalSX3");
|
plotter->Fill2D("sx3backs_gm", 100, 0, 100, 800, 0, 8192, det.stripB + 4 * id, backE, "hCalSX3");
|
||||||
plotter->Fill1D("sx3backs_calib", 800, 0, 8192, backE, "hCalSX3");
|
plotter->Fill1D("sx3backs_calib", 800, 0, 8192, backE, "hCalSX3");
|
||||||
|
|
||||||
// plotter->Fill2D("SX3CartesianPlot", 200, -100, 100, 200, -100, 100, 88.0*TMath::Cos(phi_n),88.0*TMath::Sin(phi_n), "hCalSX3");
|
plotter->Fill2D("SX3CartesianPlot", 200, -100, 100, 200, -100, 100, 88.0*TMath::Cos(phi_n),88.0*TMath::Sin(phi_n), "hCalSX3");
|
||||||
plotter->Fill2D("SX3CartesianPlot" + std::to_string(id), 200, -100, 100, 200, -100, 100, rho_at_strip * TMath::Cos(phi_n), rho_at_strip * TMath::Sin(phi_n), "hCalSX3");
|
plotter->Fill2D("SX3CartesianPlot" + std::to_string(id), 200, -100, 100, 200, -100, 100, rho_at_strip * TMath::Cos(phi_n), rho_at_strip * TMath::Sin(phi_n), "hCalSX3");
|
||||||
}
|
}
|
||||||
if (diagnostic_tplots)
|
if (diagnostic_tplots)
|
||||||
|
|
@ -4078,8 +4079,8 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
||||||
double vertex_z = r_rhoMin_fix.Z();
|
double vertex_z = r_rhoMin_fix.Z();
|
||||||
const bool axisSafe = (topo1 == "a1c2fix") || (topo1 == "a2c0") || (topo2 == "a1c1_inband");
|
const bool axisSafe = (topo1 == "a1c2fix") || (topo1 == "a2c0") || (topo2 == "a1c1_inband");
|
||||||
fillBeamProfile(plotter, r_rhoMin_fix, sievent.pos, x2f - sievent.pos,
|
fillBeamProfile(plotter, r_rhoMin_fix, sievent.pos, x2f - sievent.pos,
|
||||||
"reaction_" + globaltag + "_" + det, axisSafe);
|
"reaction_" + rx + "_" + det, axisSafe);
|
||||||
if (beamPerp(r_rhoMin_fix) > perp_cut || vertex_z < z_entrance)
|
if (beamPerp(r_rhoMin_fix) > perp_cut || vertex_z < z_entrance || vertex_z > 30.0)
|
||||||
return;
|
return;
|
||||||
|
|
||||||
double theta = (sievent.pos - r_rhoMin_fix).Theta();
|
double theta = (sievent.pos - r_rhoMin_fix).Theta();
|
||||||
|
|
@ -4217,7 +4218,7 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
||||||
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_VertexZ" + ejtag + sfx, 800, -400, 400, 800, 0, 0.6, 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_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 + "_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, 720, 0, 180, 800, 0, 0.6, theta * 180 / M_PI, anodeE_MeV, 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_phi" + ejtag + sfx, 90, -180, 180, 800, 0, 0.6, phi * 180 / M_PI, anodeE_MeV, pmlabel);
|
||||||
// if (anodeCh >= 0)
|
// 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_E" + ejtag + sfx + "_anode" + pad2(anodeCh), 400, 0, ef_max, 800, 0, 0.6, sievent.Energy1, anodeE_MeV, pmlabel);
|
||||||
|
|
@ -4245,7 +4246,11 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
||||||
{
|
{
|
||||||
if (!(pcevent.multi1 == 1 && (pcevent.multi2 == 1 || pcevent.multi2 == 2)))
|
if (!(pcevent.multi1 == 1 && (pcevent.multi2 == 1 || pcevent.multi2 == 2)))
|
||||||
continue;
|
continue;
|
||||||
if (TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) > phi_win)
|
// phi + time gate, matching protonAlphaElastic_core's equivalent a1c1/a1c2
|
||||||
|
// dispatch loop -- this branch previously gated on phi only.
|
||||||
|
bool phicut = TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) <= phi_win;
|
||||||
|
bool timecut = siPcCoincident(sievent.Time1, pcevent.Time1);
|
||||||
|
if (!(phicut && timecut))
|
||||||
continue;
|
continue;
|
||||||
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
double anodeE_MeV = (pcevent.Anodech >= 0 && pcevent.Anodech < 24)
|
||||||
? pcEnergySlope[pcevent.Anodech] * pcevent.Energy1
|
? pcEnergySlope[pcevent.Anodech] * pcevent.Energy1
|
||||||
|
|
@ -4290,13 +4295,18 @@ static void reaction_ax_core(HistPlotter *plotter, const std::vector<Event> &Si_
|
||||||
auto aPw = pwinstance.GetPseudoWire(aCl, "ANODE");
|
auto aPw = pwinstance.GetPseudoWire(aCl, "ANODE");
|
||||||
auto apwire = std::get<0>(aPw);
|
auto apwire = std::get<0>(aPw);
|
||||||
double apSumE = std::get<1>(aPw);
|
double apSumE = std::get<1>(aPw);
|
||||||
|
double apTSMaxE = std::get<3>(aPw);
|
||||||
|
|
||||||
bool isA2C0 = (aCl.size() == 2);
|
bool isA2C0 = (aCl.size() == 2);
|
||||||
const std::string a0tag = isA2C0 ? "a2c0" : "a1c0";
|
const std::string a0tag = isA2C0 ? "a2c0" : "a1c0";
|
||||||
TVector3 pc = isA2C0 ? a2c0_wirePos(apwire, sievent.pos.Phi(), isQQQ)
|
TVector3 pc = isA2C0 ? a2c0_wirePos(apwire, sievent.pos.Phi(), isQQQ)
|
||||||
: a1c0_wirePos(apwire, sievent.pos.Phi(), isQQQ);
|
: a1c0_wirePos(apwire, sievent.pos.Phi(), isQQQ);
|
||||||
|
|
||||||
if (TMath::Abs(sievent.pos.DeltaPhi(pc)) > phi_win)
|
// phi + time gate, matching protonAlphaElastic_core's equivalent a1c0/a2c0
|
||||||
|
// loop -- this branch previously gated on phi only.
|
||||||
|
bool phicut = TMath::Abs(sievent.pos.DeltaPhi(pc)) <= phi_win;
|
||||||
|
bool timecut = siPcCoincident(sievent.Time1, apTSMaxE);
|
||||||
|
if (!(phicut && timecut))
|
||||||
continue;
|
continue;
|
||||||
|
|
||||||
std::string pmlabel = folderPrefix + globaltag + "_" + rx + "+misc_" + det + "_" + a0tag;
|
std::string pmlabel = folderPrefix + globaltag + "_" + rx + "+misc_" + det + "_" + a0tag;
|
||||||
|
|
|
||||||
Binary file not shown.
|
Before Width: | Height: | Size: 201 KiB |
|
|
@ -34,7 +34,7 @@ export pressure_in_torr=250
|
||||||
export CO2percent=3
|
export CO2percent=3
|
||||||
|
|
||||||
# --- Block 1: 27Al Source Runs No Gas (1-8) ---
|
# --- Block 1: 27Al Source Runs No Gas (1-8) ---
|
||||||
if [[ 1 -eq 0 ]]; then
|
if [[ 1 -eq 1 ]]; then
|
||||||
export DATASET="27Al"
|
export DATASET="27Al"
|
||||||
export PREFIX="Run_"
|
export PREFIX="Run_"
|
||||||
export OUT_DIR="Output_av"
|
export OUT_DIR="Output_av"
|
||||||
|
|
@ -42,6 +42,7 @@ if [[ 1 -eq 0 ]]; then
|
||||||
rm -f ${OUT_DIR}/all.root
|
rm -f ${OUT_DIR}/all.root
|
||||||
|
|
||||||
parallel --bar -j 6 process_run ::: {1..8}
|
parallel --bar -j 6 process_run ::: {1..8}
|
||||||
|
exit
|
||||||
fi
|
fi
|
||||||
|
|
||||||
# --- Block 2: 17F Source Runs (5-14) ---
|
# --- Block 2: 17F Source Runs (5-14) ---
|
||||||
|
|
|
||||||
|
|
@ -5,9 +5,11 @@
|
||||||
//
|
//
|
||||||
// Usage (Single File):
|
// Usage (Single File):
|
||||||
// root -l -b -q 'make_pretty.C("myfile.root", "histName")'
|
// root -l -b -q 'make_pretty.C("myfile.root", "histName")'
|
||||||
|
// root -l -b -q 'make_pretty.C("myfile.root", "histName", "XLabel", "YLabel", minZ, minY, minX, maxX, maxY, rebinX, rebinY)'
|
||||||
//
|
//
|
||||||
// Usage (Multi-File Overlay):
|
// Usage (Multi-File Overlay):
|
||||||
// root -l -b -q 'make_pretty.C("file1.root,file2.root", "Run 1,Run 2", "histName", "X-Axis", "Y-Axis")'
|
// root -l -b -q 'make_pretty.C("file1.root,file2.root", "Run 1,Run 2", "histName", "X-Axis", "Y-Axis")'
|
||||||
|
// root -l -b -q 'make_pretty.C("file1.root,file2.root", "Run 1,Run 2", "histName", "X", "Y", yMin, yMax, xMin, xMax, canvasW, canvasH, rebinX)'
|
||||||
//
|
//
|
||||||
// =============================================================================
|
// =============================================================================
|
||||||
|
|
||||||
|
|
@ -78,7 +80,9 @@ void make_prettyplots(const char *rootFile,
|
||||||
double minY = -9999.0,
|
double minY = -9999.0,
|
||||||
double minX = -9999.0,
|
double minX = -9999.0,
|
||||||
double maxX = -9999.0,
|
double maxX = -9999.0,
|
||||||
double maxY = -9999.0)
|
double maxY = -9999.0,
|
||||||
|
int rebinX = 1,
|
||||||
|
int rebinY = 1) // rebinY only used for 2D histograms
|
||||||
{
|
{
|
||||||
|
|
||||||
SetStyle();
|
SetStyle();
|
||||||
|
|
@ -129,6 +133,15 @@ void make_prettyplots(const char *rootFile,
|
||||||
TH2 *h = (TH2 *)clone;
|
TH2 *h = (TH2 *)clone;
|
||||||
h->SetStats(0);
|
h->SetStats(0);
|
||||||
|
|
||||||
|
// Rebin() and Rebin2D() modify in place and return `this` when called
|
||||||
|
// with no newname, but reassign anyway in case that ever changes.
|
||||||
|
if (rebinX > 1 || rebinY > 1)
|
||||||
|
{
|
||||||
|
TH2 *rebinned = (TH2 *)h->Rebin2D(rebinX, rebinY);
|
||||||
|
if (rebinned)
|
||||||
|
h = rebinned;
|
||||||
|
}
|
||||||
|
|
||||||
// --- Apply Minimum Z if provided ---
|
// --- Apply Minimum Z if provided ---
|
||||||
if (minZ != -9999.0)
|
if (minZ != -9999.0)
|
||||||
{
|
{
|
||||||
|
|
@ -154,13 +167,19 @@ void make_prettyplots(const char *rootFile,
|
||||||
h->GetXaxis()->CenterTitle(true);
|
h->GetXaxis()->CenterTitle(true);
|
||||||
h->GetYaxis()->CenterTitle(true);
|
h->GetYaxis()->CenterTitle(true);
|
||||||
|
|
||||||
// gPad->SetLogz();
|
gPad->SetLogz();
|
||||||
|
|
||||||
h->Draw("colz");
|
h->Draw("colz");
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
TH1 *h = (TH1 *)clone;
|
TH1 *h = (TH1 *)clone;
|
||||||
|
if (rebinX > 1)
|
||||||
|
{
|
||||||
|
TH1 *rebinned = (TH1 *)h->Rebin(rebinX);
|
||||||
|
if (rebinned)
|
||||||
|
h = rebinned;
|
||||||
|
}
|
||||||
if (minX != -9999.0)
|
if (minX != -9999.0)
|
||||||
h->GetXaxis()->SetRangeUser(minX, h->GetXaxis()->GetXmax());
|
h->GetXaxis()->SetRangeUser(minX, h->GetXaxis()->GetXmax());
|
||||||
if (maxX != -9999.0)
|
if (maxX != -9999.0)
|
||||||
|
|
@ -200,7 +219,7 @@ void make_prettyplots(const char *rootFile,
|
||||||
// NOTE: canvasW/canvasH added as trailing optional args (default = old
|
// NOTE: canvasW/canvasH added as trailing optional args (default = old
|
||||||
// hardcoded values, so every existing call site keeps working unchanged).
|
// hardcoded values, so every existing call site keeps working unchanged).
|
||||||
// Pass larger values for a bigger/higher-resolution output image.
|
// Pass larger values for a bigger/higher-resolution output image.
|
||||||
void make_prettyplots(TString filesCSV, TString labelsCSV, TString histName, TString xAxisLabel="", TString yAxisLabel="", double yMin=-9999, double yMax=-9999, double xMin=-9999, double xMax=-9999, int canvasW=2100, int canvasH=1575)
|
void make_prettyplots(TString filesCSV, TString labelsCSV, TString histName, TString xAxisLabel="", TString yAxisLabel="", double yMin=-9999, double yMax=-9999, double xMin=-9999, double xMax=-9999, int canvasW=2100, int canvasH=1575, int rebinX=1)
|
||||||
{
|
{
|
||||||
SetStyle();
|
SetStyle();
|
||||||
|
|
||||||
|
|
@ -243,6 +262,13 @@ void make_prettyplots(TString filesCSV, TString labelsCSV, TString histName, TSt
|
||||||
clone->SetDirectory(0); // Detach from file
|
clone->SetDirectory(0); // Detach from file
|
||||||
f->Close();
|
f->Close();
|
||||||
|
|
||||||
|
if (rebinX > 1)
|
||||||
|
{
|
||||||
|
TH1 *rebinned = (TH1 *)clone->Rebin(rebinX);
|
||||||
|
if (rebinned)
|
||||||
|
clone = rebinned;
|
||||||
|
}
|
||||||
|
|
||||||
// Apply X-range BEFORE checking the max Y-height, otherwise peaks
|
// Apply X-range BEFORE checking the max Y-height, otherwise peaks
|
||||||
// outside the viewing range might scale the Y-axis unnecessarily!
|
// outside the viewing range might scale the Y-axis unnecessarily!
|
||||||
if (xMin != -9999.0 && xMax != -9999.0)
|
if (xMin != -9999.0 && xMax != -9999.0)
|
||||||
|
|
|
||||||
Loading…
Reference in New Issue
Block a user