FairRoot/PandaRoot
Functions
analysis_cluster_energyCorrection_Ntuple_suplement.C File Reference
#include <stdio.h>
#include <stdlib.h>
#include <vector>
#include "TF1.h"
#include "MyClass.h"

Go to the source code of this file.

Functions

int analysis_cluster_energyCorrection_Ntuple_suplement ()
 

Function Documentation

int analysis_cluster_energyCorrection_Ntuple_suplement ( )

Definition at line 7 of file analysis_cluster_energyCorrection_Ntuple_suplement.C.

References c1, ctime, Double_t, i, ncounts, printf(), rtime, timer, and TString.

8 {
9  TStopwatch timer;
10  timer.Start();
11 
13  cf=new TFile("/opt/fairroot/pandaroot/macro/params/gamma_en_th_corr_TGeant3.root");//the histogram hisEnergyDelta has 26 bins in energy and 61 in theta
14  TH2F *HistClusterEnergyCorrectionFactor_johan = (TH2F*)hisEnergyDelta->Clone();//making a copy of the histogram hisEnergyDelta
15  HistClusterEnergyCorrectionFactor_johan->SetName("HistClusterEnergyCorrectionFactor_johan");//renaming the new one is recommended
16  for (int i=1; i<=26; i++)//cluster energy
17  for (int j=1; j<=61; j++)//cluster theta
18  HistClusterEnergyCorrectionFactor_johan->SetBinContent(i,j,0);//emptying the histogram to be filled later in this code; this would have the same binning as of Johan map!
19 
20  TH2F *HistClusterEnergyCorrectionFactor_johan_original = (TH2F*)hisEnergyDelta->Clone();//making another copy
21  HistClusterEnergyCorrectionFactor_johan_original->SetName("HistClusterEnergyCorrectionFactor_johan_original");
22  for (int i=1; i<=26; i++)//cluster energy
23  for (int j=1; j<=61; j++)//cluster theta
24  HistClusterEnergyCorrectionFactor_johan_original->SetBinContent(i,j,0);//emptying the histogram
25 
26  TH2F *HistClusterEnergyCorrectionFactor_johan_count = (TH2F*)hisEnergyDelta->Clone();//making another copy to keep the bin contents
27  HistClusterEnergyCorrectionFactor_johan_count->SetName("HistClusterEnergyCorrectionFactor_johan_count");
28  for (int i=1; i<=26; i++)//cluster energy
29  for (int j=1; j<=61; j++)//cluster theta
30  HistClusterEnergyCorrectionFactor_johan_count->SetBinContent(i,j,0);//emptying the histogram
31 
33  Float_t xlowarray[36] = {0, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 1.,
34  1.05, 1.1, 1.15, 1.2, 1.25, 1.3, 1.35, 1.4, 1.45, 1.5, 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5.}; //energy bin edges (GeV)
35  Float_t ylowarray[54] = {0, 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 20., 21., 22., 23., 24., 25., 30., 35., 40., 45., 50., 55., 60.,
36  65., 70., 75., 80., 85., 90., 95., 100., 105., 110., 115., 120., 125., 130., 135., 140., 142., 147., 152.,
37  154., 156., 158., 161., 165., 171.};//theta bin edges
38  TH2F *HistClusterEnergyCorrectionFactor_Hossein = new TH2F("HistClusterEnergyCorrectionFactor_Hossein","", 35, xlowarray, 53, ylowarray);
39  TH2F *HistClusterEnergyCorrectionFactor_Hossein_original = new TH2F("HistClusterEnergyCorrectionFactor_Hossein_original","", 35, xlowarray, 53, ylowarray);
40  TH2F *HistClusterEnergyCorrectionFactor_Hossein_count = new TH2F("HistClusterEnergyCorrectionFactor_Hossein_count","", 35, xlowarray, 53, ylowarray);
41 
42  TH1I *HISTnumberOfClusters_thresholdPassed = new TH1I("HISTnumberOfClusters_thresholdPassed","Total number of clusters per event satisfying cluster threshold",140,0,70);
43 
44  TH2F *HistClusterEnergyCorrectionFactor_original = new TH2F("HistClusterEnergyCorrectionFactor_original","cluster energy correction factor",180, 0, 180, 360, -180, 180);//theta, phi
45  TH2F *HistClusterEnergyCorrectionFactor_count = new TH2F("HistClusterEnergyCorrectionFactor_count","cluster energy correction factor (bin content container)",180, 0, 180, 360, -180, 180);
46  TH2F *HistClusterEnergyCorrectionFactor = new TH2F("HistClusterEnergyCorrectionFactor","cluster energy correction factor (normalized)",180, 0, 180, 360, -180, 180);
47 
48  TH2F *HistClusterEnergyCorrectionFactor_theta_Egamma_original = new TH2F("HistClusterEnergyCorrectionFactor_theta_Egamma_original","cluster energy correction factor (normalized)",500, 0, 5., 80, 0, 180);//energy in GeV
49  TH2F *HistClusterEnergyCorrectionFactor_theta_Egamma_count = new TH2F("HistClusterEnergyCorrectionFactor_theta_Egamma_count","cluster energy correction factor (bin content container)",500, 0, 5., 80, 0, 180);//energy in GeV
50  TH2F *HistClusterEnergyCorrectionFactor_theta_Egamma = new TH2F("HistClusterEnergyCorrectionFactor_theta_Egamma","cluster energy correction factor (normalized)",500, 0, 5., 80, 0, 180);//energy in GeV
51 
53  TString inputfile("./energyCorrectionFactors_30MeVthreshold_cond2_Ntuple_10MillionEvt.root");
54  TString outputfile("./energyCorrectionFactors_30MeVthreshold_cond2_Ntuple_10MillionEvt_supplement.root");
56  TChain *chain = new TChain("EventParams");
57  chain->Add(inputfile);
58  chain->Print();
59 
60  //Declaration of branch types:
61  float eMc, eCluster, ClustEnCorrFactor, clustTheta, clustPhi;
62  //Set branch addresses:
63  chain->SetBranchAddress("eMc", &eMc);
64  chain->SetBranchAddress("eCluster", &eCluster);
65  chain->SetBranchAddress("ClustEnCorrFactor", &ClustEnCorrFactor);
66  chain->SetBranchAddress("clustTheta", &clustTheta);
67  chain->SetBranchAddress("clustPhi", &clustPhi);
68 
69  Int_t ncounts = chain->GetEntries();
70  cout << "Number of events = " << ncounts<< endl;
71  //ncounts=20000; // analyse specific number of events
72 
73  for (int j = 0; j < ncounts; j++)
74  {
75  chain->GetEntry(j);
76  HistClusterEnergyCorrectionFactor_original->Fill(clustTheta, clustPhi, ClustEnCorrFactor);//filling the 2D histogram with a weight of clusterEnergyCorrectionFactor
77  HistClusterEnergyCorrectionFactor_count->Fill(clustTheta, clustPhi);//this histogram is filled just to keep track of the number of events in a specific (x,y) bin
78  HistClusterEnergyCorrectionFactor_theta_Egamma_original->Fill(eCluster, clustTheta, ClustEnCorrFactor);
79  HistClusterEnergyCorrectionFactor_theta_Egamma_count->Fill(eCluster, clustTheta);
80  HistClusterEnergyCorrectionFactor_johan_original->Fill(eCluster, clustTheta, ClustEnCorrFactor);
81  HistClusterEnergyCorrectionFactor_johan_count->Fill(eCluster, clustTheta);
82  HistClusterEnergyCorrectionFactor_Hossein_original->Fill(eCluster, clustTheta, ClustEnCorrFactor);
83  HistClusterEnergyCorrectionFactor_Hossein_count->Fill(eCluster, clustTheta);
84 
85 
86  if (j%100000 ==0)
87  cout<<"Evt. No. "<< j <<endl;
88  }
89 
90 
92 for (int i=1; i<=180; i++)//theta
93  for (int j=1; j<=360; j++)//phi
94  if (HistClusterEnergyCorrectionFactor_count->GetBinContent(i,j) != 0)
95  HistClusterEnergyCorrectionFactor->SetBinContent(i,j,(HistClusterEnergyCorrectionFactor_original->GetBinContent(i,j))/(HistClusterEnergyCorrectionFactor_count->GetBinContent(i,j)));
96  else
97  HistClusterEnergyCorrectionFactor->SetBinContent(i,j,0);
98 
99 for (int i=1; i<=500; i++)//Egamma
100  for (int j=1; j<=80; j++)//theta
101  if (HistClusterEnergyCorrectionFactor_theta_Egamma_count->GetBinContent(i,j) != 0)
102  HistClusterEnergyCorrectionFactor_theta_Egamma->SetBinContent(i,j,(HistClusterEnergyCorrectionFactor_theta_Egamma_original->GetBinContent(i,j))/(HistClusterEnergyCorrectionFactor_theta_Egamma_count->GetBinContent(i,j)));
103  else
104  HistClusterEnergyCorrectionFactor_theta_Egamma->SetBinContent(i,j,0);
105 
106 
107 for (int i=1; i<=26; i++)//Egamma
108  for (int j=1; j<=61; j++)//theta
109  if (HistClusterEnergyCorrectionFactor_johan_count->GetBinContent(i,j) != 0)
110  HistClusterEnergyCorrectionFactor_johan->SetBinContent(i,j,(HistClusterEnergyCorrectionFactor_johan_original->GetBinContent(i,j))/(HistClusterEnergyCorrectionFactor_johan_count->GetBinContent(i,j)));
111  else
112  HistClusterEnergyCorrectionFactor_johan->SetBinContent(i,j,0);
113 
114 
115 for (int i=1; i<=35; i++)//Egamma
116  for (int j=1; j<=53; j++)//theta
117  if (HistClusterEnergyCorrectionFactor_Hossein_count->GetBinContent(i,j) != 0)
118  HistClusterEnergyCorrectionFactor_Hossein->SetBinContent(i,j,(HistClusterEnergyCorrectionFactor_Hossein_original->GetBinContent(i,j))/(HistClusterEnergyCorrectionFactor_Hossein_count->GetBinContent(i,j)));
119  else
120  HistClusterEnergyCorrectionFactor_Hossein->SetBinContent(i,j,0);
121 
122 
124  TFile* file_out = new TFile(outputfile,"RECREATE");
125 
126  HISTnumberOfClusters_thresholdPassed->Write();
127  HistClusterEnergyCorrectionFactor_original->Write();
128  HistClusterEnergyCorrectionFactor_count->Write();
129  HistClusterEnergyCorrectionFactor->Write();
130  HistClusterEnergyCorrectionFactor_theta_Egamma_original->Write();
131  HistClusterEnergyCorrectionFactor_theta_Egamma_count->Write();
132  HistClusterEnergyCorrectionFactor_theta_Egamma->Write();
133  HistClusterEnergyCorrectionFactor_johan_original->Write();
134  HistClusterEnergyCorrectionFactor_johan_count->Write();
135  HistClusterEnergyCorrectionFactor_johan->Write();
136  HistClusterEnergyCorrectionFactor_Hossein_original->Write();
137  HistClusterEnergyCorrectionFactor_Hossein_count->Write();
138  HistClusterEnergyCorrectionFactor_Hossein->Write();
139 
140  file_out->Close();
141 
143  gStyle->SetPalette(1);
144  TCanvas* c1 = new TCanvas("c1","", 100, 100, 1200, 800);
145  HistClusterEnergyCorrectionFactor_johan->Draw("surf4");
146  HistClusterEnergyCorrectionFactor_johan->GetZaxis()->SetRangeUser(1,1.4);
147  HistClusterEnergyCorrectionFactor_johan->GetXaxis()->SetRangeUser(0,5.1);
149 
150  timer.Stop();
151  Double_t rtime = timer.RealTime();
152  Double_t ctime = timer.CpuTime();
153  printf("RealTime=%f seconds, CpuTime=%f seconds\n",rtime,ctime);
154  return 0;
155 }
printf
printf("RealTime=%f seconds, CpuTime=%f seconds\n", rtime, ctime)
i
Int_t i
Definition: run_full.C:25
ncounts
Int_t ncounts
Definition: anal_hit_digi_cluster_fwendcap.C:155
Double_t
Double_t
Definition: SimCompleteLinkDef.h:6
timer
TStopwatch timer
Definition: hit_dirc.C:51
c1
c1
Definition: plot_dirc.C:35
ctime
Double_t ctime
Definition: hit_dirc.C:114
rtime
Double_t rtime
Definition: hit_dirc.C:113
Generated on Wed Apr 3 2019 10:03:10 for FairRoot/PandaRoot by   doxygen 1.8.5