hlit-ana  0.0.0
 All Classes Functions Variables Pages
HlitAnalysisSelectorIO.cxx
1 #include <TFile.h>
2 #include <TTree.h>
3 #include <TCanvas.h>
4 #include "HlitEvent.h"
5 
6 #include "HlitAnalysisSelectorIO.h"
7 
9 ClassImp(HlitAnalysisSelectorIO);
11 
13  : TSelector(), fChain(0), fEvent(0), fHistPx(0), fHistPy(0), fHistPz(0),
14  fHistPxPy(0), fOutputFile("hlit-ana.root") {
18 }
19 
24  delete fEvent;
25 }
26 
27 void HlitAnalysisSelectorIO::Init(TTree *tree) {
31  if (!tree)
32  return;
33  fChain = tree;
34  fChain->SetBranchAddress("HlitEvent", &fEvent);
35 }
36 
41  if (fChain->GetCurrentFile())
42  Printf("File : %s", fChain->GetCurrentFile()->GetName());
43  return kTRUE;
44 }
45 
46 void HlitAnalysisSelectorIO::Begin(TTree * /*tree*/) {
50 }
51 
52 void HlitAnalysisSelectorIO::SlaveBegin(TTree * /*tree*/) {
56 
57  TString option = GetOption();
58 
59  fHistPx = new TH1D("hPx", "p_{x}", 100, -5.0, 5.0);
60  fOutput->Add(fHistPx);
61  fHistPy = new TH1D("hPy", "p_{y}", 100, -5.0, 5.0);
62  fOutput->Add(fHistPy);
63  fHistPz = new TH1D("hPz", "p_{z}", 100, -5.0, 5.0);
64  fOutput->Add(fHistPz);
65 
66  fHistPxPy =
67  new TH2D("hPxPy", "p_{x} vs p_{y}", 100, -5.0, 5.0, 100, -5.0, 5.0);
68  fOutput->Add(fHistPxPy);
69 }
70 
71 Bool_t HlitAnalysisSelectorIO::Process(Long64_t entry) {
75 
76  GetEntry(entry);
77 
78  HlitTrack *t;
79  for (Int_t iTrack = 0; iTrack < fEvent->GetNTrack(); iTrack++) {
80  t = fEvent->GetTrack(iTrack);
81  fHistPx->Fill(t->GetPx());
82  fHistPy->Fill(t->GetPy());
83  fHistPz->Fill(t->GetPz());
84  fHistPxPy->Fill(t->GetPx(), t->GetPy());
85  }
86 
87  return kTRUE;
88 }
89 
94 }
95 
100 
101  fHistPx = dynamic_cast<TH1D *>(fOutput->FindObject("hPx"));
102  fHistPy = dynamic_cast<TH1D *>(fOutput->FindObject("hPy"));
103  fHistPz = dynamic_cast<TH1D *>(fOutput->FindObject("hPz"));
104  fHistPxPy = dynamic_cast<TH2D *>(fOutput->FindObject("hPxPy"));
105 
106  if ((fHistPx) && (fHistPy) && (fHistPz) && (fHistPxPy)) {
107  TCanvas *c = new TCanvas("cP", "Momentum distributions");
108  c->Divide(2, 2);
109  c->cd(1);
110  fHistPx->Draw();
111  c->cd(2);
112  fHistPy->Draw();
113  c->cd(3);
114  fHistPz->Draw();
115  c->cd(4);
116  fHistPxPy->Draw();
117  }
118 
119  TFile *output = TFile::Open(fOutputFile.Data(),"RECREATE");
120  fOutput->Write();
121  output->Close();
122 }
123 
124 Int_t HlitAnalysisSelectorIO::GetEntry(Long64_t entry, Int_t getall) {
128  return fChain ? fChain->GetTree()->GetEntry(entry, getall) : 0;
129 }
Track object.
Definition: HlitTrack.h:13
TH2D * fHistPxPy
px vs py distribution
Double_t GetPz() const
Momentum z component.
Definition: HlitTrack.h:46
Double_t GetPy() const
Momentum y component.
Definition: HlitTrack.h:45
TString fOutputFile
Output file name.
TTree * fChain
Pointer to the analyzed TTree or TChain.
virtual void Init(TTree *tree)
HlitTrack * GetTrack(Long64_t id)
Definition: HlitEvent.h:76
virtual void Begin(TTree *)
virtual Bool_t Process(Long64_t entry)
HlitEvent * fEvent
Current Event.
TH1D * fHistPz
pz distribution
Double_t GetPx() const
Momentum x component.
Definition: HlitTrack.h:44
TH1D * fHistPx
px distribution
TH1D * fHistPy
py distribution
Long64_t GetNTrack() const
Definition: HlitEvent.h:75
virtual Int_t GetEntry(Long64_t entry, Int_t getall=0)
virtual void SlaveBegin(TTree *tree)