Skip to content
Snippets Groups Projects
Commit c1cbe3d5 authored by Theodore Efremov's avatar Theodore Efremov :hibiscus:
Browse files

Added two macro, one to make cut on the DE_E and the other to apply a...

Added two macro, one to make cut on the DE_E and the other to apply a decreasing expo fit on each segment of the IC
parent 1e068a96
No related branches found
No related tags found
1 merge request!27Draft: [Epic] Preparation of the environement for the new GaseousDetectorScorers...
Pipeline #362526 failed
#include "IC_Dampen_Recombination.cxx"
void HistoFillerDE_E();
/**
* * Draw the DE_E plot and save the histo into a file
* */
void CutZ(bool GenDE=true){
if (GenDE==true) HistoFillerDE_E() ;
TFile* infile = TFile::Open("Output/DE_E.root");
//===========================================================================================================
// Load Histo
//===========================================================================================================
// Get number of TH2F
int HistoCount = 0 ;
TIter next(infile->GetListOfKeys());
TKey* key;
while ((key=(TKey*)next())){
if (std::string(key->GetClassName()) == "TH2F"){
HistoCount ++;
}
}
vector<TH2F*> hDE_E;
//for each histo in the file
for (int i=0 ; i< HistoCount ; i++){
if (i==0) hDE_E.push_back((TH2F*)infile->Get("DE_E_online"));
if (i==1) hDE_E.push_back((TH2F*)infile->Get("DE_Etot_online"));
}
//===========================================================================================================
// Draw histo
//===========================================================================================================
TCanvas *cDE = new TCanvas("cDE","cDE",1000,500);
cDE->Divide(2);
cDE->cd(1);
hDE_E[0]->Draw();
cDE->cd(2);
hDE_E[1]->Draw();
}
void HistoFillerDE_E(){
//===========================================================================================================
// Loading var
//===========================================================================================================
TChain* chain = new TChain("PhysicsTree");
chain->Add("../../../root/analysis/VamosCalib241.root");
TICPhysics* IC = new TICPhysics() ;
double FF_IC_X, FF_IC_Y;
chain->SetBranchStatus("FF_IC_X", true);
chain->SetBranchAddress("FF_IC_X", &FF_IC_X);
chain->SetBranchStatus("FF_IC_Y", true);
chain->SetBranchAddress("FF_IC_Y", &FF_IC_Y);
chain->SetBranchStatus("IC", true);
chain->SetBranchAddress("IC", &IC);
//===========================================================================================================
// Defining histo
//===========================================================================================================
TH2F* hDE_E_nocorr = new
TH2F("DE_E_online","DE_E_online",1000,0,20000,1000,0,25000);
TH2F* hDE_Etot = new
TH2F("DE_Etot_online","DE_Etot_online",1000,0,50000,1000,0,25000);
//===========================================================================================================
// FILL HISTO
//===========================================================================================================
int Nentries = chain->GetEntries();
for (int e = 0; e < Nentries; e++) {
// Printing status
if (e % 100000 == 0){
cout << "Please wait, we are " << e * 100. / Nentries << "% done ..\r" << flush;
}
chain->GetEntry(e);
double DE = 0 ,E = 0, Etot = 0 ;
//===========================================================================================================
// Fill De from seg 1 to 3
//===========================================================================================================
for (int seg = 0 ; seg < sizeof(IC->fIC_raw)/sizeof(IC->fIC_raw[0]) ; seg++ ){
if(seg > 1 && seg < 4) DE += IC->fIC_raw[seg] ;
else if (seg >4) E+= IC->fIC_raw[seg] ;
}
//===========================================================================================================
// Filling DE with ICO and 4
//===========================================================================================================
DE += IC->fIC_raw[0] ;
DE = DE *0.5 + IC->fIC_raw[4];
double IC1corr = (IC->fIC_raw[1]*(1-0.000686068*FF_IC_Y))*(1-4.88238e-05*FF_IC_Y+7.40395e-06*FF_IC_Y*FF_IC_Y);
double DE_Bis = 0.5*(IC1corr+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
Etot = E+DE_Bis;
//===========================================================================================================
// Fill histo
//===========================================================================================================
hDE_E_nocorr->Fill(E,DE_Bis);
hDE_Etot->Fill(Etot,DE_Bis);
}
TFile* outfile = new TFile("Output/DE_E.root","recreate");
hDE_E_nocorr->Write();
hDE_Etot->Write();
outfile->Close();
}
...@@ -51,7 +51,9 @@ void ICCorr() { ...@@ -51,7 +51,9 @@ void ICCorr() {
TH2F* hDE_E_IC0_Corr = new TH2F* hDE_E_IC0_Corr = new
TH2F("hDE_E_IC0_Corr","hDE_E_IC0_Corr",1000,0,20000,1000,0,25000); TH2F("hDE_E_IC0_Corr","hDE_E_IC0_Corr",1000,0,20000,1000,0,25000);
TH2F* hDEcut=
new TH2F("DEcut","DEcut",1000,0,20000,1000,0,25000);
hICorr.push_back(hICorr_X); hICorr.push_back(hICorr_X);
hICorr.push_back(hICorr_Y); hICorr.push_back(hICorr_Y);
...@@ -169,6 +171,7 @@ void ICCorr() { ...@@ -169,6 +171,7 @@ void ICCorr() {
// Chio online // Chio online
//=========================================================================================================== //===========================================================================================================
double IC1corr = (IC->fIC_raw[1]*(1-0.000686068*FF_IC_Y))*(1-4.88238e-05*FF_IC_Y+7.40395e-06*FF_IC_Y*FF_IC_Y); double IC1corr = (IC->fIC_raw[1]*(1-0.000686068*FF_IC_Y))*(1-4.88238e-05*FF_IC_Y+7.40395e-06*FF_IC_Y*FF_IC_Y);
double DE_Bis = 0.5*(IC1corr+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4]; double DE_Bis = 0.5*(IC1corr+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
//double DE_Bis = 0.5*(IC->fIC_raw[1]+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4]; //double DE_Bis = 0.5*(IC->fIC_raw[1]+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
...@@ -178,7 +181,11 @@ void ICCorr() { ...@@ -178,7 +181,11 @@ void ICCorr() {
hDE_E->Fill(E,DE); hDE_E->Fill(E,DE);
hDE_Ecorr01->Fill(E,DEcorr); hDE_Ecorr01->Fill(E,DEcorr);
hDE_E_IC1_Corr->Fill(E,DEIC1corr); hDE_E_IC1_Corr->Fill(E,DEIC1corr);
hDE_E_IC0_Corr->Fill(E,DE_IC0_Corr_Y);
double DEcut;
if ((FF_IC_Y>-60 && FF_IC_Y<-55)|| (FF_IC_Y>55 && FF_IC_Y<60)){ DEcut=DE_Bis;
hDEcut->Fill(E,DEcut);
}
} //End loop on event } //End loop on event
//=========================================================================================================== //===========================================================================================================
...@@ -205,6 +212,6 @@ void ICCorr() { ...@@ -205,6 +212,6 @@ void ICCorr() {
c2->cd(3); c2->cd(3);
hDEBis_Ecorr01->Draw(); hDEBis_Ecorr01->Draw();
c2->cd(4); c2->cd(4);
hDE_E_IC0_Corr->Draw(); hDEcut->Draw();
} // End spline chio XY } // End spline chio XY
#include <TFile.h>
#include <TICPhysics.h>
#include <TChain.h>
#include <TF1.h>
#include <TH2.h>
#include <Math/Functor.h>
#include <TTreeReaderValue.h>
#include <vector>
using namespace std;
void HistoFillerIC();
void IC_Dampen_Recombination(bool Fill = false){
if (Fill==true) HistoFillerIC();
TFile *infile = TFile::Open("Output/HistoDampen.root");
vector<TH2F*> hICX, hICY ;
vector<TH1F*> hICX_Profile, hICY_Profile ;
for (int seg = 0 ; seg < 11; seg++ ){
hICX.push_back((TH2F*)infile->Get(Form("hChio_X_%d",seg)));
hICY.push_back((TH2F*)infile->Get(Form("hChio_Y_%d",seg)));
hICX_Profile.push_back((TH1F*)infile->Get(Form("hChio_X_%d_pfx",seg)));
hICY_Profile.push_back((TH1F*)infile->Get(Form("hChio_Y_%d_pfx",seg)));
}
//===========================================================================================================
// Fit Histo
//===========================================================================================================
vector<TF1*> Fit_Y(11);
for (int seg=0; seg<11 ; seg++){
Fit_Y[seg] = new
TF1(Form("Fit_Y_Expo_IC%d",seg),[](double *x, double*p) {
if (x[0]>-50 && x[0]<50) {
TF1::RejectPoint();
}
//Define here the function
return p[0] * exp(p[1]* (x[0]-p[3])) + p[2] ;
},-60,63,4);
Fit_Y[seg]->SetParameters(1,-10,4800,-60);
cout << "****************" << endl << "IC"<<seg << endl;
hICY_Profile[seg]->Fit(Fit_Y[seg] ,"Rnone"); vector<TH2F*> hICX, hICY ;
vector<TH1F*> hICX_Profile, hICY_Profile ;
}//end loop on histo
//===========================================================================================================
// Draw histo
//===========================================================================================================
vector<TCanvas*> cIC(11);
for (int seg=0; seg<1; seg++){
cIC[seg] = new TCanvas
(Form("cIC%d",seg),Form("cIC%d",seg),100,100,1000,600);
cIC[seg]->Divide(2);
cIC[seg]->cd(1);
hICY[seg]->Draw();
cIC[seg]->cd(2);
hICY_Profile[seg]->Draw();
Fit_Y[seg]->Draw("same");
}
}
void HistoFillerIC(){
//===========================================================================================================
// INPUT
//===========================================================================================================
TChain* chain = new TChain("PhysicsTree");
chain->Add("../../../root/analysis/VamosCalib241.root");
double FF_IC_X, FF_IC_Y;
TICPhysics* IC = new TICPhysics();
chain->SetBranchStatus("FF_IC_X", true);
chain->SetBranchAddress("FF_IC_X", &FF_IC_X);
chain->SetBranchStatus("FF_IC_Y", true);
chain->SetBranchAddress("FF_IC_Y", &FF_IC_Y);
chain->SetBranchStatus("IC", true);
chain->SetBranchAddress("IC", &IC);
//===========================================================================================================
// Declare and fill histo
//===========================================================================================================
vector<TH2F*> hICX, hICY ;
vector<TH1F*> hICX_Profile, hICY_Profile ;
int NSegmentChio = 11 ;
for (int seg = 0 ; seg < 11 ; seg++){
TH2F* hChio_X = new TH2F(Form("hChio_X_%d",seg) ,Form("hChio_X_%d",seg) ,1000 ,-800 ,800 ,500 ,2000,9500 );
TH2F* hChio_Y = new TH2F(Form("hChio_Y_%d",seg) ,Form("hChio_Y_%d",seg) ,1000 ,-160 ,160 ,500 ,2000,9500 );
TH1F* hProfile_X = (TH1F*)hChio_X->ProfileX();
TH1F* hProfile_Y = (TH1F*)hChio_Y->ProfileX();
hICX.push_back(hChio_X);
hICY.push_back(hChio_Y);
hICX_Profile.push_back(hProfile_X);
hICY_Profile.push_back(hProfile_Y);
} //end loop on segment
// Loop on event
int Nentries = chain->GetEntries();
for (int e = 0 ; e<Nentries ; e++){ // to be parallelised
if ((e/Nentries)*100 %100 ==0){
cout << "Please wait, we are " << e * 100. / Nentries << "% done ..\r" << flush;
}
chain->GetEntry(e);
for (int seg=0 ; seg<11 ; seg++ ){
hICX[seg]->Fill(FF_IC_X,IC->fIC_raw[seg]);
hICY[seg]->Fill(FF_IC_Y,IC->fIC_raw[seg]);
}//end loop seg
} //end loop entries
//ReProfile
for (int seg=0; seg<11 ;seg++){
hICX_Profile[seg] = (TH1F*)hICX[seg]->ProfileX();
hICY_Profile[seg] = (TH1F*)hICY[seg]->ProfileX();
}
TFile* outfile = new TFile("Output/HistoDampen.root","recreate");
for (int seg=0 ; seg<11 ;seg++){
hICX[seg]->Write();
hICY[seg]->Write();
hICX_Profile[seg]->Write();
hICY_Profile[seg]->Write();
}
outfile->Close();
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment