1 /**************************************************************************
2 * Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 /////////////////////////////////////////////////////////////
20 // class used to extract and store info at event level
22 // Author: X-M. Zhang, zhang@clermont.in2p3.fr
23 // zhangxm@iopp.ccnu.edu.cn
24 /////////////////////////////////////////////////////////////
31 #include "AliInputEventHandler.h"
32 #include "AliAODMCHeader.h"
33 #include "AliVEvent.h"
34 #include "AliAODEvent.h"
35 #include "AliESDEvent.h"
36 #include "AliVVertex.h"
37 #include "AliCentrality.h"
38 #include "AliEventplane.h"
40 #include "AliMuonTrackCuts.h"
41 #include "AliMuonInfoStoreRD.h"
42 #include "AliMuonInfoStoreMC.h"
43 #include "AliDimuInfoStoreRD.h"
44 #include "AliDimuInfoStoreMC.h"
45 #include "AliMuonsHFHeader.h"
47 ClassImp(AliMuonsHFHeader)
49 const TString AliMuonsHFHeader::fgkStdBranchName("MuEvsH");
50 Int_t AliMuonsHFHeader::fgAnaMode = 0;
51 Bool_t AliMuonsHFHeader::fgIsMC = kFALSE;
52 Double_t AliMuonsHFHeader::fgCuts[5] = { -999999., 999999., 999999., -999999., 999999. };
54 //_____________________________________________________________________________
55 AliMuonsHFHeader::AliMuonsHFHeader() :
57 fSelMask(AliVEvent::kAny),
68 // default constructor
70 for (Int_t i=3; i--;) fVtx[i] = 0.;
71 for (Int_t i=3; i--;) fVMC[i] = 0.;
74 //_____________________________________________________________________________
75 AliMuonsHFHeader::AliMuonsHFHeader(const AliMuonsHFHeader &src) :
77 fSelMask(src.fSelMask),
80 fIsPileupSPD(src.fIsPileupSPD),
81 fVtxContrsN(src.fVtxContrsN),
82 fFiredTriggerClass(src.fFiredTriggerClass),
83 fCentrality(src.fCentrality),
85 fEventPlane(src.fEventPlane)
90 for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i];
91 for (Int_t i=3; i--;) fVMC[i] = src.fVMC[i];
94 //_____________________________________________________________________________
95 AliMuonsHFHeader& AliMuonsHFHeader::operator=(const AliMuonsHFHeader &src)
98 // assignment constructor
101 if(&src==this) return *this;
103 fSelMask = src.fSelMask;
106 fIsPileupSPD = src.fIsPileupSPD;
107 fVtxContrsN = src.fVtxContrsN;
108 fFiredTriggerClass = src.fFiredTriggerClass;
109 fCentrality = src.fCentrality;
110 fCentQA = src.fCentQA;
111 fEventPlane = src.fEventPlane;
112 for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i];
113 for (Int_t i=3; i--;) fVMC[i] = src.fVMC[i];
118 //_____________________________________________________________________________
119 AliMuonsHFHeader::~AliMuonsHFHeader()
122 // default destructor
126 //_____________________________________________________________________________
127 void AliMuonsHFHeader::SetEventInfo(AliInputEventHandler* const handler, AliMCEvent* const eventMC)
129 // fill info at event level
131 AliVEvent *event = handler->GetEvent();
132 AliAODEvent *aod = dynamic_cast<AliAODEvent*>(event);
133 AliESDEvent *esd = dynamic_cast<AliESDEvent*>(event);
135 fSelMask = handler->IsEventSelected();
136 if (aod) fFiredTriggerClass = aod->GetFiredTriggerClasses();
137 if (esd) fFiredTriggerClass = esd->GetFiredTriggerClasses();
138 fIsMB = fSelMask & AliVEvent::kMB;
139 fIsMU = fSelMask & AliVEvent::kMUON;
141 const AliVVertex *vertex = event->GetPrimaryVertex();
142 vertex->GetXYZ(fVtx);
143 fVtxContrsN = vertex->GetNContributors();
145 if (esd) eventMC->GetPrimaryVertex()->GetXYZ(fVMC);
146 if (aod) ((AliAODMCHeader*)aod->FindListObject(AliAODMCHeader::StdBranchName()))->GetVertex(fVMC);
147 } this->SetTitle(vertex->GetTitle());
148 fIsPileupSPD = (aod && !aod->GetTracklets()) ? event->IsPileupFromSPD(3,0.8,3.,2.,5.) : event->IsPileupFromSPDInMultBins();
150 AliCentrality *cent = event->GetCentrality();
152 fCentrality = cent->GetCentralityPercentileUnchecked("V0M");
153 fCentQA = cent->GetQuality();
156 AliEventplane *evnP = event->GetEventplane();
157 if (evnP) fEventPlane = evnP->GetEventplane("Q");
158 //if (evnP) fEventPlane = evnP->GetEventplane("V0A");
163 //_____________________________________________________________________________
164 Bool_t AliMuonsHFHeader::IsSelected()
166 // select event according to the event selection cuts
167 if (this->VtxContrsN()<fgCuts[0]) return kFALSE;
168 if (TMath::Abs(this->Vz())>fgCuts[1]) return kFALSE;
169 if (this->Vt()>fgCuts[2]) return kFALSE;
171 // centrality selection
172 Float_t centr = this->Centrality();
173 if (centr<fgCuts[3] || centr>fgCuts[4]) return kFALSE;
177 //_____________________________________________________________________________
178 void AliMuonsHFHeader::CreateHistograms(TList *list)
180 // create output histos of muon analysis according to the analysis mode & MC flag
183 this->CreateHistosEvnH(list);
185 TString sName[7] = { "Unidentified", "Hadron", "SecondaryMu", "PrimaryMu", "CharmMu", "BottomMu", "" };
186 for (Int_t i=7; i--;) this->CreateHistosMuon(list, sName[i]);
189 TString sName[7] = { "Uncorr", "Resonance", "DDsame", "DDdiff", "BBsame", "BBdiff", "" };
190 for (Int_t i=7; i--;) this->CreateHistosDimu(list, sName[i]);
195 this->CreateHistosEvnH(list,"MB"); this->CreateHistosEvnH(list,"MU");
196 if (fgAnaMode!=2) { this->CreateHistosMuon(list,"MB"); this->CreateHistosMuon(list,"MU"); }
197 if (fgAnaMode!=1) { this->CreateHistosDimu(list,"MB"); this->CreateHistosDimu(list,"MU"); }
201 //_____________________________________________________________________________
202 void AliMuonsHFHeader::CreateHistosEvnH(TList *list, TString sName)
204 // create histograms at event level
206 if (!list) list = new TList();
208 Bool_t oldStatus = TH1::AddDirectoryStatus();
209 TH1::AddDirectory(kFALSE);
212 TString tName[nhs] = { "Vz", "Vt", "VtxNcontr" };
213 Int_t nbins[nhs] = { 800 , 40 , 202 };
214 Double_t xlow[nhs] = { -40., 0., -2.5 };
215 Double_t xup[nhs] = { 40., 4., 199.5 };
218 for (Int_t i=0; i<nhs; i++) {
219 char *hName = Form("h%s_%s", sName.Data(), tName[i].Data());
220 histo = new TH1D(hName, hName, nbins[i], xlow[i], xup[i]);
221 histo->Sumw2(); list->Add(histo); histo = 0;
224 TH1::AddDirectory(oldStatus);
228 //_____________________________________________________________________________
229 void AliMuonsHFHeader::CreateHistosMuon(TList *list, TString sName)
231 // create histograms for single muon
233 if (!list) list = new TList();
235 Bool_t oldStatus = TH1::AddDirectoryStatus();
236 TH1::AddDirectory(kFALSE);
239 TString tName[nhs] = { "P", "Pt", "Eta", "DCA", "TrM", "Charge", "Rabs" };
240 Int_t nbins[nhs] = { 1500 , 300 , 15 , 1000 , 4 , 3 , 48 };
241 Double_t xlow[nhs] = { 0., 0., -4.0, 0., -0.5, -1.5, 17.6 };
242 Double_t xup[nhs] = { 150., 30., -2.5, 500., 3.5, 1.5, 80.0 };
245 for (Int_t i=0; i<nhs; i++) {
246 char *hName = Form("h%s_%s", sName.Data(), tName[i].Data());
247 histo = new TH1D(hName, hName, nbins[i], xlow[i], xup[i]);
248 histo->Sumw2(); list->Add(histo); histo = 0;
251 TH1::AddDirectory(oldStatus);
255 //_____________________________________________________________________________
256 void AliMuonsHFHeader::CreateHistosDimu(TList *list, TString sName)
258 // create histograms for dimuon
260 if (!list) list = new TList();
262 Bool_t oldStatus = TH1::AddDirectoryStatus();
263 TH1::AddDirectory(kFALSE);
267 TString tName[nhs] = { "P", "Pt", "InvM" };
268 Int_t nbins[nhs] = { 1500 , 300 , 300 };
269 Double_t xlow[nhs] = { 0., 0., 0. };
270 Double_t xup[nhs] = { 150., 30., 30. };
271 TString dimuName[3] = { "DimuNN", "DimuNP", "DimuPP" };
272 for (Int_t i=0; i<3; i++) {
273 for (Int_t j=0; j<nhs; j++) {
274 char *hName = Form("h%s_%s_%s", sName.Data(), dimuName[i].Data(), tName[j].Data());
275 histo = new TH1D(hName, hName, nbins[j], xlow[j], xup[j]);
276 histo->Sumw2(); list->Add(histo); histo = 0;
280 TH1::AddDirectory(oldStatus);
284 //_____________________________________________________________________________
285 void AliMuonsHFHeader::FillHistosEvnH(TList *list)
287 // fill histograms at event level according to event selection cuts
290 if (!this->IsSelected()) return;
293 TString tName[nhs] = { "Vz", "Vt", "VtxNcontr" };
294 Double_t dist[nhs] = { this->Vz(), this->Vt(), this->VtxContrsN() };
295 if (fgIsMC && (fSelMask & AliVEvent::kAny)) {
296 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h_%s",tName[i].Data())))->Fill(dist[i]);
298 if (fIsMB && (fSelMask & AliVEvent::kMB)) { for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s","MB",tName[i].Data())))->Fill(dist[i]); }
299 if (fIsMU && (fSelMask & AliVEvent::kMUON)) { for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s","MU",tName[i].Data())))->Fill(dist[i]); }
304 //_____________________________________________________________________________
305 void AliMuonsHFHeader::FillHistosMuon(TList *list, AliMuonInfoStoreRD* const infoStore, Int_t s)
307 // fill histograms for single muon according to event & muon track selection cuts
310 if (!this->IsSelected()) return;
311 if (!infoStore->IsSelected(0)) return;
314 TString tName[nhs] = { "P", "Pt", "Eta", "DCA", "TrM", "Charge", "Rabs" };
315 Double_t dist[nhs] = { infoStore->MomentumAtVtx().Mag(),
316 infoStore->MomentumAtVtx().Pt(),
317 infoStore->MomentumAtVtx().Eta(),
319 infoStore->MatchTrigger(),
321 infoStore->RabsEnd() };
323 if (fgIsMC && (fSelMask & AliVEvent::kAny)) {
324 TString sName[7] = { "BottomMu", "CharmMu", "PrimaryMu", "SecondaryMu", "Hadron", "Unidentified", "" };
325 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s",sName[6].Data(),tName[i].Data())))->Fill(dist[i]);
326 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s",sName[s].Data(),tName[i].Data())))->Fill(dist[i]);
328 if (fIsMB && (fSelMask & AliVEvent::kMB)) { for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s","MB",tName[i].Data())))->Fill(dist[i]); }
329 if (fIsMU && (fSelMask & AliVEvent::kMUON)) { for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s","MU",tName[i].Data())))->Fill(dist[i]); }
335 //_____________________________________________________________________________
336 void AliMuonsHFHeader::FillHistosDimu(TList *list, AliDimuInfoStoreRD* const infoStore, Int_t s)
338 // fill histograms for dimuon according to evnet & dimuon candidates selection cuts
341 if (!this->IsSelected()) return;
342 if (!infoStore->IsSelected(0)) return;
344 TString dimuName = "DimuNN";
345 if (infoStore->Charge()==0) dimuName = "DimuNP";
346 else if (infoStore->Charge()>0) dimuName = "DimuPP";
349 TString tName[nhs] = { "P", "Pt", "InvM" };
350 Double_t dist[nhs] = { infoStore->Momentum().Mag(),
351 infoStore->Momentum().Pt(),
354 if (fgIsMC && (fSelMask & AliVEvent::kAny)) {
355 TString sName[7] = { "BBdiff", "BBsame", "DDdiff", "DDsame", "Resonance", "Uncorr", "" };
356 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s_%s",sName[6].Data(),dimuName.Data(),tName[i].Data())))->Fill(dist[i]);
357 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s_%s",sName[s].Data(),dimuName.Data(),tName[i].Data())))->Fill(dist[i]);
359 if (fIsMB && (fSelMask & AliVEvent::kMB)) {
360 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s_%s","MB",dimuName.Data(),tName[i].Data())))->Fill(dist[i]);
362 if (fIsMU && (fSelMask & AliVEvent::kMUON)) {
363 for (Int_t i=nhs; i--;) ((TH1D*)list->FindObject(Form("h%s_%s_%s","MU",dimuName.Data(),tName[i].Data())))->Fill(dist[i]);