modified: Analyzer.C

new file:   PCGainMatch.C
	new file:   PCGainMatch.h
This commit is contained in:
vs19g 2024-08-28 10:49:36 -04:00
parent 48ede97992
commit a6e754b958
3 changed files with 519 additions and 4 deletions

View File

@ -279,7 +279,13 @@ Bool_t Analyzer::Process(Long64_t entry){
std::unordered_set<int> excludeQQQ = {0, 17, 109, 110, 111, 112, 113, 119, 127, 128}; std::unordered_set<int> excludeQQQ = {0, 17, 109, 110, 111, 112, 113, 119, 127, 128};
inCut=false; inCut=false;
for( int i = 0; i < pc.multi; i ++){ for( int i = 0; i < pc.multi; i ++){
if(pc.e[i]>50 && pc.multi<10){ for(int j=0; j<pc.multi;j++){
if(pc.id[j]==0){
anodeCount++;
}
}
if(pc.e[i]>50 && anodeCount==1){
// hpcIndexVE->Fill( pc.index[i], pc.e[i] ); // hpcIndexVE->Fill( pc.index[i], pc.e[i] );
// for( int j = i+1; j < pc.multi; j++){ // for( int j = i+1; j < pc.multi; j++){
// hpcCoin->Fill( pc.index[i], pc.index[j]); // hpcCoin->Fill( pc.index[i], pc.index[j]);
@ -303,8 +309,7 @@ Bool_t Analyzer::Process(Long64_t entry){
} }
} }
// hpcIndexVE->Fill( pc.index[i], pc.e[i] ); // hpcIndexVE->Fill( pc.index[i], pc.e[i] );
anodeCount=pc.multi; hAnodeMultiplicity->Fill(pc.multi);
hAnodeMultiplicity->Fill(anodeCount);
float aESum = 0; float aESum = 0;
if (anodeHits.size()>=1 && cathodeHits.size() >= 1){ if (anodeHits.size()>=1 && cathodeHits.size() >= 1){
@ -420,5 +425,5 @@ void Analyzer::Terminate(){
// hZProj->Draw(); // hZProj->Draw();
hanVScatsum->Draw("colz"); hanVScatsum->Draw("colz");
// hAnodeHits->Draw("colz"); // hAnodeHits->Draw("colz");
// hAnodeMultiplicity->Draw();? // hAnodeMultiplicity->Draw();
} }

396
PCGainMatch.C Normal file
View File

@ -0,0 +1,396 @@
#define PCGainMatch_cxx
#include "PCGainMatch.h"
#include <TH2.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 * hsx3IndexVE;
TH2F * hqqqIndexVE;
TH2F * hpcIndexVE;
TH2F * hsx3Coin;
TH2F * hqqqCoin;
TH2F * hpcCoin;
TH2F * hqqqPolar;
TH2F * hsx3VpcIndex;
TH2F * hqqqVpcIndex;
TH2F * hqqqVpcE;
TH2F * hsx3VpcE;
TH2F * hanVScatsum;
TH2F * hanVScatsum_hcut;
TH2F * hanVScatsum_lcut;
TH2F * hAnodeHits;
TH1F * hAnodeMultiplicity;
int padID = 0;
SX3 sx3_contr;
PW pw_contr;
TVector3 hitPos;
bool HitNonZero;
TH1F * hZProj;
TCutG *AnCatSum_high;
TCutG *AnCatSum_low;
bool inCuth;
bool inCut;
void PCGainMatch::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); hsx3IndexVE->SetNdivisions( -612, "x");
hqqqIndexVE = new TH2F("hqqqIndexVE", "QQQ index vs Energy; QQQ index ; Energy", 4*2*16, 0, 4*2*16, 400, 0, 5000); hqqqIndexVE->SetNdivisions( -1204, "x");
hpcIndexVE = new TH2F("hpcIndexVE", "PC index vs Energy; PC index ; Energy", 2*24, 0, 2*24, 400, 0, 4000); hpcIndexVE->SetNdivisions( -1204, "x");
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);
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);
hsx3VpcIndex->SetNdivisions( -612, "x");
hsx3VpcIndex->SetNdivisions( -12, "y");
hqqqVpcIndex = new TH2F("hqqqVpcindex", "qqq vs pc; qqq index; pc index", 4*2*16, 0, 4*2*16, 48, 0, 48);
hqqqVpcIndex->SetNdivisions( -612, "x");
hqqqVpcIndex->SetNdivisions( -12, "y");
hqqqVpcE = new TH2F("hqqqVpcEnergy", "qqq vs pc; qqq energy; pc energy", 400, 0, 5000, 400, 0, 5000);
hqqqVpcE->SetNdivisions( -612, "x");
hqqqVpcE->SetNdivisions( -12, "y");
hsx3VpcE = new TH2F("hsx3VpcEnergy", "sx3 vs pc; sx3 energy; pc energy", 400, 0, 5000, 400, 0, 5000);
hsx3VpcE->SetNdivisions( -612, "x");
hsx3VpcE->SetNdivisions( -12, "y");
hZProj = new TH1F("hZProj", "Nos of anodes", 20, 0, 19);
hAnodeHits = new TH2F("hAnodeHits", "Anode vs Anode Energy, Anode ID; Anode E", 24,0 , 23, 400, 0 , 20000);
hAnodeMultiplicity = new TH1F("hAnodeMultiplicity", "Number of Anodes/Event", 40, 0, 40);
hanVScatsum = new TH2F("hanVScatsum", "Anode vs Cathode Sum; Anode E; Cathode E", 400,0 , 10000, 400, 0 , 16000);
hanVScatsum_lcut = new TH2F("hanVScatsum_LCUT", "Anode vs Cathode Sum; Anode E; Cathode E", 400,0 , 10000, 400, 0 , 16000);
hanVScatsum_hcut = new TH2F("hanVScatsum_HCUT", "Anode vs Cathode Sum; Anode E; Cathode E", 400,0 , 10000, 400, 0 , 16000);
sx3_contr.ConstructGeo();
pw_contr.ConstructGeo();
TFile *f1 = new TFile("AnCatSum_high.root");
TFile *f2 = new TFile("AnCatSum_low.root");
AnCatSum_high= (TCutG*)f1->Get("AnCatSum_high");
AnCatSum_low= (TCutG*)f1->Get("AnCatSum_low");
}
Bool_t PCGainMatch::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
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));
}
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 && pc.index[j]<24 && pc.e[j]>50 ){
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
for( int i = 0; i < qqq.multi; i ++){
for( int j = 0; j < pc.multi; j++){
if(pc.index[j]<24 && pc.e[j]>50 ){
hqqqVpcE->Fill( qqq.e[i], pc.e[j] );
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
hqqqVpcIndex->Fill( qqq.index[i], pc.index[j] );
}
// }
}
}
// //======================= PC
std::vector<std::pair<int, double>> anodeHits;
std::vector<std::pair<int, double>> cathodeHits;
int aID = 0;
int cID = 0;
int anodeCount = 0;
float aE = 0;
float cE = 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;
inCut=false;
for( int i = 0; i < pc.multi; i ++){
for(int j=0; j<pc.multi;j++){
if(pc.id[j]==0){
anodeCount++;
}
}
if(pc.e[i]>50){
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
// for( int j = i+1; j < pc.multi; j++){
// hpcCoin->Fill( pc.index[i], pc.index[j]);
// }
// for (int j=0;j<sx3.multi;j++){
// if(excludeSX3.find(sx3.index[j]) == excludeSX3.end()){
// if(pc.e[i]>100){
if (pc.index[i] < 24 ){
anodeHits.push_back(std::pair<int, double>(pc.index[i], pc.e[i]));
// anodeCount++;
} else if (pc.index[i] >= 24){
cathodeHits.push_back(std::pair<int, double>(pc.index[i], pc.e[i]));
}
// }
// }
// }
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
float aESum = 0;
if (anodeHits.size()==1 && cathodeHits.size() >= 1){//gain matching for anodes using simple single anode events
// hpcIndexVE->Fill( pc.index[i], pc.e[i] );
// for( int j = i+1; j < pc.multi; j++){
// hpcCoin->Fill( pc.index[i], pc.index[j]);
// }
for (const auto& anode : anodeHits) {
float cESum = 0;
// for(int l=0; l<sx3.multi; l++){
// if (sx3.index[l]==80){
int aID = anode.first;
float aE = anode.second;
aESum += aE;
for (const auto& cathode : cathodeHits) {
int cID = cathode.first;
float cE = cathode.second;
cESum += cE;
}
if( AnCatSum_high->IsInside(aE, cESum)){
inCuth = true;
}
// if( AnCatSum_low->IsInside(aE, cESum)){
// inCut = true;
// }
if(inCuth){
// hanVScatsum_hcut->Fill(aE,cESum);
// }
// if(inCut){
// hanVScatsum_lcut->Fill(aE,cESum);
// }
hpcIndexVE->Fill( pc.index[i], pc.e[i] );
for( int j = i+1; j < pc.multi; j++){
hpcCoin->Fill( pc.index[i], pc.index[j]);
}
hAnodeMultiplicity->Fill(anodeHits.size());
hanVScatsum->Fill(aE,cESum);
hAnodeHits->Fill(aID, aE);
}
}
}
}
}
// hanVScatsum->Fill(aE,cE);
if( HitNonZero){
pw_contr.CalTrack( hitPos, aID, cID);
hZProj->Fill(pw_contr.GetZ0());
}
//########################################################### Track constrcution
//############################## DO THE KINEMATICS
return kTRUE;
}
void PCGainMatch::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);
hanVScatsum->Draw("colz");
// hAnodeHits->Draw("colz");
// hAnodeMultiplicity->Draw();
}

114
PCGainMatch.h Normal file
View File

@ -0,0 +1,114 @@
#ifndef PCGainMatch_h
#define PCGainMatch_h
#include <TROOT.h>
#include <TChain.h>
#include <TFile.h>
#include <TSelector.h>
#include "Armory/ClassDet.h"
class PCGainMatch : 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; //!
PCGainMatch(TTree * /*tree*/ =0) : fChain(0) { }
virtual ~PCGainMatch() { }
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(PCGainMatch,0);
};
#endif
#ifdef PCGainMatch_cxx
void PCGainMatch::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 PCGainMatch::Notify(){
return kTRUE;
}
void PCGainMatch::SlaveBegin(TTree * /*tree*/){
TString option = GetOption();
}
void PCGainMatch::SlaveTerminate(){
}
#endif // #ifdef Analyzer_cxx