1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //////////////////////////////////////////////////////////
17 // Filter the ESDCaloClusters and ESDCaloCells of EMCAL,
18 // PHOS or both, creating the corresponing AODCaloClusters
20 // Keep also the AODHeader information and the vertex.
21 // Keep tracks, v0s, VZERO if requested
22 // Select events containing a cluster or track avobe a given threshold
23 // Copy of AliAnalysisTaskESDfilter.
24 // Author: Gustavo Conesa Balbastre (INFN - Frascati)
25 //////////////////////////////////////////////////////////
28 #include "TGeoManager.h"
31 #include "TInterpreter.h"
34 #include "AliESDEvent.h"
35 #include "AliAODEvent.h"
37 #include "AliVCluster.h"
38 #include "AliVCaloCells.h"
39 #include "AliVEventHandler.h"
40 #include "AliAnalysisManager.h"
41 #include "AliInputEventHandler.h"
44 #include "AliEMCALRecoUtils.h"
45 #include "AliEMCALGeometry.h"
47 #include "AliAnalysisTaskCaloFilter.h"
49 ClassImp(AliAnalysisTaskCaloFilter)
51 ////////////////////////////////////////////////////////////////////////
53 AliAnalysisTaskCaloFilter::AliAnalysisTaskCaloFilter():
54 AliAnalysisTaskSE("CaloFilterTask"),
55 fCaloFilter(0), fEventSelection(),
56 fAcceptAllMBEvent(kFALSE), fCorrect(kFALSE),
57 fEMCALGeo(0x0), fEMCALGeoName("EMCAL_COMPLETE12SMV1"),
58 fEMCALRecoUtils(new AliEMCALRecoUtils),
59 fLoadEMCALMatrices(kFALSE), //fLoadPHOSMatrices(kFALSE),
60 fGeoMatrixSet(kFALSE),
61 fConfigName(""), fFillAODFile(kTRUE),
62 fFillTracks(kFALSE), fFillHybridTracks(kFALSE),
63 fFillAllVertices(kFALSE), fFillv0s(kFALSE),
65 fEMCALEnergyCut(0.), fEMCALNcellsCut (0),
66 fPHOSEnergyCut(0.), fPHOSNcellsCut (0),
68 fVzCut(100.), fEvent(0x0),
69 fESDEvent(0x0), fAODEvent(0x0)
71 // Default constructor
73 fEventSelection[0] = kFALSE;
74 fEventSelection[1] = kFALSE;
75 fEventSelection[2] = kFALSE;
77 for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = 0 ;
78 //for(Int_t i = 0; i < 5 ; i++) fPHOSMatrix[i] = 0 ;
82 //__________________________________________________
83 AliAnalysisTaskCaloFilter::AliAnalysisTaskCaloFilter(const char* name):
84 AliAnalysisTaskSE(name),
85 fCaloFilter(0), fEventSelection(),
86 fAcceptAllMBEvent(kFALSE), fCorrect(kFALSE),
87 fEMCALGeo(0x0), fEMCALGeoName("EMCAL_COMPLETE12SMV1"),
88 fEMCALRecoUtils(new AliEMCALRecoUtils),
89 fLoadEMCALMatrices(kFALSE), //fLoadPHOSMatrices(kFALSE),
90 fGeoMatrixSet(kFALSE),
91 fConfigName(""), fFillAODFile(kTRUE),
92 fFillTracks(kFALSE), fFillHybridTracks(kFALSE),
93 fFillAllVertices(kFALSE), fFillv0s(kFALSE),
95 fEMCALEnergyCut(0.), fEMCALNcellsCut(0),
96 fPHOSEnergyCut(0.), fPHOSNcellsCut(0),
98 fVzCut(100.), fEvent(0x0),
99 fESDEvent(0x0), fAODEvent(0x0)
103 fEventSelection[0] = kFALSE;
104 fEventSelection[1] = kFALSE;
105 fEventSelection[2] = kFALSE;
107 for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = 0 ;
108 //for(Int_t i = 0; i < 5 ; i++) fPHOSMatrix[i] = 0 ;
112 //__________________________________________________
113 AliAnalysisTaskCaloFilter::~AliAnalysisTaskCaloFilter()
117 if(fEMCALGeo) delete fEMCALGeo;
118 if(fEMCALRecoUtils) delete fEMCALRecoUtils;
122 //_____________________________________________
123 Bool_t AliAnalysisTaskCaloFilter::AcceptEvent()
125 // Define conditions to accept the event to be filtered
127 if(!AcceptEventVertex()) return kFALSE;
129 Bool_t eventSel = kFALSE;
131 Bool_t isMB = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
133 if ( isMB && fAcceptAllMBEvent ) eventSel = kTRUE; // accept any MB event
135 else if( fEventSelection[0] && AcceptEventEMCAL() ) eventSel = kTRUE; // accept event depending on EMCAL activity
137 else if( fEventSelection[1] && AcceptEventPHOS () ) eventSel = kTRUE; // accept event depending on PHOS activity
139 else if( fEventSelection[2] && AcceptEventTrack() ) eventSel = kTRUE; // accept event depending on Track activity
145 //__________________________________________________
146 Bool_t AliAnalysisTaskCaloFilter::AcceptEventEMCAL()
148 // Accept event given there is a EMCAL cluster with enough energy, and not noisy, exotic
150 if(fCaloFilter==kPHOS) return kTRUE; // accept
152 if(fEMCALEnergyCut <= 0) return kTRUE; // accept
154 Int_t nCluster = InputEvent() -> GetNumberOfCaloClusters();
155 AliVCaloCells * caloCell = InputEvent() -> GetEMCALCells();
156 Int_t bc = InputEvent() -> GetBunchCrossNumber();
158 for(Int_t icalo = 0; icalo < nCluster; icalo++)
160 AliVCluster *clus = (AliVCluster*) (InputEvent()->GetCaloCluster(icalo));
162 if( ( clus->IsEMCAL() ) && ( clus->GetNCells() > fEMCALNcellsCut ) && ( clus->E() > fEMCALEnergyCut ) &&
163 fEMCALRecoUtils->IsGoodCluster(clus,fEMCALGeo,caloCell,bc))
166 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventEMCAL() - Accept : E %2.2f > %2.2f, nCells %d > %d \n",
167 clus->E(), fEMCALEnergyCut, clus->GetNCells(), fEMCALNcellsCut);
174 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventEMCAL() - Reject \n");
180 //_________________________________________________
181 Bool_t AliAnalysisTaskCaloFilter::AcceptEventPHOS()
183 // Accept event given there is a PHOS cluster with enough energy and not noisy/exotic
185 if(fCaloFilter==kEMCAL) return kTRUE; // accept
187 if(fPHOSEnergyCut <= 0) return kTRUE; // accept
189 Int_t nCluster = InputEvent() -> GetNumberOfCaloClusters();
191 for(Int_t icalo = 0; icalo < nCluster; icalo++)
193 AliVCluster *clus = (AliVCluster*) (InputEvent()->GetCaloCluster(icalo));
195 if( ( clus->IsPHOS() ) && ( clus->GetNCells() > fPHOSNcellsCut ) && ( clus->E() > fPHOSEnergyCut ))
198 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventPHOS() - Accept : E %2.2f > %2.2f, nCells %d > %d \n",
199 clus->E(), fPHOSEnergyCut, clus->GetNCells(), fPHOSNcellsCut);
206 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventPHOS() - Reject \n");
212 //__________________________________________________
213 Bool_t AliAnalysisTaskCaloFilter::AcceptEventTrack()
215 // Accept event if there is a track avobe a certain pT
217 if(fTrackPtCut <= 0) return kTRUE; // accept
219 Double_t pTrack[3] = {0,0,0};
221 for (Int_t nTrack = 0; nTrack < fEvent->GetNumberOfTracks(); ++nTrack)
223 AliVTrack *track = (AliVTrack*) fEvent->GetTrack(nTrack);
225 // Select only hybrid tracks?
226 if(fAODEvent && fFillHybridTracks && !((AliAODTrack*)track)->IsHybridGlobalConstrainedGlobal()) continue;
228 track->GetPxPyPz(pTrack) ;
230 TLorentzVector momentum(pTrack[0],pTrack[1],pTrack[2],0);
232 if(momentum.Pt() > fTrackPtCut)
234 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventTrack() - Accept : pT %2.2f > %2.2f \n",
235 momentum.Pt(), fTrackPtCut);
242 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventTrack() - Reject \n");
248 //___________________________________________________
249 Bool_t AliAnalysisTaskCaloFilter::AcceptEventVertex()
251 // Accept event with good vertex
254 InputEvent()->GetPrimaryVertex()->GetXYZ(v) ;
257 if(TMath::Abs(v[2]) > fVzCut)
259 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::AcceptEventVertex() - Vz Reject : vz %2.2f > %2.2f\n",v[2],fVzCut);
264 return CheckForPrimaryVertex();
267 //_______________________________________________________
268 Bool_t AliAnalysisTaskCaloFilter::CheckForPrimaryVertex()
270 //Check if the vertex was well reconstructed, copy from v0Reader of conversion group
271 //It only works for ESDs
276 // Check that the vertex is not (0,0,0)
278 InputEvent()->GetPrimaryVertex()->GetXYZ(v) ;
280 if(TMath::Abs(v[2]) < 1e-6 &&
281 TMath::Abs(v[1]) < 1e-6 &&
282 TMath::Abs(v[0]) < 1e-6 )
284 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::CheckForPrimaryVertex() - Reject v(0,0,0) \n");
293 if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors() > 0)
298 if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors() < 1)
301 if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors() > 0)
303 //cout<<"spd vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
307 if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors() < 1)
309 // cout<<"bad vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
310 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::CheckForPrimaryVertex() - Reject, GetPrimaryVertexSPD()->GetNContributors() < 1 \n");
316 if (fDebug > 0) printf("AliAnalysisTaskCaloFilter::CheckForPrimaryVertex() - Reject, GetPrimaryVertexTracks()->GetNContributors() > 1 \n");
322 //__________________________________________________
323 void AliAnalysisTaskCaloFilter::CorrectionsInEMCAL()
325 //If EMCAL, and requested, correct energy, position ...
327 //Need to do this in a separate loop before filling the ESDs because of the track matching recalculations
329 if(fCorrect && (fCaloFilter==kEMCAL || fCaloFilter==kBoth) )
333 if(fLoadEMCALMatrices)
335 printf("AliAnalysisTaskCaloFilter::UserExec() - Load user defined EMCAL geometry matrices\n");
336 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
337 if(fEMCALMatrix[mod]){
339 fEMCALMatrix[mod]->Print();
340 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod) ;
345 else if(!gGeoManager)
347 printf("AliAnalysisTaskCaloFilter::UserExec() - Get geo matrices from data\n");
348 //Still not implemented in AOD, just a workaround to be able to work at least with ESDs
349 if(!strcmp(InputEvent()->GetName(),"AliAODEvent"))
352 printf("AliAnalysisTaskCaloFilter Use ideal geometry, values geometry matrix not kept in AODs.\n");
356 if(DebugLevel() > 1) printf("AliAnalysisTaskCaloFilter Load Misaligned matrices. \n");
357 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent()) ;
360 printf("AliAnalysisTaskCaloFilter::UserExec() - This event does not contain ESDs?");
363 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
365 //if(DebugLevel() > 1)
366 esd->GetEMCALMatrix(mod)->Print();
367 if(esd->GetEMCALMatrix(mod)) fEMCALGeo->SetMisalMatrix(esd->GetEMCALMatrix(mod),mod) ;
371 }//Load matrices from Data
376 Int_t nCaloClus = InputEvent()->GetNumberOfCaloClusters();
378 for (Int_t iClust=0; iClust<nCaloClus; ++iClust)
380 AliVCluster * cluster = InputEvent()->GetCaloCluster(iClust);
382 if(cluster->IsPHOS()) continue ;
384 Float_t position[]={0,0,0};
386 printf("Check cluster %d for bad channels and close to border\n",cluster->GetID());
387 if(fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo,cluster->GetCellsAbsId(), cluster->GetNCells())) continue;
391 printf("Filter, before : i %d, E %f, dispersion %f, m02 %f, m20 %f, distToBad %f\n",iClust,cluster->E(),
392 cluster->GetDispersion(),cluster->GetM02(),cluster->GetM20(), cluster->GetDistanceToBadChannel());
393 cluster->GetPosition(position);
394 printf("Filter, before : i %d, x %f, y %f, z %f\n",cluster->GetID(), position[0], position[1], position[2]);
397 //Recalculate distance to bad channels, if new list of bad channels provided
398 fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, InputEvent()->GetEMCALCells(), cluster);
400 if(fEMCALRecoUtils->IsRecalibrationOn())
402 fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, cluster, InputEvent()->GetEMCALCells());
403 fEMCALRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, InputEvent()->GetEMCALCells(),cluster);
404 fEMCALRecoUtils->RecalculateClusterPID(cluster);
407 fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, InputEvent()->GetEMCALCells(),cluster);
411 printf("Filter, after : i %d, E %f, dispersion %f, m02 %f, m20 %f, distToBad %f\n",cluster->GetID(),cluster->E(),
412 cluster->GetDispersion(),cluster->GetM02(),cluster->GetM20(), cluster->GetDistanceToBadChannel());
413 cluster->GetPosition(position);
414 printf("Filter, after : i %d, x %f, y %f, z %f\n",cluster->GetID(), position[0], position[1], position[2]);
417 cluster->SetE(fEMCALRecoUtils->CorrectClusterEnergyLinearity(cluster));
421 //Recalculate track-matching
422 fEMCALRecoUtils->FindMatches(InputEvent(),0,fEMCALGeo);
424 } // corrections in EMCAL
427 //________________________________________________
428 void AliAnalysisTaskCaloFilter::FillAODCaloCells()
430 // Fill EMCAL/PHOS cell info
433 if ((fCaloFilter==kBoth || fCaloFilter==kEMCAL) && fEvent->GetEMCALCells())
434 { // protection against missing ESD information
435 AliVCaloCells &eventEMcells = *(fEvent->GetEMCALCells());
436 Int_t nEMcell = eventEMcells.GetNumberOfCells() ;
438 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
439 aodEMcells.CreateContainer(nEMcell);
440 aodEMcells.SetType(AliVCaloCells::kEMCALCell);
441 Double_t calibFactor = 1.;
442 for (Int_t iCell = 0; iCell < nEMcell; iCell++)
444 Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
445 fEMCALGeo->GetCellIndex(eventEMcells.GetCellNumber(iCell),imod,iTower,iIphi,iIeta);
446 fEMCALGeo->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
448 if(fCorrect && fEMCALRecoUtils->IsRecalibrationOn())
450 calibFactor = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi);
453 if(!fEMCALRecoUtils->GetEMCALChannelStatus(imod, ieta, iphi))
454 { //Channel is not declared as bad
455 aodEMcells.SetCell(iCell,eventEMcells.GetCellNumber(iCell),eventEMcells.GetAmplitude(iCell)*calibFactor,
456 eventEMcells.GetTime(iCell),eventEMcells.GetMCLabel(iCell),eventEMcells.GetEFraction(iCell));
457 //printf("GOOD channel\n");
461 aodEMcells.SetCell(iCell,eventEMcells.GetCellNumber(iCell),0,-1,-1,0);
462 //printf("BAD channel\n");
469 if ((fCaloFilter==kBoth || fCaloFilter==kPHOS) && fEvent->GetPHOSCells())
470 { // protection against missing ESD information
471 AliVCaloCells &eventPHcells = *(fEvent->GetPHOSCells());
472 Int_t nPHcell = eventPHcells.GetNumberOfCells() ;
474 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
475 aodPHcells.CreateContainer(nPHcell);
476 aodPHcells.SetType(AliVCaloCells::kPHOSCell);
478 for (Int_t iCell = 0; iCell < nPHcell; iCell++)
480 aodPHcells.SetCell(iCell,eventPHcells.GetCellNumber(iCell),eventPHcells.GetAmplitude(iCell),
481 eventPHcells.GetTime(iCell),eventPHcells.GetMCLabel(iCell),eventPHcells.GetEFraction(iCell));
489 //___________________________________________________
490 void AliAnalysisTaskCaloFilter::FillAODCaloClusters()
492 // Fill the AOD with caloclusters
494 // Access to the AOD container of clusters
496 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
500 Int_t nCaloClus = fEvent->GetNumberOfCaloClusters();
501 for (Int_t iClust=0; iClust<nCaloClus; ++iClust)
504 AliVCluster * cluster = fEvent->GetCaloCluster(iClust);
506 //Check which calorimeter information we want to keep.
508 if(fCaloFilter!=kBoth)
510 if (fCaloFilter==kPHOS && cluster->IsEMCAL()) continue ;
511 else if(fCaloFilter==kEMCAL && cluster->IsPHOS()) continue ;
514 // Get original residuals, in case of previous recalculation, reset them
515 Float_t dR = cluster->GetTrackDx();
516 Float_t dZ = cluster->GetTrackDz();
519 printf("Original residuals : dZ %f, dR %f\n ",dZ, dR);
521 //--------------------------------------------------------------
522 //If EMCAL and corrections done, get the new matching parameters, do not copy noisy clusters
523 if(cluster->IsEMCAL() && fCorrect)
526 printf("Check cluster %d for bad channels and close to border\n",cluster->GetID());
528 if(fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo,cluster->GetCellsAbsId(), cluster->GetNCells())) continue;
530 if(fEMCALRecoUtils->IsExoticCluster(cluster, InputEvent()->GetEMCALCells(),InputEvent()->GetBunchCrossNumber())) continue;
532 fEMCALRecoUtils->GetMatchedResiduals(cluster->GetID(),dR,dZ);
533 cluster->SetTrackDistance(dR,dZ);
538 if(cluster->IsEMCAL()) printf("EMCAL Track-Cluster Residuals : dZ %f, dR %f\n ",dZ, dR);
539 if(cluster->IsPHOS()) printf("PHOS Track-Cluster Residuals : dZ %f, dR %f\n ",dZ, dR);
542 //--------------------------------------------------------------
546 Int_t id = cluster->GetID();
547 Float_t energy = cluster->E();
548 cluster->GetPosition(posF);
550 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++])
551 AliAODCaloCluster(id,
552 cluster->GetNLabels(),
553 cluster->GetLabels(),
559 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
560 cluster->GetDispersion(),
561 cluster->GetM20(), cluster->GetM02(),
563 cluster->GetNExMax(),cluster->GetTOF()) ;
565 caloCluster->SetPIDFromESD(cluster->GetPID());
566 caloCluster->SetNCells(cluster->GetNCells());
567 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
568 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
569 caloCluster->SetTrackDistance(dR, dZ);
573 printf("Filter, aod : i %d, E %f, dispersion %f, m02 %f, m20 %f\n",caloCluster->GetID(),caloCluster->E(),
574 caloCluster->GetDispersion(),caloCluster->GetM02(),caloCluster->GetM20());
575 caloCluster->GetPosition(posF);
576 printf("Filter, aod : i %d, x %f, y %f, z %f\n",caloCluster->GetID(), posF[0], posF[1], posF[2]);
579 //Matched tracks, just to know if there was any match, the track pointer is useless if tracks not stored
580 if(TMath::Abs(dR) < 990 && TMath::Abs(dZ) < 990)
581 { //Default value in PHOS 999, in EMCAL 1024, why?
582 caloCluster->AddTrackMatched(new AliAODTrack);
584 // TO DO, in case Tracks available, think how to put the matched track in AOD
587 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots
591 //__________________________________________________
592 void AliAnalysisTaskCaloFilter::FillAODCaloTrigger()
594 // AOD CaloTrigger copy
596 AliAODCaloTrigger* triggerEM = AODEvent()->GetCaloTrigger("EMCAL");
597 AliAODCaloTrigger* triggerPH = AODEvent()->GetCaloTrigger("PHOS");
602 if(fCaloFilter==kBoth || fCaloFilter==kPHOS) *triggerPH = *(fAODEvent->GetCaloTrigger("PHOS"));
604 if(fCaloFilter==kBoth || fCaloFilter==kEMCAL) *triggerEM = *(fAODEvent->GetCaloTrigger("EMCAL"));
611 //______________________________________________
612 void AliAnalysisTaskCaloFilter::FillAODHeader()
616 AliAODHeader* header = AODEvent()->GetHeader();
621 *header = *(fAODEvent->GetHeader());
625 if(!fESDEvent) return;
629 header->SetRunNumber(fEvent->GetRunNumber());
631 TTree* tree = fInputHandler->GetTree();
634 TFile* file = tree->GetCurrentFile();
635 if (file) header->SetESDFileName(file->GetName());
638 header->SetBunchCrossNumber(fEvent->GetBunchCrossNumber());
639 header->SetOrbitNumber(fEvent->GetOrbitNumber());
640 header->SetPeriodNumber(fEvent->GetPeriodNumber());
641 header->SetEventType(fEvent->GetEventType());
644 if(fEvent->GetCentrality())
646 header->SetCentrality(new AliCentrality(*(fEvent->GetCentrality())));
650 header->SetCentrality(0);
654 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
655 header->SetFiredTriggerClasses(fESDEvent->GetFiredTriggerClasses());
656 header->SetTriggerMask(fEvent->GetTriggerMask());
657 header->SetTriggerCluster(fEvent->GetTriggerCluster());
658 header->SetL0TriggerInputs(fESDEvent->GetHeader()->GetL0TriggerInputs());
659 header->SetL1TriggerInputs(fESDEvent->GetHeader()->GetL1TriggerInputs());
660 header->SetL2TriggerInputs(fESDEvent->GetHeader()->GetL2TriggerInputs());
662 header->SetMagneticField(fEvent->GetMagneticField());
663 //header->SetMuonMagFieldScale(fESDEvent->GetCurrentDip()/6000.);
665 header->SetZDCN1Energy(fEvent->GetZDCN1Energy());
666 header->SetZDCP1Energy(fEvent->GetZDCP1Energy());
667 header->SetZDCN2Energy(fEvent->GetZDCN2Energy());
668 header->SetZDCP2Energy(fEvent->GetZDCP2Energy());
669 header->SetZDCEMEnergy(fEvent->GetZDCEMEnergy(0),fEvent->GetZDCEMEnergy(1));
671 Float_t diamxy[2]={fEvent->GetDiamondX(),fEvent->GetDiamondY()};
673 fEvent->GetDiamondCovXY(diamcov);
674 header->SetDiamond(diamxy,diamcov);
675 header->SetDiamondZ(fESDEvent->GetDiamondZ(),fESDEvent->GetSigma2DiamondZ());
679 //_____________________________________________
680 void AliAnalysisTaskCaloFilter::FillAODTracks()
684 if(!fFillTracks) return;
686 AliAODTrack* aodTrack(0x0);
688 Double_t pos[3] = { 0. };
689 Double_t covTr[21]= { 0. };
690 Double_t pid[10] = { 0. };
691 Double_t p[3] = { 0. };
696 //TClonesArray* inTracks = fAODEvent ->GetTracks();
697 TClonesArray* ouTracks = AODEvent()->GetTracks();
698 //new (ouTracks) TClonesArray(*inTracks);
700 //printf("N tracks %d\n",fAODEvent->GetNumberOfTracks());
701 Int_t nCopyTrack = 0;
702 for (Int_t nTrack = 0; nTrack < fAODEvent->GetNumberOfTracks(); ++nTrack)
704 AliAODTrack *track = fAODEvent->GetTrack(nTrack);
706 // Select only hybrid tracks?
707 if(fFillHybridTracks && !track->IsHybridGlobalConstrainedGlobal()) continue;
709 // Remove PID object to save space
710 //track->SetDetPID(0x0);
712 //new((*ouTracks)[nCopyTrack++]) AliAODTrack(*track);
715 Bool_t isDCA = track->GetPosition(pos);
716 track->GetCovMatrix(covTr);
719 AliAODVertex* primVertex = (AliAODVertex*) AODEvent()->GetVertices()->At(0); // primary vertex, copied previously!!!
721 aodTrack = new((*ouTracks)[nCopyTrack++]) AliAODTrack(
730 track->GetITSClusterMap(),
733 track->GetUsedForVtxFit(),
734 track->GetUsedForPrimVtxFit(),
735 (AliAODTrack::AODTrk_t) track->GetType(),
736 track->GetFilterMap(),
737 track->Chi2perNDF());
740 aodTrack->SetIsHybridGlobalConstrainedGlobal(track->IsHybridGlobalConstrainedGlobal());
741 aodTrack->SetIsHybridTPCConstrainedGlobal (track->IsHybridTPCConstrainedGlobal());
742 aodTrack->SetIsGlobalConstrained (track->IsGlobalConstrained());
743 aodTrack->SetIsTPCConstrained (track->IsTPCConstrained());
745 aodTrack->SetTPCFitMap (track->GetTPCFitMap());
746 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
747 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
749 aodTrack->SetChi2MatchTrigger(track->GetChi2MatchTrigger());
751 // set the DCA values to the AOD track
753 aodTrack->SetPxPyPzAtDCA(track->PxAtDCA(),track->PyAtDCA(),track->PzAtDCA());
754 aodTrack->SetXYAtDCA (track->XAtDCA() ,track->YAtDCA());
756 aodTrack->SetFlags (track->GetFlags());
757 aodTrack->SetTPCPointsF (track->GetTPCNclsF());
761 if(track->IsEMCAL()) aodTrack->SetEMCALcluster(track->GetEMCALcluster());
762 if(track->IsPHOS()) aodTrack->SetPHOScluster (track->GetPHOScluster());
763 aodTrack->SetTrackPhiEtaOnEMCal( track->GetTrackPhiOnEMCal(), track->GetTrackPhiOnEMCal() );
767 //printf("Final N tracks %d\n",nCopyTrack);
774 //_________________________________________
775 void AliAnalysisTaskCaloFilter::FillAODv0s()
777 // Copy v0s (use if you know what you do, use quite a lot of memory)
779 if(!fFillv0s) return;
784 TClonesArray* inv0 = fAODEvent ->GetV0s();
785 TClonesArray* ouv0 = AODEvent()->GetV0s();
787 //new (ouv0s) TClonesArray(*inv0s);
789 Int_t allv0s = inv0->GetEntriesFast();
791 for (Int_t nv0s = 0; nv0s < allv0s; ++nv0s)
793 AliAODv0 *v0 = (AliAODv0*)inv0->At(nv0s);
795 new((*ouv0)[nv0s]) AliAODv0(*v0);
803 //____________________________________________
804 void AliAnalysisTaskCaloFilter::FillAODVZERO()
808 if(!fFillVZERO) return;
810 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
812 if(fESDEvent) *vzeroData = *(fESDEvent->GetVZEROData());
813 else *vzeroData = *(fAODEvent->GetVZEROData());
817 //_______________________________________________
818 void AliAnalysisTaskCaloFilter::FillAODVertices()
822 // set arrays and pointers
825 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
830 TClonesArray* inVertices = fAODEvent ->GetVertices();
831 TClonesArray* ouVertices = AODEvent()->GetVertices();
833 //new (ouVertices) TClonesArray(*inVertices);
835 //Keep only the first 3 vertices if not requested
836 Int_t allVertices = inVertices->GetEntriesFast();
838 //printf("n Vertices %d\n",allVertices);
840 if(!fFillAllVertices)
842 if(allVertices > 3) allVertices = 3;
845 //printf("Final n Vertices %d\n",allVertices);
847 for (Int_t nVertices = 0; nVertices < allVertices; ++nVertices)
849 AliAODVertex *vertex = (AliAODVertex*)inVertices->At(nVertices);
851 new((*ouVertices)[nVertices]) AliAODVertex(*vertex);
857 if(!fESDEvent) return;
861 // Access to the AOD container of vertices
863 TClonesArray &vertices = *(AODEvent()->GetVertices());
865 // Add primary vertex. The primary tracks will be defined
866 // after the loops on the composite objects (v0, cascades, kinks)
867 fEvent ->GetPrimaryVertex()->GetXYZ(pos);
868 fESDEvent->GetPrimaryVertex()->GetCovMatrix(covVtx);
869 Float_t chi = fESDEvent->GetPrimaryVertex()->GetChi2toNDF();
871 AliAODVertex * primary = new(vertices[jVertices++])
872 AliAODVertex(pos, covVtx, chi, NULL, -1, AliAODVertex::kPrimary);
873 primary->SetName(fEvent->GetPrimaryVertex()->GetName());
874 primary->SetTitle(fEvent->GetPrimaryVertex()->GetTitle());
878 //____________________________________
879 void AliAnalysisTaskCaloFilter::Init()
881 //Init analysis with configuration macro if available
883 if(gROOT->LoadMacro(fConfigName) >=0)
885 printf("Configure analysis with %s\n",fConfigName.Data());
887 AliAnalysisTaskCaloFilter *filter = (AliAnalysisTaskCaloFilter*)gInterpreter->ProcessLine("ConfigCaloFilter()");
889 fEMCALGeoName = filter->fEMCALGeoName;
890 fLoadEMCALMatrices = filter->fLoadEMCALMatrices;
891 fFillAODFile = filter->fFillAODFile;
892 fFillTracks = filter->fFillTracks;
893 fFillHybridTracks = filter->fFillHybridTracks;
894 fFillv0s = filter->fFillv0s;
895 fFillVZERO = filter->fFillVZERO;
896 fFillAllVertices = filter->fFillAllVertices;
897 fEMCALRecoUtils = filter->fEMCALRecoUtils;
898 fConfigName = filter->fConfigName;
899 fCaloFilter = filter->fCaloFilter;
900 fEventSelection[0] = filter->fEventSelection[0];
901 fEventSelection[1] = filter->fEventSelection[1];
902 fEventSelection[2] = filter->fEventSelection[2];
903 fAcceptAllMBEvent = filter->fAcceptAllMBEvent;
904 fCorrect = filter->fCorrect;
905 fEMCALEnergyCut = filter->fEMCALEnergyCut;
906 fEMCALNcellsCut = filter->fEMCALNcellsCut;
907 fPHOSEnergyCut = filter->fPHOSEnergyCut;
908 fPHOSNcellsCut = filter->fPHOSNcellsCut;
909 fTrackPtCut = filter->fTrackPtCut;
910 fVzCut = filter->fVzCut;
912 for(Int_t i = 0; i < 12; i++) fEMCALMatrix[i] = filter->fEMCALMatrix[i] ;
916 //_________________________________________
917 void AliAnalysisTaskCaloFilter::PrintInfo()
921 printf("AnalysisCaloFilter::PrintInfo() \n");
923 printf("\t Not only filter, correct Clusters? %d\n",fCorrect);
924 printf("\t Calorimeter Filtering Option ? %d\n",fCaloFilter);
926 //printf("\t Use handmade geo matrices? EMCAL %d, PHOS %d\n",fLoadEMCALMatrices, fLoadPHOSMatrices);
927 printf("\t Use handmade geo matrices? EMCAL %d, PHOS 0\n",fLoadEMCALMatrices);
929 printf("\t Fill: AOD file? %d Tracks? %d; all Vertex? %d; v0s? %d; VZERO ? %d\n",
930 fFillAODFile,fFillTracks,fFillAllVertices, fFillv0s, fFillVZERO);
932 printf("\t Event Selection based : EMCAL? %d, PHOS? %d Tracks? %d - Accept all MB? %d\n",
933 fEventSelection[0],fEventSelection[1],fEventSelection[2],fAcceptAllMBEvent);
935 printf("\t \t EMCAL E > %2.2f, EMCAL nCells >= %d, PHOS E > %2.2f, PHOS nCells >= %d, Track pT > %2.2f, |vz| < %2.2f\n",
936 fEMCALEnergyCut,fEMCALNcellsCut,fPHOSEnergyCut,fPHOSNcellsCut, fTrackPtCut,fVzCut);
939 //_______________________________________________________
940 void AliAnalysisTaskCaloFilter::UserCreateOutputObjects()
944 fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
948 //____________________________________________________________
949 void AliAnalysisTaskCaloFilter::UserExec(Option_t */*option*/)
951 // Execute analysis for current event
952 // Copy input ESD or AOD header, vertex, CaloClusters and CaloCells to output AOD
955 printf("CaloFilter: Analysing event # %d\n", (Int_t)Entry());
957 fEvent = InputEvent();
958 fAODEvent = dynamic_cast<AliAODEvent*> (fEvent);
959 fESDEvent = dynamic_cast<AliESDEvent*> (fEvent);
963 printf("AliAnalysisTaskCaloFilter::UserExec - This event does not contain Input?");
967 // printf("Start processing : %s\n",fAODEvent->GetFiredTriggerClasses().Data());
971 if(!AcceptEvent()) return ;
973 //Magic line to write events to file
975 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(fFillAODFile);
980 if(fFillv0s) nVertices = fEvent->GetNumberOfV0s();
981 Int_t nCaloClus = fEvent->GetNumberOfCaloClusters();
982 Int_t nTracks = fEvent->GetNumberOfTracks();
984 AODEvent()->ResetStd(nTracks, nVertices, 0, 0, 0, nCaloClus, 0, 0);
994 FillAODVertices(); // Do it before the track filtering to have the reference to the vertex
1003 CorrectionsInEMCAL();
1006 FillAODCaloClusters();
1012 FillAODCaloTrigger();
1014 //printf("Filtered event, end processing : %s\n",fAODEvent->GetFiredTriggerClasses().Data());