modified: Analyzer.C
new file: GainMatch.C new file: GainMatch.h
This commit is contained in:
parent
39e8f41ab1
commit
a7a765c059
233
Analyzer.C
233
Analyzer.C
|
@ -194,10 +194,13 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
|
||||
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] );
|
||||
// }
|
||||
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] );
|
||||
// }
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
@ -300,18 +303,20 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
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++)
|
||||
{
|
||||
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]);
|
||||
}
|
||||
// }
|
||||
}
|
||||
// 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
|
||||
|
@ -333,7 +338,7 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
// 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 = 10. + 40. / 16. * (chRing + 0.5);
|
||||
double rho = 50. + 40. / 16. * (chRing + 0.5);
|
||||
// if(qqq.e[i]>50){
|
||||
hqqqPolar->Fill(theta, rho);
|
||||
// }
|
||||
|
@ -385,37 +390,20 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
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][2].x has nothing to do with the geometry of the crossover and is being used to store the alpha value-
|
||||
}
|
||||
// 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*/);
|
||||
// }
|
||||
// }
|
||||
// }
|
||||
}
|
||||
}
|
||||
|
||||
// if (i == 4)
|
||||
{
|
||||
// for (int k = -4; k < 3; k++)
|
||||
// {
|
||||
// if ((i + 24 + k) % 24 == 23 - j) // the 23-j is used to accomodate for the fact that the order of the cathodes was reversed
|
||||
// {
|
||||
// Crossover[i][j][1].y = 1;
|
||||
// }
|
||||
// }
|
||||
// for (int k = -4; k < 3; k++)
|
||||
// {
|
||||
// if (Crossover[i][j][k].y != 1)
|
||||
// // if (alpha < 1 && alpha >= -1)
|
||||
// {
|
||||
// // printf("i an: %d %f %f %f \n", i, an.X(), an.Y(), an.Z());
|
||||
// Crossover[i][j][0].z = 9999999; // this is a placeholder value to indicate that the anode and cathode wires do not intersect
|
||||
// }
|
||||
// }
|
||||
}
|
||||
// printf("Anode and cathode indices, alpha, denom, andiff, cndiff : %d %d %f %f %f %f\n", i, j, alpha, denom, adiff, cdiff);
|
||||
|
||||
// anodeIntersection.Clear();
|
||||
|
@ -502,11 +490,11 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
{ return a.second > b.second; });
|
||||
|
||||
bool SiPCflag;
|
||||
corrcatMax.clear();
|
||||
|
||||
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())))
|
||||
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)
|
||||
|
@ -543,15 +531,17 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
// 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)
|
||||
for (int j = -4; j < 3; j++)
|
||||
{
|
||||
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;
|
||||
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;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
@ -591,10 +581,13 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
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());
|
||||
// hPCZProj->Fill(anodeIntersection.Z());
|
||||
|
||||
// }
|
||||
|
||||
|
@ -645,7 +638,7 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
hCat0An->Fill(cathodeHits.size());
|
||||
}
|
||||
|
||||
if (HitNonZero)
|
||||
if (HitNonZero && anodeIntersection.Z() != 0)
|
||||
{
|
||||
pw_contr.CalTrack2(hitPos, anodeIntersection);
|
||||
hZProj->Fill(pw_contr.GetZ0());
|
||||
|
@ -661,99 +654,99 @@ Bool_t Analyzer::Process(Long64_t entry)
|
|||
void Analyzer::Terminate()
|
||||
{
|
||||
|
||||
gStyle->SetOptStat("neiou");
|
||||
TCanvas *canvas = new TCanvas("cANASEN", "ANASEN", 2000, 2000);
|
||||
canvas->Divide(3, 3);
|
||||
// gStyle->SetOptStat("neiou");
|
||||
// TCanvas *canvas = new TCanvas("cANASEN", "ANASEN", 2000, 2000);
|
||||
// canvas->Divide(3, 3);
|
||||
|
||||
// hsx3VpcIndex->Draw("colz");
|
||||
// // hsx3VpcIndex->Draw("colz");
|
||||
|
||||
//=============================================== pad-1
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-1
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
hsx3IndexVE->Draw("colz");
|
||||
// hsx3IndexVE->Draw("colz");
|
||||
|
||||
//=============================================== pad-2
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-2
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
hqqqIndexVE->Draw("colz");
|
||||
// hqqqIndexVE->Draw("colz");
|
||||
|
||||
//=============================================== pad-3
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-3
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
hpcIndexVE->Draw("colz");
|
||||
// hpcIndexVE->Draw("colz");
|
||||
|
||||
//=============================================== pad-4
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-4
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
hsx3Coin->Draw("colz");
|
||||
// hsx3Coin->Draw("colz");
|
||||
|
||||
//=============================================== pad-5
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-5
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
canvas->cd(padID)->SetLogz(true);
|
||||
// canvas->cd(padID)->SetLogz(true);
|
||||
|
||||
hqqqCoin->Draw("colz");
|
||||
// hqqqCoin->Draw("colz");
|
||||
|
||||
//=============================================== pad-6
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-6
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
hpcCoin->Draw("colz");
|
||||
// hpcCoin->Draw("colz");
|
||||
|
||||
//=============================================== pad-7
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-7
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
// hsx3VpcIndex ->Draw("colz");
|
||||
hsx3VpcE->Draw("colz");
|
||||
// // hsx3VpcIndex ->Draw("colz");
|
||||
// hsx3VpcE->Draw("colz");
|
||||
|
||||
//=============================================== pad-8
|
||||
padID++;
|
||||
canvas->cd(padID);
|
||||
canvas->cd(padID)->SetGrid(1);
|
||||
// //=============================================== pad-8
|
||||
// padID++;
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
|
||||
// hqqqVpcIndex ->Draw("colz");
|
||||
// // hqqqVpcIndex ->Draw("colz");
|
||||
|
||||
hqqqVpcE->Draw("colz");
|
||||
//=============================================== pad-9
|
||||
padID++;
|
||||
// hqqqVpcE->Draw("colz");
|
||||
// //=============================================== pad-9
|
||||
// padID++;
|
||||
|
||||
// canvas->cd(padID)->DrawFrame(-50, -50, 50, 50);
|
||||
// hqqqPolar->Draw("same colz pol");
|
||||
// // 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");
|
||||
// canvas->cd(padID);
|
||||
// canvas->cd(padID)->SetGrid(1);
|
||||
// // hZProj->Draw();
|
||||
// hanVScatsum->Draw("colz");
|
||||
|
||||
// TFile *outRoot = new TFile("Histograms.root", "RECREATE");
|
||||
// // TFile *outRoot = new TFile("Histograms.root", "RECREATE");
|
||||
|
||||
// if (!outRoot->IsOpen())
|
||||
// {
|
||||
// std::cerr << "Error opening file for writing!" << std::endl;
|
||||
// return;
|
||||
// }
|
||||
// // 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
|
||||
// }
|
||||
// }
|
||||
// // // 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();
|
||||
// // outRoot->Close();
|
||||
}
|
||||
|
|
273
GainMatch.C
Normal file
273
GainMatch.C
Normal file
|
@ -0,0 +1,273 @@
|
|||
#define GainMatch_cxx
|
||||
|
||||
#include "GainMatch.h"
|
||||
#include <TH2.h>
|
||||
#include <TF1.h>
|
||||
#include <TStyle.h>
|
||||
#include <TCanvas.h>
|
||||
#include <TMath.h>
|
||||
#include <TCutG.h>
|
||||
|
||||
#include <utility>
|
||||
#include <algorithm>
|
||||
|
||||
#include "Armory/ClassSX3.h"
|
||||
#include "Armory/ClassPW.h"
|
||||
|
||||
#include "TVector3.h"
|
||||
|
||||
TH2F *hSX3FvsB;
|
||||
TH2F *hQQQFVB;
|
||||
|
||||
int padID = 0;
|
||||
|
||||
SX3 sx3_contr;
|
||||
PW pw_contr;
|
||||
TVector3 hitPos;
|
||||
bool HitNonZero;
|
||||
|
||||
void GainMatch::Begin(TTree * /*tree*/)
|
||||
{
|
||||
TString option = GetOption();
|
||||
|
||||
hSX3FvsB = new TH2F("hSX3FvsB", "SX3 Front vs Back; Front E; Back E", 400, 0, 16000, 400, 0, 16000);
|
||||
|
||||
hQQQFVB = new TH2F("hQQQFVB", "number of good QQQ vs QQQ id", 400, 0, 16000, 400, 0, 16000);
|
||||
sx3_contr.ConstructGeo();
|
||||
pw_contr.ConstructGeo();
|
||||
}
|
||||
|
||||
Bool_t GainMatch::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
|
||||
// Initialize vector to store pairs of ID and index
|
||||
std::vector<std::pair<int, int>> ID; // first = id, second = index
|
||||
|
||||
// Loop through each entry in sx3.multi
|
||||
for (int i = 0; i < sx3.multi; i++)
|
||||
{
|
||||
// Store ID and index pair in ID vector
|
||||
ID.push_back(std::pair<int, int>(sx3.id[i], i));
|
||||
}
|
||||
|
||||
// Check if the ID vector is not empty
|
||||
if (ID.size() > 0)
|
||||
{
|
||||
// Sort the ID vector by the first element (ID)
|
||||
std::sort(ID.begin(), ID.end(), [](const std::pair<int, int> &a, const std::pair<int, int> &b)
|
||||
{ return a.first < b.first; });
|
||||
|
||||
// Create a new vector to store IDs that are the same
|
||||
std::vector<std::pair<int, int>> sx3ID;
|
||||
sx3ID.push_back(ID[0]); // Start with the first ID
|
||||
bool found = false;
|
||||
|
||||
// Loop through the sorted IDs and group by same IDs
|
||||
for (size_t i = 1; i < ID.size(); i++)
|
||||
{
|
||||
if (ID[i].first == sx3ID.back().first)
|
||||
{
|
||||
// If the ID matches the last one in sx3ID, add to the group
|
||||
sx3ID.push_back(ID[i]);
|
||||
|
||||
// If 3 or more IDs are grouped, set found to true
|
||||
if (sx3ID.size() >= 3)
|
||||
{
|
||||
found = true;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// If a new ID is encountered and no group is found, reset the group
|
||||
if (!found)
|
||||
{
|
||||
sx3ID.clear();
|
||||
sx3ID.push_back(ID[i]);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
// If a group of 3 or more IDs is found, process the channels and energies
|
||||
if (found)
|
||||
{
|
||||
int sx3ChUp = -1, sx3ChDn = -1, sx3ChBk = -1;
|
||||
float sx3EUp = 0.0, sx3EDn = 0.0, sx3EBk = 0.0;
|
||||
|
||||
// Loop through the grouped IDs
|
||||
for (size_t i = 0; i < sx3ID.size(); i++)
|
||||
{
|
||||
int index = sx3ID[i].second; // Get the index from the pair
|
||||
|
||||
// Categorize channel and energy
|
||||
if (sx3.ch[index] < 8)
|
||||
{
|
||||
// If channel is less than 8, assign it to up or down
|
||||
if (sx3.ch[index] % 2 == 0)
|
||||
{
|
||||
sx3ChDn = sx3.ch[index];
|
||||
sx3EDn = sx3.e[index];
|
||||
}
|
||||
else
|
||||
{
|
||||
sx3ChUp = sx3.ch[index];
|
||||
sx3EUp = sx3.e[index];
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// If channel is >= 8, assign it to back
|
||||
sx3ChBk = sx3.ch[index];
|
||||
sx3EBk = sx3.e[index];
|
||||
}
|
||||
}
|
||||
// make a histogram for Sx3 back energy vs sum of sx3 front energies
|
||||
hSX3FvsB->Fill(sx3EUp + sx3EDn, sx3EBk);
|
||||
// Further energy matching calculations can be added here...
|
||||
}
|
||||
}
|
||||
// //======================= QQQ
|
||||
TH2F *hist1 = NULL;
|
||||
// Loop through each entry in qqq.multi
|
||||
for (int i = 0; i < qqq.multi; i++)
|
||||
{
|
||||
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
|
||||
|
||||
int chWedge = -1;
|
||||
int chRing = -1;
|
||||
float eWedge = 0.0;
|
||||
float eRing = 0.0;
|
||||
if (qqq.ch[i] < qqq.ch[j])
|
||||
{
|
||||
chRing = qqq.ch[j] - 16;
|
||||
eRing = qqq.e[j];
|
||||
chWedge = qqq.ch[i];
|
||||
eWedge = qqq.e[i];
|
||||
}
|
||||
else
|
||||
{
|
||||
chRing = qqq.ch[i];
|
||||
eRing = qqq.e[i];
|
||||
chWedge = qqq.ch[j] - 16;
|
||||
eWedge = qqq.e[j];
|
||||
}
|
||||
// printf(" ID : %d , chWedge : %d, chRing : %d \n", qqq.id[i], chWedge, chRing);
|
||||
hQQQFVB->Fill(qqq.e[i], qqq.e[j]);
|
||||
// 1. Create/get individual 2D histogram
|
||||
TString histName = Form("hQQQFVB_r%d_w%d_id%d", chRing, chWedge, qqq.id[i]);
|
||||
TH2F *hist2d = (TH2F *)gDirectory->Get(histName);
|
||||
if (!hist2d)
|
||||
{
|
||||
hist2d = new TH2F(histName, Form("QQQ GainMatch ID%d R%d W%d;Wedge E;Ring E", qqq.id[i], chRing, chWedge),
|
||||
400, 0, 16000, 400, 0, 16000);
|
||||
}
|
||||
|
||||
hist2d->Fill(eWedge, eRing);
|
||||
}
|
||||
}
|
||||
}
|
||||
return kTRUE;
|
||||
}
|
||||
void GainMatch::Terminate()
|
||||
{
|
||||
TIter next(gDirectory->GetList());
|
||||
TObject *obj;
|
||||
|
||||
while ((obj = next()))
|
||||
{
|
||||
if (!obj->InheritsFrom("TH2F"))
|
||||
continue;
|
||||
|
||||
TH2F *hist2d = (TH2F *)obj;
|
||||
TString name = hist2d->GetName();
|
||||
if (!name.BeginsWith("hQQQFVB_r"))
|
||||
continue;
|
||||
|
||||
if (hist2d->GetEntries() < 1000)
|
||||
continue;
|
||||
|
||||
std::vector<double> wedge_vals;
|
||||
std::vector<double> peak_vals;
|
||||
|
||||
for (int binX = 1; binX <= hist2d->GetNbinsX(); ++binX)
|
||||
{
|
||||
TH1D *projY = hist2d->ProjectionY("_py", binX, binX);
|
||||
if (projY->GetEntries() < 30)
|
||||
{
|
||||
delete projY;
|
||||
continue;
|
||||
}
|
||||
|
||||
double mean = projY->GetMean();
|
||||
TF1 *fit = new TF1("fit", "gaus", mean - 100, mean + 100);
|
||||
projY->Fit(fit, "QNR");
|
||||
|
||||
double wedgeE = hist2d->GetXaxis()->GetBinCenter(binX);
|
||||
double ringPeak = fit->GetParameter(1);
|
||||
|
||||
wedge_vals.push_back(wedgeE);
|
||||
peak_vals.push_back(ringPeak);
|
||||
|
||||
delete projY;
|
||||
delete fit;
|
||||
}
|
||||
|
||||
if (wedge_vals.size() >= 5)
|
||||
{
|
||||
TGraph *gr = new TGraph(wedge_vals.size(), &wedge_vals[0], &peak_vals[0]);
|
||||
gr->SetName(name + "_fit");
|
||||
TF1 *line = new TF1("line", "pol1", 0, 16000);
|
||||
gr->Fit(line, "Q");
|
||||
|
||||
double gain = line->GetParameter(1);
|
||||
double offset = line->GetParameter(0);
|
||||
printf("Gain match %s → Gain = %.4f, Offset = %.2f\n", name.Data(), gain, offset);
|
||||
|
||||
TCanvas *c1 = new TCanvas(name + "_c", name + "_c");
|
||||
gr->SetTitle(Form("Gain Match: %s", name.Data()));
|
||||
gr->GetXaxis()->SetTitle("Wedge Energy");
|
||||
gr->GetYaxis()->SetTitle("Ring Energy");
|
||||
gr->Draw("AP");
|
||||
line->Draw("same");
|
||||
c1->SaveAs(Form("%s_fit.png", name.Data()));
|
||||
delete c1;
|
||||
}
|
||||
}
|
||||
}
|
114
GainMatch.h
Normal file
114
GainMatch.h
Normal file
|
@ -0,0 +1,114 @@
|
|||
#ifndef GainMatch_h
|
||||
#define GainMatch_h
|
||||
|
||||
#include <TROOT.h>
|
||||
#include <TChain.h>
|
||||
#include <TFile.h>
|
||||
#include <TSelector.h>
|
||||
|
||||
#include "Armory/ClassDet.h"
|
||||
|
||||
class GainMatch : 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 ;
|
||||
|
||||
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; //!
|
||||
|
||||
GainMatch(TTree * /*tree*/ =0) : fChain(0) { }
|
||||
virtual ~GainMatch() { }
|
||||
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(GainMatch,0);
|
||||
};
|
||||
|
||||
#endif
|
||||
|
||||
#ifdef GainMatch_cxx
|
||||
void GainMatch::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);
|
||||
|
||||
}
|
||||
|
||||
Bool_t GainMatch::Notify(){
|
||||
|
||||
return kTRUE;
|
||||
}
|
||||
|
||||
void GainMatch::SlaveBegin(TTree * /*tree*/){
|
||||
|
||||
TString option = GetOption();
|
||||
|
||||
}
|
||||
|
||||
void GainMatch::SlaveTerminate(){
|
||||
|
||||
}
|
||||
|
||||
|
||||
#endif // #ifdef GainMatch_cxx
|
Loading…
Reference in New Issue
Block a user