]>
Commit | Line | Data |
---|---|---|
b1f673b5 | 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 | // // | |
18 | // Basic Analysis Task // | |
19 | // // | |
20 | /////////////////////////////////////////////////////////////////////////// | |
21 | ||
22 | #include <TChain.h> | |
23 | #include <TH1D.h> | |
24 | ||
25 | #include <AliCFContainer.h> | |
26 | #include <AliInputEventHandler.h> | |
27 | #include <AliESDInputHandler.h> | |
28 | #include <AliAODInputHandler.h> | |
29 | #include <AliAnalysisManager.h> | |
30 | #include <AliVEvent.h> | |
31 | #include <AliTriggerAnalysis.h> | |
32 | #include <AliPIDResponse.h> | |
33 | #include <AliTPCPIDResponse.h> | |
34 | ||
35 | #include "AliDielectron.h" | |
36 | #include "AliDielectronHistos.h" | |
37 | #include "AliDielectronCF.h" | |
38 | #include "AliDielectronMC.h" | |
39 | #include "AliDielectronMixingHandler.h" | |
40 | #include "AliAnalysisTaskMultiDielectronTG.h" | |
41 | ||
42 | #include "AliESDtrack.h" | |
43 | #include "AliESDtrackCuts.h" | |
44 | ||
45 | #include "AliLog.h" | |
46 | ||
47 | #include <vector> | |
48 | #include <deque> | |
49 | #include <cstdlib> | |
50 | #include <TRandom.h> | |
51 | ||
7329335b | 52 | using std::vector; |
53 | using std::deque; | |
54 | ||
55 | const char* kPairClassNamesTG[7] = { | |
b1f673b5 | 56 | "unlike", |
57 | "like_pp", | |
58 | "like_ee", | |
59 | "mixunlike_pe", | |
60 | "mixunlike_ep", | |
61 | "mixlike_pp", | |
62 | "mixlike_ee" | |
63 | }; | |
64 | ||
65 | ||
66 | ClassImp(AliDielectronSingleTG) | |
67 | ClassImp(AliAnalysisTaskMultiDielectronTG) | |
68 | ||
69 | //_________________________________________________________________________________ | |
70 | AliAnalysisTaskMultiDielectronTG::AliAnalysisTaskMultiDielectronTG() : | |
71 | AliAnalysisTaskSE(), | |
72 | fListDielectron(), | |
73 | fListHistos(), | |
74 | fListCF(), | |
7329335b | 75 | fQAElectron(), |
b1f673b5 | 76 | fSelectPhysics(kFALSE), |
77 | fTriggerMask(AliVEvent::kMB), | |
78 | fExcludeTriggerMask(0), | |
79 | fTriggerOnV0AND(kFALSE), | |
80 | fRejectPileup(kFALSE), | |
81 | fTriggerLogic(kAny), | |
82 | fTriggerAnalysis(0x0), | |
83 | fEventFilter(0x0), | |
84 | fCutsMother(0x0), | |
85 | fEventStat(0x0), | |
86 | fEvent(0x0), | |
87 | fdEdXvsPt(0x0), | |
88 | fdEdXnSigmaElecvsPt(0x0), | |
89 | fdEdXvsPtTOF(0x0), | |
90 | fdEdXnSigmaElecvsPtTOF(0x0), | |
91 | fTOFbetavsPt(0x0), | |
92 | fTOFnSigmaElecvsPt(0x0), | |
7329335b | 93 | fNCrossedRowsTPC(0x0), |
94 | fChi2ClusTPC(0x0), | |
95 | fRatioCrossClusTPC(0x0), | |
96 | fVem(0x0), | |
97 | fVep(0x0), | |
98 | fVemtmp(0x0), | |
99 | fVeptmp(0x0), | |
100 | fdconvphiv(acos(-1.0)), | |
805bd069 | 101 | fdop(0), |
7329335b | 102 | fbz(0), |
620b34a5 | 103 | fdv0mixing(kTRUE) |
b1f673b5 | 104 | { |
105 | // | |
106 | // Constructor | |
107 | // | |
108 | } | |
109 | ||
110 | //_________________________________________________________________________________ | |
111 | AliAnalysisTaskMultiDielectronTG::AliAnalysisTaskMultiDielectronTG(const char *name) : | |
112 | AliAnalysisTaskSE(name), | |
113 | fListDielectron(), | |
114 | fListHistos(), | |
115 | fListCF(), | |
7329335b | 116 | fQAElectron(), |
b1f673b5 | 117 | fSelectPhysics(kFALSE), |
118 | fTriggerMask(AliVEvent::kMB), | |
119 | fExcludeTriggerMask(0), | |
120 | fTriggerOnV0AND(kFALSE), | |
121 | fRejectPileup(kFALSE), | |
122 | fTriggerLogic(kAny), | |
123 | fTriggerAnalysis(0x0), | |
124 | fEventFilter(0x0), | |
125 | fCutsMother(0x0), | |
126 | fEventStat(0x0), | |
127 | fEvent(0x0), | |
128 | fdEdXvsPt(0x0), | |
129 | fdEdXnSigmaElecvsPt(0x0), | |
130 | fdEdXvsPtTOF(0x0), | |
131 | fdEdXnSigmaElecvsPtTOF(0x0), | |
132 | fTOFbetavsPt(0x0), | |
133 | fTOFnSigmaElecvsPt(0x0), | |
7329335b | 134 | fNCrossedRowsTPC(0x0), |
135 | fChi2ClusTPC(0x0), | |
136 | fRatioCrossClusTPC(0x0), | |
137 | fVem(0x0), | |
138 | fVep(0x0), | |
139 | fVemtmp(0x0), | |
140 | fVeptmp(0x0), | |
141 | fdconvphiv(acos(-1.0)), | |
805bd069 | 142 | fdop(0), |
7329335b | 143 | fbz(0), |
620b34a5 | 144 | fdv0mixing(kTRUE) |
b1f673b5 | 145 | { |
146 | // | |
147 | // Constructor | |
148 | // | |
149 | DefineInput(0,TChain::Class()); | |
150 | DefineOutput(1, TList::Class()); | |
151 | DefineOutput(2, TList::Class()); | |
152 | DefineOutput(3, TH1D::Class()); | |
153 | fListHistos.SetName("Dielectron_Histos_Multi"); | |
154 | fListCF.SetName("Dielectron_CF_Multi"); | |
155 | fListDielectron.SetOwner(); | |
156 | fListHistos.SetOwner(); | |
157 | fListCF.SetOwner(); | |
158 | ||
159 | /////////////// | |
7329335b | 160 | for(int i=0;i<fgkNDIE; i++){ |
161 | for(int j=0;j<fgkNZBIN;j++){ | |
162 | for(int k=0;k<fgkNCENT;k++){ | |
163 | for(int l=0; l<fgkNRPBIN; l++){ | |
164 | fibuf[i][j][k][l] = 0; | |
165 | fpoolp[i][j][k][l].clear(); | |
166 | fpoolm[i][j][k][l].clear(); | |
167 | for(int ib=0;ib<fgkNBUF; ib++){ | |
168 | fvep[ib][i][j][k][l].clear(); | |
169 | fvem[ib][i][j][k][l].clear(); | |
b1f673b5 | 170 | } |
171 | } | |
172 | } | |
173 | } | |
174 | } | |
175 | ||
176 | ||
177 | ||
178 | } | |
179 | ||
180 | //_________________________________________________________________________________ | |
181 | AliAnalysisTaskMultiDielectronTG::~AliAnalysisTaskMultiDielectronTG() | |
182 | { | |
183 | // | |
184 | // Destructor | |
185 | // | |
186 | ||
187 | //histograms and CF are owned by the dielectron framework. | |
188 | //however they are streamed to file, so in the first place the | |
189 | //lists need to be owner... | |
190 | fListHistos.SetOwner(kFALSE); | |
191 | fListCF.SetOwner(kFALSE); | |
192 | ||
7329335b | 193 | for(int i=0;i<fgkNDIE; i++){ |
194 | for(int j=0;j<fgkNZBIN;j++){ | |
195 | for(int k=0;k<fgkNCENT;k++){ | |
196 | for(int l=0; l<fgkNRPBIN; l++){ | |
197 | fibuf[i][j][k][l] = 0; | |
198 | fpoolp[i][j][k][l].clear(); | |
199 | fpoolm[i][j][k][l].clear(); | |
200 | for(int ib=0;ib<fgkNBUF; ib++){ | |
201 | fvep[ib][i][j][k][l].clear(); | |
202 | fvem[ib][i][j][k][l].clear(); | |
b1f673b5 | 203 | } |
204 | } | |
205 | } | |
206 | } | |
207 | } | |
208 | } | |
209 | //_________________________________________________________________________________ | |
210 | void AliAnalysisTaskMultiDielectronTG::UserCreateOutputObjects() | |
211 | { | |
212 | // | |
213 | // Add all histogram manager histogram lists to the output TList | |
214 | // | |
215 | ||
216 | if (!fListHistos.IsEmpty()||!fListCF.IsEmpty()) return; //already initialised | |
217 | ||
218 | // AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager(); | |
219 | // Bool_t isESD=man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class(); | |
220 | // Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class(); | |
221 | ||
222 | TIter nextDie(&fListDielectron); | |
223 | AliDielectron *die=0; | |
224 | while ( (die=static_cast<AliDielectron*>(nextDie())) ){ | |
225 | die->Init(); | |
226 | if (die->GetHistogramList()) fListHistos.Add(const_cast<THashList*>(die->GetHistogramList())); | |
227 | if (die->GetCFManagerPair()) fListCF.Add(const_cast<AliCFContainer*>(die->GetCFManagerPair()->GetContainer())); | |
228 | } | |
229 | ||
230 | Int_t cuts=fListDielectron.GetEntries(); | |
231 | Int_t nbins=kNbinsEvent+2*cuts; | |
232 | if (!fEventStat){ | |
233 | fEventStat=new TH1D("hEventStat","Event statistics",nbins,0,nbins); | |
234 | fEventStat->GetXaxis()->SetBinLabel(1,"Before Phys. Sel."); | |
235 | fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel."); | |
236 | ||
237 | //default names | |
238 | fEventStat->GetXaxis()->SetBinLabel(3,"Bin3 not used"); | |
239 | fEventStat->GetXaxis()->SetBinLabel(4,"Bin4 not used"); | |
240 | fEventStat->GetXaxis()->SetBinLabel(5,"Bin5 not used"); | |
241 | ||
242 | if(fTriggerOnV0AND) fEventStat->GetXaxis()->SetBinLabel(3,"V0and triggers"); | |
243 | if (fEventFilter) fEventStat->GetXaxis()->SetBinLabel(4,"After Event Filter"); | |
244 | if (fRejectPileup) fEventStat->GetXaxis()->SetBinLabel(5,"After Pileup rejection"); | |
245 | ||
246 | for (Int_t i=0; i<cuts; ++i){ | |
247 | fEventStat->GetXaxis()->SetBinLabel((kNbinsEvent+1)+2*i,Form("#splitline{1 candidate}{%s}",fListDielectron.At(i)->GetName())); | |
248 | fEventStat->GetXaxis()->SetBinLabel((kNbinsEvent+2)+2*i,Form("#splitline{With >1 candidate}{%s}",fListDielectron.At(i)->GetName())); | |
249 | } | |
250 | } | |
251 | ||
252 | if (!fTriggerAnalysis) fTriggerAnalysis=new AliTriggerAnalysis; | |
253 | fTriggerAnalysis->EnableHistograms(); | |
254 | fTriggerAnalysis->SetAnalyzeMC(AliDielectronMC::Instance()->HasMC()); | |
255 | ||
256 | ||
257 | ||
258 | int nbinx=400; | |
7329335b | 259 | float maxx=20; |
260 | float minx=0.2; | |
261 | float binw = (TMath::Log(maxx)-TMath::Log(minx))/nbinx; | |
b1f673b5 | 262 | double xbin[401]; |
263 | for(int ii=0;ii<nbinx+1;ii++){ | |
7329335b | 264 | xbin[ii] = TMath::Exp(TMath::Log(minx) + 0.5*binw+binw*ii); |
b1f673b5 | 265 | } |
266 | ||
267 | ||
7329335b | 268 | fQAElectron = new TList(); |
269 | fQAElectron->SetName("QAElectron"); | |
270 | fQAElectron->SetOwner(); | |
b1f673b5 | 271 | |
272 | ||
273 | fEvent = new TH1D("Event","centrality", 100,0,100); | |
7329335b | 274 | fQAElectron->Add(fEvent); |
b1f673b5 | 275 | fdEdXvsPt = new TH2D("dEdXvsPt","dE/dX vs. PT of TPC", nbinx, xbin, 2000,0,200); |
7329335b | 276 | fQAElectron->Add(fdEdXvsPt); |
b1f673b5 | 277 | fdEdXnSigmaElecvsPt = new TH2D("fdEdXnSigmaElecvsPt"," dE/dX normalized to electron vs. pT of TPC", |
278 | nbinx, xbin, 2000, -10, 10); | |
7329335b | 279 | fQAElectron->Add(fdEdXnSigmaElecvsPt); |
b1f673b5 | 280 | |
281 | fdEdXvsPtTOF = new TH2D("dEdXvsPtTOF","dE/dX vs. PT of TPC", nbinx, xbin, 2000,0,200); | |
7329335b | 282 | fQAElectron->Add(fdEdXvsPtTOF); |
b1f673b5 | 283 | fdEdXnSigmaElecvsPtTOF = new TH2D("fdEdXnSigmaElecvsPtTOF"," dE/dX normalized to electron vs. pT of TPC", |
284 | nbinx, xbin, 2000, -10, 10); | |
7329335b | 285 | fQAElectron->Add(fdEdXnSigmaElecvsPtTOF); |
b1f673b5 | 286 | |
287 | ||
288 | ||
289 | fTOFbetavsPt = new TH2D("fTOFbetavsPt","TOF beta vs. p", 400, 0, 20, 1200, 0, 1.2); | |
7329335b | 290 | fQAElectron->Add(fTOFbetavsPt); |
b1f673b5 | 291 | fTOFnSigmaElecvsPt = new TH2D("fTOFnSigmaElecvsPt","TOF nsigma for electron", 400, 0, 20, 2000, -10, 10); |
7329335b | 292 | fQAElectron->Add(fTOFnSigmaElecvsPt); |
b1f673b5 | 293 | |
7329335b | 294 | fNCrossedRowsTPC = new TH2F("fNCrossedRowsTPC", "TPC nCrossed Rows vs. pT", 200, 0, 20, 200, 0, 200); |
295 | fQAElectron->Add(fNCrossedRowsTPC); | |
296 | fChi2ClusTPC = new TH2F("fChi2ClusTPC", "hChi2ClusTPC vs. pT", 200, 0, 20, 200, 0, 10); | |
297 | fQAElectron->Add(fChi2ClusTPC); | |
b1f673b5 | 298 | |
7329335b | 299 | fRatioCrossClusTPC = new TH2F("fRatioCrossClusTPC", "hRatioCrossClusTPC vs. pT", 200, 0, 20, 200, 0, 10); |
300 | fQAElectron->Add(fRatioCrossClusTPC); | |
b1f673b5 | 301 | |
7329335b | 302 | fListHistos.Add(fQAElectron); |
b1f673b5 | 303 | |
304 | fListHistos.SetOwner(); | |
305 | ||
306 | PostData(1, &fListHistos); | |
307 | PostData(2, &fListCF); | |
308 | PostData(3, fEventStat); | |
309 | ||
310 | fCutsMother = new AliESDtrackCuts; | |
311 | fCutsMother->SetDCAToVertex2D(kTRUE); | |
312 | fCutsMother->SetMaxDCAToVertexZ(3.0); | |
313 | fCutsMother->SetMaxDCAToVertexXY(1.0); | |
314 | fCutsMother->SetPtRange( 0.05 , 200.0); | |
315 | fCutsMother->SetEtaRange( -0.84 , 0.84 ); | |
316 | fCutsMother->SetAcceptKinkDaughters(kFALSE); | |
317 | fCutsMother->SetRequireITSRefit(kTRUE); | |
318 | fCutsMother->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst); | |
319 | fCutsMother->SetMinNClustersITS(3); | |
320 | fCutsMother->SetRequireTPCRefit(kTRUE); | |
321 | ||
322 | ||
805bd069 | 323 | AliInfo(Form("PairCutType = %d %d %d %d %d", |
324 | fRejectPairFlag[0], | |
325 | fRejectPairFlag[1], | |
326 | fRejectPairFlag[2], | |
327 | fRejectPairFlag[3], | |
328 | fRejectPairFlag[4])); | |
329 | ||
b1f673b5 | 330 | } |
331 | ||
332 | //_________________________________________________________________________________ | |
333 | void AliAnalysisTaskMultiDielectronTG::UserExec(Option_t *) | |
334 | { | |
335 | // | |
336 | // Main loop. Called for every event | |
337 | // | |
338 | ||
339 | if (fListHistos.IsEmpty()&&fListCF.IsEmpty()) return; | |
340 | ||
341 | AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager(); | |
342 | Bool_t isESD=man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class(); | |
343 | Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class(); | |
344 | ||
345 | AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler()); | |
346 | if (!inputHandler) return; | |
347 | ||
348 | // AliPIDResponse *pidRes=inputHandler->GetPIDResponse(); | |
349 | if ( inputHandler->GetPIDResponse() ){ | |
350 | // for the 2.76 pass2 MC private train. Together with a sigma shift of -0.169 | |
351 | // pidRes->GetTPCResponse().SetSigma(4.637e-3,2.41332105409873257e+04); | |
352 | AliDielectronVarManager::SetPIDResponse( inputHandler->GetPIDResponse() ); | |
353 | } else { | |
354 | AliFatal("This task needs the PID response attached to the input event handler!"); | |
355 | } | |
356 | ||
357 | // Was event selected ? | |
358 | ULong64_t isSelected = AliVEvent::kAny; | |
359 | Bool_t isRejected = kFALSE; | |
360 | if( fSelectPhysics && inputHandler){ | |
361 | if((isESD && inputHandler->GetEventSelection()) || isAOD){ | |
362 | isSelected = inputHandler->IsEventSelected(); | |
363 | if (fExcludeTriggerMask && (isSelected&fExcludeTriggerMask)) isRejected=kTRUE; | |
364 | if (fTriggerLogic==kAny) isSelected&=fTriggerMask; | |
365 | else if (fTriggerLogic==kExact) isSelected=((isSelected&fTriggerMask)==fTriggerMask); | |
366 | } | |
367 | } | |
368 | ||
369 | ||
370 | //Before physics selection | |
371 | fEventStat->Fill(kAllEvents); | |
372 | if (isSelected==0||isRejected) { | |
373 | PostData(3,fEventStat); | |
374 | return; | |
375 | } | |
376 | //after physics selection | |
377 | fEventStat->Fill(kSelectedEvents); | |
378 | ||
379 | //V0and | |
380 | if(fTriggerOnV0AND){ | |
805bd069 | 381 | if(isESD){if (!fTriggerAnalysis->IsOfflineTriggerFired(static_cast<AliESDEvent*>(InputEvent()), AliTriggerAnalysis::kV0AND)) |
b1f673b5 | 382 | return;} |
805bd069 | 383 | if(isAOD){if(!((static_cast<AliAODEvent*>(InputEvent()))->GetVZEROData()->GetV0ADecision() == AliVVZERO::kV0BB && |
384 | (static_cast<AliAODEvent*>(InputEvent()))->GetVZEROData()->GetV0CDecision() == AliVVZERO::kV0BB) ) | |
385 | return;} | |
386 | } | |
b1f673b5 | 387 | |
388 | ||
389 | fEventStat->Fill(kV0andEvents); | |
390 | ||
391 | //Fill Event histograms before the event filter | |
392 | TIter nextDie(&fListDielectron); | |
393 | AliDielectron *die=0; | |
394 | Double_t values[AliDielectronVarManager::kNMaxValues]={0}; | |
395 | Double_t valuesMC[AliDielectronVarManager::kNMaxValues]={0}; | |
396 | AliDielectronVarManager::SetEvent(InputEvent()); | |
397 | AliDielectronVarManager::Fill(InputEvent(),values); | |
398 | AliDielectronVarManager::Fill(InputEvent(),valuesMC); | |
399 | Bool_t hasMC=AliDielectronMC::Instance()->HasMC(); | |
400 | if (hasMC) { | |
401 | if (AliDielectronMC::Instance()->ConnectMCEvent()) | |
402 | AliDielectronVarManager::Fill(AliDielectronMC::Instance()->GetMCEvent(),valuesMC); | |
403 | } | |
404 | ||
405 | while ( (die=static_cast<AliDielectron*>(nextDie())) ){ | |
406 | AliDielectronHistos *h=die->GetHistoManager(); | |
407 | if (h){ | |
408 | if (h->GetHistogramList()->FindObject("Event_noCuts")) | |
409 | h->FillClass("Event_noCuts",AliDielectronVarManager::kNMaxValues,values); | |
410 | if (hasMC && h->GetHistogramList()->FindObject("MCEvent_noCuts")) | |
411 | h->FillClass("Event_noCuts",AliDielectronVarManager::kNMaxValues,valuesMC); | |
412 | } | |
413 | } | |
414 | nextDie.Reset(); | |
415 | ||
416 | //event filter | |
417 | if (fEventFilter) { | |
418 | if (!fEventFilter->IsSelected(InputEvent())) return; | |
419 | } | |
420 | fEventStat->Fill(kFilteredEvents); | |
421 | ||
422 | //pileup | |
423 | if (fRejectPileup){ | |
424 | if (InputEvent()->IsPileupFromSPD(3,0.8,3.,2.,5.)) return; | |
425 | } | |
426 | fEventStat->Fill(kPileupEvents); | |
427 | ||
428 | //bz for AliKF | |
7329335b | 429 | fbz = InputEvent()->GetMagneticField(); |
430 | AliKFParticle::SetField( fbz ); | |
b1f673b5 | 431 | AliDielectronPID::SetCorrVal((Double_t)InputEvent()->GetRunNumber()); |
432 | ||
433 | //Process event in all AliDielectron instances | |
434 | // TIter nextDie(&fListDielectron); | |
435 | // AliDielectron *die=0; | |
436 | Int_t idie=0; | |
437 | while ( (die=static_cast<AliDielectron*>(nextDie())) ){ | |
805bd069 | 438 | //AliInfo(Form(" **************** die->Process(InputEvent()) : %d **************************", idie)); |
b1f673b5 | 439 | die->SetDontClearArrays(kTRUE); |
440 | die->Process(InputEvent()); | |
441 | if (die->HasCandidates()){ | |
442 | Int_t ncandidates=die->GetPairArray(1)->GetEntriesFast(); | |
443 | if (ncandidates==1) fEventStat->Fill((kNbinsEvent)+2*idie); | |
444 | else if (ncandidates>1) fEventStat->Fill((kNbinsEvent+1)+2*idie); | |
445 | } | |
446 | ||
447 | ||
448 | AliDielectronVarManager::Fill(InputEvent(), fgValues); | |
449 | for (Int_t ii=0; ii<2; ++ii){ | |
450 | TObjArray *obj = (TObjArray*)die->GetTrackArray(ii); | |
451 | Int_t ntracks=obj->GetEntriesFast(); | |
452 | //AliInfo(Form(" ************** # of tracks = %d", ntracks)); | |
453 | for (Int_t itrack=0; itrack<ntracks; ++itrack){ | |
454 | ||
455 | //////////////////////////////////////////////////////////////////// | |
456 | AliDielectronVarManager::Fill(obj->UncheckedAt(itrack), fgValues); | |
457 | //////////////////////////////////////////////////////////////////// | |
458 | ||
459 | if(fgValues[AliDielectronVarManager::kCharge]>0){ | |
7329335b | 460 | fVeptmp.push_back(new AliDielectronSingleTG(1, |
b1f673b5 | 461 | fgValues[AliDielectronVarManager::kCentrality], |
462 | fgValues[AliDielectronVarManager::kXv], | |
463 | fgValues[AliDielectronVarManager::kYv], | |
464 | fgValues[AliDielectronVarManager::kZv], | |
465 | fgValues[AliDielectronVarManager::kPx], | |
466 | fgValues[AliDielectronVarManager::kPy], | |
467 | fgValues[AliDielectronVarManager::kPz], | |
468 | fgValues[AliDielectronVarManager::kPt], | |
469 | fgValues[AliDielectronVarManager::kEta], | |
470 | fgValues[AliDielectronVarManager::kPhi], | |
471 | fgValues[AliDielectronVarManager::kTheta], | |
472 | 1, 1, static_cast<AliVTrack*>(obj->UncheckedAt(itrack))) | |
473 | ); | |
474 | }else if(fgValues[AliDielectronVarManager::kCharge]<0){ | |
7329335b | 475 | fVemtmp.push_back(new AliDielectronSingleTG(-1, |
b1f673b5 | 476 | fgValues[AliDielectronVarManager::kCentrality], |
477 | fgValues[AliDielectronVarManager::kXv], | |
478 | fgValues[AliDielectronVarManager::kYv], | |
479 | fgValues[AliDielectronVarManager::kZv], | |
480 | fgValues[AliDielectronVarManager::kPx], | |
481 | fgValues[AliDielectronVarManager::kPy], | |
482 | fgValues[AliDielectronVarManager::kPz], | |
483 | fgValues[AliDielectronVarManager::kPt], | |
484 | fgValues[AliDielectronVarManager::kEta], | |
485 | fgValues[AliDielectronVarManager::kPhi], | |
486 | fgValues[AliDielectronVarManager::kTheta], | |
487 | 1, 1, static_cast<AliVTrack*>(obj->UncheckedAt(itrack))) | |
488 | ); | |
489 | } | |
490 | } | |
491 | } | |
7329335b | 492 | //AliInfo(Form("size of e and p = %d %d", (int)fVeptmp.size(), (int)fVemtmp.size())); |
493 | ||
b1f673b5 | 494 | |
7329335b | 495 | CheckGhostPairs(fVeptmp); |
496 | CheckGhostPairs(fVemtmp); | |
805bd069 | 497 | if(fRejectPairFlag[idie]==1 || fRejectPairFlag[idie]==2){ |
498 | RejectPairs(fVeptmp, fVemtmp, idie); | |
499 | } | |
7329335b | 500 | RandomizePool(fVeptmp, fVemtmp); |
501 | CalcPair(fVep, fVem, die, idie); | |
b1f673b5 | 502 | |
7329335b | 503 | // AliInfo(Form("size of e and p (after) = %d %d", (int)fVep.size(), (int)fVem.size())); |
b1f673b5 | 504 | |
7329335b | 505 | double dwcent = 100.0/fgkNCENT; |
506 | double dwiz = 20.0/fgkNZBIN; | |
507 | double dwrp = acos(-1.0)/fgkNRPBIN; | |
b1f673b5 | 508 | |
7329335b | 509 | int icent = (int)(fgValues[AliDielectronVarManager::kCentrality]/dwcent); |
510 | int izbin = (int)((fgValues[AliDielectronVarManager::kZvPrim]+10)/dwiz); | |
511 | int irp = (int)((fgValues[AliDielectronVarManager::kV0ACrpH2])/dwrp); | |
b1f673b5 | 512 | |
513 | if(icent<0) icent=0; | |
7329335b | 514 | if(icent>=fgkNCENT) icent=fgkNCENT-1; |
b1f673b5 | 515 | if(izbin<0) izbin=0; |
7329335b | 516 | if(izbin>=fgkNZBIN) izbin=fgkNZBIN-1; |
b1f673b5 | 517 | if(irp<0) irp=0; |
7329335b | 518 | if(irp>=fgkNRPBIN) irp=fgkNRPBIN-1; |
b1f673b5 | 519 | |
7329335b | 520 | fvep[fibuf[idie][izbin][icent][irp]][idie][izbin][icent][irp].clear(); |
521 | for(int iep = 0; iep<(int)fVep.size();iep++) { | |
522 | fvep[fibuf[idie][izbin][icent][irp]][idie][izbin][icent][irp].push_back(fVep[iep]); | |
523 | fpoolp[idie][izbin][icent][irp].push_back(fVep[iep]); | |
524 | if(fpoolp[idie][izbin][icent][irp].size()>fgkMAXPOOL) { | |
525 | fpoolp[idie][izbin][icent][irp].pop_front(); | |
b1f673b5 | 526 | } |
527 | } | |
7329335b | 528 | fvem[fibuf[idie][izbin][icent][irp]][idie][izbin][icent][irp].clear(); |
529 | for(int iem = 0; iem<(int)fVem.size();iem++) { | |
530 | fvem[fibuf[idie][izbin][icent][irp]][idie][izbin][icent][irp].push_back(fVem[iem]); | |
531 | fpoolm[idie][izbin][icent][irp].push_back(fVem[iem]); | |
532 | if(fpoolm[idie][izbin][icent][irp].size()>fgkMAXPOOL) { | |
533 | fpoolm[idie][izbin][icent][irp].pop_front(); | |
b1f673b5 | 534 | } |
535 | } | |
536 | ||
537 | ||
7329335b | 538 | fibuf[idie][izbin][icent][irp]++; |
539 | if(fibuf[idie][izbin][icent][irp]>= fgkNBUF) fibuf[idie][izbin][icent][irp]=0; | |
b1f673b5 | 540 | |
541 | ||
7329335b | 542 | fVeptmp.clear(); |
543 | fVemtmp.clear(); | |
544 | fVep.clear(); | |
545 | fVem.clear(); | |
b1f673b5 | 546 | |
547 | ||
548 | ++idie; | |
549 | } | |
550 | ||
551 | ||
552 | AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent()); | |
553 | fEvent->Fill(values[AliDielectronVarManager::kCentrality]); | |
554 | for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) { | |
555 | AliESDtrack* track = fESD->GetTrack(iTracks); | |
556 | if (!track) { | |
557 | Printf("ERROR: Could not receive track %d", iTracks); | |
558 | continue; | |
559 | } | |
560 | if(!fCutsMother->AcceptTrack(track)) continue; | |
561 | fdEdXvsPt->Fill(track->GetTPCmomentum(), track->GetTPCsignal()); | |
562 | fdEdXnSigmaElecvsPt->Fill(track->GetTPCmomentum(), | |
563 | AliDielectronVarManager::GetPIDResponse()->NumberOfSigmasTPC(track, | |
564 | AliPID::kElectron) | |
565 | -AliDielectronPID::GetCorrVal()); | |
566 | /// for beta caliculaton | |
567 | Double_t l = track->GetIntegratedLength(); // cm | |
568 | Double_t t = track->GetTOFsignal(); | |
569 | Double_t t0 = AliDielectronVarManager::GetPIDResponse()->GetTOFResponse().GetTimeZero(); // ps | |
570 | Double_t beta = 0; | |
571 | if( (l < 360. || l > 800.) || (t <= 0.) || (t0 >999990.0) ) { | |
572 | beta=-9999; | |
573 | } | |
574 | else { | |
575 | t -= t0; // subtract the T0 | |
576 | l *= 0.01; // cm ->m | |
577 | t *= 1e-12; //ps -> s | |
578 | ||
579 | Double_t v = l / t; | |
580 | beta = v / TMath::C(); | |
581 | } | |
582 | ||
583 | fTOFbetavsPt->Fill(track->GetTPCmomentum(), beta); | |
584 | fTOFnSigmaElecvsPt->Fill(track->GetTPCmomentum(), | |
585 | AliDielectronVarManager::GetPIDResponse()->NumberOfSigmasTOF(track, | |
586 | AliPID::kElectron)); | |
587 | ////rough PID is required | |
588 | if( fabs(AliDielectronVarManager::GetPIDResponse()->NumberOfSigmasTOF(track, AliPID::kElectron))<3){ | |
589 | fdEdXvsPtTOF->Fill(track->GetTPCmomentum(), track->GetTPCsignal()); | |
590 | fdEdXnSigmaElecvsPtTOF->Fill(track->GetTPCmomentum(), | |
591 | AliDielectronVarManager::GetPIDResponse()->NumberOfSigmasTPC(track, | |
592 | AliPID::kElectron) | |
593 | -AliDielectronPID::GetCorrVal()); | |
594 | ||
595 | ||
596 | if(track->GetTPCsignal()>70 && track->GetTPCsignal()<90){ | |
7329335b | 597 | fNCrossedRowsTPC->Fill(track->GetTPCmomentum(), track->GetTPCCrossedRows()); |
598 | //fChi2ClusTPC->Fill(track->GetTPCmomentum(), track->GetTPCchi2()/track->GetTPCNcls()); | |
599 | fChi2ClusTPC->Fill(track->GetTPCmomentum(), track->GetTPCchi2()); | |
600 | fRatioCrossClusTPC->Fill(track->GetTPCmomentum(), track->GetTPCCrossedRows()/track->GetTPCNclsF()); | |
b1f673b5 | 601 | } |
602 | } | |
603 | } | |
604 | ||
605 | PostData(1, &fListHistos); | |
606 | PostData(2, &fListCF); | |
607 | PostData(3,fEventStat); | |
608 | } | |
609 | ||
610 | //_________________________________________________________________________________ | |
611 | void AliAnalysisTaskMultiDielectronTG::FinishTaskOutput() | |
612 | { | |
613 | // | |
614 | // Write debug tree | |
615 | // | |
616 | TIter nextDie(&fListDielectron); | |
617 | AliDielectron *die=0; | |
618 | while ( (die=static_cast<AliDielectron*>(nextDie())) ){ | |
619 | die->SaveDebugTree(); | |
620 | AliDielectronMixingHandler *mix=die->GetMixingHandler(); | |
621 | // printf("\n\n\n===============\ncall mix in Terminate: %p (%p)\n=================\n\n",mix,die); | |
622 | if (mix) mix->MixRemaining(die); | |
623 | } | |
624 | PostData(1, &fListHistos); | |
625 | PostData(2, &fListCF); | |
626 | } | |
627 | ||
628 | //_________________________________________________________________________________ | |
7329335b | 629 | void AliAnalysisTaskMultiDielectronTG::CheckGhostPairs(vector<AliDielectronSingleTG*> e1) |
630 | { | |
631 | ////// Check whether the pairs are in ghost pairs | |
632 | ///// ghost means that one of the pairs is highly associated but reconstructed as a true particle | |
633 | ||
634 | ||
b1f673b5 | 635 | bool reject = false; |
636 | if(e1.size()>1){ | |
637 | for(int i1=0; i1<(int)e1.size(); i1++){ | |
638 | reject = false; | |
639 | for(int i2=i1+1; i2<(int)e1.size(); i2++){ | |
620b34a5 | 640 | if( fabs(e1[i1]->Phi() - e1[i2]->Phi())<0.01 ){ |
b1f673b5 | 641 | reject = true; |
642 | e1[i2]->SetGstFlag(0); | |
643 | } | |
644 | } | |
645 | if(reject==true)e1[i1]->SetGstFlag(0); | |
646 | } | |
647 | } | |
648 | } | |
649 | ||
650 | //_________________________________________________________________________________ | |
7329335b | 651 | void AliAnalysisTaskMultiDielectronTG::RandomizePool(vector<AliDielectronSingleTG*> e1, vector<AliDielectronSingleTG*> e2) |
652 | { | |
b1f673b5 | 653 | |
7329335b | 654 | ///// Randomize the pool constent to cancel the filling scheme of the single tracks |
655 | //// | |
656 | ||
657 | ||
b1f673b5 | 658 | int size1 = e1.size(); |
7329335b | 659 | int usedindex[1000]; |
b1f673b5 | 660 | for(int i=0;i<1000;i++){ |
7329335b | 661 | usedindex[i] = -1; |
b1f673b5 | 662 | } |
663 | for(int i=0;i<size1;i++){ | |
7329335b | 664 | usedindex[i] = 0; |
b1f673b5 | 665 | } |
666 | ||
667 | for(int i=0;i<size1;i++){ | |
668 | int j = (int)(gRandom->Uniform(0,size1)); | |
7329335b | 669 | while(usedindex[j]==1){ |
b1f673b5 | 670 | j = (int)(gRandom->Uniform(0,size1)); |
671 | } | |
672 | if( (e1[j]->GetGstFlag()==1) && | |
673 | (e1[j]->GetConvFlag()==1) | |
674 | ){ | |
7329335b | 675 | fVep.push_back(e1[j]); |
b1f673b5 | 676 | } |
7329335b | 677 | usedindex[j] = 1; |
b1f673b5 | 678 | } |
679 | ||
680 | ||
681 | int size2 = e2.size(); | |
682 | for(int i=0;i<1000;i++){ | |
7329335b | 683 | usedindex[i] = -1; |
b1f673b5 | 684 | } |
685 | for(int i=0;i<size2;i++){ | |
7329335b | 686 | usedindex[i] = 0; |
b1f673b5 | 687 | } |
688 | ||
689 | for(int i=0;i<size2;i++){ | |
690 | int j = (int)(gRandom->Uniform(0,size2)); | |
7329335b | 691 | while(usedindex[j]==1){ |
b1f673b5 | 692 | j = (int)(gRandom->Uniform(0,size2)); |
693 | } | |
694 | if( (e2[j]->GetGstFlag()==1) && | |
695 | (e2[j]->GetConvFlag()==1) | |
696 | ){ | |
7329335b | 697 | fVem.push_back(e2[j]); |
b1f673b5 | 698 | } |
7329335b | 699 | usedindex[j] = 1; |
b1f673b5 | 700 | } |
701 | } | |
702 | ||
703 | ||
704 | //_________________________________________________________________________________ | |
7329335b | 705 | void AliAnalysisTaskMultiDielectronTG::CalcPair(vector<AliDielectronSingleTG*> ve1, |
706 | vector<AliDielectronSingleTG*> ve2, AliDielectron *die, Int_t idie) | |
707 | { | |
708 | ||
709 | // | |
710 | // main routine for the pair procesing | |
711 | // | |
712 | ||
b1f673b5 | 713 | |
714 | for(int i1=0; i1<(int)ve1.size(); i1++){ | |
715 | for(int i2=0; i2<(int)ve2.size(); i2++){ | |
805bd069 | 716 | FillPair(ve1[i1], ve2[i2], 0, die, idie); |
b1f673b5 | 717 | } |
718 | } | |
719 | ||
720 | for(int i1=0; i1<(int)ve1.size(); i1++){ | |
721 | for(int i2=i1+1; i2<(int)ve1.size(); i2++){ | |
805bd069 | 722 | FillPair(ve1[i1], ve1[i2], 1, die, idie); |
b1f673b5 | 723 | } |
724 | } | |
725 | ||
726 | for(int i1=0; i1<(int)ve2.size(); i1++){ | |
727 | for(int i2=i1+1; i2<(int)ve2.size(); i2++){ | |
805bd069 | 728 | FillPair(ve2[i1], ve2[i2], 2, die, idie ); |
b1f673b5 | 729 | } |
730 | } | |
731 | ||
732 | ||
7329335b | 733 | double dwcent = 100.0/fgkNCENT; |
734 | double dwiz = 20.0/fgkNZBIN; | |
735 | double dwrp = acos(-1.0)/fgkNRPBIN; | |
b1f673b5 | 736 | |
7329335b | 737 | int icent = (int)(fgValues[AliDielectronVarManager::kCentrality]/dwcent); |
738 | int izbin = (int)((fgValues[AliDielectronVarManager::kZvPrim]+10)/dwiz); | |
739 | int irp = (int)((fgValues[AliDielectronVarManager::kV0ACrpH2])/dwrp); | |
b1f673b5 | 740 | |
741 | if(icent<0) icent=0; | |
7329335b | 742 | if(icent>=fgkNCENT) icent=fgkNCENT-1; |
b1f673b5 | 743 | if(izbin<0) izbin=0; |
7329335b | 744 | if(izbin>=fgkNZBIN) izbin=fgkNZBIN-1; |
b1f673b5 | 745 | if(irp<0) irp=0; |
7329335b | 746 | if(irp>=fgkNRPBIN) irp=fgkNRPBIN-1; |
b1f673b5 | 747 | |
748 | int nmixed; | |
749 | if(ve1.size()>0) { | |
750 | // | |
751 | // Now mixed event for +- pairs | |
752 | // | |
753 | nmixed = 0; | |
7329335b | 754 | for(int ibuf=0;(nmixed<fgkNMix);ibuf++) { |
b1f673b5 | 755 | int ntry = 0; |
620b34a5 | 756 | while(ntry<fgkMAXTRY) { |
7329335b | 757 | ReshuffleBuffer(fvem[ibuf][idie][izbin][icent][irp],fpoolm[idie][izbin][icent][irp]); |
b1f673b5 | 758 | ntry++; |
759 | } | |
760 | for(int i1=0; i1<(int)ve1.size(); i1++){ | |
7329335b | 761 | for(int i2=0; i2<(int)fvem[ibuf][idie][izbin][icent][irp].size(); i2++){ |
805bd069 | 762 | FillPair(ve1[i1],fvem[ibuf][idie][izbin][icent][irp][i2], 3, die, idie); |
b1f673b5 | 763 | } |
764 | } | |
765 | ++nmixed; | |
766 | }//for(ibuf) | |
767 | } | |
768 | if(ve2.size()>0) { | |
769 | // | |
805bd069 | 770 | // Now mixed event for -+ pairs |
b1f673b5 | 771 | // |
772 | nmixed = 0; | |
7329335b | 773 | for(int ibuf=0;(nmixed<fgkNMix);ibuf++) { |
b1f673b5 | 774 | int ntry = 0; |
620b34a5 | 775 | while(ntry<fgkMAXTRY) { |
7329335b | 776 | ReshuffleBuffer(fvep[ibuf][idie][izbin][icent][irp],fpoolp[idie][izbin][icent][irp]); |
b1f673b5 | 777 | ntry++; |
778 | } | |
779 | for(int i1=0; i1<(int)ve2.size(); i1++){ | |
7329335b | 780 | for(int i2=0; i2<(int)fvep[ibuf][idie][izbin][icent][irp].size(); i2++){ |
805bd069 | 781 | FillPair(fvep[ibuf][idie][izbin][icent][irp][i2],ve2[i1],4, die, idie); |
b1f673b5 | 782 | } |
783 | } | |
784 | ++nmixed; | |
785 | }//for(ibuf) | |
786 | } | |
787 | ||
788 | if(ve1.size()>0) { | |
789 | // | |
790 | // Now mixed event for ++ pairs | |
791 | // | |
792 | nmixed = 0; | |
7329335b | 793 | for(int ibuf=0;(nmixed<fgkNMix);ibuf++) { |
b1f673b5 | 794 | int ntry = 0; |
620b34a5 | 795 | while(ntry<fgkMAXTRY) { |
7329335b | 796 | ReshuffleBuffer(fvep[ibuf][idie][izbin][icent][irp],fpoolp[idie][izbin][icent][irp]); |
b1f673b5 | 797 | ntry++; |
798 | } | |
799 | for(int i1=0; i1<(int)ve1.size(); i1++){ | |
7329335b | 800 | for(int i2=0;i2<(int)fvep[ibuf][idie][izbin][icent][irp].size(); i2++){ |
805bd069 | 801 | FillPair(ve1[i1],fvep[ibuf][idie][izbin][icent][irp][i2], 5, die, idie); |
b1f673b5 | 802 | } |
803 | } | |
804 | ++nmixed; | |
805 | }//for(ibuf) | |
806 | } | |
807 | ||
808 | if(ve2.size()>0) { | |
809 | // | |
810 | // Now mixed event for +- pairs | |
811 | // | |
812 | nmixed = 0; | |
7329335b | 813 | for(int ibuf=0;(nmixed<fgkNMix);ibuf++) { |
b1f673b5 | 814 | int ntry = 0; |
620b34a5 | 815 | while(ntry<fgkMAXTRY) { |
7329335b | 816 | ReshuffleBuffer(fvem[ibuf][idie][izbin][icent][irp],fpoolm[idie][izbin][icent][irp]); |
b1f673b5 | 817 | ntry++; |
818 | } | |
819 | for(int i1=0; i1<(int)ve2.size(); i1++){ | |
7329335b | 820 | for(int i2=0; i2<(int)fvem[ibuf][idie][izbin][icent][irp].size(); i2++){ |
805bd069 | 821 | FillPair(ve2[i1],fvem[ibuf][idie][izbin][icent][irp][i2],6, die, idie); |
b1f673b5 | 822 | } |
823 | } | |
824 | ++nmixed; | |
825 | }//for(ibuf) | |
826 | } | |
827 | ||
828 | } | |
829 | ||
830 | ||
831 | //_________________________________________________________________________________ | |
7329335b | 832 | void AliAnalysisTaskMultiDielectronTG::FillPair(AliDielectronSingleTG *iep, |
805bd069 | 833 | AliDielectronSingleTG *iem, int type, AliDielectron *die, Int_t idie) |
7329335b | 834 | { |
835 | ||
836 | // | |
837 | // main routine for filling kinematics of pairs | |
838 | // | |
839 | ||
840 | ||
b1f673b5 | 841 | |
7329335b | 842 | double dmass, dphiv, dpxpair, dpypair, dpzpair; |
843 | double dptpair, depair, dphipair, detapair, dcos, dpsi; | |
844 | ||
845 | CalcVars(iep, iem, dmass, dphiv, dpxpair, dpypair, dpzpair, | |
846 | dptpair, depair, dphipair, detapair, dcos, dpsi); | |
847 | ||
848 | ||
849 | double dopeningangle = -9999; | |
850 | double dv0mass = -9999; | |
851 | double dv0pxpair = -9999; | |
852 | double dv0pypair = -9999; | |
853 | double dv0pzpair = -9999; | |
854 | double dv0ptpair = -9999; | |
855 | double dv0epair = -9999; | |
856 | double dv0xvpair = -9999; | |
857 | double dv0yvpair = -9999; | |
858 | double dv0zvpair = -9999; | |
859 | double dv0phipair = -9999; | |
860 | double dv0etapair = -9999; | |
861 | double dv0rpair = -9999; | |
862 | double dpsipair = -9999; | |
805bd069 | 863 | double dphivpair = -9999; |
b1f673b5 | 864 | |
865 | //////////////////////////// | |
866 | ///// calculate v0 //////// | |
867 | /////////////////////////// | |
7329335b | 868 | Bool_t v0OFF=kFALSE; |
b1f673b5 | 869 | /// for the moment, this doesn't work for MC |
7329335b | 870 | if(fdv0mixing == kFALSE && (type==3 || type==4 || type==5 || type==6)){ |
871 | v0OFF = kTRUE; | |
b1f673b5 | 872 | } |
873 | if(die->GetHasMC()==kTRUE && (type==3 || type==4 || type==5 || type==6)){ | |
7329335b | 874 | v0OFF = kTRUE; |
b1f673b5 | 875 | } |
876 | if(type==0 || type==1 || type==2){ | |
7329335b | 877 | v0OFF = kFALSE; |
b1f673b5 | 878 | } |
879 | ||
880 | ||
7329335b | 881 | if(v0OFF==kFALSE){ |
b1f673b5 | 882 | AliVTrack* trackob1= iep->GetTrack(); |
883 | AliVTrack* trackob2= iem->GetTrack(); | |
884 | ||
885 | if(!trackob1 || !trackob2){ | |
886 | return; | |
887 | } | |
888 | ||
889 | AliDielectronPair candidate; | |
890 | candidate.SetTracks(trackob1, (int)(11*iep->Charge()), | |
891 | trackob2, (int)(11*iem->Charge())); | |
892 | ||
893 | if(trackob1==trackob2){ | |
894 | AliInfo("Same Objects!!"); | |
895 | return; | |
896 | } | |
805bd069 | 897 | |
898 | //const AliKFParticle &kfPairLeg1 = candidate.GetKFFirstDaughter(); | |
899 | //const AliKFParticle &kfPairLeg2 = candidate.GetKFSecondDaughter(); | |
900 | ||
b1f673b5 | 901 | const AliKFParticle &kfPair = candidate.GetKFParticle(); |
805bd069 | 902 | |
903 | /* | |
7329335b | 904 | dv0mass = candidate.M(); |
905 | dv0pxpair = candidate.Px(); | |
906 | dv0pypair = candidate.Py(); | |
907 | dv0pzpair = candidate.Pz(); | |
908 | dv0ptpair = candidate.Pt(); | |
909 | dv0epair = candidate.E(); | |
910 | dv0xvpair = candidate.Xv(); | |
911 | dv0yvpair = candidate.Yv(); | |
912 | dv0zvpair = candidate.Zv(); | |
913 | dv0phipair = candidate.Phi(); | |
7329335b | 914 | dv0etapair = candidate.Eta(); |
805bd069 | 915 | */ |
916 | dv0mass = kfPair.GetMass(); | |
917 | dv0pxpair = kfPair.GetPx(); | |
918 | dv0pypair = kfPair.GetPy(); | |
919 | dv0pzpair = kfPair.GetPz(); | |
920 | dv0ptpair = kfPair.GetPt(); | |
921 | dv0epair = kfPair.GetE(); | |
922 | dv0xvpair = kfPair.GetX(); | |
923 | dv0yvpair = kfPair.GetY(); | |
924 | dv0zvpair = kfPair.GetZ(); | |
925 | dv0phipair = kfPair.GetPhi(); | |
926 | dv0etapair = kfPair.GetEta(); | |
7329335b | 927 | dv0rpair = kfPair.GetR(); |
805bd069 | 928 | |
929 | dopeningangle = candidate.OpeningAngle(); | |
7329335b | 930 | dpsipair = candidate.PsiPair(fbz); |
805bd069 | 931 | dphivpair = candidate.PhivPair(fbz); |
932 | ||
933 | ||
b1f673b5 | 934 | } |
935 | ||
936 | Double_t values[AliDielectronVarManager::kNMaxValues]; | |
937 | TString className1; | |
938 | TString className2; | |
7329335b | 939 | className1.Form("MyPair_%s",kPairClassNamesTG[type]); |
940 | className2.Form("MyPairV0_%s",kPairClassNamesTG[type]); | |
b1f673b5 | 941 | |
942 | AliDielectronHistos *fHistos = die->GetHistoManager(); | |
943 | Bool_t pairClass1=fHistos->GetHistogramList()->FindObject(className1.Data())!=0x0; | |
944 | Bool_t pairClass2=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0; | |
945 | ||
805bd069 | 946 | if (pairClass1 && PairTrackcut(dphiv, dcos, idie)==true){ |
b1f673b5 | 947 | ///import pair variables to values!! |
7329335b | 948 | values[AliDielectronVarManager::kPx] = dpxpair; |
949 | values[AliDielectronVarManager::kPy] = dpypair; | |
950 | values[AliDielectronVarManager::kPz] = dpzpair; | |
951 | values[AliDielectronVarManager::kPt] = dptpair; | |
952 | values[AliDielectronVarManager::kXv] = dv0xvpair; | |
953 | values[AliDielectronVarManager::kYv] = dv0yvpair; | |
954 | values[AliDielectronVarManager::kZv] = dv0zvpair; | |
955 | values[AliDielectronVarManager::kR] = dv0rpair; | |
956 | values[AliDielectronVarManager::kE] = depair; | |
957 | values[AliDielectronVarManager::kEta] = detapair; | |
958 | values[AliDielectronVarManager::kM] = dmass; | |
959 | values[AliDielectronVarManager::kPsiPair] = dphiv; | |
805bd069 | 960 | values[AliDielectronVarManager::kPhivPair] = dphiv; |
7329335b | 961 | values[AliDielectronVarManager::kPhi] = dphipair; |
962 | values[AliDielectronVarManager::kOpeningAngle] = dcos; | |
b1f673b5 | 963 | fHistos->FillClass(className1, AliDielectronVarManager::kNMaxValues, values); |
964 | } | |
965 | ||
966 | ||
805bd069 | 967 | if (pairClass2 && PairTrackcut(dphiv, dopeningangle, idie)==true){ |
7329335b | 968 | values[AliDielectronVarManager::kPx] = dv0pxpair; |
969 | values[AliDielectronVarManager::kPy] = dv0pypair; | |
970 | values[AliDielectronVarManager::kPz] = dv0pzpair; | |
971 | values[AliDielectronVarManager::kPt] = dv0ptpair; | |
972 | values[AliDielectronVarManager::kXv] = dv0xvpair; | |
973 | values[AliDielectronVarManager::kYv] = dv0yvpair; | |
974 | values[AliDielectronVarManager::kZv] = dv0zvpair; | |
975 | values[AliDielectronVarManager::kR] = dv0rpair; | |
976 | values[AliDielectronVarManager::kE] = dv0epair; | |
977 | values[AliDielectronVarManager::kEta] = dv0etapair; | |
978 | values[AliDielectronVarManager::kM] = dv0mass; | |
979 | values[AliDielectronVarManager::kPsiPair] = dpsipair; | |
805bd069 | 980 | values[AliDielectronVarManager::kPhivPair] = dphivpair; |
7329335b | 981 | values[AliDielectronVarManager::kPhi] = dv0phipair; |
982 | values[AliDielectronVarManager::kOpeningAngle] = dopeningangle; | |
b1f673b5 | 983 | fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values); |
984 | } | |
985 | ||
986 | ||
987 | ||
988 | } | |
989 | ||
990 | //_________________________________________________________________________________ | |
805bd069 | 991 | bool AliAnalysisTaskMultiDielectronTG::PairTrackcut(double phiv, double op, int idie) |
7329335b | 992 | { |
993 | ||
994 | // | |
995 | // pair-by-pair cuts | |
996 | // | |
b1f673b5 | 997 | |
998 | ||
999 | bool pairOK = true; | |
1000 | ||
805bd069 | 1001 | if(fRejectPairFlag[idie] == 1 || fRejectPairFlag[idie] == 2 || |
1002 | fRejectPairFlag[idie] == 3 || fRejectPairFlag[idie] == 4 ){ | |
1003 | if(fRejectPairFlag[idie] == 2 || fRejectPairFlag[idie] == 4 ){ | |
1004 | if(fbz>0 && phiv>fdconvphiv){ | |
1005 | pairOK = false; | |
1006 | }else if(fbz<0 && phiv<acos(-1.0)-fdconvphiv){ | |
1007 | pairOK = false; | |
1008 | } | |
1009 | }else if(fRejectPairFlag[idie] == 1 || fRejectPairFlag[idie] == 3){ | |
1010 | if(op<fdop){ | |
1011 | pairOK = false; | |
1012 | } | |
1013 | } | |
b1f673b5 | 1014 | } |
1015 | ||
1016 | return pairOK; | |
1017 | ||
1018 | } | |
1019 | ||
1020 | ||
1021 | //_________________________________________________________________________________ | |
7329335b | 1022 | void AliAnalysisTaskMultiDielectronTG::CalcVars(AliDielectronSingleTG *iep, AliDielectronSingleTG *iem, |
b1f673b5 | 1023 | double &mass, double &phiv, double &px, double &py, double&pz, |
1024 | double &pt, double &e, double &phi, | |
7329335b | 1025 | double &eta, double &cos, double &psi) |
1026 | { | |
1027 | ||
1028 | ||
1029 | // | |
1030 | // standalone calculator for the pair variables | |
1031 | // | |
b1f673b5 | 1032 | |
1033 | px = iep->Px()+iem->Px(); | |
1034 | py = iep->Py()+iem->Py(); | |
1035 | pz = iep->Pz()+iem->Pz(); | |
1036 | pt = sqrt(px*px+py*py); | |
7329335b | 1037 | double dppair = sqrt(pt*pt+pz*pz); |
b1f673b5 | 1038 | static const double me=0.0005109989; |
1039 | e = sqrt(me*me+iep->Px()*iep->Px()+iep->Py()*iep->Py()+iep->Pz()*iep->Pz()) | |
1040 | + sqrt(me*me+iem->Px()*iem->Px()+iem->Py()*iem->Py()+iem->Pz()*iem->Pz()); | |
1041 | ||
1042 | mass = e*e-px*px-py*py-pz*pz; | |
1043 | if(mass<0){ | |
1044 | mass = mass; | |
1045 | }else{ | |
1046 | mass = sqrt(mass); | |
1047 | } | |
1048 | ||
1049 | ||
1050 | phi = atan2(py, px); | |
7329335b | 1051 | eta = -0.5*TMath::Log((dppair+pz)/(dppair-pz)); |
b1f673b5 | 1052 | double p1 = sqrt(pow(iep->Px(),2)+pow(iep->Py(),2)+pow(iep->Pz(),2)); |
1053 | double p2 = sqrt(pow(iem->Px(),2)+pow(iem->Py(),2)+pow(iem->Pz(),2)); | |
1054 | cos = acos((iep->Px()*iem->Px()+iep->Py()*iem->Py()+iep->Pz()*iem->Pz())/(p1*p2)); | |
1055 | ||
1056 | double dtheta = iep->Theta()-iem->Theta(); | |
1057 | psi = asin(dtheta/cos); | |
1058 | ||
1059 | ||
1060 | //unit vector of (pep+pem) | |
7329335b | 1061 | double pl = dppair; |
b1f673b5 | 1062 | double ux = px/pl; |
1063 | double uy = py/pl; | |
1064 | double uz = pz/pl; | |
1065 | double ax = uy/sqrt(ux*ux+uy*uy); | |
1066 | double ay = -ux/sqrt(ux*ux+uy*uy); | |
1067 | ||
1068 | //momentum of e+ and e- in (ax,ay,az) axis. Note that az=0 by | |
1069 | //definition. | |
1070 | //double ptep = iep->Px()*ax + iep->Py()*ay; | |
1071 | //double ptem = iem->Px()*ax + iem->Py()*ay; | |
1072 | ||
1073 | double pxep = iep->Px(); | |
1074 | double pyep = iep->Py(); | |
1075 | double pzep = iep->Pz(); | |
1076 | double pxem = iem->Px(); | |
1077 | double pyem = iem->Py(); | |
1078 | double pzem = iem->Pz(); | |
1079 | ||
1080 | ||
1081 | //vector product of pep X pem | |
1082 | double vpx = pyep*pzem - pzep*pyem; | |
1083 | double vpy = pzep*pxem - pxep*pzem; | |
1084 | double vpz = pxep*pyem - pyep*pxem; | |
1085 | double vp = sqrt(vpx*vpx+vpy*vpy+vpz*vpz); | |
1086 | //double thev = acos(vpz/vp); | |
1087 | ||
1088 | //unit vector of pep X pem | |
1089 | double vx = vpx/vp; | |
1090 | double vy = vpy/vp; | |
1091 | double vz = vpz/vp; | |
1092 | ||
1093 | //The third axis defined by vector product (ux,uy,uz)X(vx,vy,vz) | |
1094 | double wx = uy*vz - uz*vy; | |
1095 | double wy = uz*vx - ux*vz; | |
1096 | double wz = ux*vy - uy*vx; | |
1097 | double wl = sqrt(wx*wx+wy*wy+wz*wz); | |
1098 | // by construction, (wx,wy,wz) must be a unit vector. | |
1099 | if(fabs(wl - 1.0) > 0.00001) std::cout << "Calculation error in W vector"<<std::endl; | |
1100 | // measure angle between (wx,wy,wz) and (ax,ay,0). The angle between them | |
1101 | // should be small if the pair is conversion | |
1102 | // | |
1103 | double cosPhiV = wx*ax + wy*ay; | |
1104 | phiv = acos(cosPhiV); | |
1105 | ||
1106 | } | |
1107 | ||
1108 | //_________________________________________________________________________________ | |
7329335b | 1109 | void AliAnalysisTaskMultiDielectronTG::ReshuffleBuffer(vector<AliDielectronSingleTG*> ve, deque<AliDielectronSingleTG*> pool) |
1110 | { | |
b1f673b5 | 1111 | //If there is not enough electron in the pool, give up |
7329335b | 1112 | // |
1113 | // ReshuffleBuffer for th event mixing | |
1114 | // | |
b1f673b5 | 1115 | |
1116 | unsigned int ne = ve.size(); | |
1117 | unsigned int poolsize = pool.size(); | |
7329335b | 1118 | int used[fgkMAXPOOL]; |
1119 | for(int i=0;i<(int)fgkMAXPOOL;i++){ | |
b1f673b5 | 1120 | used[i]=0; |
1121 | } | |
1122 | ||
1123 | if(poolsize < ne) { | |
1124 | std::cout <<" pool size="<<poolsize<<" ne"<<ne<<std::endl; | |
1125 | return; | |
1126 | } | |
1127 | for(unsigned int ie=0; ie < ne; ie++) { | |
1128 | int j = rand()%poolsize; | |
1129 | while(used[j]==1){ | |
1130 | j = rand()%poolsize; | |
1131 | } | |
1132 | ve[ie] = pool[j]; | |
1133 | used[j]=1; | |
1134 | } | |
1135 | ||
1136 | } | |
805bd069 | 1137 | |
1138 | //_________________________________________________________________________________ | |
1139 | void AliAnalysisTaskMultiDielectronTG::RejectPairs(vector<AliDielectronSingleTG*> e1, vector<AliDielectronSingleTG*> e2, Int_t idie){ | |
1140 | ||
1141 | //////////////////////////////////// | |
1142 | ///// to reject pairs from track arrays | |
1143 | ///// conversions, ghost .. | |
1144 | /////////////////////////////////// | |
1145 | if(e1.size()>0 && e2.size()>0){ | |
1146 | for(int i1=0; i1<(int)e1.size(); i1++){ | |
1147 | for(int i2=0; i2<(int)e2.size(); i2++){ | |
1148 | if(fRejectPairFlag[idie]==1){ | |
1149 | Double_t cosine = GetOpeningAngle(e1[i1], e2[i2]); | |
1150 | if(cosine<fdop){ | |
1151 | e1[i1]->SetConvFlag(0); | |
1152 | e2[i2]->SetConvFlag(0); | |
1153 | } | |
1154 | }else if(fRejectPairFlag[idie]==2){ | |
1155 | Double_t phiv = GetPhiv(e1[i1], e2[i2]); | |
1156 | if(fbz>0 && phiv>fdconvphiv){ | |
1157 | e1[i1]->SetConvFlag(0); | |
1158 | e2[i2]->SetConvFlag(0); | |
1159 | }else if(fbz<0 && phiv<acos(-1.0)-fdconvphiv){ | |
1160 | e1[i1]->SetConvFlag(0); | |
1161 | e2[i2]->SetConvFlag(0); | |
1162 | } | |
1163 | } | |
1164 | } | |
1165 | } | |
1166 | } | |
620b34a5 | 1167 | |
805bd069 | 1168 | if(e1.size()>0){ |
1169 | for(int i1=0; i1<(int)e1.size(); i1++){ | |
1170 | for(int i2=i1+1; i2<(int)e1.size(); i2++){ | |
1171 | if(fRejectPairFlag[idie]==1){ | |
1172 | Double_t cosine = GetOpeningAngle(e1[i1], e1[i2]); | |
1173 | if(cosine<fdop){ | |
1174 | e1[i1]->SetConvFlag(0); | |
1175 | e1[i2]->SetConvFlag(0); | |
1176 | } | |
1177 | }else if(fRejectPairFlag[idie]==2){ | |
1178 | Double_t phiv = GetPhiv(e1[i1], e1[i2]); | |
1179 | if(fbz>0 && phiv>fdconvphiv){ | |
1180 | e1[i1]->SetConvFlag(0); | |
1181 | e1[i2]->SetConvFlag(0); | |
1182 | }else if(fbz<0 && phiv<acos(-1.0)-fdconvphiv){ | |
1183 | e1[i1]->SetConvFlag(0); | |
1184 | e1[i2]->SetConvFlag(0); | |
1185 | } | |
1186 | } | |
1187 | } | |
1188 | } | |
1189 | } | |
1190 | ||
1191 | if(e2.size()>0){ | |
1192 | for(int i1=0; i1<(int)e2.size(); i1++){ | |
1193 | for(int i2=i1+1; i2<(int)e2.size(); i2++){ | |
1194 | if(fRejectPairFlag[idie]==1){ | |
1195 | Double_t cosine = GetOpeningAngle(e2[i1], e2[i2]); | |
1196 | if(cosine<fdop){ | |
1197 | e2[i1]->SetConvFlag(0); | |
1198 | e2[i2]->SetConvFlag(0); | |
1199 | } | |
1200 | }else if(fRejectPairFlag[idie]==2){ | |
1201 | Double_t phiv = GetPhiv(e2[i1], e2[i2]); | |
1202 | if(fbz>0 && phiv>fdconvphiv){ | |
1203 | e2[i1]->SetConvFlag(0); | |
1204 | e2[i2]->SetConvFlag(0); | |
1205 | }else if(fbz<0 && phiv<acos(-1.0)-fdconvphiv){ | |
1206 | e2[i1]->SetConvFlag(0); | |
1207 | e2[i2]->SetConvFlag(0); | |
1208 | } | |
1209 | } | |
1210 | } | |
1211 | } | |
1212 | } | |
1213 | } | |
1214 | ||
1215 | ||
1216 | //_________________________________________________________________________________ | |
1217 | Double_t AliAnalysisTaskMultiDielectronTG::GetOpeningAngle(AliDielectronSingleTG* e1, AliDielectronSingleTG* e2){ | |
1218 | ||
1219 | ////////////////////// | |
1220 | //////// calculate pairs and get opening angle | |
1221 | ////////////////////// | |
1222 | ||
1223 | double dmass, dphiv, dpxpair, dpypair, dpzpair; | |
1224 | double dptpair, depair, dphipair, detapair, dcos, dpsi; | |
1225 | ||
1226 | CalcVars(e1, e2, dmass, dphiv, dpxpair, dpypair, dpzpair, | |
1227 | dptpair, depair, dphipair, detapair, dcos, dpsi); | |
1228 | ||
1229 | return dcos; | |
1230 | } | |
1231 | ||
1232 | //_________________________________________________________________________________ | |
1233 | Double_t AliAnalysisTaskMultiDielectronTG::GetPhiv(AliDielectronSingleTG* e1, AliDielectronSingleTG* e2){ | |
1234 | ||
1235 | ////////////////////// | |
1236 | //////// calculate pairs and get opening angle | |
1237 | ////////////////////// | |
1238 | ||
1239 | double dmass, dphiv, dpxpair, dpypair, dpzpair; | |
1240 | double dptpair, depair, dphipair, detapair, dcos, dpsi; | |
1241 | ||
1242 | CalcVars(e1, e2, dmass, dphiv, dpxpair, dpypair, dpzpair, | |
1243 | dptpair, depair, dphipair, detapair, dcos, dpsi); | |
1244 | ||
1245 | return dphiv; | |
1246 | } |