ddca290912907b746576a5b57c576646d907cedd
[u/mrichter/AliRoot.git] / PWGHF / correlationHF / AliDxHFECorrelation.cxx
1 // $Id$
2
3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE Project            * 
5 //* ALICE Experiment at CERN, All rights reserved.                         *
6 //*                                                                        *
7 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
8 //*                  Sedat Altinpinar <Sedat.Altinpinar@cern.ch>           *
9 //*                  Hege Erdal       <hege.erdal@gmail.com>               *
10 //*                                                                        *
11 //* Permission to use, copy, modify and distribute this software and its   *
12 //* documentation strictly for non-commercial purposes is hereby granted   *
13 //* without fee, provided that the above copyright notice appears in all   *
14 //* copies and that both the copyright notice and this permission notice   *
15 //* appear in the supporting documentation. The authors make no claims     *
16 //* about the suitability of this software for any purpose. It is          *
17 //* provided "as is" without express or implied warranty.                  *
18 //**************************************************************************
19
20 /// @file   AliDxHFECorrelation.cxx
21 /// @author Sedat Altinpinar, Hege Erdal, Matthias Richter
22 /// @date   2012-04-25
23 /// @brief  Worker class for D0-HF electron correlation
24 ///
25
26 #include "AliDxHFECorrelation.h"
27 #include "AliVParticle.h"
28 #include "AliLog.h"
29 #include "TObjArray.h"
30 #include "TH1D.h"
31 #include "TH2D.h"
32 #include "TMath.h"
33 #include "TFile.h"
34 #include "TCanvas.h"
35 #include <iostream>
36 #include <cerrno>
37 #include <memory>
38
39 using namespace std;
40
41 ClassImp(AliDxHFECorrelation)
42
43 AliDxHFECorrelation::AliDxHFECorrelation(const char* name)
44   : TNamed(name?name:"AliDxHFECorrelation", "")
45   , fHistograms(NULL)
46 {
47   // default constructor
48   // 
49   //
50
51 }
52
53 AliDxHFECorrelation::~AliDxHFECorrelation()
54 {
55   // destructor
56   //
57   //
58   if (fHistograms) {
59     delete fHistograms;
60     fHistograms=NULL;
61   }
62 }
63
64 int AliDxHFECorrelation::Init()
65 {
66   /// init class and create histograms
67   if (fHistograms) delete fHistograms;
68   fHistograms=new TObjArray;
69   if (!fHistograms) return -ENOMEM;
70   fHistograms->SetOwner(kTRUE);
71   AliInfo(Form("initializing %s ", GetName()));
72
73   // TODO: This is just a mockup to illustrate the functionality
74   // the specific histograms need to be defined
75
76   // avoid the objects to be added to the global directory 
77   bool statusAddDirectory=TH1::AddDirectoryStatus();
78   TH1::AddDirectory(false);
79   const double Pii=TMath::Pi();
80   TObjArray* a=fHistograms;
81   a->AddAt(new TH1D("hD0pT"       , "D0 pT"              ,100,0.0,10.0)   , khD0pT        );
82   a->AddAt(new TH1D("hD0Phi"      , "D0 phi"             ,180,0.0,2*Pii)  , khD0Phi       );
83   a->AddAt(new TH1D("hD0Eta"      , "D0 eta"             ,100,-10.0,10.0) , khD0Eta       );
84   a->AddAt(new TH1D("hElectronpT" , "Electron pT"        ,100,0.0,10.0)   , khElectronpT  );
85   a->AddAt(new TH1D("hElectronPhi", "Electron phi"       ,180,0.0,2*Pii)  , khElectronPhi );
86   a->AddAt(new TH1D("hElectronEta", "Electron eta"       ,100,-10.0,10.0) , khElectronEta );
87   a->AddAt(new TH1D("hDeltaPhi"   , "#Delta#Phi D0 - HFE",180,0.0,2*Pii)  , khDeltaPhi    );
88
89   TH1::AddDirectory(statusAddDirectory);
90   return 0;
91 }
92
93 int AliDxHFECorrelation::Fill(const TObjArray* candidatesD0, const TObjArray* candidatesElectron)
94 {
95   /// fill histograms from array of AliVParticle objects
96   if (!candidatesD0 || !candidatesElectron) return -EINVAL;
97   if (!fHistograms) {
98     Init();
99   }
100   if (!fHistograms) {
101     AliError("Initialisation failed, can not fill histograms");
102   }
103
104   const double Pii=TMath::Pi();
105
106   TIter itrigger(candidatesD0);
107   TObject* otrigger=NULL;
108   int ctrigger=-1;
109   while ((otrigger=itrigger())!=NULL) {
110     // loop over trigger D0 particle
111     ctrigger++;
112     AliVParticle* ptrigger=reinterpret_cast<AliVParticle*>(otrigger);
113     //    if (!ptrigger) continue;
114     ((TH1D*)fHistograms->At(khD0pT)) ->Fill(ptrigger->Pt());
115     ((TH1D*)fHistograms->At(khD0Phi))->Fill(ptrigger->Phi());
116     ((TH1D*)fHistograms->At(khD0Eta))->Fill(ptrigger->Eta());
117     // TODO: add further correlation specific cuts here, e.g acceptance
118     // which are no primarily part of the particle selection
119
120     TIter iElectron(candidatesElectron);
121     TObject* oElectron=NULL;
122     int cElectron=-1;
123     while ((oElectron=iElectron())!=NULL) {
124       // loop over electrons
125       cElectron++;
126       AliVParticle* pElectron=reinterpret_cast<AliVParticle*>(oElectron);
127       //      if (!pElectron) continue;
128       ((TH1D*)fHistograms->At(khElectronpT)) ->Fill(pElectron->Pt());
129       ((TH1D*)fHistograms->At(khElectronPhi))->Fill(pElectron->Phi());
130       ((TH1D*)fHistograms->At(khElectronEta))->Fill(pElectron->Eta());
131
132       // phi difference to trigger particle
133       Double_t DeltaPhi = ptrigger->Phi() - pElectron->Phi();
134       if (DeltaPhi<-0.5*Pii) DeltaPhi += 2*Pii;
135       if (DeltaPhi>1.5*Pii)  DeltaPhi -= 2*Pii;
136       ((TH1D*)fHistograms->At(khDeltaPhi))->Fill(DeltaPhi);
137
138     } // loop over electrons
139   } // loop over D0 trigger particle
140
141   return 0;
142 }
143
144 void AliDxHFECorrelation::Clear(Option_t * /*option*/)
145 {
146   /// overloaded from TObject: cleanup
147
148   // nothing to be done so far
149   return TObject::Clear();
150 }
151
152 void AliDxHFECorrelation::Print(Option_t */*option*/) const
153 {
154   /// overloaded from TObject: print info
155   cout << "====================================================================" << endl;
156   TObject::Print();
157   if (fHistograms) {
158     fHistograms->Print();
159   }
160 }
161
162 void AliDxHFECorrelation::Draw(Option_t */*option*/)
163 {
164   /// overloaded from TObject: draw histograms
165   if (fHistograms) {
166     TString name;
167     int canvasno=1;
168     int padno=1;
169     const char* drawoption="";
170     name.Form("%s_%d", GetName(), canvasno++);
171     TCanvas* c=new TCanvas(name);
172     c->SetWindowSize(1600,800);
173     c->SetTitle(Form("%s: particle properties", GetName()));
174     c->Divide(3,2);
175     padno=1;
176     c->cd(padno++); fHistograms->At(khD0pT)       ->Draw(drawoption);
177     c->cd(padno++); fHistograms->At(khD0Phi)      ->Draw(drawoption);
178     c->cd(padno++); fHistograms->At(khD0Eta)      ->Draw(drawoption);
179     c->cd(padno++); fHistograms->At(khDeltaPhi)   ->Draw(drawoption);
180     c->cd(padno++); fHistograms->At(khElectronpT) ->Draw(drawoption);
181     c->cd(padno++); fHistograms->At(khElectronPhi)->Draw(drawoption);
182     c->cd(padno++); fHistograms->At(khElectronEta)->Draw(drawoption);
183     c->Print(".png");
184   }
185 }
186
187 TObject* AliDxHFECorrelation::FindObject(const char *name) const
188 {
189   /// overloaded from TObject: find object by name
190   if (fHistograms) {
191     return fHistograms->FindObject(name);
192   }
193   return NULL;
194 }
195
196 TObject* AliDxHFECorrelation::FindObject(const TObject *obj) const
197 {
198   /// overloaded from TObject: find object by pointer
199   if (fHistograms) {
200     return fHistograms->FindObject(obj);
201   }
202   return NULL;
203 }
204
205 void     AliDxHFECorrelation::SaveAs(const char *filename,Option_t */*option*/) const
206 {
207   /// overloaded from TObject: save to file
208   std::auto_ptr<TFile> output(TFile::Open(filename, "RECREATE"));
209   if (!output.get() || output->IsZombie()) {
210     AliError(Form("can not open file %s from writing", filename));
211     return;
212   }
213   output->cd();
214   if (fHistograms) fHistograms->Write();
215   output->Close();
216 }
217
218 AliDxHFECorrelation& AliDxHFECorrelation::operator+=(const AliDxHFECorrelation& other)
219 {
220   /// add histograms from another instance
221   if (!fHistograms || !other.fHistograms) return *this;
222   
223   for (int i=0; i<kNofHistograms; i++) {
224     if (fHistograms->At(i)==NULL || other.fHistograms->At(i)==NULL) continue;
225     TH1* target=reinterpret_cast<TH1*>(fHistograms->At(i));
226     TH1* source=reinterpret_cast<TH1*>(other.fHistograms->At(i));
227     if (!target || !source) continue;
228     TString name(fHistograms->At(i)->GetName());
229     if (name.CompareTo(target->GetName())!=0) {
230       AliWarning(Form("skipping incompatible objects at position %d: %s vs %s", i, source->GetName(), target->GetName()));
231       continue;
232     }
233     if (source->IsA()!=target->IsA()) {
234       AliWarning(Form("skipping incompatible classes at position %d: %s vs %s", i, source->ClassName(), target->ClassName()));
235       continue;
236     }
237     target->Add(source);
238   }
239   return *this;
240 }