Fixes for CAF
[u/mrichter/AliRoot.git] / MFT / AliMuonForwardTrackAnalysis.cxx
CommitLineData
e806e863 1//====================================================================================================================================================
2//
3// Class for the analysis of the ALICE muon forward tracks (MUON + MFT)
4//
5// Contact author: antonio.uras@cern.ch
6//
7//====================================================================================================================================================
8
9#include "TObject.h"
10#include "TClonesArray.h"
11#include "AliMuonForwardTrack.h"
12#include "AliMuonForwardTrackPair.h"
13#include "TMatrixD.h"
14#include "TTree.h"
15#include "TH1D.h"
16#include "AliLog.h"
17#include "TFile.h"
18#include "TParticle.h"
19#include "AliMUONTrackParam.h"
20#include "AliMUONTrackExtrap.h"
21#include "TGeoManager.h"
22#include "TRandom.h"
23#include "TLorentzVector.h"
24#include "TDatabasePDG.h"
b5ab1ac4 25#include "TGraph.h"
e806e863 26#include "AliMuonForwardTrackAnalysis.h"
27
28ClassImp(AliMuonForwardTrackAnalysis)
29
30//====================================================================================================================================================
31
32AliMuonForwardTrackAnalysis::AliMuonForwardTrackAnalysis():
33 TObject(),
34 fInputDir(0),
35 fOutputDir(0),
7e3dd1af 36 fInputTreeWithBranson(0x0),
37 fInputTreeWithoutBranson(0x0),
38 fMuonForwardTracksWithBranson(0),
39 fMuonForwardTrackPairsWithBranson(0),
40 fMuonForwardTracksWithoutBranson(0),
41 fMuonForwardTrackPairsWithoutBranson(0),
42 fMFTTrackWithBranson(0),
43 fMFTTrackWithoutBranson(0),
e806e863 44 fMFTTrack(0),
7e3dd1af 45 fMFTTrackPairWithBranson(0),
46 fMFTTrackPairWithoutBranson(0),
e806e863 47 fMFTTrackPair(0),
48 fMCRefTrack(0),
49 fEv(0),
50 fFirstEvent(-1),
51 fLastEvent(-1),
52 fNTracksOfEvent(0),
53 fNTracksAnalyzedOfEvent(0),
54 fNTracksAnalyzed(0),
55 fNPairsOfEvent(0),
56 fNPairsAnalyzedOfEvent(0),
d8c2cc3e 57 fNTracksAnalyzedOfEventAfterCut(0),
58 fNPairsAnalyzedOfEventAfterCut(0),
e806e863 59 fHistOffsetSingleMuonsX(0x0),
60 fHistOffsetSingleMuonsY(0x0),
61 fHistOffsetSingleMuons(0x0),
62 fHistWOffsetSingleMuons(0x0),
63 fHistErrorSingleMuonsX(0x0),
64 fHistErrorSingleMuonsY(0x0),
d8c2cc3e 65 fHistZOriginSingleMuonsMC(0x0),
cabfe25b 66 fHistZROriginSingleMuonsMC(0x0),
67 fHistSingleMuonsPtRapidity(0x0),
b5ab1ac4 68 fHistSingleMuonsOffsetChi2(0x0),
d8c2cc3e 69 fHistRapidityPtMuonPairs(0x0),
cabfe25b 70 fHistMassMuonPairsMCVsPt(0x0),
71 fHistMassMuonPairsVsPt(0x0),
72 fHistMassMuonPairsWithoutMFTVsPt(0x0),
73 fHistMassMuonPairsVsPtLSp(0x0),
74 fHistMassMuonPairsWithoutMFTVsPtLSp(0x0),
75 fHistMassMuonPairsVsPtLSm(0x0),
76 fHistMassMuonPairsWithoutMFTVsPtLSm(0x0),
77 fHistWOffsetMuonPairsAtPrimaryVtxVsPt(0x0),
78 fHistWOffsetMuonPairsAtPCAVsPt(0x0),
79 fHistDistancePrimaryVtxPCAVsPt(0x0),
80 fHistPCAQualityVsPt(0x0),
81 fHistPseudoProperDecayLengthVsPt(0x0),
82 fHistDimuonVtxResolutionXVsPt(0x0),
83 fHistDimuonVtxResolutionYVsPt(0x0),
84 fHistDimuonVtxResolutionZVsPt(0x0),
85 fTrueMass(0.),
e806e863 86 fMassMin(0),
cabfe25b 87 fMassMax(9999),
e806e863 88 fSingleMuonAnalysis(1),
89 fMuonPairAnalysis(1),
e806e863 90 fOption(0),
cabfe25b 91 fTriggerLevel(0),
b5ab1ac4 92 fXVertResMC(50.e-4),
93 fYVertResMC(50.e-4),
94 fZVertResMC(50.e-4),
7e3dd1af 95 fPrimaryVtxX(0.),
96 fPrimaryVtxY(0.),
97 fPrimaryVtxZ(0.),
e806e863 98 fMaxNWrongClustersMC(999),
7e3dd1af 99 fPtMinSingleMuons(0),
100 fUseBransonForCut(kFALSE),
101 fUseBransonForKinematics(kFALSE),
cabfe25b 102 fCorrelateCutOnOffsetChi2(kFALSE),
103 fMaxChi2SingleMuons(1.e9),
104 fMaxOffsetSingleMuons(1.e9),
105 fMaxWOffsetMuonPairsAtPrimaryVtx(1.e9),
106 fMaxWOffsetMuonPairsAtPCA(1.e9),
107 fMaxDistancePrimaryVtxPCA(1.e9),
108 fMinPCAQuality(0.)
e806e863 109{
110
111 // default constructor
112
e806e863 113}
114
115//====================================================================================================================================================
116
117Bool_t AliMuonForwardTrackAnalysis::Init(Char_t *inputFileName) {
118
119 BookHistos();
120
121 TFile *inputFile = new TFile(Form("%s/%s",fInputDir.Data(),inputFileName));
122 if (!inputFile || !inputFile->IsOpen()) {
123 AliError(Form("Error opening file %s", inputFileName));
124 return kFALSE;
125 }
7e3dd1af 126 fInputTreeWithBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithBranson");
127 if (!fInputTreeWithBranson) {
128 AliError("Error reading input tree");
129 return kFALSE;
130 }
131 fInputTreeWithoutBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithoutBranson");
132 if (!fInputTreeWithoutBranson) {
e806e863 133 AliError("Error reading input tree");
134 return kFALSE;
135 }
136
789c2e43 137 if (fFirstEvent<0 || fLastEvent<0 || fFirstEvent>fLastEvent || fFirstEvent>=fInputTreeWithBranson->GetEntriesFast()) {
e806e863 138 fFirstEvent = 0;
789c2e43 139 fLastEvent = fInputTreeWithBranson->GetEntriesFast()-1;
e806e863 140 }
141 else {
789c2e43 142 fLastEvent = TMath::Min(fLastEvent, Int_t(fInputTreeWithBranson->GetEntriesFast()-1));
e806e863 143 }
144
145 AliInfo(Form("Analysing events %d to %d", fFirstEvent, fLastEvent));
146
789c2e43 147 fMuonForwardTracksWithBranson = new TClonesArray("AliMuonForwardTrack",30);
7e3dd1af 148 fInputTreeWithBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithBranson);
149
789c2e43 150 fMuonForwardTracksWithoutBranson = new TClonesArray("AliMuonForwardTrack",30);
7e3dd1af 151 fInputTreeWithoutBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithoutBranson);
e806e863 152
153 TGeoManager::Import(Form("%s/geometry.root",fInputDir.Data()));
154
155 AliMUONTrackExtrap::SetField();
156
789c2e43 157 fMuonForwardTrackPairsWithBranson = new TClonesArray("AliMuonForwardTrackPair",10);
158 fMuonForwardTrackPairsWithoutBranson = new TClonesArray("AliMuonForwardTrackPair",10);
cabfe25b 159 fMuonForwardTrackPairsWithBranson -> SetOwner(kTRUE);
160 fMuonForwardTrackPairsWithoutBranson -> SetOwner(kTRUE);
161
e806e863 162 return kTRUE;
163
164}
165
166//====================================================================================================================================================
167
168Bool_t AliMuonForwardTrackAnalysis::LoadNextEvent() {
169
170 if (fEv>fLastEvent) return kFALSE;
171 if (fEv<fFirstEvent) { fEv++; return kTRUE; }
789c2e43 172 fMuonForwardTracksWithBranson -> Clear("");
173 fMuonForwardTracksWithoutBranson -> Clear("");
7e3dd1af 174 fInputTreeWithBranson->GetEvent(fEv);
175 fInputTreeWithoutBranson->GetEvent(fEv);
789c2e43 176 AliDebug(2,Form("**** analyzing event # %4d (%3d tracks) ****", fEv, fMuonForwardTracksWithBranson->GetEntriesFast()));
7e3dd1af 177
cabfe25b 178 AliInfo(Form("**** analyzing event # %6d of %6d ****", fEv, fLastEvent));
179
7e3dd1af 180 fPrimaryVtxX = gRandom->Gaus(0., fXVertResMC);
181 fPrimaryVtxY = gRandom->Gaus(0., fYVertResMC);
182 fPrimaryVtxZ = gRandom->Gaus(0., fZVertResMC);
e806e863 183
184 if (fSingleMuonAnalysis) {
185 fNTracksAnalyzedOfEvent = 0;
d8c2cc3e 186 fNTracksAnalyzedOfEventAfterCut = 0;
789c2e43 187 fNTracksOfEvent = fMuonForwardTracksWithBranson->GetEntriesFast();
e806e863 188 while (AnalyzeSingleMuon()) continue;
189 }
190
191 if (fMuonPairAnalysis) {
7e3dd1af 192 if (fMuonForwardTrackPairsWithBranson) {
789c2e43 193 fMuonForwardTrackPairsWithBranson->Clear("C");
194 fMuonForwardTrackPairsWithoutBranson->Clear("C");
e806e863 195 }
196 BuildMuonPairs();
197 fNPairsAnalyzedOfEvent = 0;
d8c2cc3e 198 fNPairsAnalyzedOfEventAfterCut = 0;
cabfe25b 199 fNPairsOfEvent = fMuonForwardTrackPairsWithBranson->GetEntriesFast();
e806e863 200 while (AnalyzeMuonPair()) continue;
201 }
202
789c2e43 203 AliDebug(2,Form("**** analyzed event # %4d (%3d tracks and %3d pairs analyzed) ****", fEv, fNTracksAnalyzedOfEventAfterCut, fNPairsAnalyzedOfEventAfterCut));
d8c2cc3e 204
e806e863 205 fEv++;
206
207 return kTRUE;
208
209}
210
211//====================================================================================================================================================
212
213Bool_t AliMuonForwardTrackAnalysis::AnalyzeSingleMuon() {
214
215 if (fNTracksAnalyzedOfEvent>=fNTracksOfEvent) return kFALSE;
216
7e3dd1af 217 fMFTTrackWithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(fNTracksAnalyzedOfEvent);
218 fMFTTrackWithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(fNTracksAnalyzedOfEvent);
219
e806e863 220 fNTracksAnalyzedOfEvent++;
7e3dd1af 221
222 Bool_t passedCut = kFALSE;
223 if (fUseBransonForCut) passedCut = PassedCutSingleMuon(fMFTTrackWithBranson);
224 else passedCut = PassedCutSingleMuon(fMFTTrackWithoutBranson);
225 if (!passedCut) return kTRUE;
226
227 if (fUseBransonForKinematics) fMFTTrack = fMFTTrackWithBranson;
228 else fMFTTrack = fMFTTrackWithoutBranson;
cabfe25b 229 if (fMFTTrack->GetMatchTrigger() < fTriggerLevel) return kTRUE;
e806e863 230 fMCRefTrack = fMFTTrack->GetMCTrackRef();
231
cabfe25b 232 if (fMaxNWrongClustersMC<6) {
233 if (!fMCRefTrack) return kTRUE;
234 if (fMFTTrack->GetNWrongClustersMC()>fMaxNWrongClustersMC) return kTRUE;
235 }
d8c2cc3e 236
237 if (fMCRefTrack) {
238 fHistZOriginSingleMuonsMC -> Fill(-1.*fMCRefTrack->Vz());
239 fHistZROriginSingleMuonsMC -> Fill(-1.*fMCRefTrack->Vz(), TMath::Sqrt(fMCRefTrack->Vx()*fMCRefTrack->Vx()+fMCRefTrack->Vy()*fMCRefTrack->Vy()));
240 }
e806e863 241
e806e863 242 AliMUONTrackParam *param = fMFTTrack->GetTrackParamAtMFTCluster(0);
7e3dd1af 243 AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
e806e863 244
245 TLorentzVector pMu;
246 Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
247 Double_t energy = TMath::Sqrt(param->P()*param->P() + mMu*mMu);
248 pMu.SetPxPyPzE(param->Px(), param->Py(), param->Pz(), energy);
249
e806e863 250 TMatrixD cov(5,5);
251 cov = param->GetCovariances();
252
253 fHistErrorSingleMuonsX -> Fill(1.e4*TMath::Sqrt(cov(0,0)));
254 fHistErrorSingleMuonsY -> Fill(1.e4*TMath::Sqrt(cov(2,2)));
255
7e3dd1af 256 Double_t dX = fMFTTrack->GetOffsetX(fPrimaryVtxX, fPrimaryVtxZ);
257 Double_t dY = fMFTTrack->GetOffsetY(fPrimaryVtxY, fPrimaryVtxZ);
e806e863 258
7e3dd1af 259 Double_t offset = fMFTTrack->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
260 Double_t weightedOffset = fMFTTrack->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
e806e863 261
262 // AliDebug(2, Form("pdg code = %d\n", fMCRefTrack->GetPdgCode()));
263
264 fHistOffsetSingleMuonsX -> Fill(1.e4*dX);
265 fHistOffsetSingleMuonsY -> Fill(1.e4*dY);
cd2f51d2 266
cabfe25b 267 fHistSingleMuonsPtRapidity -> Fill(pMu.Rapidity(), pMu.Pt());
268 fHistOffsetSingleMuons -> Fill(1.e4*offset);
269 fHistWOffsetSingleMuons -> Fill(weightedOffset);
b5ab1ac4 270 Double_t chi2OverNdf = fMFTTrack->GetGlobalChi2()/Double_t(fMFTTrack->GetNMFTClusters()+fMFTTrack->GetNMUONClusters());
271 fHistSingleMuonsOffsetChi2 -> Fill(1.e4*offset, chi2OverNdf);
e806e863 272
273 fNTracksAnalyzed++;
d8c2cc3e 274 fNTracksAnalyzedOfEventAfterCut++;
e806e863 275
276 return kTRUE;
277
278}
279
280//====================================================================================================================================================
281
282Bool_t AliMuonForwardTrackAnalysis::AnalyzeMuonPair() {
283
284 if (fNPairsAnalyzedOfEvent>=fNPairsOfEvent) return kFALSE;
285
7e3dd1af 286 fMFTTrackPairWithBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithBranson->At(fNPairsAnalyzedOfEvent);
287 fMFTTrackPairWithoutBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithoutBranson->At(fNPairsAnalyzedOfEvent);
e806e863 288
7e3dd1af 289 fNPairsAnalyzedOfEvent++;
290
7e3dd1af 291 if (fUseBransonForKinematics) fMFTTrackPair = fMFTTrackPairWithBranson;
292 else fMFTTrackPair = fMFTTrackPairWithoutBranson;
e806e863 293
7e3dd1af 294 if ( fMFTTrackPair->GetTrack(0)->GetNWrongClustersMC()>fMaxNWrongClustersMC ||
295 fMFTTrackPair->GetTrack(1)->GetNWrongClustersMC()>fMaxNWrongClustersMC ) return kTRUE;
e806e863 296
cabfe25b 297 if (fOption==kResonanceOnly && !fMFTTrackPair->IsResonance()) return kTRUE;
789c2e43 298 if (fOption==kResonanceOnly && fMFTTrackPair->GetCharge() != 0) return kTRUE;
299
300 Double_t pca[3] = {0};
301 fMFTTrackPair -> GetPointOfClosestApproach(pca);
302 Double_t distancePrimaryVtxPCA = TMath::Sqrt(TMath::Power(fPrimaryVtxX-pca[0],2)+
303 TMath::Power(fPrimaryVtxY-pca[1],2)+
304 TMath::Power(fPrimaryVtxZ-pca[2],2));
305
d8c2cc3e 306 fMFTTrackPair -> SetKinem(fPrimaryVtxZ);
307
cabfe25b 308 Bool_t passedCut = kFALSE;
309 if (fUseBransonForCut) passedCut = PassedCutMuonPair(fMFTTrackPairWithBranson);
310 else passedCut = PassedCutMuonPair(fMFTTrackPairWithoutBranson);
311
312 if (!passedCut) return kTRUE;
cd2f51d2 313
cabfe25b 314 // --------------- Filling dimuon histograms --------------------------------
cd2f51d2 315
789c2e43 316 if (fMFTTrackPair->GetCharge() == 0) {
cabfe25b 317 if (fOption==kResonanceOnly) fHistMassMuonPairsMCVsPt -> Fill(fMFTTrackPair->GetMassMC(), fMFTTrackPair->GetPt());
318 fHistMassMuonPairsVsPt -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
319 fHistMassMuonPairsWithoutMFTVsPt -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
320 fHistWOffsetMuonPairsAtPrimaryVtxVsPt -> Fill(fMFTTrackPair->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
321 fHistWOffsetMuonPairsAtPCAVsPt -> Fill(fMFTTrackPair->GetWeightedOffsetAtPCA(), fMFTTrackPair->GetPt());
322 fHistDistancePrimaryVtxPCAVsPt -> Fill(distancePrimaryVtxPCA*1.e4, fMFTTrackPair->GetPt());
323 fHistPCAQualityVsPt -> Fill(fMFTTrackPair->GetPCAQuality(), fMFTTrackPair->GetPt());
324 if (fOption==kResonanceOnly) fHistPseudoProperDecayLengthVsPt->Fill(GetPseudoProperDecayLength(fMFTTrackPair, fTrueMass), fMFTTrackPair->GetPt());
325 if (fEvalDimuonVtxResolution) {
326 fHistDimuonVtxResolutionXVsPt->Fill(pca[0]*1.e4, fMFTTrackPair->GetPt());
327 fHistDimuonVtxResolutionYVsPt->Fill(pca[1]*1.e4, fMFTTrackPair->GetPt());
328 fHistDimuonVtxResolutionZVsPt->Fill(pca[2]*1.e4, fMFTTrackPair->GetPt());
329 }
330 fHistRapidityPtMuonPairs -> Fill(fMFTTrackPair->GetRapidityMC(), fMFTTrackPair->GetPt());
789c2e43 331 }
332 else if (fMFTTrackPair->GetCharge() == -2) {
333 fHistMassMuonPairsVsPtLSm -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
334 fHistMassMuonPairsWithoutMFTVsPtLSm -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
335 }
336 else if (fMFTTrackPair->GetCharge() == 2) {
337 fHistMassMuonPairsVsPtLSp -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
338 fHistMassMuonPairsWithoutMFTVsPtLSp -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
339 }
340
d8c2cc3e 341 AliDebug(1, Form("mass = %f MC = %f", fMFTTrackPair->GetMass(), fMFTTrackPair->GetMassMC()));
342
343 fNPairsAnalyzedOfEventAfterCut++;
e806e863 344
345 return kTRUE;
346
347}
348
349//====================================================================================================================================================
350
351void AliMuonForwardTrackAnalysis::BuildMuonPairs() {
352
789c2e43 353 Int_t nMuonPairs = 0;
354
355 for (Int_t iTrack=0; iTrack<fMuonForwardTracksWithBranson->GetEntriesFast(); iTrack++) {
e806e863 356 for (Int_t jTrack=0; jTrack<iTrack; jTrack++) {
789c2e43 357
7e3dd1af 358 AliMuonForwardTrack *track0_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(iTrack);
359 AliMuonForwardTrack *track1_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(jTrack);
789c2e43 360
7e3dd1af 361 AliMuonForwardTrack *track0_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(iTrack);
362 AliMuonForwardTrack *track1_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(jTrack);
e806e863 363
cabfe25b 364 if (track0_WithBranson->GetMatchTrigger()<fTriggerLevel || track1_WithBranson->GetMatchTrigger()<fTriggerLevel) continue;
e806e863 365
789c2e43 366 new ((*fMuonForwardTrackPairsWithBranson)[nMuonPairs]) AliMuonForwardTrackPair(track0_WithBranson, track1_WithBranson);
367 new ((*fMuonForwardTrackPairsWithoutBranson)[nMuonPairs]) AliMuonForwardTrackPair(track0_WithoutBranson, track1_WithoutBranson);
368
369 AliMuonForwardTrackPair *trackPairWithBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithBranson->At(nMuonPairs);
370 if (!(fOption==kResonanceOnly && !trackPairWithBranson->IsResonance())) nMuonPairs++;
371
e806e863 372 }
373 }
374
375}
376
377//====================================================================================================================================================
378
7e3dd1af 379Bool_t AliMuonForwardTrackAnalysis::PassedCutSingleMuon(AliMuonForwardTrack *track) {
380
381 AliMUONTrackParam *param = track->GetTrackParamAtMFTCluster(0);
382 AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
383
384 if (track->Pt()<fPtMinSingleMuons) return kFALSE;
385
cabfe25b 386 Double_t offset = 1.e4*track->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
387 Double_t chi2OverNdf = track->GetGlobalChi2() / Double_t(track->GetNMFTClusters()+track->GetNMUONClusters());
388 offset /= fMaxOffsetSingleMuons;
389 chi2OverNdf /= fMaxChi2SingleMuons;
390
391 if (fCorrelateCutOnOffsetChi2) {
392 if (TMath::Sqrt(offset*offset + chi2OverNdf*chi2OverNdf) > 1) return kFALSE;
393 }
394 else {
395 if (offset>1 || chi2OverNdf>1) return kFALSE;
7e3dd1af 396 }
397
398 return kTRUE;
399
400}
401
402//====================================================================================================================================================
403
404Bool_t AliMuonForwardTrackAnalysis::PassedCutMuonPair(AliMuonForwardTrackPair *pair) {
405
cabfe25b 406 if (!PassedCutSingleMuon(pair->GetTrack(0)) || !PassedCutSingleMuon(pair->GetTrack(1))) return kFALSE;
407
408 if (pair->GetMass()>fMassMax || pair->GetMass()<fMassMin) return kFALSE;
7e3dd1af 409
cabfe25b 410 if (pair->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ) > fMaxWOffsetMuonPairsAtPrimaryVtx) return kFALSE;
411 if (pair->GetWeightedOffsetAtPCA() > fMaxWOffsetMuonPairsAtPCA) return kFALSE;
412
413 Double_t pca[3] = {0};
414 pair -> GetPointOfClosestApproach(pca);
415 Double_t distancePrimaryVtxPCA = TMath::Sqrt(TMath::Power(fPrimaryVtxX-pca[0],2)+
416 TMath::Power(fPrimaryVtxY-pca[1],2)+
417 TMath::Power(fPrimaryVtxZ-pca[2],2));
418
419 if (distancePrimaryVtxPCA > fMaxDistancePrimaryVtxPCA) return kFALSE;
420 if (pair->GetPCAQuality() < fMinPCAQuality) return kFALSE;
421
422 return kTRUE;
423
424}
425
426//====================================================================================================================================================
427
428Double_t AliMuonForwardTrackAnalysis::GetPseudoProperDecayLength(AliMuonForwardTrackPair *pair, Double_t trueMass) {
429
430 Double_t pca[3] = {0};
431 pair -> GetPointOfClosestApproach(pca);
432
433 TVector2 vecVertexPCA(pca[0]-fPrimaryVtxX, pca[1]-fPrimaryVtxY);
434 TVector2 dimuonPt(pair->GetPx(), pair->GetPy());
435 dimuonPt.Unit();
436
437 Double_t l_xy = vecVertexPCA*dimuonPt;
438 Double_t pseudoProperDecayLength = (l_xy * trueMass / pair->GetPt()) * 1.e4 ; // in micron
439
440 return pseudoProperDecayLength;
441
7e3dd1af 442}
443
444//====================================================================================================================================================
445
e806e863 446void AliMuonForwardTrackAnalysis::Terminate(Char_t *outputFileName) {
447
e806e863 448 TFile *fileOut = new TFile(Form("%s/%s",fOutputDir.Data(),outputFileName), "recreate");
449
450 printf("Writing output objects to file %s\n", fileOut->GetName());
451
cabfe25b 452 fHistOffsetSingleMuonsX -> Write();
453 fHistOffsetSingleMuonsY -> Write();
454 fHistErrorSingleMuonsX -> Write();
455 fHistErrorSingleMuonsY -> Write();
456 fHistOffsetSingleMuons -> Write();
457 fHistWOffsetSingleMuons -> Write();
e806e863 458
cabfe25b 459 fHistSingleMuonsPtRapidity -> Write();
b5ab1ac4 460 fHistSingleMuonsOffsetChi2 -> Write();
d8c2cc3e 461 fHistZOriginSingleMuonsMC -> Write();
462 fHistZROriginSingleMuonsMC -> Write();
463
cabfe25b 464 fHistMassMuonPairsVsPt -> Write();
465 fHistMassMuonPairsWithoutMFTVsPt -> Write();
466 fHistMassMuonPairsVsPtLSp -> Write();
467 fHistMassMuonPairsWithoutMFTVsPtLSp -> Write();
468 fHistMassMuonPairsVsPtLSm -> Write();
469 fHistMassMuonPairsWithoutMFTVsPtLSm -> Write();
470 fHistMassMuonPairsMCVsPt -> Write();
471 fHistWOffsetMuonPairsAtPrimaryVtxVsPt -> Write();
472 fHistWOffsetMuonPairsAtPCAVsPt -> Write();
473 fHistDistancePrimaryVtxPCAVsPt -> Write();
474 fHistPCAQualityVsPt -> Write();
475 if (fOption==kResonanceOnly) fHistPseudoProperDecayLengthVsPt->Write();
476 if (fEvalDimuonVtxResolution) {
477 fHistDimuonVtxResolutionXVsPt -> Write();
478 fHistDimuonVtxResolutionYVsPt -> Write();
479 fHistDimuonVtxResolutionZVsPt -> Write();
cd2f51d2 480 }
d8c2cc3e 481 fHistRapidityPtMuonPairs -> Write();
e806e863 482
483 fileOut -> Close();
484
485}
486
487//====================================================================================================================================================
488
489void AliMuonForwardTrackAnalysis::BookHistos() {
490
cabfe25b 491 // -------------------- single muons
492
e806e863 493 fHistOffsetSingleMuonsX = new TH1D("fHistOffsetSingleMuonsX", "Offset for single muons along X", 200, -1000, 1000);
494 fHistOffsetSingleMuonsY = new TH1D("fHistOffsetSingleMuonsY", "Offset for single muons along Y", 200, -1000, 1000);
495 fHistErrorSingleMuonsX = new TH1D("fHistErrorSingleMuonsX", "Coordinate Error for single muons along X", 200, 0, 1000);
496 fHistErrorSingleMuonsY = new TH1D("fHistErrorSingleMuonsY", "Coordinate Error for single muons along Y", 200, 0, 1000);
b5ab1ac4 497 fHistOffsetSingleMuons = new TH1D("fHistOffsetSingleMuons", "Offset for single muons", 200, 0, 2000);
e806e863 498 fHistWOffsetSingleMuons = new TH1D("fHistWOffsetSingleMuons", "Weighted Offset for single muons", 300, 0, 15);
499
cabfe25b 500 fHistSingleMuonsPtRapidity = new TH2D("fHistSingleMuonsPtRapidity", "Phase Space for single muons", 100, -4.5, -2., 100, 0., 10.);
789c2e43 501 fHistSingleMuonsOffsetChi2 = new TH2D("fHistSingleMuonsOffsetChi2", "Offset vs #chi^{2}/ndf for single muons", 400, 0, 4000, 200, 0, 10);
d8c2cc3e 502 fHistZOriginSingleMuonsMC = new TH1D("fHistZOriginSingleMuonsMC", "Z origin for single muons (from MC)", 1000, 0., 500.);
503 fHistZROriginSingleMuonsMC = new TH2D("fHistZROriginSingleMuonsMC", "Z-R origin for single muons (from MC)", 1000, 0., 500., 1000, 0., 100.);
e806e863 504
505 fHistOffsetSingleMuonsX -> SetXTitle("Offset(X) [#mum]");
506 fHistOffsetSingleMuonsY -> SetXTitle("Offset(Y) [#mum]");
507 fHistErrorSingleMuonsX -> SetXTitle("Err. on track position at z_{vtx} (X) [#mum]");
508 fHistErrorSingleMuonsY -> SetXTitle("Err. on track position at z_{vtx} (Y) [#mum]");
509 fHistOffsetSingleMuons -> SetXTitle("Offset [#mum]");
510 fHistWOffsetSingleMuons -> SetXTitle("Weighted Offset");
511
cabfe25b 512 fHistSingleMuonsPtRapidity -> SetXTitle("y^{#mu}");
513 fHistSingleMuonsPtRapidity -> SetYTitle("p_{T}^{#mu} [GeV/c]");
b5ab1ac4 514 fHistSingleMuonsOffsetChi2 -> SetXTitle("Offset [#mum]");
515 fHistSingleMuonsOffsetChi2 -> SetYTitle("#chi^{2}/ndf");
d8c2cc3e 516
517 fHistZOriginSingleMuonsMC -> SetXTitle("Z [cm]");
518 fHistZROriginSingleMuonsMC -> SetXTitle("Z [cm]");
519 fHistZROriginSingleMuonsMC -> SetXTitle("R [cm]");
e806e863 520
521 fHistOffsetSingleMuonsX -> Sumw2();
522 fHistOffsetSingleMuonsY -> Sumw2();
523 fHistErrorSingleMuonsX -> Sumw2();
524 fHistErrorSingleMuonsY -> Sumw2();
525 fHistOffsetSingleMuons -> Sumw2();
526 fHistWOffsetSingleMuons -> Sumw2();
527
cabfe25b 528 fHistSingleMuonsPtRapidity -> Sumw2();
529 fHistSingleMuonsOffsetChi2 -> Sumw2();
d8c2cc3e 530 fHistZOriginSingleMuonsMC -> Sumw2();
531 fHistZROriginSingleMuonsMC -> Sumw2();
532
cabfe25b 533 // -------------------- muon pairs
534
535 Int_t nBinsPt = 20; Double_t ptMin = 0., ptMax = 10.; // dimuon pt
536
537 fHistMassMuonPairsVsPt = new TH2D("fHistMassMuonPairsVsPt", "Dimuon Mass (MUON+MFT) vs p_{T}^{#mu#mu}", 1000, 0, 10, nBinsPt, ptMin, ptMax);
538 fHistMassMuonPairsWithoutMFTVsPt = new TH2D("fHistMassMuonPairsWithoutMFTVsPt", "Dimuon Mass (MUON only) vs p_{T}^{#mu#mu}", 1000, 0, 10, nBinsPt, ptMin, ptMax);
539 fHistMassMuonPairsVsPtLSp = new TH2D("fHistMassMuonPairsVsPtLSp", "Dimuon Mass (MUON+MFT) vs p_{T}^{#mu#mu} Like-Sign ++", 1000, 0, 10, nBinsPt, ptMin, ptMax);
540 fHistMassMuonPairsWithoutMFTVsPtLSp = new TH2D("fHistMassMuonPairsWithoutMFTVsPtLSp", "Dimuon Mass (MUON only) vs p_{T}^{#mu#mu} Like-Sign ++", 1000, 0, 10, nBinsPt, ptMin, ptMax);
541 fHistMassMuonPairsVsPtLSm = new TH2D("fHistMassMuonPairsVsPtLSm", "Dimuon Mass (MUON+MFT) vs p_{T}^{#mu#mu} Like-Sign --", 1000, 0, 10, nBinsPt, ptMin, ptMax);
542 fHistMassMuonPairsWithoutMFTVsPtLSm = new TH2D("fHistMassMuonPairsWithoutMFTVsPtLSm", "Dimuon Mass (MUON only) vs p_{T}^{#mu#mu} Like-Sign --", 1000, 0, 10, nBinsPt, ptMin, ptMax);
543 fHistMassMuonPairsMCVsPt = new TH2D("fHistMassMuonPairsMCVsPt", "Dimuon Mass (MC) vs p_{T}^{#mu#mu}", 1000, 0, 10, nBinsPt, ptMin, ptMax);
544 fHistWOffsetMuonPairsAtPrimaryVtxVsPt = new TH2D("fHistWOffsetMuonPairsAtPrimaryVtxVsPt", "Weighted Offset for Muon Pairs at Primary Vertex vs p_{T}^{#mu#mu}", 300, 0, 60, nBinsPt, ptMin, ptMax);
545 fHistWOffsetMuonPairsAtPCAVsPt = new TH2D("fHistWOffsetMuonPairsAtPCAVsPt", "Weighted Offset for Muon Pairs at PCA vs p_{T}^{#mu#mu}", 300, 0, 60, nBinsPt, ptMin, ptMax);
546 fHistDistancePrimaryVtxPCAVsPt = new TH2D("fHistDistancePrimaryVtxPCA_%d", "Distance between PCA and primary vertex vs p_{T}^{#mu#mu}", 1000, 0, 50000, nBinsPt, ptMin, ptMax);
547 fHistPCAQualityVsPt = new TH2D("fHistPCAQualityVsPt", "PCA Quality vs p_{T}^{#mu#mu}", 200, 0, 1, nBinsPt, ptMin, ptMax);
548 fHistPseudoProperDecayLengthVsPt = new TH2D("fHistPseudoProperDecayLengthVsPt", "Pseudo proper decay length vs p_{T}^{#mu#mu}", 1000, -5000, 5000, nBinsPt, ptMin, ptMax);
549 fHistDimuonVtxResolutionXVsPt = new TH2D("fHistDimuonVtxResolutionXVsPt", "Dimuon Vtx offset along X w.r.t. true Vtx vs p_{T}^{#mu#mu}", 2000, -1000, 1000, nBinsPt, ptMin, ptMax);
550 fHistDimuonVtxResolutionYVsPt = new TH2D("fHistDimuonVtxResolutionYVsPt", "Dimuon Vtx offset along Y w.r.t. true Vtx vs p_{T}^{#mu#mu}", 2000, -1000, 1000, nBinsPt, ptMin, ptMax);
551 fHistDimuonVtxResolutionZVsPt = new TH2D("fHistDimuonVtxResolutionZVsPt", "Dimuon Vtx offset along Z w.r.t. true Vtx vs p_{T}^{#mu#mu}", 2000, -1000, 1000, nBinsPt, ptMin, ptMax);
cd2f51d2 552
cabfe25b 553 fHistMassMuonPairsVsPt -> SetXTitle("Mass [GeV/c^{2}]");
554 fHistMassMuonPairsWithoutMFTVsPt -> SetXTitle("Mass [GeV/c^{2}]");
555 fHistMassMuonPairsVsPtLSp -> SetXTitle("Mass [GeV/c^{2}]");
556 fHistMassMuonPairsWithoutMFTVsPtLSp -> SetXTitle("Mass [GeV/c^{2}]");
557 fHistMassMuonPairsVsPtLSm -> SetXTitle("Mass [GeV/c^{2}]");
558 fHistMassMuonPairsWithoutMFTVsPtLSm -> SetXTitle("Mass [GeV/c^{2}]");
559 fHistMassMuonPairsMCVsPt -> SetXTitle("Mass [GeV/c^{2}]");
560 fHistWOffsetMuonPairsAtPrimaryVtxVsPt -> SetXTitle("Weighted Offset at Primary Vtx");
561 fHistWOffsetMuonPairsAtPCAVsPt -> SetXTitle("Weighted Offset at PCA");
562 fHistDistancePrimaryVtxPCAVsPt -> SetXTitle("PCA - Primary Vtx [#mum]");
563 fHistPCAQualityVsPt -> SetXTitle("PCA Quality");
564 fHistPseudoProperDecayLengthVsPt -> SetXTitle("l_{xy} [#mum]");
565 fHistDimuonVtxResolutionXVsPt -> SetXTitle("Offset [#mum]");
566 fHistDimuonVtxResolutionYVsPt -> SetXTitle("Offset [#mum]");
567 fHistDimuonVtxResolutionZVsPt -> SetXTitle("Offset [#mum]");
568
569 fHistMassMuonPairsVsPt -> SetYTitle("p_{T} [GeV/c]");
570 fHistMassMuonPairsWithoutMFTVsPt -> SetYTitle("p_{T} [GeV/c]");
571 fHistMassMuonPairsVsPtLSp -> SetYTitle("p_{T} [GeV/c]");
572 fHistMassMuonPairsWithoutMFTVsPtLSp -> SetYTitle("p_{T} [GeV/c]");
573 fHistMassMuonPairsVsPtLSm -> SetYTitle("p_{T} [GeV/c]");
574 fHistMassMuonPairsWithoutMFTVsPtLSm -> SetYTitle("p_{T} [GeV/c]");
575 fHistMassMuonPairsMCVsPt -> SetYTitle("p_{T} [GeV/c]");
576 fHistWOffsetMuonPairsAtPrimaryVtxVsPt -> SetYTitle("p_{T} [GeV/c]");
577 fHistWOffsetMuonPairsAtPCAVsPt -> SetYTitle("p_{T} [GeV/c]");
578 fHistDistancePrimaryVtxPCAVsPt -> SetYTitle("p_{T} [GeV/c]");
579 fHistPCAQualityVsPt -> SetYTitle("p_{T} [GeV/c]");
580 fHistPseudoProperDecayLengthVsPt -> SetYTitle("p_{T} [GeV/c]");
581 fHistDimuonVtxResolutionXVsPt -> SetYTitle("p_{T} [GeV/c]");
582 fHistDimuonVtxResolutionYVsPt -> SetYTitle("p_{T} [GeV/c]");
583 fHistDimuonVtxResolutionZVsPt -> SetYTitle("p_{T} [GeV/c]");
584
585 fHistMassMuonPairsVsPt -> Sumw2();
586 fHistMassMuonPairsWithoutMFTVsPt -> Sumw2();
587 fHistMassMuonPairsVsPtLSp -> Sumw2();
588 fHistMassMuonPairsWithoutMFTVsPtLSp -> Sumw2();
589 fHistMassMuonPairsVsPtLSm -> Sumw2();
590 fHistMassMuonPairsWithoutMFTVsPtLSm -> Sumw2();
591 fHistMassMuonPairsMCVsPt -> Sumw2();
592 fHistWOffsetMuonPairsAtPrimaryVtxVsPt -> Sumw2();
593 fHistWOffsetMuonPairsAtPCAVsPt -> Sumw2();
594 fHistDistancePrimaryVtxPCAVsPt -> Sumw2();
595 fHistPCAQualityVsPt -> Sumw2();
596 fHistPseudoProperDecayLengthVsPt -> Sumw2();
597 fHistDimuonVtxResolutionXVsPt -> Sumw2();
598 fHistDimuonVtxResolutionYVsPt -> Sumw2();
599 fHistDimuonVtxResolutionZVsPt -> Sumw2();
600
601 fHistRapidityPtMuonPairs = new TH2D("fHistRapidityPtMuonPairs", "Dimuon Phase Space (rec)", 20, -4.5, -2., 20, 0., 10.);
602 fHistRapidityPtMuonPairs -> SetXTitle("Rapidity");
603 fHistRapidityPtMuonPairs -> SetYTitle("p_{T} [GeV/c]");
604 fHistRapidityPtMuonPairs -> Sumw2();
e806e863 605
606}
607
608//====================================================================================================================================================
609