From fd1d0cb9f78e1bd1ff8597a3eba41ee2a842ba6b Mon Sep 17 00:00:00 2001 From: martinez Date: Wed, 10 Mar 2010 20:44:35 +0000 Subject: [PATCH] Update of the Xiaoming code for pp900 first muon analysis. Fixing wanirngs and violtions (Xiaoming) --- PWG3/PWG3muonLinkDef.h | 8 +- PWG3/libPWG3muon.pkg | 12 +- PWG3/muon/AddTaskMuonsHF.C | 46 +- PWG3/muon/AliAODMuonPair.cxx | 68 --- PWG3/muon/AliAODMuonPair.h | 42 -- PWG3/muon/AliAODMuonTrack.cxx | 115 ----- PWG3/muon/AliAODMuonTrack.h | 43 -- PWG3/muon/AliAnalysisTaskSEMuonsHF.cxx | 572 ++++++++----------------- PWG3/muon/AliAnalysisTaskSEMuonsHF.h | 128 ++---- PWG3/muon/AliDimuInfoStoreMC.cxx | 150 +++++++ PWG3/muon/AliDimuInfoStoreMC.h | 56 +++ PWG3/muon/AliDimuInfoStoreRD.cxx | 142 ++++++ PWG3/muon/AliDimuInfoStoreRD.h | 58 +++ PWG3/muon/AliMCMuonPair.cxx | 96 ----- PWG3/muon/AliMCMuonPair.h | 38 -- PWG3/muon/AliMCMuonTrack.cxx | 407 ------------------ PWG3/muon/AliMCMuonTrack.h | 83 ---- PWG3/muon/AliMuonInfoStoreMC.cxx | 490 +++++++++++++++++++++ PWG3/muon/AliMuonInfoStoreMC.h | 101 +++++ PWG3/muon/AliMuonInfoStoreRD.cxx | 218 ++++++++++ PWG3/muon/AliMuonInfoStoreRD.h | 92 ++++ PWG3/muon/AliMuonsHFHeader.cxx | 451 ++++++++++++++++++- PWG3/muon/AliMuonsHFHeader.h | 136 ++++-- 23 files changed, 2096 insertions(+), 1456 deletions(-) delete mode 100644 PWG3/muon/AliAODMuonPair.cxx delete mode 100644 PWG3/muon/AliAODMuonPair.h delete mode 100644 PWG3/muon/AliAODMuonTrack.cxx delete mode 100644 PWG3/muon/AliAODMuonTrack.h create mode 100644 PWG3/muon/AliDimuInfoStoreMC.cxx create mode 100644 PWG3/muon/AliDimuInfoStoreMC.h create mode 100644 PWG3/muon/AliDimuInfoStoreRD.cxx create mode 100644 PWG3/muon/AliDimuInfoStoreRD.h delete mode 100644 PWG3/muon/AliMCMuonPair.cxx delete mode 100644 PWG3/muon/AliMCMuonPair.h delete mode 100644 PWG3/muon/AliMCMuonTrack.cxx delete mode 100644 PWG3/muon/AliMCMuonTrack.h create mode 100644 PWG3/muon/AliMuonInfoStoreMC.cxx create mode 100644 PWG3/muon/AliMuonInfoStoreMC.h create mode 100644 PWG3/muon/AliMuonInfoStoreRD.cxx create mode 100644 PWG3/muon/AliMuonInfoStoreRD.h diff --git a/PWG3/PWG3muonLinkDef.h b/PWG3/PWG3muonLinkDef.h index 2dbac690274..0394fd7cfc3 100644 --- a/PWG3/PWG3muonLinkDef.h +++ b/PWG3/PWG3muonLinkDef.h @@ -20,11 +20,11 @@ #pragma link C++ class AliAnalysisTaskCreateMixedDimuons+; #pragma link C++ class AliAnalysisTaskMuonAODCreation+; #pragma link C++ class AliAnalysisTaskMuonDistributions+; +#pragma link C++ class AliMuonInfoStoreRD+; +#pragma link C++ class AliDimuInfoStoreRD+; +#pragma link C++ class AliMuonInfoStoreMC+; +#pragma link C++ class AliDimuInfoStoreMC+; #pragma link C++ class AliMuonsHFHeader+; -#pragma link C++ class AliMCMuonTrack+; -#pragma link C++ class AliMCMuonPair+; -#pragma link C++ class AliAODMuonTrack+; -#pragma link C++ class AliAODMuonPair+; #pragma link C++ class AliAnalysisTaskSEMuonsHF+; #pragma link C++ class AliAnalysisTaskDimuonCFContainerBuilder+; #endif diff --git a/PWG3/libPWG3muon.pkg b/PWG3/libPWG3muon.pkg index 639a22ef73a..867e28da26d 100644 --- a/PWG3/libPWG3muon.pkg +++ b/PWG3/libPWG3muon.pkg @@ -16,12 +16,12 @@ SRCS:= muon/AliAnalysisTaskESDMuonFilter.cxx \ muon/AliAnalysisTaskCreateMixedDimuons.cxx \ muon/AliAnalysisTaskMuonAODCreation.cxx \ muon/AliAnalysisTaskMuonDistributions.cxx \ - muon/AliMuonsHFHeader.cxx \ - muon/AliMCMuonTrack.cxx \ - muon/AliMCMuonPair.cxx \ - muon/AliAODMuonTrack.cxx \ - muon/AliAODMuonPair.cxx \ - muon/AliAnalysisTaskSEMuonsHF.cxx \ + muon/AliMuonInfoStoreRD.cxx \ + muon/AliDimuInfoStoreRD.cxx \ + muon/AliMuonInfoStoreMC.cxx \ + muon/AliDimuInfoStoreMC.cxx \ + muon/AliMuonsHFHeader.cxx \ + muon/AliAnalysisTaskSEMuonsHF.cxx \ muon/AliAnalysisTaskDimuonCFContainerBuilder.cxx diff --git a/PWG3/muon/AddTaskMuonsHF.C b/PWG3/muon/AddTaskMuonsHF.C index 9f2acf89d46..9c56da8a4b6 100644 --- a/PWG3/muon/AddTaskMuonsHF.C +++ b/PWG3/muon/AddTaskMuonsHF.C @@ -28,38 +28,36 @@ AliAnalysisTaskSEMuonsHF* AddTaskMuonsHF(Int_t mode=0, Bool_t isMC=kFALSE, Bool_ } } - // set cuts for single muon track selection - Double_t cuts[10]={-1., // 0, min of 3-momentum - 999999., // 1, max of 3-momnentum - 1., // 2, PtMin - 999999., // 3, PtMax - -4., // 4, EtaMin - -2.5, // 5, EtaMax - -1., // 6, DCAmin - 10., // 7, DCAmax - 0.5, // 8, for trigger - 3.5. // 9, for trigger - }; - + // set cuts for events or muons selection + Double_t cutsEvsH[3] ={-999999.0, // low limit of Ncontrs + 999999.0, // up limit of |vz| + 999999.0}; // up limit of vt + Double_t cutsMuon[10]={-999999.0, // 0, min of 3-momentum + 999999.0, // 1, max of 3-momnentum + -999999.0, // 2, PtMin + 999999.0, // 3, PtMax + -999999.0, // 4, EtaMin + 999999.0, // 5, EtaMax + -999999.0, // 6, DCAmin + 999999.0, // 7, DCAmax + -999999.0, // 8, for trigger + 999999.0}; // 9, for trigger + Double_t cutsDimu[10]={-999999.0, 999999.0, // single muon cuts used for dimuon selection + -999999.0, 999999.0, + -999999.0, 999999.0, + -999999.0, 999999.0, + -999999.0, 999999.0}; AliAnalysisTaskSEMuonsHF *taskMuonsHF = new AliAnalysisTaskSEMuonsHF("MuonsHF Analysis Task"); taskMuonsHF->SetAnaMode(mode); taskMuonsHF->SetIsUseMC(isMC); taskMuonsHF->SetIsOutputTree(isTree); - taskMuonsHF->SetSingleMuonCuts(cuts); + taskMuonsHF->SetEvsHCuts(cutsEvsH); + taskMuonsHF->SetMuonCuts(cutsMuon); + taskMuonsHF->SetDimuCuts(cutsDimu); mgr->AddTask(taskMuonsHF); - - TString outputfile = AliAnalysisManager::GetCommonFileName(); - outputfile += ":PWG3Muon_MuonHF"; - - AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("listHisEventH",TList::Class(),AliAnalysisManager::kOutputContainer,outputfile); - AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("listHisSingleMuon",TList::Class(),AliAnalysisManager::kOutputContainer,outputfile); - AliAnalysisDataContainer *coutput4 = mgr->CreateContainer("listHisDimuon",TList::Class(),AliAnalysisManager::kOutputContainer,outputfile); mgr->ConnectInput(taskMuonsHF,0,mgr->GetCommonInputContainer()); if (isTree) mgr->ConnectOutput(taskMuonsHF,0,mgr->GetCommonOutputContainer()); - mgr->ConnectOutput(taskMuonsHF,1,coutput2); - mgr->ConnectOutput(taskMuonsHF,2,coutput3); - mgr->ConnectOutput(taskMuonsHF,3,coutput4); return taskMuonsHF; } diff --git a/PWG3/muon/AliAODMuonPair.cxx b/PWG3/muon/AliAODMuonPair.cxx deleted file mode 100644 index 7b294dc042b..00000000000 --- a/PWG3/muon/AliAODMuonPair.cxx +++ /dev/null @@ -1,68 +0,0 @@ -#include "AliAODMuonTrack.h" -#include "AliAODMuonPair.h" - -ClassImp(AliAODMuonPair) - -//----------------------------------------------------------------------------- -AliAODMuonPair::AliAODMuonPair() : -TObject(), -fP(), -fCharge(0) -{ - // - // default constructor - // - for (Int_t i=0; i<2; i++) { - fTrk[i] = 0; - fTrigger[i] = 0; - fDca[i] = 0.; - fChi2[i] = 0.; - } -} - -//----------------------------------------------------------------------------- -AliAODMuonPair::AliAODMuonPair(AliAODMuonTrack *trk0, AliAODMuonTrack *trk1) : -TObject(), -fP(), -fCharge(0) -{ - // - // default constructor - // - for (Int_t i=0; i<2; i++) { - fTrigger[i]=0; fDca[i]=0.; fChi2[i]=0.; - } - - fTrk[0] = trk0; - fTrk[1] = trk1; - FillPairInfo(); -} - -//----------------------------------------------------------------------------- -AliAODMuonPair::~AliAODMuonPair() -{ - // - // destructor - // -} - -//----------------------------------------------------------------------------- -void AliAODMuonPair::FillPairInfo() -{ - AliAODMuonTrack *trk0 = (AliAODMuonTrack*)fTrk[0].GetObject(); - AliAODMuonTrack *trk1 = (AliAODMuonTrack*)fTrk[1].GetObject(); - - fP = trk0->GetP() + trk1->GetP(); - fCharge = trk0->GetCharge() + trk1->GetCharge(); - - fTrigger[0] = trk0->GetTrigger(); - fTrigger[1] = trk1->GetTrigger(); - - fDca[0] = trk0->GetDCA(); - fDca[1] = trk1->GetDCA(); - - fChi2[0] = trk0->GetChi2(); - fChi2[1] = trk1->GetChi2(); - - return; -} diff --git a/PWG3/muon/AliAODMuonPair.h b/PWG3/muon/AliAODMuonPair.h deleted file mode 100644 index 85c1fe7ab38..00000000000 --- a/PWG3/muon/AliAODMuonPair.h +++ /dev/null @@ -1,42 +0,0 @@ -#ifndef ALIAODMUONPAIR_H -#define ALIAODMUONPAIR_H - -#include -#include -#include - -#include "AliAODMuonTrack.h" - -class AliAODMuonPair : public TObject { - public: - - AliAODMuonPair(); - AliAODMuonPair(AliAODMuonTrack *trk0, AliAODMuonTrack *trk1); - virtual ~AliAODMuonPair(); - - AliAODMuonTrack* GetTrack(Int_t i) const { return (i<2 ? (AliAODMuonTrack*)(fTrk[i].GetObject()) : 0x0); } - - TLorentzVector GetP() const { return fP; } - Int_t GetCharge() const { return fCharge; } - Int_t GetMuMatchTrigger(Int_t i) const { return (i<2 ? fTrigger[i] : -9999 ); } - Double_t GetMuDCA(Int_t i) const { return (i<2 ? fDca[i] : -9999.); } - Double_t GetMuChi2(Int_t i) const { return (i<2 ? fChi2[i] : -9999.); } - - protected: - - TRef fTrk[2]; - void FillPairInfo(); - - private: - - - TLorentzVector fP; - Int_t fCharge; - Int_t fTrigger[2]; - Double_t fDca[2]; - Double_t fChi2[2]; - - ClassDef(AliAODMuonPair, 3); -}; - -#endif diff --git a/PWG3/muon/AliAODMuonTrack.cxx b/PWG3/muon/AliAODMuonTrack.cxx deleted file mode 100644 index dd3be40b587..00000000000 --- a/PWG3/muon/AliAODMuonTrack.cxx +++ /dev/null @@ -1,115 +0,0 @@ -#include -#include -#include -#include - -#include "AliAODTrack.h" -#include "AliESDMuonTrack.h" -#include "AliMCMuonTrack.h" -#include "AliAODMuonTrack.h" - -ClassImp(AliAODMuonTrack) - -//----------------------------------------------------------------------------- -AliAODMuonTrack::AliAODMuonTrack() : -TObject(), -fP(), -fCharge(0), -fTrigger(0), -fDca(0.), -fChi2(0.), -fCentr(0.) -{ - // - // default constructor - // -} - -//----------------------------------------------------------------------------- -AliAODMuonTrack::AliAODMuonTrack(AliAODTrack *trk) : -TObject(), -fP(), -fCharge(0), -fTrigger(0), -fDca(0.), -fChi2(0.), -fCentr(0.) -{ - // - // default constructor - // - this->FillTrackInfo(trk); -} - -//----------------------------------------------------------------------------- -AliAODMuonTrack::AliAODMuonTrack(AliESDMuonTrack *trk) : -TObject(), -fP(), -fCharge(0), -fTrigger(0), -fDca(0.), -fChi2(0.), -fCentr(0.) -{ - // - // default constructor - // - this->FillTrackInfo(trk); -} - -//----------------------------------------------------------------------------- -AliAODMuonTrack::~AliAODMuonTrack() -{ - // - // destructor - // -} - -//----------------------------------------------------------------------------- -void AliAODMuonTrack::FillTrackInfo(AliAODTrack *trk) -{ - Double_t mMu = TDatabasePDG::Instance()->GetParticle(13)->Mass(); - Double_t px = trk->Px(); - Double_t py = trk->Py(); - Double_t pz = trk->Pz(); - Double_t energy = TMath::Sqrt(mMu*mMu + px*px + py*py + pz*pz); - fP.SetPxPyPzE(px,py,pz,energy); - fCharge = trk->Charge(); - fTrigger = trk->GetMatchTrigger(); - fDca = trk->DCA(); - fChi2 = trk->Chi2perNDF(); - return; -} - -//----------------------------------------------------------------------------- -void AliAODMuonTrack::FillTrackInfo(AliESDMuonTrack *trk) -{ - trk->LorentzP(fP); - fCharge = trk->Charge(); - fTrigger = trk->GetMatchTrigger(); - fDca = trk->GetDCA(); - fChi2 = trk->GetChi2()/(2.*trk->GetNHit()-5.); - return; -} - -//----------------------------------------------------------------------------- -Bool_t AliAODMuonTrack::SelectSingleMuon(Double_t cuts[10]) -{ - TLorentzVector lorentzP = this->GetP(); - Double_t p = lorentzP.P(); - if (pcuts[1])) return kFALSE; - - Double_t pt = lorentzP.Pt(); - if (ptcuts[3]) return kFALSE; - - Double_t eta = lorentzP.Eta(); - if (etacuts[5]) return kFALSE; - - Double_t dca = this->GetDCA(); - if (dcacuts[7]) return kFALSE; - - Int_t trigger = this->GetTrigger(); - if (triggercuts[9]) return kFALSE; - - return kTRUE; -} diff --git a/PWG3/muon/AliAODMuonTrack.h b/PWG3/muon/AliAODMuonTrack.h deleted file mode 100644 index 371b84ac4cb..00000000000 --- a/PWG3/muon/AliAODMuonTrack.h +++ /dev/null @@ -1,43 +0,0 @@ -#ifndef ALIAODMUONTRACK_H -#define ALIAODMUONTRACK_H - -#include -#include -#include - -#include "AliAODTrack.h" -#include "AliESDMuonTrack.h" - -class AliAODMuonTrack : public TObject { - public: - - AliAODMuonTrack(); - AliAODMuonTrack(AliAODTrack *trk); - AliAODMuonTrack(AliESDMuonTrack *trk); - virtual ~AliAODMuonTrack(); - - Bool_t SelectSingleMuon(Double_t cuts[10]); - - TLorentzVector GetP() const { return fP; } - Int_t GetCharge() const { return fCharge; } - Int_t GetTrigger() const { return fTrigger; } - Double_t GetDCA() const { return fDca; } - Double_t GetChi2() const { return fChi2; } - Double_t GetCentr() const { return fCentr; } - - private: - - void FillTrackInfo(AliAODTrack *trk); - void FillTrackInfo(AliESDMuonTrack *trk); - - TLorentzVector fP; - Short_t fCharge; - Int_t fTrigger; - Double_t fDca; - Double_t fChi2; - Double_t fCentr; // used for PbPb conllisions - - ClassDef(AliAODMuonTrack, 3); -}; - -#endif diff --git a/PWG3/muon/AliAnalysisTaskSEMuonsHF.cxx b/PWG3/muon/AliAnalysisTaskSEMuonsHF.cxx index 09a0a3e194a..9063a6bbb0e 100644 --- a/PWG3/muon/AliAnalysisTaskSEMuonsHF.cxx +++ b/PWG3/muon/AliAnalysisTaskSEMuonsHF.cxx @@ -1,26 +1,52 @@ -#include -#include -#include +/************************************************************************** + * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +///////////////////////////////////////////////////////////// +// +// AliAnalysisTaskSE for the single muon and dimuon from HF analysis, +// using the classes AliMuonsHFHeader, +// AliMuonInfoStoreRD, +// AliDimuInfoStoreRD, +// AliMuonInfoStoreMC, +// AliDimuInfoStoreMC. +// +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +///////////////////////////////////////////////////////////// + +#include #include +#include #include "AliAnalysisManager.h" #include "AliMCEventHandler.h" #include "AliAODMCParticle.h" -#include "AliVVertex.h" -#include "AliVHeader.h" #include "AliVEvent.h" #include "AliESDEvent.h" #include "AliAODEvent.h" #include "AliESDMuonTrack.h" #include "AliAODTrack.h" #include "AliMuonsHFHeader.h" -#include "AliAODMuonTrack.h" -#include "AliAODMuonPair.h" -#include "AliMCMuonTrack.h" -#include "AliMCMuonPair.h" -#include "AliAnalysisTaskSE.h" +#include "AliMuonInfoStoreRD.h" +#include "AliMuonInfoStoreMC.h" +#include "AliDimuInfoStoreRD.h" +#include "AliDimuInfoStoreMC.h" #include "AliAnalysisTaskSEMuonsHF.h" +class AliAnalysisTaskSE; + ClassImp(AliAnalysisTaskSEMuonsHF) //_____________________________________________________________________________ @@ -30,25 +56,15 @@ fAnaMode(0), fIsOutputTree(kFALSE), fIsUseMC(kFALSE), fHeader(0), -fMuTrkClArr(0), -fMuPairClArr(0), -fListHisAtEvLevel(0), -fListHisSingleMuon(0), -fListHisDimuon(0) +fMuonClArr(0), +fDimuClArr(0), +fListHisHeader(0), +fListHisMuon(), +fListHisDimu(0) { // // Default constructor // - fSingleMuonCuts[0] = -1.; - fSingleMuonCuts[1] = 999999.; - fSingleMuonCuts[2] = -1.; - fSingleMuonCuts[3] = 999999.; - fSingleMuonCuts[4] = -1.; - fSingleMuonCuts[5] = 999999.; - fSingleMuonCuts[6] = -1.; - fSingleMuonCuts[7] = 999999.; - fSingleMuonCuts[8] = -1.; - fSingleMuonCuts[9] = 4.; } //_____________________________________________________________________________ @@ -58,29 +74,15 @@ fAnaMode(0), fIsOutputTree(kFALSE), fIsUseMC(kFALSE), fHeader(0), -fMuTrkClArr(0), -fMuPairClArr(0), -fListHisAtEvLevel(0), -fListHisSingleMuon(0), -fListHisDimuon(0) +fMuonClArr(0), +fDimuClArr(0), +fListHisHeader(0), +fListHisMuon(0), +fListHisDimu(0) { // // Constructor // - fSingleMuonCuts[0] = -1.; - fSingleMuonCuts[1] = 999999.; - fSingleMuonCuts[2] = -1.; - fSingleMuonCuts[3] = 999999.; - fSingleMuonCuts[4] = -1.; - fSingleMuonCuts[5] = 999999.; - fSingleMuonCuts[6] = -1.; - fSingleMuonCuts[7] = 999999.; - fSingleMuonCuts[8] = -1.; - fSingleMuonCuts[9] = 4.; - - DefineOutput(1, TList::Class()); - DefineOutput(2, TList::Class()); - DefineOutput(3, TList::Class()); } //_____________________________________________________________________________ @@ -89,41 +91,71 @@ AliAnalysisTaskSEMuonsHF::~AliAnalysisTaskSEMuonsHF() // // Default destructor // - if (fHeader) { delete fHeader; fHeader=NULL; } - if (fMuTrkClArr) { delete fMuTrkClArr; fMuTrkClArr =NULL; } - if (fMuPairClArr) { delete fMuPairClArr; fMuPairClArr=NULL; } + if (fHeader) { delete fHeader; fHeader =NULL; } + if (fMuonClArr) { delete fMuonClArr; fMuonClArr =NULL; } + if (fDimuClArr) { delete fDimuClArr; fDimuClArr =NULL; } - if (fListHisSingleMuon) { delete fListHisSingleMuon; fListHisSingleMuon=NULL; } - if (fListHisDimuon) { delete fListHisDimuon; fListHisDimuon=NULL; } + if (fListHisHeader) { delete fListHisHeader; fListHisHeader=NULL; } + if (fListHisMuon) { delete [] fListHisMuon; fListHisMuon =NULL; } + if (fListHisDimu) { delete [] fListHisDimu; fListHisDimu =NULL; } } //_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::UserCreateOutputObjects() +void AliAnalysisTaskSEMuonsHF::Init() { - CreateOutoutHistosAtEvnetLevel(); - if (fAnaMode!=2) CreateOutputHistosSingleMuon(); - if (fAnaMode!=1) CreateOutputHistosDimuon(); + // Initialization + // Setting and initializing the running mode and status + AliMuonsHFHeader::SetAnaMode(fAnaMode); + AliMuonsHFHeader::SetIsMC(fIsUseMC); if (!fHeader) { fHeader = new AliMuonsHFHeader(); - fHeader->SetName("MuonsHFHeader"); + fHeader->SetName(AliMuonsHFHeader::StdBranchName()); } - if (!fMuTrkClArr) { - if (fIsUseMC) fMuTrkClArr = new TClonesArray("AliMCMuonTrack", 0); - else fMuTrkClArr = new TClonesArray("AliAODMuonTrack", 0); - fMuTrkClArr->SetName("MuonTrack"); + + if (!fMuonClArr) { + if (fIsUseMC) { + fMuonClArr = new TClonesArray("AliMuonInfoStoreMC", 0); + fMuonClArr->SetName(AliMuonInfoStoreMC::StdBranchName()); + } else { + fMuonClArr = new TClonesArray("AliMuonInfoStoreRD", 0); + fMuonClArr->SetName(AliMuonInfoStoreRD::StdBranchName()); + } } - if (fAnaMode!=1 && !fMuPairClArr) { - if (fIsUseMC) fMuPairClArr = new TClonesArray("AliMCMuonPair", 0); - else fMuPairClArr = new TClonesArray("AliAODMuonPair", 0); - fMuPairClArr->SetName("MuonPair"); + if (fAnaMode!=1 && !fDimuClArr) { + if (fIsUseMC) { + fDimuClArr = new TClonesArray("AliDimuInfoStoreMC", 0); + fDimuClArr->SetName(AliDimuInfoStoreMC::StdBranchName()); + } else { + fDimuClArr = new TClonesArray("AliDimuInfoStoreRD", 0); + fDimuClArr->SetName(AliDimuInfoStoreRD::StdBranchName()); + } + } + + return; +} + +//_____________________________________________________________________________ +void AliAnalysisTaskSEMuonsHF::UserCreateOutputObjects() +{ + // Create the output container + + if (!fListHisHeader) fListHisHeader = new TList(); + if (fAnaMode!=2 && !fListHisMuon) { + if (fIsUseMC) fListHisMuon = new TList[AliMuonInfoStoreMC::NSources()]; + else fListHisMuon = new TList(); } + if (fAnaMode!=1 && !fListHisDimu) { + if (fIsUseMC) fListHisDimu = new TList[AliDimuInfoStoreMC::NSources()]; + else fListHisDimu = new TList(); + } + fHeader->CreateHistograms(fListHisHeader, fListHisMuon, fListHisDimu); if (fIsOutputTree) { AddAODBranch("AliMuonsHFHeader", &fHeader); - AddAODBranch("TClonesArray", &fMuTrkClArr); - if (fAnaMode!=1) AddAODBranch("TClonesArray", &fMuPairClArr); + AddAODBranch("TClonesArray", &fMuonClArr); + if (fAnaMode!=1) AddAODBranch("TClonesArray", &fDimuClArr); } return; @@ -132,12 +164,13 @@ void AliAnalysisTaskSEMuonsHF::UserCreateOutputObjects() //_____________________________________________________________________________ void AliAnalysisTaskSEMuonsHF::UserExec(Option_t *) { + // Execute analysis for current event: + // muon event header & (di)muon info store + AliVEvent *event = dynamic_cast(InputEvent()); //if (AODEvent() && IsStandardAOD()) event = dynamic_cast(AODEvent()); TString evName = event->IsA()->GetName(); Bool_t isAOD = ((evName=="AliAODEvent") ? kTRUE : kFALSE); - fHeader->SetHeader(((AliVHeader*)event->GetHeader())); - fHeader->SetVertex(((AliVVertex*)event->GetPrimaryVertex())); event = 0x0; Int_t ntrks = 0; @@ -148,373 +181,136 @@ void AliAnalysisTaskSEMuonsHF::UserExec(Option_t *) if (isAOD) { aod = dynamic_cast(InputEvent()); //if (AODEvent() && IsStandardAOD()) aod = dynamic_cast(AODEvent()); - if (!aod) { - AliError("AOD event not found. Nothing done!"); - return; - } + if (!aod) { AliError("AOD event not found. Nothing done!"); return; } if (fIsUseMC) { mcClArr = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()); - if (!mcClArr) { - AliError("MC Array not found. Nothing done!"); - return; - } + if (!mcClArr) { AliError("MC Array not found. Nothing done!"); return; } } + fHeader->SetEvent(aod); ntrks = aod->GetNTracks(); } else { esd = dynamic_cast(InputEvent()); - if (!esd) { - AliError("ESD event not found. Nothing done!"); - return; - } + if (!esd) { AliError("ESD event not found. Nothing done!"); return; } if (fIsUseMC) { - mcH = (AliMCEventHandler*) - ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); - if (!mcH) { - AliError("MC Handler not found. Nothing done!"); - return; - } + mcH = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); + if (!mcH) { AliError("MC Handler not found. Nothing done!"); return; } } + fHeader->SetEvent(esd); ntrks = esd->GetNumberOfMuonTracks(); } - fMuTrkClArr->Delete(); - TClonesArray &muTrkRef = *fMuTrkClArr; - Int_t countN = fMuTrkClArr->GetEntriesFast(); + fHeader->FillHistosEventH(fListHisHeader); + + fMuonClArr->Delete(); + TClonesArray &muonRef = *fMuonClArr; + Int_t countN = fMuonClArr->GetEntriesFast(); - Int_t theMult = 0; - AliAODTrack *trkAOD = 0; - AliESDMuonTrack *trkESD = 0; - AliAODMuonTrack *track = 0; - AliMCMuonTrack *trkMC = 0; + AliAODTrack *trkAOD = 0; + AliESDMuonTrack *trkESD = 0; + AliMuonInfoStoreRD *trkRD = 0; + AliMuonInfoStoreMC *trkMC = 0; for (Int_t itrk=0; itrkGetTrack(itrk); - if (!trkAOD->IsMuonTrack()) { trkAOD=0; continue; } - if (fIsUseMC) trkMC = new AliMCMuonTrack(trkAOD, mcClArr); - else track = new AliAODMuonTrack(trkAOD); + if (!trkAOD->IsMuonTrack()) { trkAOD=0; trkRD=0; trkMC=0; continue; } + if (fIsUseMC) trkMC = new AliMuonInfoStoreMC(trkAOD, mcClArr); + else trkRD = new AliMuonInfoStoreRD(trkAOD); + trkAOD = 0; } else { trkESD = (AliESDMuonTrack*)esd->GetMuonTrack(itrk); - if (!trkESD->ContainTrackerData()) { trkESD=0; continue; } - if (fIsUseMC) trkMC = new AliMCMuonTrack(trkESD, esd, mcH); - else track = new AliAODMuonTrack(trkESD); + if (!trkESD->ContainTrackerData()) { trkESD=0; trkRD=0; trkMC=0; continue; } + if (fIsUseMC) trkMC = new AliMuonInfoStoreMC(trkESD, esd, mcH); + else trkRD = new AliMuonInfoStoreRD(trkESD); + trkESD = 0; } - if (fAnaMode!=2) { - if (track) FillDistributionsSingleMuon(track); - else FillDistributionsSingleMuon((AliAODMuonTrack*)trkMC, trkMC->GetSource()); + if (trkRD) { + if (fAnaMode!=2) fHeader->FillHistosMuonRD(fListHisMuon, trkRD); + new(muonRef[countN++]) AliMuonInfoStoreRD(*trkRD); } - - if (fIsOutputTree || fAnaMode!=1) { - if (track) new(muTrkRef[countN++]) AliAODMuonTrack(*track); - else new(muTrkRef[countN++]) AliMCMuonTrack(*trkMC); + if (trkMC) { + if (fAnaMode!=2) fHeader->FillHistosMuonMC(fListHisMuon, trkMC); + new(muonRef[countN++]) AliMuonInfoStoreMC(*trkMC); } - theMult++; - trkESD = 0; - trkAOD = 0; - if (track) { delete track; track=0; } + if (trkRD) { delete trkRD; trkRD=0; } if (trkMC) { delete trkMC; trkMC=0; } } // end loop of all tracks - fHeader->SetMultSingleMuon(theMult); - FillDistributionsAtEventLeavel(); - PostData(1, fListHisAtEvLevel); - if (fAnaMode!=2) PostData(2, fListHisSingleMuon); + if (fAnaMode==1) return; + fDimuClArr->Delete(); + countN = fDimuClArr->GetEntriesFast(); + TClonesArray &dimuRef = *fDimuClArr; - theMult = 0; - fMuPairClArr->Delete(); - TClonesArray &muPairRef = *fMuPairClArr; - countN = fMuPairClArr->GetEntriesFast(); - AliAODMuonPair *pair = 0; - AliMCMuonPair *pairMC = 0; - ntrks = fMuTrkClArr->GetEntriesFast(); + AliDimuInfoStoreRD *dimuRD = 0; + AliDimuInfoStoreMC *dimuMC = 0; + ntrks = fMuonClArr->GetEntriesFast(); for (Int_t itrk=0; itrkAt(itrk), (AliMCMuonTrack*)fMuTrkClArr->At(jtrk)); + dimuMC = new AliDimuInfoStoreMC((AliMuonInfoStoreMC*)fMuonClArr->At(itrk), (AliMuonInfoStoreMC*)fMuonClArr->At(jtrk)); else - pair = new AliAODMuonPair((AliAODMuonTrack*)fMuTrkClArr->At(itrk), (AliAODMuonTrack*)fMuTrkClArr->At(jtrk)); - - if (pair) FillDistributionsDimuon(pair); - else FillDistributionsDimuon(pairMC, pairMC->GetSource()); + dimuRD = new AliDimuInfoStoreRD((AliMuonInfoStoreRD*)fMuonClArr->At(itrk), (AliMuonInfoStoreRD*)fMuonClArr->At(jtrk)); - if (fIsOutputTree) { - if (pair) new(muPairRef[countN++]) AliAODMuonPair(*pair); - else new(muPairRef[countN++]) AliMCMuonPair(*pairMC); + if (dimuRD) { + fHeader->FillHistosDimuRD(fListHisDimu, dimuRD); + if (fIsOutputTree) new(dimuRef[countN++]) AliDimuInfoStoreRD(*dimuRD); + } + if (dimuMC) { + fHeader->FillHistosDimuMC(fListHisDimu, dimuMC); + if (fIsOutputTree) new(dimuRef[countN++]) AliDimuInfoStoreMC(*dimuMC); } - theMult++; - if (pair) { delete pair; pair =0; } - if (pairMC) { delete pairMC; pairMC=0; } - } // end 2nd loop of muon traks + if (dimuRD) { delete dimuRD; dimuRD=0; } + if (dimuMC) { delete dimuMC; dimuMC=0; } + } // end 2nd loop of muon tracks } // end 1st loop of muon tracks - fHeader->SetMultDimuon(theMult); - PostData(3, fListHisDimuon); return; } //_____________________________________________________________________________ void AliAnalysisTaskSEMuonsHF::Terminate(Option_t *) { - printf("This is Terminate!\n"); - // adding the fitting - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::CreateOutoutHistosAtEvnetLevel() -{ - if (!fListHisAtEvLevel) fListHisAtEvLevel = new TList(); - Int_t nbins[kNHistEv]; - Double_t xlow[kNHistEv], xup[kNHistEv]; - TString name[kNHistEv], axis[kNHistEv], unit[kNHistEv]; - nbins[kHistVx]=200; xlow[kHistVx] =-50.; xup[kHistVx] =50.; name[kHistVx] = "vx"; axis[kHistVx] ="v_{x}"; unit[kHistVx]="cm"; - nbins[kHistVy]=200; xlow[kHistVy] =-50.; xup[kHistVy] =50.; name[kHistVy] = "vy"; axis[kHistVy] ="v_{y}"; unit[kHistVy]="cm"; - nbins[kHistVz]=200; xlow[kHistVz] =-50.; xup[kHistVz] =50.; name[kHistVz] = "vz"; axis[kHistVz] ="v_{z}"; unit[kHistVz]="cm"; - nbins[kHistMult]=9; xlow[kHistMult]=-0.5; xup[kHistMult]=8.5; name[kHistMult]="mult"; axis[kHistMult]="mult"; unit[kHistMult]=""; - - TH1F *histo = 0; - TString histName, histTitle; - for (Int_t ihist=0; ihistGetXaxis()->SetTitle(Form("%s %s", axis[ihist].Data(), unit[ihist].Data())); - histo->Sumw2(); - fListHisAtEvLevel->AddAt(histo, ihist); - } // end loop over histos - - return; -} - - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::CreateOutputHistosSingleMuon() -{ - if (!fListHisSingleMuon) fListHisSingleMuon = new TList(); - - Int_t nbins[kNHistMu]; - Double_t xlow[kNHistMu], xup[kNHistMu]; - TString name[kNHistMu], axis[kNHistMu], unit[kNHistMu]; - nbins[kHistP] =300 ; xlow[kHistP]= 0. ; xup[kHistP] =150.; name[kHistP] ="P" ; axis[kHistP] ="P" ; unit[kHistP] ="[GeV/c]"; - nbins[kHistPt] =120 ; xlow[kHistPt]= 0. ; xup[kHistPt] = 30.; name[kHistPt] ="Pt" ; axis[kHistPt] ="p_{t}" ; unit[kHistPt] ="[GeV/c]"; - nbins[kHistEta]=50 ; xlow[kHistEta]=-4. ; xup[kHistEta]=-2.5; name[kHistEta]="Eta"; axis[kHistEta]="#eta" ; unit[kHistEta]=""; - nbins[kHistDca]=1000; xlow[kHistDca]= 0. ; xup[kHistDca]=500.; name[kHistDca]="Dca"; axis[kHistDca]="DCA" ; unit[kHistDca]="cm"; - nbins[kHistTrg]= 4; xlow[kHistTrg]=-0.5; xup[kHistTrg]= 3.5; name[kHistTrg]="trg"; axis[kHistTrg]="Trigger"; unit[kHistTrg]=""; - - TH1F *histo = 0; - TString histName, histTitle; - Int_t index=-1, lowBound=0; - for (Int_t ihist=0; ihistGetXaxis()->SetTitle(Form("%s %s", axis[ihist].Data(), unit[ihist].Data())); - histo->Sumw2(); - fListHisSingleMuon->AddAt(histo, index); - } // end loop over histos - - if (fIsUseMC) CreateOutputHistosSingleMuonMC(nbins, xlow, xup, name, axis, unit); - TH2F *hVzPt = new TH2F("hVzPt", "v_{z} vs. p_{t}", 300, -30., 30., 120, 0., 30); - hVzPt->GetXaxis()->SetTitle("v_{z} [cm]"); - hVzPt->GetYaxis()->SetTitle("p_{t} [GeV/c]"); - fListHisSingleMuon->AddLast(hVzPt); - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::CreateOutputHistosDimuon() -{ - if (!fListHisDimuon) fListHisDimuon = new TList(); - - TString dimuName[kNMuMus], dimuTitle[kNMuMus]; - dimuName[kMuNMuP]="MuNMuP"; dimuTitle[kMuNMuP]="#mu^{-}#mu^{+}"; - dimuName[kMuNMuN]="MuNMuN"; dimuTitle[kMuNMuN]="#mu^{-}#mu^{-}"; - dimuName[kMuPMuP]="MuPMuP"; dimuTitle[kMuPMuP]="#mu^{+}#mu^{+}"; - - Int_t nbins[kNHistDimu]; - Double_t xlow[kNHistDimu], xup[kNHistDimu]; - TString name[kNHistDimu], axis[kNHistDimu], unit[kNHistDimu]; - nbins[kInvM] =300; xlow[kInvM] =0.; xup[kInvM] =30.; name[kInvM] ="InvM"; axis[kInvM] ="M"; unit[kInvM] ="[GeV/c^{2}]"; - nbins[kPtPair]=120; xlow[kPtPair]=0.; xup[kPtPair]=60.; name[kPtPair]="Pt"; axis[kPtPair]="p_{t}"; unit[kPtPair]="[GeV/c]"; - - TH1F *histo = 0; - TString histName, histTitle; - Int_t index=-1, lowBound=0; - for (Int_t idmu=0; idmuGetXaxis()->SetTitle(Form("%s(%s) %s", axis[ihist].Data(), dimuTitle[idmu].Data(), unit[ihist].Data())); - histo->Sumw2(); - fListHisDimuon->AddAt(histo, index); - } // end loop over histos - } // end loop of different kinds of dimuon - - if (fIsUseMC) CreateOutputHistosDimuonMC(nbins, xlow, xup, name, axis, unit, dimuName, dimuTitle); - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::CreateOutputHistosSingleMuonMC(Int_t *nbins, Double_t *xlow, Double_t *xup, - TString *name, TString *axis, TString *unit) -{ - // define the histos of single muons - TString srcName[kNSingleMuSrcs]; - srcName[kBeautyMu] = "Bottom"; - srcName[kCharmMu] = "Charm"; - srcName[kPrimaryMu] = "Primary"; - srcName[kSecondaryMu] = "Secondary"; - srcName[kNotMu] = "NotMu"; - - TH1F *histo = 0; - TString histName, histTitle; - Int_t index=-1, lowBound=kNHistMu; - for (Int_t src=0; srcGetXaxis()->SetTitle(Form("%s %s", axis[ihist].Data(), unit[ihist].Data())); - histo->Sumw2(); - fListHisSingleMuon->AddAt(histo, index); - } // end loop over histos - } // // end loop of single muon src - - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::CreateOutputHistosDimuonMC(Int_t *nbins, Double_t *xlow, Double_t *xup, - TString *name, TString *axis, TString *unit, - TString *dimuName, TString *dimuTitle) -{ - TString srcName[kNDimuSrcs]; - srcName[kBBdiff] = "BBdiff"; - srcName[kBchain] = "Bchain"; - srcName[kDDdiff] = "DDdiff"; - srcName[kDchain] = "Dchain"; - srcName[kResonance] = "Resonance"; - srcName[kUncorr] = "Uncorr"; - - TH1F *histo = 0; - TString histName, histTitle; - Int_t index=-1, lowBound=kNMuMus*kNHistDimu; - for (Int_t idmu=0; idmuGetXaxis()->SetTitle(Form("%s(%s) %s", axis[ihist].Data(), dimuTitle[idmu].Data(), unit[ihist].Data())); - fListHisDimuon->AddAt(histo, index); - } // end loop of dimuon src - } // end loop over histos - } // end loop of different kinds of dimuon - - return; -} + // Terminate analysis -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::FillDistributionsAtEventLeavel() -{ - Double_t disEv[kNHistEv]; - disEv[kHistVx] = fHeader->GetXv(); - disEv[kHistVy] = fHeader->GetYv(); - disEv[kHistVz] = fHeader->GetZv(); - disEv[kHistMult] = (Double_t)fHeader->GetMultSingleMuon(); - - Int_t index=-1, lowBound=0; - for (Int_t ihist=0; ihistAt(index))->Fill(disEv[ihist]); - } - - return; -} + TFile *fout = new TFile("muonsHF.root", "RECREATE"); + gDirectory->mkdir("header"); + gDirectory->cd("header"); + fListHisHeader->Write(); + gDirectory->cd(".."); -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::FillDistributionsSingleMuon(AliAODMuonTrack *track, Int_t src) -{ - if (fHeader->IsUnrecoVertex()) return; - if (!(track->SelectSingleMuon(fSingleMuonCuts))) return; - - TLorentzVector lorentzP = track->GetP(); - Double_t disMu[kNHistMu]; - disMu[kHistP] = lorentzP.P(); - disMu[kHistPt] = lorentzP.Pt(); - disMu[kHistEta] = lorentzP.Eta(); - disMu[kHistDca] = track->GetDCA(); - disMu[kHistTrg] = (Double_t)track->GetTrigger(); - - Int_t index=-1, lowBound=0; - for (Int_t ihist=0; ihistAt(index))->Fill(disMu[ihist]); - } - - if (fIsUseMC && src>=0) FillDistributionsSingleMuonMC(disMu, src); - - Double_t vtx[3]; - fHeader->GetXYZ(vtx); - ((TH2F*)fListHisSingleMuon->FindObject("hVzPt"))->Fill(vtx[2], disMu[kHistPt]); - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::FillDistributionsSingleMuonMC(Double_t *disMu, Int_t src) -{ - if (src<0) return; - Int_t index=-1, lowBound=kNHistMu; - for (Int_t ihist=0; ihistAt(index))->Fill(disMu[ihist]); + if (fAnaMode!=2) { + gDirectory->mkdir("muon"); + gDirectory->cd("muon"); + if (fIsUseMC) { + char *muonName[] = {"BottomMu", "CharmMu", "PrimaryMu", "SecondaryMu", "NotMu", "Undentified", "All"}; + for (Int_t i=AliMuonInfoStoreMC::NSources(); i--;) { + gDirectory->mkdir(muonName[i]); + gDirectory->cd(muonName[i]); + fListHisMuon[i].Write(); + gDirectory->cd(".."); + } + } else fListHisMuon->Write(); + gDirectory->cd(".."); } - return; -} -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::FillDistributionsDimuon(AliAODMuonPair *pair, Int_t src) -{ - if (!(pair->GetTrack(0)->SelectSingleMuon(fSingleMuonCuts)) || - !(pair->GetTrack(1)->SelectSingleMuon(fSingleMuonCuts))) return; - - TLorentzVector lorentzP = pair->GetP(); - Double_t disDimu[kNHistDimu]; - disDimu[kInvM] = lorentzP.Mag(); - disDimu[kPtPair] = lorentzP.Pt(); - - Int_t dimuK=-1, charge=pair->GetCharge(); - if (charge==0) dimuK = kMuNMuP; - if (charge<0) dimuK = kMuNMuN; - if (charge>0) dimuK = kMuPMuP; - - Int_t index=-1, lowBound=0; - for (Int_t ihist=0; ihistAt(index))->Fill(disDimu[ihist]); + if (fAnaMode!=1) { + gDirectory->mkdir("dimu"); + gDirectory->cd("dimu"); + if (fIsUseMC) { + char *dimuName[] = {"BBdiff", "Bchain", "DDdiff", "Dchain", "Resonance", "Background", "All"}; + for (Int_t i=AliDimuInfoStoreMC::NSources(); i--;) { + gDirectory->mkdir(dimuName[i]); + gDirectory->cd(dimuName[i]); + fListHisDimu[i].Write(); + gDirectory->cd(".."); + } + } else fListHisDimu->Write(); + gDirectory->cd(".."); } - if (fIsUseMC && src>=0) FillDistributionsDimuonMC(disDimu, dimuK, src); - return; -} - -//_____________________________________________________________________________ -void AliAnalysisTaskSEMuonsHF::FillDistributionsDimuonMC(Double_t *disDimu, Int_t dimuK, Int_t src) -{ - if (src<0) return; - Int_t index=-1, lowBound=kNMuMus*kNHistDimu; - for (Int_t ihist=0; ihistAt(index))->Fill(disDimu[ihist]); - } + fout->Write(); + fout->Close(); return; } diff --git a/PWG3/muon/AliAnalysisTaskSEMuonsHF.h b/PWG3/muon/AliAnalysisTaskSEMuonsHF.h index b8095b4aaf8..c48c8aae221 100644 --- a/PWG3/muon/AliAnalysisTaskSEMuonsHF.h +++ b/PWG3/muon/AliAnalysisTaskSEMuonsHF.h @@ -1,15 +1,23 @@ #ifndef ALIANALYSISTASKSEMUONSHF_H #define ALIANALYSISTASKSEMUONSHF_H -#include -#include -#include +/* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//************************************************************************* +// Class AliAnalysisTaskSEMuonsHF +// AliAnalysisTaskSE for the single muon and dimuon from HF analysis +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//************************************************************************* -#include "AliMuonsHFHeader.h" -#include "AliAODMuonTrack.h" -#include "AliAODMuonPair.h" #include "AliAnalysisTaskSE.h" +class TString; +class TList; +class TClonesArray; +class AliMuonsHFHeader; + class AliAnalysisTaskSEMuonsHF : public AliAnalysisTaskSE { public: @@ -17,106 +25,40 @@ class AliAnalysisTaskSEMuonsHF : public AliAnalysisTaskSE { AliAnalysisTaskSEMuonsHF(const char *name); virtual ~AliAnalysisTaskSEMuonsHF(); + virtual void Init(); + virtual void LocalInit() { Init(); } virtual void UserCreateOutputObjects(); virtual void UserExec(Option_t *opt); virtual void Terminate(Option_t *opt); - void SetAnaMode(Int_t mode) { fAnaMode = (mode<3 ? mode: 0); } + void SetAnaMode(Int_t mode) { fAnaMode = (mode<3 ? mode : 0); } void SetIsOutputTree(Bool_t ist) { fIsOutputTree = ist; } void SetIsUseMC(Bool_t isMC) { fIsUseMC = isMC; } - void SetSingleMuonCuts(Double_t cuts[10]) { - for (Int_t i=0; i<10; i++) fSingleMuonCuts[i]=cuts[i]; - } + void SetEvsHCuts(Double_t cuts[3]) const { AliMuonsHFHeader::SetSelectionCuts(cuts); } + void SetMuonCuts(Double_t cuts[10]) const { AliMuonInfoStoreRD::SetSelectionCuts(cuts); } + void SetDimuCuts(Double_t cuts[10]) const { AliDimuInfoStoreRD::SetSelectionCuts(cuts); } private: - void CreateOutoutHistosAtEvnetLevel(); - void CreateOutputHistosSingleMuon(); - void CreateOutputHistosDimuon(); - void CreateOutputHistosSingleMuonMC(Int_t *nbins, Double_t *xlow, Double_t *xup, - TString *name, TString *axis, TString *unit); - void CreateOutputHistosDimuonMC(Int_t *nbins, Double_t *xlow, Double_t *xup, - TString *name, TString *axis, TString *unit, - TString *dimuName, TString *dimuTitle); - - void FillDistributionsAtEventLeavel(); - void FillDistributionsSingleMuon(AliAODMuonTrack *track, Int_t src=-1); - void FillDistributionsDimuon(AliAODMuonPair *pair, Int_t src=-1); - void FillDistributionsSingleMuonMC(Double_t *disMu, Int_t src); - void FillDistributionsDimuonMC(Double_t *disDimu, Int_t dimuK, Int_t src); - - enum { // histos at event level - kHistVx, // x position of vtx - kHistVy, // y position of vtx - kHistVz, // z position of vtx - kHistMult, // multiplicity of single muon track - kNHistEv // number of histos at event level - }; - enum { // distributions of single mu - kHistP, // histo of single mu p - kHistPt, // histo of single mu pt - kHistEta, // histo of single mu eta - kHistDca, // histo of single mu DCA - kHistTrg, // histo of single mu pass trigger - kNHistMu // number of histos of single mu - }; - enum { // distribution of dimuon - kInvM, // invariance mass if dimuon - kPtPair, // pt of dimu pair - kNHistDimu // number of histos of dimuon - }; - - enum { // MC sources of single muon - kBeautyMu, // mu<-b - kCharmMu, // mu<-c - kPrimaryMu, // primary muon - kSecondaryMu, // mu<-secondary paritcle - kNotMu, // not muon - kNSingleMuSrcs // number of single muon sources - }; - enum { // MC soureces of dimuon - kBBdiff, // dimuon<-BB_diff - kBchain, // dimuon<-B-chain - kDDdiff, // dimuon<-DD_diff - kDchain, // dimuon<-D-chain - kResonance, // dimuon<-resonances - kUncorr, // uncorr dimuon - kNDimuSrcs // number of dimuon sources - }; - enum { // different kinds of correlated dimu - kMuNMuP, // mu-mu+ - kMuNMuN, // mu-mu- - kMuPMuP, // mu+mu+ - kNMuMus // number of corr dimu kinds - }; - - Int_t fAnaMode; // = 0, ana both single muon and dimuon - // = 1, ana single muon - // = 2, ana dimuon + AliAnalysisTaskSEMuonsHF(const AliAnalysisTaskSEMuonsHF&); + AliAnalysisTaskSEMuonsHF& operator=(const AliAnalysisTaskSEMuonsHF&); + + Int_t fAnaMode; // = 0, ana both single muon and dimuon + // = 1, ana single muon + // = 2, ana dimuon Bool_t fIsOutputTree; // flag used to switch on/off tree output Bool_t fIsUseMC; // flag used to switch on/off MC ana - Double_t fSingleMuonCuts[10]; // 0, max of 3-momentum - // 1, min of 3-momentum - // 2, pt_Min - // 3, pt_Max - // 4, eta_Min - // 5, eta_Max - // 6, dca_Min - // 7, dca_Max - // 8, about trigger - // 9, about trigger - - AliMuonsHFHeader *fHeader; // output clones array for info at ev level - TClonesArray *fMuTrkClArr; // output clones array for single mu - TClonesArray *fMuPairClArr; // output clones array for dimu - - TList *fListHisAtEvLevel; // output list of histos at event level - TList *fListHisSingleMuon; // output list of histos for single mu - TList *fListHisDimuon; // output list of histos for dimuon - - ClassDef(AliAnalysisTaskSEMuonsHF, 5); + AliMuonsHFHeader *fHeader; // output for info at ev level + TClonesArray *fMuonClArr; // output clones array for single mu + TClonesArray *fDimuClArr; // output clones array for dimu + + TList *fListHisHeader; // output list of histos at event level + TList *fListHisMuon; // output list of histos for single mu + TList *fListHisDimu; // output list of histos for dimuon + + ClassDef(AliAnalysisTaskSEMuonsHF, 6); }; #endif diff --git a/PWG3/muon/AliDimuInfoStoreMC.cxx b/PWG3/muon/AliDimuInfoStoreMC.cxx new file mode 100644 index 00000000000..5ba0e2c89ca --- /dev/null +++ b/PWG3/muon/AliDimuInfoStoreMC.cxx @@ -0,0 +1,150 @@ +/************************************************************************** + * Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +///////////////////////////////////////////////////////////// +// +// class used to extract and store info of MC particles +// +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +///////////////////////////////////////////////////////////// + +#include "AliMuonInfoStoreMC.h" +#include "AliDimuInfoStoreRD.h" +#include "AliDimuInfoStoreMC.h" + +ClassImp(AliDimuInfoStoreMC) + +const TString AliDimuInfoStoreMC::fgkStdBranchName("DimuMC"); +const Int_t AliDimuInfoStoreMC::fgkNSources = 7; + +//----------------------------------------------------------------------------- +AliDimuInfoStoreMC::AliDimuInfoStoreMC() : +AliDimuInfoStoreRD(), +fIsFull(kFALSE), +fLorentzP(), +fSource(-1) +{ + // + // default constructor + // +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreMC::AliDimuInfoStoreMC(AliMuonInfoStoreMC* const trk0, AliMuonInfoStoreMC* const trk1, Bool_t full) : +AliDimuInfoStoreRD(), +fIsFull(full), +fLorentzP(), +fSource(-1) +{ + // + // default constructor + // + fMuonRef[0] = trk0; + fMuonRef[1] = trk1; + fLorentzP = trk0->LorentzP() + trk1->LorentzP(); + AliDimuInfoStoreRD::FillDimuInfo(); + if (fIsFull) this->FindDimuonSourceFull(); + else this->FindDimuonSourceFast(); +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreMC::AliDimuInfoStoreMC(const AliDimuInfoStoreMC &src) : +AliDimuInfoStoreRD(src), +fIsFull(src.fIsFull), +fLorentzP(src.fLorentzP), +fSource(src.fSource) +{ + // + // copy constructor + // +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreMC& AliDimuInfoStoreMC::operator=(const AliDimuInfoStoreMC &src) +{ + // + // assignment constructor + // + if(&src==this) return *this; + + fIsFull = src.fIsFull; + fLorentzP = src.fLorentzP; + fSource = src.fSource; + + return *this; +} + + +//----------------------------------------------------------------------------- +AliDimuInfoStoreMC::~AliDimuInfoStoreMC() +{ + // + // destructor + // +} + +//----------------------------------------------------------------------------- +void AliDimuInfoStoreMC::FindDimuonSourceFast() +{ + // find corr relation of two particles (fast for p-p) + + AliMuonInfoStoreMC *trk0 = (AliMuonInfoStoreMC*)fMuonRef[0].GetObject(); + Int_t src0 = trk0->MuonSource(); + if (src0<0 || src0==4 || src0==3) { + fSource=5; return; + } + + AliMuonInfoStoreMC *trk1 = (AliMuonInfoStoreMC*)fMuonRef[1].GetObject(); + Int_t src1 = trk1->MuonSource(); + if (src1<0 || src1==4 || src1==3) { + fSource=5; return; + } + + // Drell-Yan is expected very small at LHC, we ingore it + Int_t np0 = trk0->NParents() - 1; + if (np0<0) { + fSource=5; return; + } + Int_t np1 = trk1->NParents() - 1; + if (np1<0) { + fSource=5; return; + } + + if (trk0->IsMotherAResonance(np0) && + trk1->IsMotherAResonance(np1) && + trk0->ParentIndex(np0)==trk1->ParentIndex(np1)) { + fSource=4; return; + } + + if (src0==0 && src1==0) { + if ((trk0->ParentIndex(0))==(trk1->ParentIndex(0))) + fSource = 1; + else + fSource = 0; + return; + } + + if (src0==1 && src1==1) { + if ((trk0->ParentIndex(0))==(trk1->ParentIndex(0))) + fSource = 3; + else + fSource = 2; + return; + } + + fSource = 5; + return; +} diff --git a/PWG3/muon/AliDimuInfoStoreMC.h b/PWG3/muon/AliDimuInfoStoreMC.h new file mode 100644 index 00000000000..2a3cd2f2e72 --- /dev/null +++ b/PWG3/muon/AliDimuInfoStoreMC.h @@ -0,0 +1,56 @@ +#ifndef ALIDIMUINFOSTOREMC_H +#define ALIDIMUINFOSTOREMC_H + +/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//*********************************************************** +// Class AliDimuInfoStoreMC +// class used to extract and store info of MC particles +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//*********************************************************** + +#include +#include "AliDimuInfoStoreRD.h" + +class AliMuonInfoStoreMC; +class AliDimuInfoStoreMC : public AliDimuInfoStoreRD { + public: + + AliDimuInfoStoreMC(); + AliDimuInfoStoreMC(AliMuonInfoStoreMC* const trk0, AliMuonInfoStoreMC* const trk1, Bool_t full=kFALSE); + AliDimuInfoStoreMC(const AliDimuInfoStoreMC &src); + AliDimuInfoStoreMC& operator=(const AliDimuInfoStoreMC &src); + virtual ~AliDimuInfoStoreMC(); + + AliMuonInfoStoreMC* Muon(Int_t i) const { return (i<2 ? (AliMuonInfoStoreMC*)(fMuonRef[i].GetObject()) : 0x0); } + + TLorentzVector LorentzP() const { return fLorentzP; } + Int_t DimuSource() const { return fSource; } + + static const char* StdBranchName() { return fgkStdBranchName.Data(); } + static const Int_t NSources() { return fgkNSources; } + + + private: + + void FindDimuonSourceFast(); + void FindDimuonSourceFull(); + + static const TString fgkStdBranchName; // Standard branch name + static const Int_t fgkNSources; // num. of dimuon sources + + Bool_t fIsFull; // flag of using full analysis (for Pb-Pb) + TLorentzVector fLorentzP; // lorentz momentum of MC particle + Int_t fSource; // = 0, BBdiff + // = 1, Bchain + // = 2, DDdiff + // = 3, Dchain + // = 4, Resonance + // = 5, UnCorr bkg + + ClassDef(AliDimuInfoStoreMC, 1); +}; + +#endif diff --git a/PWG3/muon/AliDimuInfoStoreRD.cxx b/PWG3/muon/AliDimuInfoStoreRD.cxx new file mode 100644 index 00000000000..8dbdaa93314 --- /dev/null +++ b/PWG3/muon/AliDimuInfoStoreRD.cxx @@ -0,0 +1,142 @@ +/************************************************************************** + * Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +///////////////////////////////////////////////////////////// +// +// class used to extract and store reco info of dimu candidate +// +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +///////////////////////////////////////////////////////////// + +#include +#include + +#include "AliMuonInfoStoreRD.h" +#include "AliDimuInfoStoreRD.h" + +ClassImp(AliDimuInfoStoreRD) + +const TString AliDimuInfoStoreRD::fgkStdBranchName("DimuRD"); +Double_t AliDimuInfoStoreRD::fgCutd[10] = {-999999., 999999., + -999999., 999999., + -999999., 999999., + -999999., 999999., + -999999., 999999.}; + +//----------------------------------------------------------------------------- +AliDimuInfoStoreRD::AliDimuInfoStoreRD() : +TObject(), +fMomentum(), +fCharge(0), +fInvM(0.) +{ + // + // default constructor + // + for (Int_t i=2; i--;) fMuonRef[i] = 0; +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreRD::AliDimuInfoStoreRD(AliMuonInfoStoreRD* const trk0, AliMuonInfoStoreRD* const trk1) : +TObject(), +fMomentum(), +fCharge(0), +fInvM(0.) +{ + // + // default constructor + // + fMuonRef[0] = trk0; + fMuonRef[1] = trk1; + FillDimuInfo(); +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreRD::AliDimuInfoStoreRD(const AliDimuInfoStoreRD &src) : +TObject(src), +fMomentum(src.fMomentum), +fCharge(src.fCharge), +fInvM(src.fInvM) +{ + // + // copy constructor + // + for (Int_t i=2; i--;) fMuonRef[i] = src.fMuonRef[i]; +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreRD& AliDimuInfoStoreRD::operator=(const AliDimuInfoStoreRD &src) +{ + // + // assignment constructor + // + if(&src==this) return *this; + + fMomentum = src.fMomentum; + fCharge = src.fCharge; + fInvM = src.fInvM; + for (Int_t i=2; i--;) fMuonRef[i] = src.fMuonRef[i]; + + return *this; +} + +//----------------------------------------------------------------------------- +AliDimuInfoStoreRD::~AliDimuInfoStoreRD() +{ + // + // destructor + // +} + +//----------------------------------------------------------------------------- +void AliDimuInfoStoreRD::FillDimuInfo() +{ + // fill dimuon candidate info from the corresponding two muon tracks + + AliMuonInfoStoreRD *trk0 = (AliMuonInfoStoreRD*)fMuonRef[0].GetObject(); + AliMuonInfoStoreRD *trk1 = (AliMuonInfoStoreRD*)fMuonRef[1].GetObject(); + + fMomentum = trk0->Momentum() + trk1->Momentum(); + fCharge = trk0->Charge() + trk1->Charge(); + + Double_t mMu = TDatabasePDG::Instance()->GetParticle(13)->Mass(); + TLorentzVector lorentzP0, lorentzP1, lorentzP; + lorentzP0.SetVectM(trk0->Momentum(), mMu); + lorentzP1.SetVectM(trk1->Momentum(), mMu); + lorentzP = lorentzP0 + lorentzP1; + fInvM = lorentzP.Mag(); + + trk0 = 0; + trk1 = 0; + return; +} + +//----------------------------------------------------------------------------- +Bool_t AliDimuInfoStoreRD::DimuSelection() +{ + // select dimuon candidates according to the corresponding two muon tracks cuts + + Double_t cutsOld[10]; + AliMuonInfoStoreRD::SelectionCust(cutsOld); + AliMuonInfoStoreRD::SetSelectionCuts(AliDimuInfoStoreRD::fgCutd); + AliMuonInfoStoreRD *trk0 = (AliMuonInfoStoreRD*)fMuonRef[0].GetObject(); + AliMuonInfoStoreRD *trk1 = (AliMuonInfoStoreRD*)fMuonRef[1].GetObject(); + if (!trk0->MuonSelection()) { AliMuonInfoStoreRD::SetSelectionCuts(cutsOld); return kFALSE; } + if (!trk1->MuonSelection()) { AliMuonInfoStoreRD::SetSelectionCuts(cutsOld); return kFALSE; } + + AliMuonInfoStoreRD::SetSelectionCuts(cutsOld); + return kTRUE; +} diff --git a/PWG3/muon/AliDimuInfoStoreRD.h b/PWG3/muon/AliDimuInfoStoreRD.h new file mode 100644 index 00000000000..fec5df8475f --- /dev/null +++ b/PWG3/muon/AliDimuInfoStoreRD.h @@ -0,0 +1,58 @@ +#ifndef ALIDIMUINFOSTORERD_H +#define ALIDIMUINFOSTORERD_H + +/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//*********************************************************** +// Class AliDimuInfoStoreRD +// class used to extract and store reco info of dimu candidate +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//*********************************************************** + +#include +#include +#include +#include + +#include "AliMuonInfoStoreRD.h" + +class AliDimuInfoStoreRD : public TObject { + public: + + AliDimuInfoStoreRD(); + AliDimuInfoStoreRD(AliMuonInfoStoreRD* const trk0, AliMuonInfoStoreRD* const trk1); + AliDimuInfoStoreRD(const AliDimuInfoStoreRD &src); + AliDimuInfoStoreRD& operator=(const AliDimuInfoStoreRD &src); + virtual ~AliDimuInfoStoreRD(); + + AliMuonInfoStoreRD* Muon(Int_t i) const { return (i<2 ? (AliMuonInfoStoreRD*)(fMuonRef[i].GetObject()) : 0x0); } + + TVector3 Momentum() const { return fMomentum; } + Short_t Charge() const { return fCharge; } + Double_t InvM() const { return fInvM; } + + Bool_t DimuSelection(); + + static const char* StdBranchName() { return fgkStdBranchName.Data(); } + static void SetSelectionCuts(Double_t cuts[10]) { for (Int_t i=10; i--;) fgCutd[i]=cuts[i]; } + + protected: + + void FillDimuInfo(); + TRef fMuonRef[2]; // ref to the two corresponding muon tracks + + private: + + static const TString fgkStdBranchName; // Standard branch name + static Double_t fgCutd[10]; // single muon cuts for dimuon selection + + TVector3 fMomentum; // 3-momentum of dimuon + Short_t fCharge; // charge of dimuon + Double_t fInvM; // invariance mass of dimuon + + ClassDef(AliDimuInfoStoreRD, 3); +}; + +#endif diff --git a/PWG3/muon/AliMCMuonPair.cxx b/PWG3/muon/AliMCMuonPair.cxx deleted file mode 100644 index 8fb7ba899c6..00000000000 --- a/PWG3/muon/AliMCMuonPair.cxx +++ /dev/null @@ -1,96 +0,0 @@ -#include "AliAODMuonPair.h" -#include "AliMCMuonTrack.h" -#include "AliMCMuonPair.h" - -ClassImp(AliAODMuonPair) - -//----------------------------------------------------------------------------- -AliMCMuonPair::AliMCMuonPair() : -AliAODMuonPair(), -fIsFull(kFALSE), -fPGen(), -fSource(-1) -{ - // - // default constructor - // - for (Int_t i=0; i<2; i++) fTrk[i] = 0; - -} - -//----------------------------------------------------------------------------- -AliMCMuonPair::AliMCMuonPair(AliMCMuonTrack *trk0, AliMCMuonTrack *trk1, Bool_t full) : -AliAODMuonPair(), -fIsFull(full), -fPGen(), -fSource(-1) -{ - // - // default constructor - // - fTrk[0] = trk0; - fTrk[1] = trk1; - fPGen = trk0->GetPGen() + trk1->GetPGen(); - AliAODMuonPair::FillPairInfo(); - if (fIsFull) this->FindDimuonSourceFull(); - else this->FindDimuonSourceFast(); -} - -//----------------------------------------------------------------------------- -AliMCMuonPair::~AliMCMuonPair() -{ - // - // destructor - // -} - -//----------------------------------------------------------------------------- -void AliMCMuonPair::FindDimuonSourceFast() -{ - AliMCMuonTrack *trk0 = (AliMCMuonTrack*)fTrk[0].GetObject(); - Int_t src0 = trk0->GetSource(); - if (src0<0 || src0==4 || src0==3) { - fSource=5; return; - } - - AliMCMuonTrack *trk1 = (AliMCMuonTrack*)fTrk[1].GetObject(); - Int_t src1 = trk1->GetSource(); - if (src1<0 || src1==4 || src1==3) { - fSource=5; return; - } - - // Drell-Yan is expected very small at LHC, we ingore it - Int_t np0 = trk0->GetNParents() - 1; - if (np0<0) { - fSource=5; return; - } - - Int_t np1 = trk1->GetNParents() - 1; - if (np1<0) { - fSource=5; return; - } - - if (trk0->IsMotherAResonance(np0) && trk1->IsMotherAResonance(np1) && - (trk0->GetParentIndex(np0))==(trk1->GetParentIndex(np1))) { - fSource=4; return; - } - - if (src0==0 && src1==0) { - if ((trk0->GetParentIndex(0))==(trk1->GetParentIndex(0))) - fSource = 1; - else - fSource = 0; - return; - } - - if (src0==1 && src1==1) { - if ((trk0->GetParentIndex(0))==(trk1->GetParentIndex(0))) - fSource = 3; - else - fSource = 2; - return; - } - - fSource = 5; - return; -} diff --git a/PWG3/muon/AliMCMuonPair.h b/PWG3/muon/AliMCMuonPair.h deleted file mode 100644 index 89003e5eb10..00000000000 --- a/PWG3/muon/AliMCMuonPair.h +++ /dev/null @@ -1,38 +0,0 @@ -#ifndef ALIMCMUONPAIR_H -#define ALIMCMUONPAIR_H - -#include - -#include "AliMCMuonTrack.h" -#include "AliAODMuonPair.h" - -class AliMCMuonPair : public AliAODMuonPair { - public: - - AliMCMuonPair(); - AliMCMuonPair(AliMCMuonTrack *trk0, AliMCMuonTrack *trk1, Bool_t full=kFALSE); - virtual ~AliMCMuonPair(); - - AliMCMuonTrack* GetTrack(Int_t i) const { return (i<2 ? (AliMCMuonTrack*)(fTrk[i].GetObject()) : 0x0); } - - TLorentzVector GetPGen() const { return fPGen; } - Int_t GetSource() const { return fSource; } - - private: - - void FindDimuonSourceFast(); - void FindDimuonSourceFull(); - - Bool_t fIsFull; - TLorentzVector fPGen; - Int_t fSource; // = 0, BBdiff - // = 1, Bchain - // = 2, DDdiff - // = 3, Dchain - // = 4, Resonance - // = 5, UnCorr bkg - - ClassDef(AliMCMuonPair, 1); -}; - -#endif diff --git a/PWG3/muon/AliMCMuonTrack.cxx b/PWG3/muon/AliMCMuonTrack.cxx deleted file mode 100644 index 6f646b5373b..00000000000 --- a/PWG3/muon/AliMCMuonTrack.cxx +++ /dev/null @@ -1,407 +0,0 @@ -#include -#include -#include - -#include "AliMCEventHandler.h" -#include "AliMCEvent.h" -#include "AliStack.h" -#include "AliAODMCParticle.h" -#include "AliESDEvent.h" -#include "AliESDMuonTrack.h" -#include "AliESDMuonCluster.h" -#include "AliAODTrack.h" -#include "AliAODMuonTrack.h" -#include "AliMCMuonTrack.h" - -ClassImp(AliMCMuonTrack) - -const Double_t AliMCMuonTrack::fgkSigmaCut = 10.; - -//----------------------------------------------------------------------------- -AliMCMuonTrack::AliMCMuonTrack() : -AliAODMuonTrack(), -fIsFull(kFALSE), -fPGen(), -fTrackIndex(-1), -fTrackPDGCode(0), -fSource(-1), -fNParents(0), -fOscillation(kFALSE), -fWeight(0.) -{ - // - // default constructor - // - for (Int_t i=0; iFindTrackRef(trkAOD, mcClArr); - if (pMC) this->SetMCInfo(pMC, mcClArr); -} - -//----------------------------------------------------------------------------- -AliMCMuonTrack::AliMCMuonTrack(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH, Bool_t full) : -AliAODMuonTrack(trkESD), -fIsFull(full), -fPGen(), -fTrackIndex(-1), -fTrackPDGCode(0), -fSource(-1), -fNParents(0), -fOscillation(kFALSE), -fWeight(0.) -{ - // - // default constructor - // - TParticle *pMC = this->FindTrackRef(trkESD, esd, mcH); // do nothing when running on official train - if (pMC) this->SetMCInfo(pMC, mcH); - // MC infomation is not implement when running on the official train for they depend on the MUON module. - // In the case of private running, one can get the information from ESD by uncommenting the definitions - // of FindTrackRef() and CovESDtoMuonTrack() - // and including the following head files in both AliMCMuonTrack.h and this file, - // #include "AliMUONRecoCheck.h" - // #include "AliMUONVTrackStore.h" - // #include "AliMUONVClusterStore.h" - // #include "AliMUONTrackParam.h" - // #include "AliMUONTrack.h" - // #include "AliMUONVCluster.h" - // #include "AliMUONESDInterface.h" -} - - - -//----------------------------------------------------------------------------- -AliMCMuonTrack::~AliMCMuonTrack() -{ - // - // destructor - // -} - -//----------------------------------------------------------------------------- -AliAODMCParticle* AliMCMuonTrack::FindTrackRef(AliAODTrack *trkAOD, TClonesArray *mcClArr) -{ - AliAODMCParticle *pMC = 0; - fTrackIndex = trkAOD->GetLabel(); - if (fTrackIndex>=0) - pMC = (AliAODMCParticle*)mcClArr->At(fTrackIndex); - return pMC; -} - -//----------------------------------------------------------------------------- -TParticle* AliMCMuonTrack::FindTrackRef(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH) -{ - TParticle *pMCRef = 0; - - /*AliMUONTrack *trkMuon = CovESDtoMuonTrack(*trkESD); - AliMUONRecoCheck rc(esd,mcH); - AliMUONVTrackStore* trkRefArr = rc.TrackRefs(-1); - TIter next(trkRefArr->CreateIterator()); - AliMUONTrack* trkRef = 0; - while ((trkRef=static_cast(next()))) { - Bool_t trkCompArr[10]; - Int_t nMatchClusters = trkMuon->CompatibleTrack(trkRef, fgkSigmaCut, trkCompArr); - Double_t matchClusterFrac = ((Double_t)nMatchClusters) / ((Double_t)trkMuon->GetNClusters()); - if ((trkCompArr[0] || trkCompArr[1] || trkCompArr[2] || trkCompArr[3]) && - (trkCompArr[6] || trkCompArr[7] || trkCompArr[8] || trkCompArr[9]) && matchClusterFrac>0.5) { - fTrackIndex = trkRef->GetUniqueID(); - trkRefArr->Remove(*trkRef); - break; - } - } - if (fTrackIndex>=0) pMCRef = mcH->MCEvent()->Stack()->Particle(fTrackIndex);*/ - - return pMCRef; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::SetMCInfo(AliAODMCParticle *pMC, TClonesArray *mcClArr) -{ - fPGen.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->E()); - fTrackPDGCode = pMC->GetPdgCode(); - if (TMath::Abs(fTrackPDGCode)!=13) { - fSource = 4; - return; - } - - Int_t lineM = pMC->GetMother(); - if (lineM<0) { - fSource = 2; - return; - } - - Bool_t isPrimary = - ((AliAODMCParticle*)mcClArr->At(lineM))->IsPrimary(); - if (!isPrimary) { - fSource = 3; - return; - } - - this->FillHistoryParents(lineM, mcClArr); - fSource = this->SelectHFMuon(); - return; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::SetMCInfo(TParticle *pMC, AliMCEventHandler *mcH) -{ - fPGen.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->Energy()); - fTrackPDGCode = pMC->GetPdgCode(); - if (TMath::Abs(fTrackPDGCode)!=13) { - fSource = 4; - return; - } - - Int_t lineM = pMC->GetFirstMother(); - if (lineM<0) { - fSource = 2; - return; - } - - AliStack *stack = mcH->MCEvent()->Stack(); - if (lineM>=stack->GetNprimary()) { - fSource = 3; - return; - } - - this->FillHistoryParents(lineM, stack); - fSource = this->SelectHFMuon(); - return; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::FillHistoryParents(Int_t lineM, TClonesArray *mcClArr) -{ - Int_t countP=-1, pdg=0; - Int_t parents[10], parLine[10]; - AliAODMCParticle *mother = 0; - while(lineM>=0){ - mother = (AliAODMCParticle*)mcClArr->At(lineM); - pdg = mother->GetPdgCode(); - if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; - parents[++countP] = pdg; - parLine[countP] = lineM; - lineM = mother->GetMother(); - } - for(Int_t i=0; i<=countP; i++){ - fParentIndex[i] = parLine[countP-i]; - fParentPDGCode[i] = parents[countP-i]; - } - fNParents = countP + 1; - - if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, mcClArr); - return; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::FillHistoryParents(Int_t lineM, AliStack *stack) -{ - Int_t countP=-1, pdg=0; - Int_t parents[10], parLine[10]; - TParticle *mother = 0; - while(lineM>=0){ - mother = stack->Particle(lineM); - pdg = mother->GetPdgCode(); - if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; - parents[++countP] = pdg; - parLine[countP] = lineM; - lineM = mother->GetFirstMother(); - } - for(Int_t i=0; i<=countP; i++){ - fParentIndex[i] = parLine[countP-i]; - fParentPDGCode[i] = parents[countP-i]; - } - fNParents = countP + 1; - - if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, stack); - return; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::FillHistoryQuarks(Int_t lineM, TClonesArray *mcClArr) -{ - // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx - if (lineM<0) return; - Int_t countP=-1, pdg=0; - AliAODMCParticle *mother = 0; - while(lineM>=0){ - mother = (AliAODMCParticle*)mcClArr->At(lineM); - pdg = mother->GetPdgCode(); - fQuarkIndex[++countP] = lineM; - fQuarkPDGCode[countP] = pdg; - lineM = mother->GetMother(); - } - - // for PYTHIA checking - countP = 1; - for(Int_t par=0; par<4; par++) { - if(TMath::Abs(this->GetQuarkPDGCode(par))<6) { countP=par; break; } - } - if(this->GetQuarkIndex(countP)>-1 && (this->GetParentFlavour(0)==4 || this->GetParentFlavour(0)==5)) { - if(this->GetParentFlavour(0)!=TMath::Abs(this->GetQuarkPDGCode(countP))) { - AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", - this->GetParentFlavour(0), TMath::Abs(this->GetQuarkPDGCode(countP)))); - - pdg = this->GetQuarkPDGCode(countP); - Int_t line = this->GetQuarkIndex(countP); - this->ResetQuarkInfo(); - while(TMath::Abs(pdg)!=this->GetParentFlavour(0)) { - pdg = ((AliAODMCParticle*)mcClArr->At(++line))->GetPdgCode(); - } - while(line>=0){ - mother = (AliAODMCParticle*)mcClArr->At(line); - pdg = mother->GetPdgCode(); - fQuarkIndex[countP] = line; - fQuarkPDGCode[countP++] = pdg; - line = mother->GetMother(); - } - } - } - return; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::FillHistoryQuarks(Int_t lineM, AliStack *stack) -{ - // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx - if (lineM<0) return; - Int_t countP=-1, pdg=0; - TParticle *mother = 0; - while(lineM>=0){ - mother = stack->Particle(lineM); - pdg = mother->GetPdgCode(); - fQuarkIndex[++countP] = lineM; - fQuarkPDGCode[countP] = pdg; - lineM = mother->GetFirstMother(); - } - - // for PYTHIA checking - countP = 1; - for(Int_t par=0; par<4; par++) { - if(TMath::Abs(this->GetQuarkPDGCode(par))<6) { countP=par; break; } - } - if(this->GetQuarkIndex(countP)>-1 && (this->GetParentFlavour(0)==4 || this->GetParentFlavour(0)==5)) { - if(this->GetParentFlavour(0)!=TMath::Abs(this->GetQuarkPDGCode(countP))) { - AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", - this->GetParentFlavour(0), TMath::Abs(this->GetQuarkPDGCode(countP)))); - - pdg = this->GetQuarkPDGCode(countP); - Int_t line = this->GetQuarkIndex(countP); - this->ResetQuarkInfo(); - while(TMath::Abs(pdg)!=this->GetParentFlavour(0)) { - pdg = stack->Particle(++line)->GetPdgCode(); - } - while(line>=0){ - mother = stack->Particle(line); - pdg = mother->GetPdgCode(); - fQuarkIndex[countP] = line; - fQuarkPDGCode[countP++] = pdg; - line = mother->GetFirstMother(); - } - } - } - return; -} - -//----------------------------------------------------------------------------- -/*AliMUONTrack* AliMCMuonTrack::CovESDtoMuonTrack(AliESDMuonTrack &trkESD) -{ - // method in $ALICE_ROOT/PWG3/muondep/AliAnalysisTaskESDMCLabelAddition.cxx - AliMUONTrack *trkMuon = new AliMUONTrack(); - if(!trkESD.ClustersStored()) return trkMuon; - - AliMUONTrackParam param; - AliMUONESDInterface::GetParamAtFirstCluster(trkESD, param); - AliMUONESDInterface::GetParamCov(trkESD, param); - AliMUONVClusterStore *cStore = AliMUONESDInterface::NewClusterStore(); - AliMUONVCluster *cluster = cStore->CreateCluster(0,0,0); - - AliESDMuonCluster *esdCluster = (AliESDMuonCluster*)trkESD.GetClusters().First(); - while (esdCluster) { - AliMUONESDInterface::ESDToMUON(*esdCluster, *cluster); - param.SetZ(cluster->GetZ()); - trkMuon->AddTrackParamAtCluster(param, *cluster, kTRUE); - esdCluster = (AliESDMuonCluster*)trkESD.GetClusters().After(esdCluster); - } - - delete cluster; cluster = 0; - delete cStore; cStore = 0; - return trkMuon; -}*/ - -//----------------------------------------------------------------------------- -Int_t AliMCMuonTrack::SelectHFMuon() -{ - Int_t flv = GetParentFlavour(0); - if (flv!=4 && flv!=5) return 2; - - Bool_t isRes = kFALSE; - Int_t i=0, nparents=this->GetNParents(); - while (i1000 && (pdg%100)<10) return kTRUE; - else return kFALSE; -} - -//----------------------------------------------------------------------------- -void AliMCMuonTrack::ResetQuarkInfo() -{ - for(Int_t pos=1; pos<4; pos++) { - fQuarkIndex[pos] = -1; - fQuarkPDGCode[pos] = 0; - } - return; -} - -//----------------------------------------------------------------------------- -Int_t AliMCMuonTrack::GetParentFlavour(Int_t i) const -{ - Int_t pdg = GetParentPDGCode(i); - pdg = TMath::Abs(pdg/100); - if(pdg>9) pdg /= 10; - return pdg; -} - -//----------------------------------------------------------------------------- -Bool_t AliMCMuonTrack::IsMotherAResonance(Int_t i) const -{ - // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx - Int_t pdg = GetParentPDGCode(i); - Int_t id=pdg%100000; - return (!((id-id%10)%110)); -} diff --git a/PWG3/muon/AliMCMuonTrack.h b/PWG3/muon/AliMCMuonTrack.h deleted file mode 100644 index 6e2a97fbaea..00000000000 --- a/PWG3/muon/AliMCMuonTrack.h +++ /dev/null @@ -1,83 +0,0 @@ -#ifndef ALIMCMUONTRACK_H -#define ALIMCMUONTRACK_H - -#include -#include -#include - -#include "AliMCEventHandler.h" -#include "AliStack.h" -#include "AliAODMCParticle.h" -#include "AliESDEvent.h" -#include "AliESDMuonTrack.h" -#include "AliAODTrack.h" -#include "AliAODMuonTrack.h" - -class AliMCMuonTrack : public AliAODMuonTrack { - public: - - AliMCMuonTrack(); - AliMCMuonTrack(AliAODTrack *trkAOD, TClonesArray *mcClArr, Bool_t full=kFALSE); - AliMCMuonTrack(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH, Bool_t full=kFALSE); - AliMCMuonTrack(AliESDMuonTrack *trkESD, AliStack *stack); - virtual ~AliMCMuonTrack(); - - Int_t GetParentFlavour(Int_t i=0) const; - Bool_t IsMotherAResonance(Int_t i) const; - Bool_t IsAMuon() const { return (fSource>=0 && fSource!=4); } - - TLorentzVector GetPGen() const { return fPGen; } - Int_t GetSource() const { return fSource; } - Int_t GetTrackIndex() const { return fTrackIndex; } - Int_t GetTrackPDGCode() const { return fTrackPDGCode; } - Int_t GetNParents() const { return fNParents; } - Int_t GetParentIndex(Int_t i=0) const { return (i +#include + +#include "AliMCEvent.h" +#include "AliMCEventHandler.h" +#include "AliStack.h" +#include "AliAODMCParticle.h" +#include "AliESDMuonTrack.h" +#include "AliAODTrack.h" +#include "AliMuonInfoStoreRD.h" +#include "AliMuonInfoStoreMC.h" + +class AliESDEvent; + +ClassImp(AliMuonInfoStoreMC) + +const TString AliMuonInfoStoreMC::fgkStdBranchName("MuonMC"); +const Int_t AliMuonInfoStoreMC::fgkNSources = 7; + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC::AliMuonInfoStoreMC() : +AliMuonInfoStoreRD(), +fIsFull(kFALSE), +fLorentzP(), +fTrackIndex(-1), +fTrackPDGCode(0), +fSource(-1), +fNParents(0), +fOscillation(kFALSE), +fWeight(0.) +{ + // + // default constructor + // + for (Int_t i=5; i--;) { fParentIndex[i] = -1; fParentPDGCode[i] = 0; } + for (Int_t i=4; i--;) { fQuarkIndex[i] = -1; fQuarkPDGCode[i] = 0; } +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC::AliMuonInfoStoreMC(AliAODTrack *trkAOD, TClonesArray *mcClArr, Bool_t full) : +AliMuonInfoStoreRD(trkAOD), +fIsFull(full), +fLorentzP(), +fTrackIndex(-1), +fTrackPDGCode(0), +fSource(-1), +fNParents(0), +fOscillation(kFALSE), +fWeight(0.) +{ + // + // default constructor + // + for (Int_t i=5; i--;) { fParentIndex[i] = -1; fParentPDGCode[i] = 0; } + for (Int_t i=4; i--;) { fQuarkIndex[i] = -1; fQuarkPDGCode[i] = 0; } + + AliAODMCParticle *pMC = this->FindTrackRef(trkAOD, mcClArr); + if (pMC) this->SetMCInfo(pMC, mcClArr); +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC::AliMuonInfoStoreMC(AliESDMuonTrack *trkESD, AliMCEventHandler *mcH, Bool_t full) : +AliMuonInfoStoreRD(trkESD), +fIsFull(full), +fLorentzP(), +fTrackIndex(-1), +fTrackPDGCode(0), +fSource(-1), +fNParents(0), +fOscillation(kFALSE), +fWeight(0.) +{ + // + // default constructor + // + for (Int_t i=5; i--;) { fParentIndex[i] = -1; fParentPDGCode[i] = 0; } + for (Int_t i=4; i--;) { fQuarkIndex[i] = -1; fQuarkPDGCode[i] = 0; } + + TParticle *pMC = this->FindTrackRef(trkESD, mcH); + if (pMC) this->SetMCInfo(pMC, mcH); +} + +//----------------------------------------------------------------------------- +/*AliMuonInfoStoreMC::AliMuonInfoStoreMC(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH, Bool_t full) : +AliMuonInfoStoreRD(trkESD), +fIsFull(full), +fLorentzP(), +fTrackIndex(-1), +fTrackPDGCode(0), +fSource(-1), +fNParents(0), +fOscillation(kFALSE), +fWeight(0.) +{ +#include "AliMUONRecoCheck.h" +#include "AliMUONVTrackStore.h" +#include "AliMUONTrack.h" +#include "AliMUONESDInterface.h" + // + // default constructor + // + for (Int_t i=5; i--;) { fParentIndex[i] = -1; fParentPDGCode[i] = 0; } + for (Int_t i=4; i--;) { fQuarkIndex[i] = -1; fQuarkPDGCode[i] = 0; } + + TParticle *pMC = this->FindTrackRef(trkESD, esd, mcH); + if (pMC) this->SetMCInfo(pMC, mcH); +}*/ + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC::AliMuonInfoStoreMC(const AliMuonInfoStoreMC &src) : +AliMuonInfoStoreRD(src), +fIsFull(src.fIsFull), +fLorentzP(src.fLorentzP), +fTrackIndex(src.fTrackIndex), +fTrackPDGCode(src.fTrackPDGCode), +fSource(src.fSource), +fNParents(src.fNParents), +fOscillation(src.fOscillation), +fWeight(src.fWeight) +{ + // + // copy constructor + // + for (Int_t i=5; i--;) { + fParentIndex[i] = src.fParentIndex[i]; + fParentPDGCode[i] = src.fParentPDGCode[i]; + } + for (Int_t i=4; i--;) { + fQuarkIndex[i] = src.fQuarkIndex[i]; + fQuarkPDGCode[i] = src.fQuarkPDGCode[i]; + } +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC& AliMuonInfoStoreMC::operator=(const AliMuonInfoStoreMC &src) +{ + // + // assignment constructor + // + if(&src==this) return *this; + AliMuonInfoStoreRD::operator=(src); + + fIsFull = src.fIsFull; + fLorentzP = src.fLorentzP; + fTrackIndex = src.fTrackIndex; + fTrackPDGCode = src.fTrackPDGCode; + fSource = src.fSource; + fNParents = src.fNParents; + fOscillation = src.fOscillation; + fWeight = src.fWeight; + for (Int_t i=5; i--;) { + fParentIndex[i] = src.fParentIndex[i]; + fParentPDGCode[i] = src.fParentPDGCode[i]; + } + for (Int_t i=4; i--;) { + fQuarkIndex[i] = src.fQuarkIndex[i]; + fQuarkPDGCode[i] = src.fQuarkPDGCode[i]; + } + + return *this; +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreMC::~AliMuonInfoStoreMC() +{ + // + // destructor + // +} + +//----------------------------------------------------------------------------- +AliAODMCParticle* AliMuonInfoStoreMC::FindTrackRef(AliAODTrack* const trkAOD, TClonesArray* const mcClArr) +{ + // find MC track ref with AOD base + + AliAODMCParticle *pMC = 0; + fTrackIndex = trkAOD->GetLabel(); + if (fTrackIndex>=0) + pMC = (AliAODMCParticle*)mcClArr->At(fTrackIndex); + return pMC; +} + +//----------------------------------------------------------------------------- +TParticle* AliMuonInfoStoreMC::FindTrackRef(AliESDMuonTrack* const trkESD, AliMCEventHandler* const mcH) +{ + // find MC track ref with ESD base + + TParticle *pMCRef = 0; + fTrackIndex = trkESD->GetLabel(); + if (fTrackIndex>=0) pMCRef = mcH->MCEvent()->Stack()->Particle(fTrackIndex); + return pMCRef; +} + +//----------------------------------------------------------------------------- +/*TParticle* AliMuonInfoStoreMC::FindTrackRef(AliESDMuonTrack* const trkESD, AliESDEvent* const esd, AliMCEventHandler* const mcH) +{ + // find MC track ref with ESD trackRef base + + TParticle *pMCRef = 0; + AliMUONRecoCheck rc(esd,mcH); + AliMUONVTrackStore *trkRefArr = rc.TrackRefs(-1); + + AliMUONTrack trkMuon; + AliMUONESDInterface::ESDToMUON(*trkESD, trkMuon, kFALSE); + + Int_t nMatchClusters = 0; + AliMUONTrack *trkRef = rc.FindCompatibleTrack(trkMuon, *trkRefArr, nMatchClusters, kFALSE, 10.); + if (trkRef) fTrackIndex = trkRef->GetUniqueID(); + if (fTrackIndex>=0) pMCRef = mcH->MCEvent()->Stack()->Particle(fTrackIndex); + return pMCRef; +}*/ + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::SetMCInfo(AliAODMCParticle *pMC, TClonesArray *mcClArr) +{ + // fill track MC info with AOD base + + fLorentzP.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->E()); + fTrackPDGCode = pMC->GetPdgCode(); + if (TMath::Abs(fTrackPDGCode)!=13) { + fSource = 4; + return; + } + + Int_t lineM = pMC->GetMother(); + if (lineM<0) { + fSource = 2; + return; + } + + Bool_t isPrimary = ((AliAODMCParticle*)mcClArr->At(lineM))->IsPrimary(); + if (!isPrimary) { + fSource = 3; + return; + } + + this->FillHistoryParents(lineM, mcClArr); + fSource = this->SelectHFMuon(); + return; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::SetMCInfo(TParticle *pMC, AliMCEventHandler* const mcH) +{ + // fill track MC info with ESD base + + fLorentzP.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->Energy()); + fTrackPDGCode = pMC->GetPdgCode(); + if (TMath::Abs(fTrackPDGCode)!=13) { + fSource = 4; + return; + } + + Int_t lineM = pMC->GetFirstMother(); + if (lineM<0) { + fSource = 2; + return; + } + + AliStack *stack = mcH->MCEvent()->Stack(); + if (lineM>=stack->GetNprimary()) { + fSource = 3; + return; + } + + this->FillHistoryParents(lineM, stack); + fSource = this->SelectHFMuon(); + return; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::FillHistoryParents(Int_t lineM, TClonesArray *mcClArr) +{ + // find track hadron parents with AOD base + + Int_t countP=-1, pdg=0; + Int_t parents[10], parLine[10]; + AliAODMCParticle *mother = 0; + while(lineM>=0){ + mother = (AliAODMCParticle*)mcClArr->At(lineM); + pdg = mother->GetPdgCode(); + if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; + parents[++countP] = pdg; + parLine[countP] = lineM; + lineM = mother->GetMother(); + } + for(Int_t i=0; i<=countP; i++){ + fParentIndex[i] = parLine[countP-i]; + fParentPDGCode[i] = parents[countP-i]; + } + fNParents = countP + 1; + + if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, mcClArr); + return; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::FillHistoryParents(Int_t lineM, AliStack *stack) +{ + // find track hadron parents with ESD base + + Int_t countP=-1, pdg=0; + Int_t parents[10], parLine[10]; + TParticle *mother = 0; + while(lineM>=0){ + mother = stack->Particle(lineM); + pdg = mother->GetPdgCode(); + if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; + parents[++countP] = pdg; + parLine[countP] = lineM; + lineM = mother->GetFirstMother(); + } + for(Int_t i=0; i<=countP; i++){ + fParentIndex[i] = parLine[countP-i]; + fParentPDGCode[i] = parents[countP-i]; + } + fNParents = countP + 1; + + if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, stack); + return; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::FillHistoryQuarks(Int_t lineM, TClonesArray* const mcClArr) +{ + // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx + + if (lineM<0) return; + Int_t countP=-1, pdg=0; + AliAODMCParticle *mother = 0; + while(lineM>=0){ + mother = (AliAODMCParticle*)mcClArr->At(lineM); + pdg = mother->GetPdgCode(); + fQuarkIndex[++countP] = lineM; + fQuarkPDGCode[countP] = pdg; + lineM = mother->GetMother(); + } + + // for PYTHIA checking + countP = 1; + for(Int_t par=0; par<4; par++) { + if(TMath::Abs(this->QuarkPDGCode(par))<6) { countP=par; break; } + } + if(this->QuarkIndex(countP)>-1 && (this->ParentFlavour(0)==4 || this->ParentFlavour(0)==5)) { + if(this->ParentFlavour(0)!=TMath::Abs(this->QuarkPDGCode(countP))) { + AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", + this->ParentFlavour(0), TMath::Abs(this->QuarkPDGCode(countP)))); + + pdg = this->QuarkPDGCode(countP); + Int_t line = this->QuarkIndex(countP); + this->ResetQuarkInfo(); + while(TMath::Abs(pdg)!=this->ParentFlavour(0)) { + pdg = ((AliAODMCParticle*)mcClArr->At(++line))->GetPdgCode(); + } + while(line>=0){ + mother = (AliAODMCParticle*)mcClArr->At(line); + pdg = mother->GetPdgCode(); + fQuarkIndex[countP] = line; + fQuarkPDGCode[countP++] = pdg; + line = mother->GetMother(); + } + } + } + return; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::FillHistoryQuarks(Int_t lineM, AliStack* const stack) +{ + // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx + + if (lineM<0) return; + Int_t countP=-1, pdg=0; + TParticle *mother = 0; + while(lineM>=0){ + mother = stack->Particle(lineM); + pdg = mother->GetPdgCode(); + fQuarkIndex[++countP] = lineM; + fQuarkPDGCode[countP] = pdg; + lineM = mother->GetFirstMother(); + } + + // for PYTHIA checking + countP = 1; + for(Int_t par=0; par<4; par++) { + if(TMath::Abs(this->QuarkPDGCode(par))<6) { countP=par; break; } + } + if(this->QuarkIndex(countP)>-1 && (this->ParentFlavour(0)==4 || this->ParentFlavour(0)==5)) { + if(this->ParentFlavour(0)!=TMath::Abs(this->QuarkPDGCode(countP))) { + AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", + this->ParentFlavour(0), TMath::Abs(this->QuarkPDGCode(countP)))); + + pdg = this->QuarkPDGCode(countP); + Int_t line = this->QuarkIndex(countP); + this->ResetQuarkInfo(); + while(TMath::Abs(pdg)!=this->ParentFlavour(0)) { + pdg = stack->Particle(++line)->GetPdgCode(); + } + while(line>=0){ + mother = stack->Particle(line); + pdg = mother->GetPdgCode(); + fQuarkIndex[countP] = line; + fQuarkPDGCode[countP++] = pdg; + line = mother->GetFirstMother(); + } + } + } + return; +} + +//----------------------------------------------------------------------------- +Int_t AliMuonInfoStoreMC::SelectHFMuon() +{ + // set info of muon from HF + + Int_t flv = ParentFlavour(0); + if (flv!=4 && flv!=5) return 2; + + Bool_t isRes = kFALSE; + Int_t i=0, nparents=this->NParents(); + while (i1000 && (pdg%100)<10) return kTRUE; + else return kFALSE; +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreMC::ResetQuarkInfo() +{ + // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx + for(Int_t pos=1; pos<4; pos++) { + fQuarkIndex[pos] = -1; + fQuarkPDGCode[pos] = 0; + } + return; +} + +//----------------------------------------------------------------------------- +Int_t AliMuonInfoStoreMC::ParentFlavour(Int_t i) const +{ + // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx + Int_t pdg = ParentPDGCode(i); + pdg = TMath::Abs(pdg/100); + if(pdg>9) pdg /= 10; + return pdg; +} + +//----------------------------------------------------------------------------- +Bool_t AliMuonInfoStoreMC::IsMotherAResonance(Int_t i) const +{ + // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx + Int_t pdg = ParentPDGCode(i); + Int_t id=pdg%100000; + return (!((id-id%10)%110)); +} diff --git a/PWG3/muon/AliMuonInfoStoreMC.h b/PWG3/muon/AliMuonInfoStoreMC.h new file mode 100644 index 00000000000..b8a2ddd946d --- /dev/null +++ b/PWG3/muon/AliMuonInfoStoreMC.h @@ -0,0 +1,101 @@ +#ifndef ALIMUONINFOSTOREMC_H +#define ALIMUONINFOSTOREMC_H + +/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//*********************************************************** +// Class AliMuonInfoStoreRD +// class used to extract and store info of MC particle +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//*********************************************************** + +#include +#include +#include +#include + +#include "AliMCEventHandler.h" +#include "AliStack.h" +#include "AliAODMCParticle.h" +#include "AliESDEvent.h" +#include "AliESDMuonTrack.h" +#include "AliAODTrack.h" +#include "AliMuonInfoStoreRD.h" + +class AliMuonInfoStoreMC : public AliMuonInfoStoreRD { + public: + + AliMuonInfoStoreMC(); + AliMuonInfoStoreMC(AliAODTrack *trkAOD, TClonesArray *mcClArr, Bool_t full=kFALSE); + AliMuonInfoStoreMC(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH, Bool_t full=kFALSE); + AliMuonInfoStoreMC(AliESDMuonTrack *trkESD, AliMCEventHandler *mcH, Bool_t full=kFALSE); + AliMuonInfoStoreMC(AliESDMuonTrack *trkESD, AliStack *stack); + AliMuonInfoStoreMC(const AliMuonInfoStoreMC &src); + AliMuonInfoStoreMC& operator=(const AliMuonInfoStoreMC &src); + virtual ~AliMuonInfoStoreMC(); + + Int_t ParentFlavour(Int_t i=0) const; + Bool_t IsMotherAResonance(Int_t i) const; + Bool_t IsAMuon() const { return (fSource>=0 && fSource!=4); } + + TLorentzVector LorentzP() const { return fLorentzP; } + Int_t MuonSource() const { return fSource; } + Int_t TrackIndex() const { return fTrackIndex; } + Int_t TrackPDGCode() const { return fTrackPDGCode; } + Int_t NParents() const { return fNParents; } + Int_t ParentIndex(Int_t i=0) const { return (iFillMuonInfo(trk); +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreRD::AliMuonInfoStoreRD(AliESDMuonTrack *trk) : +TObject(), +fMomentum(), +fMomentumAtDCA(), +fCharge(0), +fMatchTrigger(-1), +fNClusters(0), +fMUONClusterMap(0), +fChi2FitMomentum(0.), +fChi2MatchTrigger(0.) +{ + // + // ESD-base constructor + // + for (Int_t i=3; i--;) fDCA[i] = 0.; + this->FillMuonInfo(trk); +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreRD::AliMuonInfoStoreRD(const AliMuonInfoStoreRD &src) : +TObject(src), +fMomentum(src.fMomentum), +fMomentumAtDCA(src.fMomentumAtDCA), +fCharge(src.fCharge), +fMatchTrigger(src.fMatchTrigger), +fNClusters(src.fNClusters), +fMUONClusterMap(src.fMUONClusterMap), +fChi2FitMomentum(src.fChi2FitMomentum), +fChi2MatchTrigger(src.fChi2MatchTrigger) +{ + // + // copy constructor + // + for (Int_t i=3; i--;) fDCA[i] = src.fDCA[i]; +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreRD& AliMuonInfoStoreRD::operator=(const AliMuonInfoStoreRD &src) +{ + // + // assignment constructor + // + if(&src==this) return *this; + + fMomentum = src.fMomentum; + fMomentumAtDCA = src.fMomentumAtDCA; + + fCharge = src.fCharge; + fMatchTrigger = src.fMatchTrigger; + fNClusters = src.fNClusters; + fMUONClusterMap = src.fMUONClusterMap; + fChi2FitMomentum = src.fChi2FitMomentum; + fChi2MatchTrigger = src.fChi2MatchTrigger; + + for (Int_t i=3; i--;) fDCA[i] = src.fDCA[i]; + + return *this; +} + +//----------------------------------------------------------------------------- +AliMuonInfoStoreRD::~AliMuonInfoStoreRD() +{ + // + // destructor + // +} + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreRD::FillMuonInfo(AliAODTrack *trk) +{ + // extract reco info of muon track from AOD + + Double_t arr[3]; + trk->PxPyPz(arr); + this->SetMomentum(arr); + + trk->PxPyPzAtDCA(arr); + this->SetMomentumAtDCA(arr); + + trk->XYZAtDCA(arr); + this->SetDCA(arr); + + this->SetCharge(trk->Charge()); + this->SetMatchTrigger(trk->GetMatchTrigger()); + this->SetMUONClusterMap(trk->GetMUONClusterMap()); + this->SetChi2FitMomentum(trk->Chi2perNDF()); + this->SetChi2MatchTrigger(trk->GetChi2MatchTrigger()); + + return; +} + + +//----------------------------------------------------------------------------- +void AliMuonInfoStoreRD::FillMuonInfo(AliESDMuonTrack *trk) +{ + // extract reco info of muon track from ESD + // tract params before related to vertex are extracted + + Double_t arr[3]; + trk->PxPyPz(arr); + this->SetMomentum(arr); + + trk->PxPyPzAtDCA(arr); + this->SetMomentumAtDCA(arr); + + arr[0] = trk->GetNonBendingCoorAtDCA(); + arr[1] = trk->GetBendingCoorAtDCA(); + arr[2] = trk->GetZ(); + this->SetDCA(arr); + + this->SetCharge(trk->Charge()); + this->SetMatchTrigger(trk->GetMatchTrigger()); + this->SetMUONClusterMap(trk->GetMuonClusterMap()); + this->SetChi2FitMomentum(trk->GetChi2()/(2.*trk->GetNHit()-5.)); + this->SetChi2MatchTrigger(trk->GetChi2MatchTrigger()); + + this->SetNClusters(trk->GetNClusters()); + return; +} + +//----------------------------------------------------------------------------- +Bool_t AliMuonInfoStoreRD::MuonSelection() +{ + // select muon tracks according to the selection cuts + + Double_t p = Momentum().Mag(); + if (pAliMuonInfoStoreRD::fgCuts[1]) return kFALSE; + + Double_t pt = Momentum().Pt(); + if (ptAliMuonInfoStoreRD::fgCuts[3]) return kFALSE; + + Double_t eta = Momentum().Eta(); + if (etaAliMuonInfoStoreRD::fgCuts[5]) return kFALSE; + + Double_t dca = this->DCA(); + if (dcaAliMuonInfoStoreRD::fgCuts[7]) return kFALSE; + + Int_t trigger = this->MatchTrigger(); + if (triggerAliMuonInfoStoreRD::fgCuts[9]) return kFALSE; + + return kTRUE; +} diff --git a/PWG3/muon/AliMuonInfoStoreRD.h b/PWG3/muon/AliMuonInfoStoreRD.h new file mode 100644 index 00000000000..365b15c40ef --- /dev/null +++ b/PWG3/muon/AliMuonInfoStoreRD.h @@ -0,0 +1,92 @@ +#ifndef ALIMUONINFOSTORERD_H +#define ALIMUONINFOSTORERD_H + +/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//*********************************************************** +// Class AliMuonInfoStoreRD +// class used to extract and store reco info of muon track +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//*********************************************************** + +#include +#include +#include + +class AliAODTrack; +class AliESDMuonTrack; + +class AliMuonInfoStoreRD : public TObject { + public: + + AliMuonInfoStoreRD(); + AliMuonInfoStoreRD(AliAODTrack *trk); + AliMuonInfoStoreRD(AliESDMuonTrack *trk); + AliMuonInfoStoreRD(const AliMuonInfoStoreRD &src); + AliMuonInfoStoreRD& operator=(const AliMuonInfoStoreRD &src); + virtual ~AliMuonInfoStoreRD(); + + TVector3 Momentum() const { return fMomentum; } + TVector3 MomentumAtDCA() const { return fMomentumAtDCA; } + + void XYZAtDCA(Double_t dca[3]) const { for (Int_t i=3; i--;) dca[i]=fDCA[i]; } + Double_t DCA() const { return TMath::Sqrt(fDCA[0]*fDCA[0]+fDCA[1]*fDCA[1]); } + + Short_t Charge() const { return fCharge; } + Int_t MatchTrigger() const { return fMatchTrigger; } + Int_t NClusters() const { return fNClusters; } + UInt_t MUONClusterMap() const { return fMUONClusterMap; } + Double_t Chi2FitMomentum() const { return fChi2FitMomentum; } + Double_t Chi2MatchTrigger() const { return fChi2MatchTrigger; } + + Bool_t MuonSelection(); + + static const char* StdBranchName() { return fgkStdBranchName.Data(); } + static const void SelectionCust(Double_t cuts[10]) { for (Int_t i=10; i--;) cuts[i]=fgCuts[i]; } + static void SetSelectionCuts(Double_t cuts[10]) { for (Int_t i=10; i--;) fgCuts[i]=cuts[i]; } + + private: + + void FillMuonInfo(AliAODTrack *trk); + void FillMuonInfo(AliESDMuonTrack *trk); + + void SetMomentum(Double_t p[3]) { fMomentum.SetXYZ(p[0],p[1],p[2]); } + void SetMomentumAtDCA(Double_t p[3]) { fMomentumAtDCA.SetXYZ(p[0],p[1],p[2]); } + + void SetDCA(Double_t dca[3]) { for (Int_t i=3; i--;) fDCA[i]=dca[i]; } + void SetCharge(Short_t charge) { fCharge = charge; } + void SetNClusters(Int_t ncls) { fNClusters = ncls; } + void SetMatchTrigger(Int_t trigger) { fMatchTrigger = trigger; } + void SetMUONClusterMap(UInt_t clMap) { fMUONClusterMap = clMap; } + void SetChi2FitMomentum(Double_t chi2) { fChi2FitMomentum = chi2; } + void SetChi2MatchTrigger(Double_t chi2) { fChi2MatchTrigger = chi2; } + + static const TString fgkStdBranchName; // Standard branch name + static Double_t fgCuts[10]; // 0, min of 3-momentum + // 1, max of 3-momentum + // 2, pt_Min + // 3, pt_Max + // 4, eta_Min + // 5, eta_Max + // 6, dca_Min + // 7, dca_Max + // 8, about trigger matching + // 9, about trigger matching + + TVector3 fMomentum; // momentum corrected w vtx + TVector3 fMomentumAtDCA; // momentum at DCA in vtx plane + + Double_t fDCA[3]; // distance of closet approach + Short_t fCharge; // track charge + Int_t fMatchTrigger; // type of match trigger + Int_t fNClusters; // number of clusters in the track + UInt_t fMUONClusterMap; // map of MUON clusters + Double_t fChi2FitMomentum; // chi2/NDF of momentum fit + Double_t fChi2MatchTrigger; // chi2 of trigger matching + + ClassDef(AliMuonInfoStoreRD, 4); +}; + +#endif diff --git a/PWG3/muon/AliMuonsHFHeader.cxx b/PWG3/muon/AliMuonsHFHeader.cxx index df87c52fcb0..515f81b1975 100644 --- a/PWG3/muon/AliMuonsHFHeader.cxx +++ b/PWG3/muon/AliMuonsHFHeader.cxx @@ -1,26 +1,126 @@ -#include +/************************************************************************** + * Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +///////////////////////////////////////////////////////////// +// +// class used to extract and store info at event level +// +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +///////////////////////////////////////////////////////////// + #include +#include +#include +#include +#include +#include +#include -#include "AliVHeader.h" -#include "AliVVertex.h" +#include "AliTriggerAnalysis.h" +#include "AliBackgroundSelection.h" +#include "AliAODVertex.h" +#include "AliAODEvent.h" +#include "AliESDEvent.h" +#include "AliMultiplicity.h" +#include "AliMuonInfoStoreRD.h" +#include "AliMuonInfoStoreMC.h" +#include "AliDimuInfoStoreRD.h" +#include "AliDimuInfoStoreMC.h" #include "AliMuonsHFHeader.h" +class TNamed; +class AliESDVertex; + ClassImp(AliMuonsHFHeader) +const TString AliMuonsHFHeader::fgkStdBranchName("MuEvsH"); +Int_t AliMuonsHFHeader::fgAnaMode = 0; +Bool_t AliMuonsHFHeader::fgIsMC = kFALSE; +Bool_t AliMuonsHFHeader::fgIsEventSelected = kFALSE; +Double_t AliMuonsHFHeader::fgCuts[3] = { -999999., 999999., 999999.}; + //_____________________________________________________________________________ AliMuonsHFHeader::AliMuonsHFHeader() : TNamed(), fTriggerMask(0), +fFiredTrigger(0), +fNFiredTrigger(0), +fIsPhysicsTriggered(kFALSE), +fIsPhysicsAccepted(kFALSE), +fEventType(0), +fUnrecoVertex(kFALSE), fNContributors(0), -fMultMuon(0), -fMultDimuon(0), -fCentrality(0.), -fUnrecoVertex(kFALSE) +fUnrecoVtxSPD(kFALSE), +fNContributorsSPD(0), +fNTrackletsSPD(0), +fCentrality(0.) { // // default constructor // - for (Int_t i=0; i<3; i++) fPosition[i]=0.; + for (Int_t i=3; i--;) fVtx[i] = 0.; + for (Int_t i=3; i--;) fVtxSPD[i] = 0.; +} + +//_____________________________________________________________________________ +AliMuonsHFHeader::AliMuonsHFHeader(const AliMuonsHFHeader &src) : +TNamed(), +fTriggerMask(src.fTriggerMask), +fFiredTrigger(src.fFiredTrigger), +fNFiredTrigger(src.fNFiredTrigger), +fIsPhysicsTriggered(src.fIsPhysicsTriggered), +fIsPhysicsAccepted(src.fIsPhysicsAccepted), +fEventType(src.fEventType), +fUnrecoVertex(src.fUnrecoVertex), +fNContributors(src.fNContributors), +fUnrecoVtxSPD(src.fUnrecoVtxSPD), +fNContributorsSPD(src.fNContributorsSPD), +fNTrackletsSPD(src.fNTrackletsSPD), +fCentrality(src.fCentrality) +{ + // + // copy constructor + // + for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i]; + for (Int_t i=3; i--;) fVtxSPD[i] = src.fVtxSPD[i]; +} + +//_____________________________________________________________________________ +AliMuonsHFHeader& AliMuonsHFHeader::operator=(const AliMuonsHFHeader &src) +{ + // + // assignment constructor + // + fTriggerMask = src.fTriggerMask; + fFiredTrigger = src.fFiredTrigger; + fNFiredTrigger = src.fNFiredTrigger; + fIsPhysicsTriggered = src.fIsPhysicsTriggered; + fIsPhysicsAccepted = src.fIsPhysicsAccepted; + fEventType = src.fEventType; + fUnrecoVertex = src.fUnrecoVertex; + fNContributors = src.fNContributors; + fUnrecoVtxSPD = src.fUnrecoVtxSPD; + fNContributorsSPD = src.fNContributorsSPD; + fNTrackletsSPD = src.fNTrackletsSPD; + fCentrality = src.fCentrality; + + for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i]; + for (Int_t i=3; i--;) fVtxSPD[i] = src.fVtxSPD[i]; + + return *this; } //_____________________________________________________________________________ @@ -32,18 +132,341 @@ AliMuonsHFHeader::~AliMuonsHFHeader() } //_____________________________________________________________________________ -void AliMuonsHFHeader::SetHeader(AliVHeader *header) +void AliMuonsHFHeader::SetEvent(AliAODEvent *event) { - fTriggerMask = header->GetTriggerMask(); + // extract event info from AOD event + + fTriggerMask = event->GetTriggerMask(); + + AliAODVertex *vertex = event->GetPrimaryVertex(); + vertex->GetXYZ(fVtx); + fNContributors = vertex->GetNContributors(); + fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 && + TMath::Abs(fVtx[2])<1e-6); + this->SetTitle(vertex->GetTitle()); + + this->SetFiredTrigger(event->GetFiredTriggerClasses()); + this->EventSelection(); return; } //_____________________________________________________________________________ -void AliMuonsHFHeader::SetVertex(AliVVertex *vertex) +void AliMuonsHFHeader::SetEvent(AliESDEvent *event) { - vertex->GetXYZ(fPosition); + // extract event info from ESD event + // SPD vertex and Physics event selection are implimented + + fTriggerMask = event->GetTriggerMask(); + + const AliESDVertex *vertex = event->GetPrimaryVertex(); + vertex->GetXYZ(fVtx); fNContributors = vertex->GetNContributors(); - fUnrecoVertex = (TMath::Abs(fPosition[0])<1e-6 && TMath::Abs(fPosition[1])<1e-6 && - TMath::Abs(fPosition[2])<1e-6); + fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 && + TMath::Abs(fVtx[2])<1e-6); + this->SetTitle(vertex->GetTitle()); + + const AliESDVertex *vtxSPD = event->GetPrimaryVertexSPD(); + vtxSPD->GetXYZ(fVtxSPD); + fNContributorsSPD = vtxSPD->GetNContributors(); + fUnrecoVtxSPD = (TMath::Abs(fVtxSPD[0])<1e-6 && TMath::Abs(fVtxSPD[1])<1e-6 && + TMath::Abs(fVtxSPD[2])<1e-6); + fNTrackletsSPD = event->GetMultiplicity()->GetNumberOfTracklets(); + + this->SetFiredTrigger(event->GetFiredTriggerClasses()); + this->PhysicsTriggerAna(event); + this->EventSelection(); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::PhysicsTriggerAna(const AliESDEvent *esd) +{ + // ESD event trigger condition analysis + // according to the method in $ALICE_ROOT/ANALYSIS/AliPhysicsSelection.cxx + + fEventType = esd->GetHeader()->GetEventType(); + fIsPhysicsTriggered = kFALSE; + fIsPhysicsAccepted = kFALSE; + + AliTriggerAnalysis *triggerAna = new AliTriggerAnalysis(); + triggerAna->SetAnalyzeMC(fgIsMC); + triggerAna->SetSPDGFOThreshhold(1); + + Int_t triggerHW = triggerAna->SPDFiredChips(esd, 1); + Bool_t isFiredv0A = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A); + Bool_t isFiredv0C = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C); + if (triggerHW==0 && !isFiredv0A && !isFiredv0C) { + delete triggerAna; + triggerAna = 0; + return; + } + + Bool_t triggerBG = (triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG) || + triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG)); + Int_t isFiredSPD = triggerAna->SPDFiredChips(esd, 0); + Bool_t triggerFD = ((isFiredSPD>1) || (isFiredSPD>0 && (isFiredv0A || isFiredv0C)) || (isFiredv0A && isFiredv0C)); + if ((!triggerBG) && triggerFD) fIsPhysicsTriggered = kTRUE; + delete triggerAna; + triggerAna = 0; + + if (fIsPhysicsTriggered) { + AliBackgroundSelection *bkgId = new AliBackgroundSelection(); + fIsPhysicsAccepted = bkgId->IsSelected(const_cast(esd)); + delete bkgId; bkgId=0; + } + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::SetFiredTrigger(TString trigger) +{ + // get info of fired trigger classes of event + + fFiredTrigger = trigger; + TObjArray *tokens = trigger.Tokenize(" "); + fNFiredTrigger = tokens->GetEntries(); + return; +} + +//_____________________________________________________________________________ +Bool_t AliMuonsHFHeader::IsTriggerFired(TString trigger) +{ + // check whether the trigger class "trigger" is fired in this event + + if (fNFiredTrigger<=0) return kFALSE; + + TObjArray *tokens = fFiredTrigger.Tokenize(" "); + for (Int_t i=fNFiredTrigger; i--;) { + TString fired = ((TObjString*)tokens->At(i))->String(); + if (fired.CompareTo(trigger)==0) return kTRUE; + } + return kFALSE; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::EventSelection(TString triggerName) +{ + // select event according to the "triggerName" & event selection cuts + + fgIsEventSelected = kFALSE; + if (!this->IsPhysicsAccepted()) return; + if (!fgIsMC && !this->IsTriggerFired(triggerName)) return; + this->EventSelection(); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::EventSelection() +{ + // select event according to the event selection cuts + + fgIsEventSelected = kFALSE; + if (this->NVtxContributorsSPD()VzSPD())>fgCuts[1]) return; + if (this->VtSPD()>fgCuts[2]) return; + fgIsEventSelected = kTRUE; + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistograms(TList *listEvent, TList *listMuon, TList *listDimu) +{ + // create output histos of muon analysis according to the analysis mode & MC flag + + this->CreateHistosEventH(listEvent); + if (fgIsMC) { + if (fgAnaMode!=2) this->CreateHistosMuonMC(listMuon); + if (fgAnaMode!=1) this->CreateHistosDimuMC(listDimu); + } else { + if (fgAnaMode!=2) this->CreateHistosMuonRD(listMuon); + if (fgAnaMode!=1) this->CreateHistosDimuRD(listDimu); + } + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistosEventH(TList *list) +{ + // create histograms at event level + + if (!list) list = new TList(); + list->SetOwner(); + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + const Int_t nHistos = 5; + char *name[nHistos] = { "hVztSPD", "hVzxSPD", "hVzySPD", "hVzNcontr", "hVtNcontr" }; + Int_t nbinsX[nHistos] = { 800 , 800 , 800 , 800 , 400 }; + Double_t xlow[nHistos] = { -40., -40., -40., -40. , 0. }; + Double_t xup[nHistos] = { 40., 40., 40., 40. , 4. }; + Int_t nbinsY[nHistos] = { 400 , 600 , 600 , 202 , 202 }; + Double_t ylow[nHistos] = { 0., -3., -3., -2.5, -2.5 }; + Double_t yup[nHistos] = { 4., 3., 3., 199.5, 199.5 }; + + TH2F *histo = 0; + for (Int_t i=0; iSumw2(); list->AddAt(histo, i); histo = 0; + } + + TH1::AddDirectory(oldStatus); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistosMuonRD(TList *list) +{ + // create histograms for single muon + + if (!list) list = new TList(); + list->SetOwner(); + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + const Int_t nHistos = 8; + char *name[nHistos] = { "hP", "hPt", "hEta", "hDCA", "hTrg", "hCharge", "hUnfVtx" , "hEtaTimesDCA"}; + Int_t nbins[nHistos] = { 1500 , 300 , 100 , 500 , 4 , 3 , 80 , 10000 }; + Double_t xlow[nHistos] = { 0., 0., -10., 0., -0.5, -1.5, -40. , 0.}; + Double_t xup[nHistos] = { 150., 30., 0., 500., 3.5, 1.5, 40. , 1000.}; + + TH1F *histo = 0; + for (Int_t i=0; iSumw2(); list->AddAt(histo, i); histo = 0; + } + + TH1::AddDirectory(oldStatus); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistosDimuRD(TList *list) +{ + // create histograms for dimuon + + if (!list) list = new TList(); + list->SetOwner(); + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + TH1F *histo = 0; + const Int_t nHistos = 3; + char *name[nHistos] = { "hP", "hPt", "hInvM" }; + Int_t nbins[nHistos] = { 1500 , 300 , 300 }; + Double_t xlow[nHistos] = { 0., 0., 0. }; + Double_t xup[nHistos] = { 150., 30., 30. }; + char *dimuName[3] = {"DimuNN", "DimuNP", "DimuPP"}; + for (Int_t i=0; i<3; i++) { + for (Int_t j=0; jSumw2(); list->AddAt(histo,i*nHistos+j); histo = 0; + } + } + + TH1::AddDirectory(oldStatus); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistosMuonMC(TList *list) +{ + // create histograms for single muon with MC info + + if (!list) list = new TList[AliMuonInfoStoreMC::NSources()]; + for (Int_t i=AliMuonInfoStoreMC::NSources(); i--;) this->CreateHistosMuonRD(&list[i]); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::CreateHistosDimuMC(TList *list) +{ + // create histograms for dimuon with MC info + + if (!list) list = new TList[AliDimuInfoStoreMC::NSources()]; + for (Int_t i=AliDimuInfoStoreMC::NSources(); i--;) this->CreateHistosDimuRD(&list[i]); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::FillHistosEventH(TList *list) +{ + // fill histograms at event level according to event selection cuts + + if (!list) return; + if (!AliMuonsHFHeader::IsEventSelected()) return; + + const Int_t nHistos = 5; + Double_t vz = this->VzSPD(); + Double_t vt = this->VtSPD(); + Int_t nc = this->NVtxContributorsSPD(); + Double_t distX[nHistos] = { vz, vz, vz, vz, vt }; + Double_t distY[nHistos] = { vt, this->VxSPD(), this->VySPD(), nc, nc }; + for (Int_t i=nHistos; i--;) ((TH2F*)list->At(i))->Fill(distX[i], distY[i]); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::FillHistosMuonRD(TList *list, AliMuonInfoStoreRD* const muonStoreRD) +{ + // fill histograms for single muon according to event & muon track selection cuts + + if (!list) return; + if (!AliMuonsHFHeader::IsEventSelected()) return; + if (!muonStoreRD->MuonSelection()) return; + + const Int_t nHistos = 8; + Double_t dist[nHistos] = {muonStoreRD->Momentum().Mag(), + muonStoreRD->Momentum().Pt(), + muonStoreRD->Momentum().Eta(), + muonStoreRD->DCA(), + muonStoreRD->MatchTrigger(), + muonStoreRD->Charge(), + this->VzSPD(), + muonStoreRD->DCA()*TMath::Exp(-2.*muonStoreRD->Momentum().Eta())/1000.}; + for (Int_t i=nHistos; i--;) ((TH1F*)list->At(i))->Fill(dist[i]); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::FillHistosDimuRD(TList *list, AliDimuInfoStoreRD* const dimuStoreRD) +{ + // fill histograms for dimuon according to evnet & dimuon candidates selection cuts + + if (!list) return; + if (!AliMuonsHFHeader::IsEventSelected()) return; + if (!dimuStoreRD->DimuSelection()) return; + + Int_t theDimu = 0; + if (dimuStoreRD->Charge()==0) theDimu = 1; + else if (dimuStoreRD->Charge()>0) theDimu = 2; + + const Int_t nHistos = 3; + Double_t dist[nHistos] = {dimuStoreRD->Momentum().Mag(), + dimuStoreRD->Momentum().Pt(), + dimuStoreRD->InvM()}; + for (Int_t i=nHistos; i--;) ((TH1F*)list->At(theDimu*nHistos+i))->Fill(dist[i]); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::FillHistosMuonMC(TList *list, AliMuonInfoStoreMC* const muonStoreMC) +{ + // fill histograms for single muon with MC info + + Int_t srcMuon = muonStoreMC->MuonSource(); + this->FillHistosMuonRD(&list[6], (AliMuonInfoStoreRD*)muonStoreMC); + if (srcMuon<0) this->FillHistosMuonRD(&list[5], (AliMuonInfoStoreRD*)muonStoreMC); + else this->FillHistosMuonRD(&list[srcMuon], (AliMuonInfoStoreRD*)muonStoreMC); + return; +} + +//_____________________________________________________________________________ +void AliMuonsHFHeader::FillHistosDimuMC(TList *list, AliDimuInfoStoreMC* const dimuStoreMC) +{ + // fill histograms for dimuon with MC info + + Int_t srcDimu = dimuStoreMC->DimuSource(); + this->FillHistosDimuRD(&list[6], (AliDimuInfoStoreRD*)dimuStoreMC); + this->FillHistosDimuRD(&list[srcDimu], (AliDimuInfoStoreRD*)dimuStoreMC); return; } diff --git a/PWG3/muon/AliMuonsHFHeader.h b/PWG3/muon/AliMuonsHFHeader.h index 56120e598d2..10c66a7ba45 100644 --- a/PWG3/muon/AliMuonsHFHeader.h +++ b/PWG3/muon/AliMuonsHFHeader.h @@ -1,53 +1,119 @@ #ifndef ALIMUONSHFHEADER_H #define ALIMUONSHFHEADER_H +/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//*********************************************************** +// Class AliMuonsHFHeader +// class used to extract and store info at event level +// Author: X-M. Zhang, zhang@clermont.in2p3.fr +// zhangxm@iopp.ccnu.edu.cn +//*********************************************************** + #include +#include + +class TList; +class AliAODEvent; +class AliESDEvent; -#include "AliVHeader.h" -#include "AliVVertex.h" +class AliMuonInfoStoreRD; +class AliDimuInfoStoreRD; +class AliMuonInfoStoreMC; +class AliDimuInfoStoreMC; class AliMuonsHFHeader : public TNamed { public : AliMuonsHFHeader(); + AliMuonsHFHeader(const AliMuonsHFHeader &src); + AliMuonsHFHeader& operator=(const AliMuonsHFHeader &src); ~AliMuonsHFHeader(); - void SetHeader(AliVHeader *header); - void SetVertex(AliVVertex *vertex); - - void SetMultSingleMuon(Int_t mul) { fMultMuon = mul; } - void SetMultDimuon(Int_t mul) { fMultDimuon = mul; } - void SetCentrality(Double_t cen) { fCentrality = cen; } - - void GetXYZ(Double_t *pos) const - { for (Int_t i=0; i<3; i++) pos[i]=fPosition[i]; } - Double_t GetXv() const { return fPosition[0]; } - Double_t GetYv() const { return fPosition[1]; } - Double_t GetZv() const { return fPosition[2]; } - Int_t GetNContributors() const { return fNContributors; } - - ULong64_t GetTriggerMask() const { return fTriggerMask; } - - Int_t GetMultSingleMuon() const { return fMultMuon; } - Int_t GetMultDimuon() const { return fMultDimuon; } - Double_t GetCentrality() const { return fCentrality; } - - Bool_t IsUnrecoVertex() const { return fUnrecoVertex; } + ULong64_t TriggerMask() const { return fTriggerMask; } + TString FiredTrigger() const { return fFiredTrigger; } + Int_t NFiredTrigger() const { return fNFiredTrigger; } + Int_t EventType() const { return fEventType; } + Bool_t IsPhysicsTriggered() const { return fIsPhysicsTriggered; } + Bool_t IsPhysicsAccepted() const { return fIsPhysicsAccepted; } + Bool_t IsTriggerFired(TString trigger); + + void GetXYZ(Double_t *vtx) const { for (Int_t i=3; i--;) vtx[i]=fVtx[i]; } + Double_t Vx() const { return fVtx[0]; } + Double_t Vy() const { return fVtx[1]; } + Double_t Vz() const { return fVtx[2]; } + Double_t Vt() const { return TMath::Sqrt(fVtx[0]*fVtx[0] + fVtx[1]*fVtx[1]); } + Bool_t IsUnrecoVertex() const { return fUnrecoVertex; } + Int_t NVtxContributors() const { return fNContributors; } + + void GetXYZSPD(Double_t *vtx) const { for (Int_t i=3; i--;) vtx[i]=fVtxSPD[i]; } + Double_t VxSPD() const { return fVtxSPD[0]; } + Double_t VySPD() const { return fVtxSPD[1]; } + Double_t VzSPD() const { return fVtxSPD[2]; } + Double_t VtSPD() const { return TMath::Sqrt(fVtxSPD[0]*fVtxSPD[0] + fVtxSPD[1]*fVtxSPD[1]); } + Bool_t IsUnrecoVtxSPD() const { return fUnrecoVtxSPD; } + Int_t NVtxContributorsSPD() const { return fNContributorsSPD; } + Int_t NTrackletsSPD() const { return fNTrackletsSPD; } + + Double_t Centrality() const { return fCentrality; } + + void SetEvent(AliAODEvent *event); + void SetEvent(AliESDEvent *event); + + void EventSelection(TString triggerName); + void CreateHistograms(TList *listEvent=0, TList *listMuon=0, TList *listDimu=0); + void CreateHistosEventH(TList *list); + void CreateHistosMuonRD(TList *list); + void CreateHistosDimuRD(TList *list); + void CreateHistosMuonMC(TList *list); + void CreateHistosDimuMC(TList *list); + + void FillHistosEventH(TList *list); + void FillHistosMuonRD(TList *list, AliMuonInfoStoreRD* const muonStoreRD); + void FillHistosDimuRD(TList *list, AliDimuInfoStoreRD* const dimuStoreRD); + void FillHistosMuonMC(TList *list, AliMuonInfoStoreMC* const muonStoreMC); + void FillHistosDimuMC(TList *list, AliDimuInfoStoreMC* const dimuStoreMC); + + static const char* StdBranchName() { return fgkStdBranchName.Data(); } + static const Bool_t IsEventSelected() { return fgIsEventSelected; } + static void SetAnaMode(Int_t anaMode=0) { fgAnaMode=anaMode; } + static void SetIsMC(Int_t isMC=kFALSE) { fgIsMC =isMC; } + static void SetSelectionCuts(Double_t cuts[3]) { for (Int_t i=3; i--;) fgCuts[i]=cuts[i]; } private : - ULong64_t fTriggerMask; // trigger mask - - Double32_t fPosition[3]; // position of vtx - Int_t fNContributors; // number of contributor for vtx - - Int_t fMultMuon; // event multiplicity - Int_t fMultDimuon; - Double_t fCentrality; // event centrality class - - Bool_t fUnrecoVertex; - - ClassDef(AliMuonsHFHeader, 1) + void SetFiredTrigger(TString str); + void PhysicsTriggerAna(const AliESDEvent *esd); + void EventSelection(); + + static const TString fgkStdBranchName; // Standard branch name + static Bool_t fgIsEventSelected; // flag for event selection + static Int_t fgAnaMode; // analysis mode + static Bool_t fgIsMC; // flag to use MC + static Double_t fgCuts[3]; // 0, low limit of num. of vtx contributors + // 1, up limit of vz + // 2, up limit of vt + + ULong64_t fTriggerMask; // trigger mask + TString fFiredTrigger; // fired of trigger class of event + Int_t fNFiredTrigger; // num. of fired trigger class + Bool_t fIsPhysicsTriggered; // flag of final physics trigger from AliPhysicsSelection + Bool_t fIsPhysicsAccepted; // flag of physiscs selection w/ BKG Id + Int_t fEventType; // event type + + Double32_t fVtx[3]; // position of vtx + Bool_t fUnrecoVertex; // flag for unreco vtx + Int_t fNContributors; // num. of contributors of vtx rec + + Double32_t fVtxSPD[3]; // position of vtx + Bool_t fUnrecoVtxSPD; // flag for unreco vtx + Int_t fNContributorsSPD; // num. of contributors of vtx rec + Int_t fNTrackletsSPD; // num. of SPD tracklets + + Double_t fCentrality; // event centrality class + + ClassDef(AliMuonsHFHeader, 2) }; #endif -- 2.43.0