Add a draw class for the CORRFW (produces a warning, will be fixed
[u/mrichter/AliRoot.git] / PWG3 / dielectron / AliDielectron.cxx
1 /*************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3 *                                                                        *
4 * Author: The ALICE Off-line Project.                                    *
5 * Contributors are mentioned in the code where appropriate.              *
6 *                                                                        *
7 * Permission to use, copy, modify and distribute this software and its   *
8 * documentation strictly for non-commercial purposes is hereby granted   *
9 * without fee, provided that the above copyright notice appears in all   *
10 * copies and that both the copyright notice and this permission notice   *
11 * appear in the supporting documentation. The authors make no claims     *
12 * about the suitability of this software for any purpose. It is          *
13 * provided "as is" without express or implied warranty.                  *
14 **************************************************************************/
15
16 ///////////////////////////////////////////////////////////////////////////
17 //                Dielectron Analysis Main class                         //
18 //                                                                       //
19 /*
20 Framework to perform event selectoin, single track selection and track pair
21 selection.
22
23 Convention for the signs of the pair in fPairCandidates:
24 The names are available via the function PairClassName(Int_t i)
25
26 0: ev1+ ev1+  (same event like sign +)
27 1: ev1+ ev1-  (same event unlike sign)
28 2: ev1- ev1-  (same event like sign -)
29
30 3: ev1+ ev2+  (mixed event like sign +)
31 4: ev1- ev2+  (mixed event unlike sign -+)
32 6: ev1+ ev2-  (mixed event unlike sign +-)
33 7: ev1- ev2-  (mixed event like sign -)
34
35 5: ev2+ ev2+  (same event like sign +)
36 8: ev2+ ev2-  (same event unlike sign)
37 9: ev2- ev2-  (same event like sign -)
38
39
40
41 */
42 //                                                                       //
43 ///////////////////////////////////////////////////////////////////////////
44
45 #include <TString.h>
46 #include <TList.h>
47 #include <TMath.h>
48
49 #include <AliESDEvent.h>
50 #include <AliESDtrack.h>
51
52 #include <AliVEvent.h>
53 #include <AliVParticle.h>
54 #include <AliVTrack.h>
55 #include "AliDielectronPair.h"
56 #include "AliDielectronHistos.h"
57 #include "AliDielectronCF.h"
58 #include "AliDielectronMC.h"
59 #include "AliDielectronVarManager.h"
60 #include "AliDielectron.h"
61
62 ClassImp(AliDielectron)
63
64 const char* AliDielectron::fgkTrackClassNames[4] = {
65   "ev1+",
66   "ev1-",
67   "ev2+",
68   "ev2-"
69 };
70
71 const char* AliDielectron::fgkPairClassNames[10] = {
72   "ev1+_ev1+",
73   "ev1+_ev1-",
74   "ev1-_ev1-",
75   "ev1+_ev2+",
76   "ev1-_ev2+",
77   "ev2+_ev2+",
78   "ev1+_ev2-",
79   "ev1-_ev2-",
80   "ev2+_ev2-",
81   "ev2-_ev2-"
82 };
83
84 //________________________________________________________________
85 AliDielectron::AliDielectron() :
86   TNamed("AliDielectron","AliDielectron"),
87   fEventFilter("EventFilter"),
88   fTrackFilter("TrackFilter"),
89   fPairFilter("PairFilter"),
90   fPdgMother(443),
91   fHistos(0x0),
92   fPairCandidates(new TObjArray(10)),
93   fCfManagerPair(0x0)
94 {
95   //
96   // Default constructor
97   //
98
99 }
100
101 //________________________________________________________________
102 AliDielectron::AliDielectron(const char* name, const char* title) :
103   TNamed(name,title),
104   fEventFilter("EventFilter"),
105   fTrackFilter("TrackFilter"),
106   fPairFilter("PairFilter"),
107   fPdgMother(443),
108   fHistos(0x0),
109   fPairCandidates(new TObjArray(10)),
110   fCfManagerPair(0x0)
111 {
112   //
113   // Named constructor
114   //
115   
116 }
117
118 //________________________________________________________________
119 AliDielectron::~AliDielectron()
120 {
121   //
122   // Default destructor
123   //
124   if (fHistos) delete fHistos;
125   if (fPairCandidates) delete fPairCandidates;
126 }
127
128 //________________________________________________________________
129 void AliDielectron::Init()
130 {
131   //
132   // Initialise objects
133   //
134   if (fCfManagerPair) fCfManagerPair->InitialiseContainer(fPairFilter);
135   AliDielectronVarManager::InitESDpid(); //currently this will take the values for MC
136 }
137
138 //________________________________________________________________
139 void AliDielectron::Process(AliVEvent *ev1, AliVEvent *ev2)
140 {
141   //
142   // Process the events
143   //
144
145   AliDielectronVarManager::SetEvent(0x0);
146    
147   //in case we have MC load the MC event and process the MC particles
148   if (AliDielectronMC::Instance()->ConnectMCEvent()) ProcessMC();
149   
150   //if candidate array doesn't exist, create it
151   if (!fPairCandidates->UncheckedAt(0)) {
152     InitPairCandidateArrays();
153   } else {
154     ClearArrays();
155   }
156
157   //mask used to require that all cuts are fulfilled
158   UInt_t selectedMask=(1<<fEventFilter.GetCuts()->GetEntries())-1;
159
160   //apply event cuts
161     if ((ev1&&fEventFilter.IsSelected(ev1)!=selectedMask) ||
162         (ev2&&fEventFilter.IsSelected(ev2)!=selectedMask)) return;
163   
164   AliDielectronVarManager::SetEvent(ev1);
165   
166   //fill track arrays for the first event
167   if (ev1) FillTrackArrays(ev1);
168
169   //fill track arrays for the second event
170   if (ev2) FillTrackArrays(ev2,1);
171
172   // create pairs and fill pair candidate arrays
173   for (Int_t itrackArr1=0; itrackArr1<4; ++itrackArr1){
174     for (Int_t itrackArr2=itrackArr1; itrackArr2<4; ++itrackArr2){
175       FillPairArrays(itrackArr1, itrackArr2);
176     }
177   }
178
179   //in case there is a histogram manager, fill the QA histograms
180   if (fHistos) FillHistograms(ev1);
181   
182 }
183
184 //________________________________________________________________
185 void AliDielectron::ProcessMC()
186 {
187   //
188   // Process the MC data
189   //
190
191   //loop over all MC data and Fill the CF container if it exist
192   if (!fCfManagerPair) return;
193   fCfManagerPair->SetPdgMother(fPdgMother);
194   
195   AliDielectronMC *dieMC=AliDielectronMC::Instance();
196   for (Int_t ipart=0; ipart<dieMC->GetNMCTracks();++ipart){
197     //TODO: MC truth cut properly!!!
198     AliVParticle *mcPart=dieMC->GetMCTrackFromMCEvent(ipart);
199     if (!dieMC->IsMCMotherToEE(mcPart, fPdgMother)) continue;
200     fCfManagerPair->FillMC(mcPart);
201   }
202 }
203
204 //________________________________________________________________
205 void AliDielectron::FillHistograms(const AliVEvent *ev)
206 {
207   //
208   // Fill Histogram information for tracks and pairs
209   //
210   
211   TString  className;
212   Double_t values[AliDielectronVarManager::kNMaxValues];
213   //Fill event information
214   AliDielectronVarManager::Fill(ev, values);
215   fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, values);
216   
217   //Fill track information, separately for the track array candidates
218   for (Int_t i=0; i<4; ++i){
219     className.Form("Track_%s",fgkTrackClassNames[i]);
220     Int_t ntracks=fTracks[i].GetEntriesFast();
221     for (Int_t itrack=0; itrack<ntracks; ++itrack){
222       AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
223       fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
224     }
225   }
226
227   //Fill Pair information, separately for all pair candidate arrays
228   for (Int_t i=0; i<10; ++i){
229     className.Form("Pair_%s",fgkPairClassNames[i]);
230     Int_t ntracks=PairArray(i)->GetEntriesFast();
231     for (Int_t ipair=0; ipair<ntracks; ++ipair){
232       AliDielectronVarManager::Fill(PairArray(i)->UncheckedAt(ipair), values);
233       fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
234     }
235   }
236   
237 }
238
239 //________________________________________________________________
240 void AliDielectron::FillTrackArrays(AliVEvent * const ev, Int_t eventNr)
241 {
242   //
243   // select tracks and fill track candidate arrays
244   // eventNr = 0: First  event, use track arrays 0 and 1
245   // eventNr = 1: Second event, use track arrays 2 and 3
246   //
247   
248   Int_t ntracks=ev->GetNumberOfTracks();
249   UInt_t selectedMask=(1<<fTrackFilter.GetCuts()->GetEntries())-1;
250   for (Int_t itrack=0; itrack<ntracks; ++itrack){
251     //get particle
252     AliVParticle *particle=ev->GetTrack(itrack);
253     //TODO: temporary solution, perhaps think about a better implementation
254     //      This is needed to use AliESDpidCuts, which relies on the ESD event
255     //      is set as a AliESDtrack attribute... somehow ugly!
256     if (ev->IsA()==AliESDEvent::Class()){
257       AliESDtrack *track=static_cast<AliESDtrack*>(particle);
258       track->SetESDEvent(static_cast<AliESDEvent*>(ev)); //only in trunk...
259     }
260     
261     //apply track cuts
262     if (fTrackFilter.IsSelected(particle)!=selectedMask) continue;
263     
264     //fill selected particle into the corresponding track arrays
265     Short_t charge=particle->Charge();
266     if (charge>0)      fTracks[eventNr*2].Add(particle);
267     else if (charge<0) fTracks[eventNr*2+1].Add(particle);
268   }
269 }
270
271 //________________________________________________________________
272 void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2) {
273   //
274   // select pairs and fill pair candidate arrays
275   //
276   Int_t pairIndex=GetPairIndex(arr1,arr2);
277
278   Int_t ntrack1=fTracks[arr1].GetEntriesFast();
279   Int_t ntrack2=fTracks[arr2].GetEntriesFast();
280
281   AliDielectronPair *candidate=new AliDielectronPair;
282
283   UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
284   
285   for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
286     Int_t end=ntrack2;
287     if (arr1==arr2) end=itrack1;
288     for (Int_t itrack2=0; itrack2<end; ++itrack2){
289       //create the pair TODO: change hardcoded PID of tracks
290       candidate->SetTracks(static_cast<AliVTrack*>(fTracks[arr1].UncheckedAt(itrack1)), 11,
291                            static_cast<AliVTrack*>(fTracks[arr2].UncheckedAt(itrack2)), 11);
292       candidate->SetType(pairIndex);
293
294       //pair cuts
295       UInt_t cutMask=fPairFilter.IsSelected(candidate);
296       
297       //CF manager for the pair
298       if (fCfManagerPair) fCfManagerPair->Fill(cutMask,candidate);
299
300       //apply cut
301       if (cutMask!=selectedMask) continue;
302
303       //add the candidate to the candidate array 
304       PairArray(pairIndex)->Add(candidate);
305       //get a new candidate
306       candidate=new AliDielectronPair;
307     }
308   }
309   //delete the surplus candidate
310   delete candidate;
311 }
312
313
314
315