]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/hfe/AliAnalysisTaskDCA.cxx
Updates + addition of EMCal
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliAnalysisTaskDCA.cxx
CommitLineData
70da6c5a 1/**************************************************************************
2* Copyright(c) 1998-1999, 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// The analysis task:
17// impact parameter resolution and pull study
18// for tracks which survivied the particle cuts
19//
20//
21// Authors:
22// Hongyan Yang <hongyan@physi.uni-heidelberg.de>
23// Carlo Bombonati <carlo.bombonati@cern.ch>
24//
6555e2ad 25#include <Riostream.h>
70da6c5a 26#include <TChain.h>
27#include <TFile.h>
28#include <TH1F.h>
29#include <TH1I.h>
30#include <TList.h>
31#include <TMath.h>
32#include <TObjArray.h>
33#include <TParticle.h>
34#include <TString.h>
35
36#include <TCanvas.h>
37
faee3b18 38#include "AliAnalysisManager.h"
39
70da6c5a 40#include "AliCFManager.h"
faee3b18 41
70da6c5a 42#include "AliESDInputHandler.h"
43#include "AliESDtrack.h"
faee3b18 44#include "AliVertexerTracks.h"
45#include "AliESDVertex.h"
46
70da6c5a 47#include "AliMCEventHandler.h"
faee3b18 48#include "AliMCEvent.h"
70da6c5a 49#include "AliMCParticle.h"
50
faee3b18 51#include "AliESDpid.h"
52#include "AliHFEpid.h"
70da6c5a 53#include "AliHFEcuts.h"
54#include "AliHFEdca.h"
6555e2ad 55#include "AliHFEtools.h"
70da6c5a 56
57#include "AliAnalysisTaskDCA.h"
58
59
60//____________________________________________________________
61AliAnalysisTaskDCA::AliAnalysisTaskDCA():
faee3b18 62 AliAnalysisTaskSE("Impact Parameter Resolution and Pull Analysis")
70da6c5a 63 , fPlugins(0)
70da6c5a 64 , fCuts(0x0)
faee3b18 65 , fDefaultPID(0x0)
66 , fHFEpid(0x0)
67 , fPIDdetectors("")
68 , fPIDstrategy(0)
70da6c5a 69 , fCFM(0x0)
70 , fDCA(0x0)
70da6c5a 71 , fNclustersITS(0x0)
faee3b18 72 , fMinNprimVtxContrbutor(0x0)
70da6c5a 73 , fNEvents(0x0)
74 , fResidualList(0x0)
75 , fPullList(0x0)
faee3b18 76 , fDcaList(0x0)
77 , fKfDcaList(0x0)
78 , fMcVertexList(0x0)
79 , fDataDcaList(0x0)
80 , fDataVertexList(0x0)
81 , fDataPullList(0x0)
82 , fMcPidList(0x0)
83 , fDataPidList(0x0)
84 , fHfeDcaList(0x0)
85 , fHfeDataDcaList(0x0)
70da6c5a 86 , fOutput(0x0)
87{
88 //
89 // Dummy constructor
90 //
70da6c5a 91}
92
93//____________________________________________________________
94AliAnalysisTaskDCA::AliAnalysisTaskDCA(const char * name):
faee3b18 95 AliAnalysisTaskSE(name)
70da6c5a 96 , fPlugins(0)
70da6c5a 97 , fCuts(0x0)
faee3b18 98 , fDefaultPID(0x0)
99 , fHFEpid(0x0)
100 , fPIDdetectors("")
101 , fPIDstrategy(0)
70da6c5a 102 , fCFM(0x0)
103 , fDCA(0x0)
70da6c5a 104 , fNclustersITS(0x0)
faee3b18 105 , fMinNprimVtxContrbutor(0x0)
70da6c5a 106 , fNEvents(0x0)
107 , fResidualList(0x0)
108 , fPullList(0x0)
faee3b18 109 , fDcaList(0x0)
110 , fKfDcaList(0x0)
111 , fMcVertexList(0x0)
112 , fDataDcaList(0x0)
113 , fDataVertexList(0x0)
114 , fDataPullList(0x0)
115 , fMcPidList(0x0)
116 , fDataPidList(0x0)
117 , fHfeDcaList(0x0)
118 , fHfeDataDcaList(0x0)
70da6c5a 119 , fOutput(0x0)
120{
121 //
122 // Default constructor
123 //
124 DefineInput(0, TChain::Class());
faee3b18 125 DefineOutput(1, TH1I::Class());
126 DefineOutput(2, TList::Class());
6555e2ad 127
128
129 //-CUTS SETTING-//
130 Int_t nMinTPCcluster = 100;
131 Float_t maxDcaXY = 0.5;
132 Float_t maxDcaZ = 1.0;
133 //--------------//
134 AliHFEcuts *hfecuts = new AliHFEcuts;
135 hfecuts->CreateStandardCuts();
136 hfecuts->SetMinNClustersTPC(nMinTPCcluster);
137 hfecuts->SetCutITSpixel(AliHFEextraCuts::kFirst);
138 hfecuts->SetCheckITSLayerStatus(kFALSE);
139 hfecuts->SetMaxImpactParam(maxDcaXY, maxDcaZ);
140 SetHFECuts(hfecuts);
bf892a6a 141
6555e2ad 142 fDefaultPID = new AliESDpid();
e3fc062d 143 fHFEpid = new AliHFEpid("PIDforDCAanalysis");
70da6c5a 144
70da6c5a 145}
146
147//____________________________________________________________
148AliAnalysisTaskDCA::AliAnalysisTaskDCA(const AliAnalysisTaskDCA &ref):
faee3b18 149 AliAnalysisTaskSE(ref)
70da6c5a 150 , fPlugins(ref.fPlugins)
faee3b18 151 , fCuts(ref.fCuts)
152 , fDefaultPID(ref.fDefaultPID)
153 , fHFEpid(ref.fHFEpid)
154 , fPIDdetectors(ref.fPIDdetectors)
155 , fPIDstrategy(ref.fPIDstrategy)
70da6c5a 156 , fCFM(ref.fCFM)
157 , fDCA(ref.fDCA)
70da6c5a 158 , fNclustersITS(ref.fNclustersITS)
faee3b18 159 , fMinNprimVtxContrbutor(ref.fMinNprimVtxContrbutor)
70da6c5a 160 , fNEvents(ref.fNEvents)
161 , fResidualList(ref.fResidualList)
162 , fPullList(ref.fPullList)
faee3b18 163 , fDcaList(ref.fDcaList)
164 , fKfDcaList(ref.fKfDcaList)
165 , fMcVertexList(ref.fMcVertexList)
166 , fDataDcaList(ref.fDataDcaList)
167 , fDataVertexList(ref.fDataVertexList)
168 , fDataPullList(ref.fDataPullList)
169 , fMcPidList(ref.fMcPidList)
170 , fDataPidList(ref.fDataPidList)
171 , fHfeDcaList(ref.fHfeDcaList)
172 , fHfeDataDcaList(ref.fHfeDataDcaList)
70da6c5a 173 , fOutput(ref.fOutput)
174{
175 //
176 // Copy Constructor
177 //
6555e2ad 178 AliInfo("Copy Constructor");
faee3b18 179 ref.Copy(*this);
70da6c5a 180}
181
182//____________________________________________________________
183AliAnalysisTaskDCA &AliAnalysisTaskDCA::operator=(const AliAnalysisTaskDCA &ref){
184 //
185 // Assignment operator
186 //
187 if(this == &ref) return *this;
faee3b18 188 AliAnalysisTaskSE::operator=(ref);
70da6c5a 189 fPlugins = ref.fPlugins;
70da6c5a 190 fCuts = ref.fCuts;
faee3b18 191 fDefaultPID = ref.fDefaultPID;
192 fHFEpid = ref.fHFEpid;
193 fPIDdetectors = ref.fPIDdetectors;
194 fPIDstrategy = ref.fPIDstrategy;
70da6c5a 195 fCFM = ref.fCFM;
196 fDCA = ref.fDCA;
70da6c5a 197 fNclustersITS = ref.fNclustersITS;
faee3b18 198 fMinNprimVtxContrbutor = ref.fMinNprimVtxContrbutor;
70da6c5a 199 fNEvents = ref.fNEvents;
70da6c5a 200 fResidualList = ref.fResidualList;
201 fPullList = ref.fPullList;
faee3b18 202 fDcaList = ref.fDcaList;
203 fKfDcaList = ref.fKfDcaList;
204 fMcVertexList = ref.fMcVertexList;
205 fDataDcaList = ref.fDataDcaList;
206 fDataVertexList = ref.fDataVertexList;
207 fDataPullList = ref.fDataPullList;
208 fMcPidList = ref.fMcPidList;
209 fDataPidList = ref.fDataPidList;
210 fHfeDcaList = ref.fHfeDcaList;
211 fHfeDataDcaList = ref.fHfeDataDcaList;
212 fOutput = ref.fOutput;
70da6c5a 213
214 return *this;
215}
216
217//____________________________________________________________
218AliAnalysisTaskDCA::~AliAnalysisTaskDCA(){
219 //
220 // Destructor
221 //
222
faee3b18 223 if(fDefaultPID) delete fDefaultPID;
224 if(fHFEpid) delete fHFEpid;
70da6c5a 225 if(fCFM) delete fCFM;
70da6c5a 226 if(fDCA) delete fDCA;
faee3b18 227 if(fNEvents) delete fNEvents;
70da6c5a 228 if(fResidualList){
229 fResidualList->Clear();
230 delete fResidualList;
231 }
232
233 if(fPullList){
234 fPullList->Clear();
235 delete fPullList;
236 }
237
faee3b18 238 if(fDcaList){
239 fDcaList->Clear();
240 delete fDcaList;
241 }
242 if(fKfDcaList){
243 fKfDcaList->Clear();
244 delete fKfDcaList;
245 }
70da6c5a 246
faee3b18 247 if(fMcVertexList){
248 fMcVertexList->Clear();
249 delete fMcVertexList;
70da6c5a 250 }
faee3b18 251
252 if(fDataDcaList){
253 fDataDcaList->Clear();
254 delete fDataDcaList;
255 }
256
257 if(fDataVertexList){
258 fDataVertexList->Clear();
259 delete fDataVertexList;
260 }
261 if(fDataPullList){
262 fDataPullList->Clear();
263 delete fDataPullList;
264 }
265
266 if(fMcPidList){
267 fMcPidList -> Clear();
268 delete fMcPidList;
269 }
270 if(fDataPidList){
271 fDataPidList -> Clear();
272 delete fDataPidList;
273 }
274
275 if(fHfeDcaList) {
276 fHfeDcaList->Clear();
277 delete fHfeDcaList;
278 }
279
280 if(fHfeDataDcaList) {
281 fHfeDataDcaList->Clear();
282 delete fHfeDataDcaList;
283 }
284
285 if(fOutput){
286 fOutput->Clear();
287 delete fOutput;
70da6c5a 288 }
faee3b18 289
70da6c5a 290}
291
292//____________________________________________________________
faee3b18 293void AliAnalysisTaskDCA::UserCreateOutputObjects(){
70da6c5a 294 // create output objects
295 // fNEvents
296 // residual and pull
bf892a6a 297 //printf("\n=====UserCreateOutputObjects=====\n");
faee3b18 298
299 // Automatic determination of the analysis mode
300 AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
69ac0e6f 301 if(!inputHandler){
302 AliError("NoEvent Handler available");
303 return;
304 }
305
faee3b18 306 if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){
307 SetAODAnalysis();
308 } else {
309 SetESDAnalysis();
310 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
311 SetHasMCData();
312 }
313
70da6c5a 314
faee3b18 315 fNEvents = new TH1I("nEvents", "Number of Events in the Analysis", 5, -0.5, 4.5); // Number of Events neccessary for the analysis and not a QA histogram
70da6c5a 316 if(!fOutput) fOutput = new TList;
317 // Initialize correction Framework and Cuts
318 fCFM = new AliCFManager;
319 MakeParticleContainer();
320 // Temporary fix: Initialize particle cuts with 0x0
321 for(Int_t istep = 0; istep < fCFM->GetParticleContainer()->GetNStep(); istep++)
322 fCFM->SetParticleCutsList(istep, 0x0);
323 if(!fCuts){
324 AliWarning("Cuts not available. Default cuts will be used");
325 fCuts = new AliHFEcuts;
326 fCuts->CreateStandardCuts();
327 }
328
70da6c5a 329 fCuts->Initialize(fCFM);
330
6555e2ad 331 if(!fHFEpid) AliWarning("Hello, fHFEpid is not available");
332 cout<<" Hello this is a cout "<<endl<<endl;
333
bf892a6a 334 if(GetPlugin(kHFEpid)) {
faee3b18 335 fHFEpid->SetHasMCData(HasMCData());
6555e2ad 336 fHFEpid->AddDetector("TOF", 0);
337 fHFEpid->AddDetector("TPC", 1);
338 cout<<endl<<" ---> TPC and TOF added to the PID"<<endl;
339 fHFEpid->ConfigureTPCrejection();
3a72645a 340 fHFEpid->InitializePID();
faee3b18 341 }
70da6c5a 342
343 // dca study----------------------------------
faee3b18 344
345
70da6c5a 346 if(!fDCA) fDCA = new AliHFEdca;
347 if(!fResidualList) fResidualList = new TList();
348 if(!fPullList) fPullList = new TList();
faee3b18 349 if(!fDcaList) fDcaList = new TList();
350 if(!fKfDcaList) fKfDcaList = new TList();
351 if(!fMcVertexList) fMcVertexList = new TList();
352 if(!fDataDcaList) fDataDcaList = new TList();
353 if(!fDataVertexList) fDataVertexList = new TList();
354 if(!fDataPullList) fDataPullList = new TList();
355 if(!fMcPidList) fMcPidList = new TList();
356 if(!fDataPidList) fDataPidList = new TList();
70da6c5a 357
faee3b18 358 if(!fHfeDcaList) fHfeDcaList = new TList();
359 if(!fHfeDataDcaList) fHfeDataDcaList = new TList();
360
361 if(HasMCData()) {
362 if(GetPlugin(kImpactPar) ) {
363 fDCA->CreateHistogramsResidual(fResidualList);
364 fDCA->CreateHistogramsPull(fPullList);
365 fDCA->CreateHistogramsDca(fDcaList);
366 fOutput->AddAt(fResidualList,0);
367 fOutput->AddAt(fPullList,1);
368 fOutput->AddAt(fDcaList,2);
369 }
370 if(GetPlugin(kKFdca)){
371 fDCA->CreateHistogramsKfDca(fKfDcaList);
372 fOutput->AddAt(fDcaList,3);
373 }
6555e2ad 374 if(GetPlugin(kPrimVtx)){//<---
faee3b18 375 fDCA->CreateHistogramsVertex(fMcVertexList);
376 fOutput->AddAt(fMcVertexList,4);
377 }
6555e2ad 378 if(GetPlugin(kCombinedPid)){//<---
faee3b18 379 fDCA->CreateHistogramsPid(fMcPidList);
380 fOutput->AddAt(fMcPidList, 5);
381 }
6555e2ad 382 if(GetPlugin(kHFEpid)){//<---
faee3b18 383 fDCA->CreateHistogramsHfeDca(fHfeDcaList);
384 fOutput->AddAt(fHfeDcaList, 6);
385 }
386 } // mc case
387
388 if(!HasMCData()) {
389
390 if(GetPlugin(kPrimVtx)){
391 fDCA->CreateHistogramsDataVertex(fDataVertexList);
392 fOutput->AddAt(fDataVertexList,0);
393 }
394
395 if(GetPlugin(kCombinedPid)){
396 fDCA->CreateHistogramsDataDca(fDataDcaList);
397 fDCA->CreateHistogramsDataPull(fDataPullList);
398 fDCA->CreateHistogramsDataPid(fDataPidList);
399 fOutput->AddAt(fDataDcaList,1);
400 fOutput->AddAt(fDataPullList,2);
401 fOutput->AddAt(fDataPidList, 3);
402 }
403 if(GetPlugin(kHFEpid)){
404 fDCA->CreateHistogramsHfeDataDca(fHfeDataDcaList);
405 fOutput->AddAt(fHfeDataDcaList, 4);
406 }
407
70da6c5a 408
faee3b18 409
410 } // data case
411
70da6c5a 412}
413
414//____________________________________________________________
faee3b18 415void AliAnalysisTaskDCA::UserExec(Option_t *){
70da6c5a 416 //
417 // Run the analysis
418 //
bf892a6a 419 //printf("\n=====UserExec=====\n");
6555e2ad 420 if(HasMCData()) printf("WITH MC!\n");
70da6c5a 421
422 AliDebug(3, "Processing ESD events");
423
faee3b18 424 if(!fInputEvent){
425 AliError("Reconstructed Event not available");
70da6c5a 426 return;
427 }
faee3b18 428 if(HasMCData()){
429 AliDebug(4, Form("MC Event: %p", fMCEvent));
430 if(!fMCEvent){
431 AliError("No MC Event, but MC Data required");
432 return;
433 }
70da6c5a 434 }
faee3b18 435
70da6c5a 436 if(!fCuts){
437 AliError("HFE cuts not available");
438 return;
439 }
faee3b18 440
441
442 // protection
443 if(IsESDanalysis() && HasMCData()){
444 // Protect against missing MC trees
445 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
69ac0e6f 446 if(!mcH){
447 AliError("No MC Event Handler available");
448 return;
449 }
faee3b18 450 if(!mcH->InitOk()) return;
451 if(!mcH->TreeK()) return;
452 if(!mcH->TreeTR()) return;
453 }
454
455 if(!IsAODanalysis()) {
456 AliESDInputHandler *inH = dynamic_cast<AliESDInputHandler *>(fInputHandler);
69ac0e6f 457 if(!inH){
458 AliError("No ESD Event Handler available");
459 return;
460 }
faee3b18 461 AliESDpid *workingPID = inH->GetESDpid();
462 if(workingPID){
463 AliDebug(1, "Using ESD PID from the input handler");
464 fHFEpid->SetESDpid(workingPID);
465 } else {
466 AliDebug(1, "Using default ESD PID");
6555e2ad 467 fHFEpid->SetESDpid(AliHFEtools::GetDefaultPID(HasMCData()));
faee3b18 468 }
469 ProcessDcaAnalysis();
470 }
471
472
473 PostData(1, fNEvents);
474 PostData(2, fOutput);
475}
476//____________________________________________________________
477void AliAnalysisTaskDCA::ProcessDcaAnalysis(){
70da6c5a 478
bf892a6a 479 //printf("\n=====ProcessDcaAnalysis=====\n");
6555e2ad 480
70da6c5a 481 //
482 // Loop ESD
483 //
69ac0e6f 484 AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
faee3b18 485 if(!fESD){
69ac0e6f 486 AliError("ESD Event required for ESD Analysis");
487 return;
488 }
489
bf892a6a 490 AliMCEvent *fMC = 0x0;
69ac0e6f 491 if(HasMCData()){
492 fMC = dynamic_cast<AliMCEvent*>(fMCEvent);
493 if(!fMC){
494 AliError("MC Event required for Analysis");
faee3b18 495 return;
69ac0e6f 496 }
faee3b18 497 }
498
499 fNEvents->Fill(1); // original event number before cut
500 fDCA->ApplyExtraCuts(fESD,fMinNprimVtxContrbutor); // cut on primVtx contributors
501 fNEvents->Fill(3); // events number after cut
70da6c5a 502 fCFM->SetRecEventInfo(fESD);
faee3b18 503
70da6c5a 504 // event cut level
505 if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
506
bf892a6a 507 AliESDtrack *track = 0x0;
508 AliMCParticle *mctrack = 0x0;
509 AliESDVertex *vtxESDSkip = 0x0;
510 AliHFEpidObject hfetrack;
511
70da6c5a 512 for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
faee3b18 513
70da6c5a 514 track = fESD->GetTrack(itrack);
bf892a6a 515 if(HasMCData()) mctrack = dynamic_cast<AliMCParticle *>(fMC->GetTrack(TMath::Abs(track->GetLabel())));
70da6c5a 516
517 // RecPrim: primary cuts
3a72645a 518 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
70da6c5a 519 // RecKine: ITSTPC cuts
3a72645a 520 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
70da6c5a 521 // HFEcuts: ITS layers cuts
3a72645a 522 if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
70da6c5a 523
faee3b18 524 if(track->GetITSclusters(0)<=fNclustersITS) continue; // require number of ITS clusters
70da6c5a 525
faee3b18 526 // track accepted, do PID
3a72645a 527 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
528 hfetrack.SetRecTrack(track);
529 if(HasMCData()) hfetrack.SetMCTrack(mctrack);
faee3b18 530
bf892a6a 531 //printf("Track %d passed all the cuts!\n",itrack);
532
faee3b18 533 if(HasMCData()){
534 if(GetPlugin(kPrimVtx))
535 fDCA->FillHistogramsVtx(fESD, fMC);
536 if(GetPlugin(kImpactPar))
537 fDCA->FillHistogramsDca(fESD, track, fMC);
538 if(GetPlugin(kKFdca))
539 fDCA->FillHistogramsKfDca(fESD, track, fMC);
540 if(GetPlugin(kCombinedPid))
541 fDCA->FillHistogramsPid(track, fMC);
6555e2ad 542 if(GetPlugin(kHFEpid)) { // data-like
bf892a6a 543 if(fHFEpid->IsSelected(&hfetrack)){
544
6555e2ad 545 // printf("Found an electron in p+p collision! from HFE pid \n");
546 if(!vtxESDSkip){
547 // method from Andrea D 28.05.2010
548 AliVertexerTracks *vertexer = new AliVertexerTracks(fESD->GetMagneticField());
549 vertexer->SetITSMode();
550 vertexer->SetMinClusters(fNclustersITS);
551 Int_t skipped[2];
552 skipped[0] = (Int_t)track->GetID();
553 vertexer->SetSkipTracks(1,skipped);
554 vtxESDSkip = (AliESDVertex*)vertexer->FindPrimaryVertex(fESD);
555 delete vertexer; vertexer = NULL;
556 if(vtxESDSkip->GetNContributors()<fMinNprimVtxContrbutor) continue;
557 }
bf892a6a 558 //printf("\n[ABOUT TO FILL HFE DCA: MC!]\n");
559 fDCA->FillHistogramsHfeDataDca(fESD, track, vtxESDSkip);
560 }
faee3b18 561 } // plugin for hfepid
562 } // MC
563
564 if(!HasMCData()){
565 if(GetPlugin(kPrimVtx))
566 fDCA->FillHistogramsDataVtx(fESD);
567 if(GetPlugin(kCombinedPid)) {
568
569 // method from Andrea D 28.05.2010
570 AliVertexerTracks *vertexer = new AliVertexerTracks(fESD->GetMagneticField());
571 vertexer->SetITSMode();
572 vertexer->SetMinClusters(fNclustersITS);
573 Int_t skipped[2];
574 skipped[0] = (Int_t)track->GetID();
575 vertexer->SetSkipTracks(1,skipped);
6555e2ad 576 vtxESDSkip = (AliESDVertex*)vertexer->FindPrimaryVertex(fESD);
faee3b18 577 delete vertexer; vertexer = NULL;
578 if(vtxESDSkip->GetNContributors()<fMinNprimVtxContrbutor) continue;
579
580 fDCA->FillHistogramsDataDca(fESD, track, vtxESDSkip);
581 fDCA->FillHistogramsDataPid(track);
582 }
583 if(GetPlugin(kHFEpid)) {
584 if(fHFEpid->IsSelected(&hfetrack)) {
585 // printf("Found an electron in p+p collision! from HFE pid \n");
6555e2ad 586 if(!vtxESDSkip){
587 // method from Andrea D 28.05.2010
588 AliVertexerTracks *vertexer = new AliVertexerTracks(fESD->GetMagneticField());
589 vertexer->SetITSMode();
590 vertexer->SetMinClusters(fNclustersITS);
591 Int_t skipped[2];
592 skipped[0] = (Int_t)track->GetID();
593 vertexer->SetSkipTracks(1,skipped);
594 vtxESDSkip = (AliESDVertex*)vertexer->FindPrimaryVertex(fESD);
595 delete vertexer; vertexer = NULL;
596 if(vtxESDSkip->GetNContributors()<fMinNprimVtxContrbutor) continue;
597 }
bf892a6a 598 printf("\n[ABOUT TO FILL HFE DCA: DATA!]\n");
6555e2ad 599 fDCA->FillHistogramsHfeDataDca(fESD, track,vtxESDSkip);
faee3b18 600 }
601 } // plugin for hfepid
602 } // data case
603
604 } // track loop
605
70da6c5a 606}
607
faee3b18 608
70da6c5a 609//____________________________________________________________
610void AliAnalysisTaskDCA::Terminate(Option_t *){
611 //
612 // Terminate not implemented at the moment
613 //
bf892a6a 614 //printf("\n=====Terminate=====\n");
70da6c5a 615
616 if(GetPlugin(kPostProcess)){
617 fOutput = dynamic_cast<TList *>(GetOutputData(1));
618 if(!fOutput){
619 AliError("Results not available");
620 return;
621 }
622 PostProcess();
623 }
624
625}
626
627
628//____________________________________________________________
629void AliAnalysisTaskDCA::Load(TString filename){
6555e2ad 630
bf892a6a 631 //printf("\n=====Load=====\n");
6555e2ad 632
70da6c5a 633 // no need for postprocessing for the moment
634 TFile *input = TFile::Open(filename.Data());
635 if(!input || input->IsZombie()){
636 AliError("Cannot read file");
637 return;
638 }
639
70da6c5a 640 input->Close();
641 delete input;
642
643
644}
645
646//____________________________________________________________
647void AliAnalysisTaskDCA::PostProcess(){
648 // do post processing
649 // should do fitting here for dca resolution
650 // moved to an external macro to do the job
651
bf892a6a 652 //printf("\n=====PostProcess=====\n");
faee3b18 653 Load("HFEdca.root");
70da6c5a 654 TCanvas *c1 = new TCanvas("c1", "number of analyzed events", 300, 400);
655 fNEvents->Draw();
656 c1->SaveAs("temp.png");
657
658}
659
660
661
662
663//____________________________________________________________
664void AliAnalysisTaskDCA::PrintStatus() const {
665
666 //
667 // Print Analysis status
668 //
669 printf("\n\tAnalysis Settings\n\t========================================\n");
faee3b18 670 printf("\t Running on %s\n", !HasMCData()?"p+p collision data":"MC sample");
671 printf("\t Cuts: %s\n", (fCuts != NULL) ? "YES" : "NO");
672 printf("\t Impact parameter analysis is %s\n", GetPlugin(kImpactPar)?"ON":"OFF");
673 printf("\t Using AliKFParticle for analysis? %s\n", GetPlugin(kKFdca)?"ON":"OFF");
674 printf("\t Primary vertex analysis is %s\n", GetPlugin(kPrimVtx)?"ON":"OFF");
675 printf("\t Combined pid analysis is %s\n", GetPlugin(kCombinedPid)?"ON":"OFF");
676 printf("\t HFE pid analysis is %s\n", GetPlugin(kHFEpid)?"ON":"OFF");
677 printf("\t Post process analysis is %s\n", GetPlugin(kPostProcess)?"ON":"OFF");
70da6c5a 678 printf("\t ");
679 printf("\n");
680}
681
682//__________________________________________
683void AliAnalysisTaskDCA::SwitchOnPlugin(Int_t plug){
684 //
685 // Switch on Plugin
686 // Available:
687 // - analyze impact parameter
6555e2ad 688 // - Post Processing
689
690 AliDebug(2,Form("SwitchOnPlugin %d",plug));
691
70da6c5a 692 switch(plug){
faee3b18 693 case kPostProcess:
694 SETBIT(fPlugins, plug);
695 break;
70da6c5a 696 case kImpactPar:
697 SETBIT(fPlugins, plug);
698 break;
faee3b18 699 case kPrimVtx:
700 SETBIT(fPlugins, plug);
701 break;
702 case kCombinedPid:
703 SETBIT(fPlugins, plug);
704 break;
705 case kHFEpid:
706 SETBIT(fPlugins, plug);
707 break;
708 case kKFdca:
70da6c5a 709 SETBIT(fPlugins, plug);
710 break;
711 default:
712 AliError("Unknown Plugin");
713 };
714}
715
716
717//____________________________________________________________
718void AliAnalysisTaskDCA::MakeParticleContainer(){
6555e2ad 719
bf892a6a 720 //printf("\n=====MakeParticleContainer=====\n");
70da6c5a 721 //
722 // Create the particle container (borrowed from AliAnalysisTaskHFE)
723 //
724 const Int_t kNvar = 3 ; //number of variables on the grid:pt,eta, phi
725 const Double_t kPtmin = 0.1, kPtmax = 10.;
726 const Double_t kEtamin = -0.9, kEtamax = 0.9;
727 const Double_t kPhimin = 0., kPhimax = 2. * TMath::Pi();
728
729 //arrays for the number of bins in each dimension
730 Int_t iBin[kNvar];
731 iBin[0] = 40; //bins in pt
732 iBin[1] = 8; //bins in eta
733 iBin[2] = 18; // bins in phi
734
735 //arrays for lower bounds :
736 Double_t* binEdges[kNvar];
737 for(Int_t ivar = 0; ivar < kNvar; ivar++)
738 binEdges[ivar] = new Double_t[iBin[ivar] + 1];
739
740 //values for bin lower bounds
741 for(Int_t i=0; i<=iBin[0]; i++) binEdges[0][i]=(Double_t)TMath::Power(10,TMath::Log10(kPtmin) + (TMath::Log10(kPtmax)-TMath::Log10(kPtmin))/iBin[0]*(Double_t)i);
742 for(Int_t i=0; i<=iBin[1]; i++) binEdges[1][i]=(Double_t)kEtamin + (kEtamax-kEtamin)/iBin[1]*(Double_t)i;
743 for(Int_t i=0; i<=iBin[2]; i++) binEdges[2][i]=(Double_t)kPhimin + (kPhimax-kPhimin)/iBin[2]*(Double_t)i;
744
745 //one "container" for MC
3a72645a 746 const Int_t kNcutStepsESDtrack = AliHFEcuts::kNcutStepsRecTrack + 1;
747 const Int_t kNcutStepsTrack = AliHFEcuts::kNcutStepsMCTrack + kNcutStepsESDtrack;
748 AliCFContainer* container = new AliCFContainer("container","container for tracks", (kNcutStepsTrack + 1 + 2*(kNcutStepsESDtrack + 1)), kNvar, iBin);
70da6c5a 749
750 //setting the bin limits
751 for(Int_t ivar = 0; ivar < kNvar; ivar++)
752 container -> SetBinLimits(ivar, binEdges[ivar]);
753 fCFM->SetParticleContainer(container);
754
755 //create correlation matrix for unfolding
756 Int_t thnDim[2*kNvar];
757 for (int k=0; k<kNvar; k++) {
758 //first half : reconstructed
759 //second half : MC
760 thnDim[k] = iBin[k];
761 thnDim[k+kNvar] = iBin[k];
762 }
763
764
faee3b18 765}
70da6c5a 766
faee3b18 767//____________________________________________________________
768void AliAnalysisTaskDCA::AddPIDdetector(TString detector){
769
770 //
771 // Adding PID detector to the task
772 //
bf892a6a 773 //printf("\n=====AddPIDdetector=====\n");
faee3b18 774
775 if(!fPIDdetectors.Length())
776 fPIDdetectors = detector;
777 else
778 fPIDdetectors += ":" + detector;
70da6c5a 779}
faee3b18 780