]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MFT/AliMuonForwardTrackAnalysis.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
1) Mods in AddTaskPHOSPi0Flow to comply with LEGO train
[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),
66 fHistZROriginSingleMuonsMC(0x0),
67 fHistSingleMuonsPtRapidityMC(0x0),
b5ab1ac4 68 fHistSingleMuonsOffsetChi2(0x0),
d8c2cc3e 69 fHistRapidityPtMuonPairs(0x0),
789c2e43 70 fHistMassMuonPairsVsPt(0),
71 fHistMassMuonPairsWithoutMFTVsPt(0),
72 fHistMassMuonPairsVsPtLSp(0),
73 fHistMassMuonPairsWithoutMFTVsPtLSp(0),
74 fHistMassMuonPairsVsPtLSm(0),
75 fHistMassMuonPairsWithoutMFTVsPtLSm(0),
76 fEvalDimuonVtxResolution(kFALSE),
cd2f51d2 77 fNMassBins(10),
78 fNPtDimuBins(1000),
e806e863 79 fMassMin(0),
80 fMassMax(10),
cd2f51d2 81 fPtDimuMin(0),
82 fPtDimuMax(5),
83 fPtAxisDimuons(0),
e806e863 84 fSingleMuonAnalysis(1),
85 fMuonPairAnalysis(1),
86 fMatchTrigger(0),
87 fOption(0),
b5ab1ac4 88 fXVertResMC(50.e-4),
89 fYVertResMC(50.e-4),
90 fZVertResMC(50.e-4),
7e3dd1af 91 fPrimaryVtxX(0.),
92 fPrimaryVtxY(0.),
93 fPrimaryVtxZ(0.),
e806e863 94 fMaxNWrongClustersMC(999),
7e3dd1af 95 fPtMinSingleMuons(0),
96 fUseBransonForCut(kFALSE),
97 fUseBransonForKinematics(kFALSE),
98 fCutOnOffsetChi2(kFALSE),
99 fCenterOffset(0.),
53b30119 100 fCenterChi2(0.),
7e3dd1af 101 fScaleOffset(250.),
789c2e43 102 fScaleChi2(1.5),
7e3dd1af 103 fRadiusCut(1.)
e806e863 104{
105
106 // default constructor
107
cd2f51d2 108 for (Int_t iPtBin=0; iPtBin<fNMaxPtBinsDimuons+1; iPtBin++) {
789c2e43 109 fHistMassMuonPairs[iPtBin] = NULL;
110 fHistMassMuonPairsWithoutMFT[iPtBin] = NULL;
111 fHistMassMuonPairsMC[iPtBin] = NULL;
112 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] = NULL;
113 fHistWOffsetMuonPairsAtPCA[iPtBin] = NULL;
114 fHistDistancePrimaryVtxPCA[iPtBin] = NULL;
115 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] = NULL;
116 fHistDimuonVtxResolutionX[iPtBin] = NULL;
117 fHistDimuonVtxResolutionY[iPtBin] = NULL;
118 fHistDimuonVtxResolutionZ[iPtBin] = NULL;
e806e863 119 }
120
121}
122
123//====================================================================================================================================================
124
125Bool_t AliMuonForwardTrackAnalysis::Init(Char_t *inputFileName) {
126
cd2f51d2 127 fPtAxisDimuons = new TAxis(fNPtDimuBins, fPtDimuMin, fPtDimuMax);
128
e806e863 129 BookHistos();
130
131 TFile *inputFile = new TFile(Form("%s/%s",fInputDir.Data(),inputFileName));
132 if (!inputFile || !inputFile->IsOpen()) {
133 AliError(Form("Error opening file %s", inputFileName));
134 return kFALSE;
135 }
7e3dd1af 136 fInputTreeWithBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithBranson");
137 if (!fInputTreeWithBranson) {
138 AliError("Error reading input tree");
139 return kFALSE;
140 }
141 fInputTreeWithoutBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithoutBranson");
142 if (!fInputTreeWithoutBranson) {
e806e863 143 AliError("Error reading input tree");
144 return kFALSE;
145 }
146
789c2e43 147 if (fFirstEvent<0 || fLastEvent<0 || fFirstEvent>fLastEvent || fFirstEvent>=fInputTreeWithBranson->GetEntriesFast()) {
e806e863 148 fFirstEvent = 0;
789c2e43 149 fLastEvent = fInputTreeWithBranson->GetEntriesFast()-1;
e806e863 150 }
151 else {
789c2e43 152 fLastEvent = TMath::Min(fLastEvent, Int_t(fInputTreeWithBranson->GetEntriesFast()-1));
e806e863 153 }
154
155 AliInfo(Form("Analysing events %d to %d", fFirstEvent, fLastEvent));
156
789c2e43 157 fMuonForwardTracksWithBranson = new TClonesArray("AliMuonForwardTrack",30);
7e3dd1af 158 fInputTreeWithBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithBranson);
159
789c2e43 160 fMuonForwardTracksWithoutBranson = new TClonesArray("AliMuonForwardTrack",30);
7e3dd1af 161 fInputTreeWithoutBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithoutBranson);
e806e863 162
163 TGeoManager::Import(Form("%s/geometry.root",fInputDir.Data()));
164
165 AliMUONTrackExtrap::SetField();
166
789c2e43 167 fMuonForwardTrackPairsWithBranson = new TClonesArray("AliMuonForwardTrackPair",10);
168 fMuonForwardTrackPairsWithoutBranson = new TClonesArray("AliMuonForwardTrackPair",10);
169 fMuonForwardTrackPairsWithBranson->SetOwner(kTRUE);
170 fMuonForwardTrackPairsWithoutBranson->SetOwner(kTRUE);
171
e806e863 172
173 return kTRUE;
174
175}
176
177//====================================================================================================================================================
178
179Bool_t AliMuonForwardTrackAnalysis::LoadNextEvent() {
180
181 if (fEv>fLastEvent) return kFALSE;
182 if (fEv<fFirstEvent) { fEv++; return kTRUE; }
789c2e43 183 fMuonForwardTracksWithBranson -> Clear("");
184 fMuonForwardTracksWithoutBranson -> Clear("");
7e3dd1af 185 fInputTreeWithBranson->GetEvent(fEv);
186 fInputTreeWithoutBranson->GetEvent(fEv);
789c2e43 187 AliDebug(2,Form("**** analyzing event # %4d (%3d tracks) ****", fEv, fMuonForwardTracksWithBranson->GetEntriesFast()));
7e3dd1af 188
189 fPrimaryVtxX = gRandom->Gaus(0., fXVertResMC);
190 fPrimaryVtxY = gRandom->Gaus(0., fYVertResMC);
191 fPrimaryVtxZ = gRandom->Gaus(0., fZVertResMC);
e806e863 192
193 if (fSingleMuonAnalysis) {
194 fNTracksAnalyzedOfEvent = 0;
d8c2cc3e 195 fNTracksAnalyzedOfEventAfterCut = 0;
789c2e43 196 fNTracksOfEvent = fMuonForwardTracksWithBranson->GetEntriesFast();
e806e863 197 while (AnalyzeSingleMuon()) continue;
198 }
199
200 if (fMuonPairAnalysis) {
7e3dd1af 201 if (fMuonForwardTrackPairsWithBranson) {
789c2e43 202 fMuonForwardTrackPairsWithBranson->Clear("C");
203 fMuonForwardTrackPairsWithoutBranson->Clear("C");
e806e863 204 }
205 BuildMuonPairs();
206 fNPairsAnalyzedOfEvent = 0;
d8c2cc3e 207 fNPairsAnalyzedOfEventAfterCut = 0;
789c2e43 208 // fNPairsOfEvent = fMuonForwardTrackPairsWithBranson->GetEntriesFast();
e806e863 209 while (AnalyzeMuonPair()) continue;
210 }
211
789c2e43 212 AliDebug(2,Form("**** analyzed event # %4d (%3d tracks and %3d pairs analyzed) ****", fEv, fNTracksAnalyzedOfEventAfterCut, fNPairsAnalyzedOfEventAfterCut));
d8c2cc3e 213
e806e863 214 fEv++;
215
216 return kTRUE;
217
218}
219
220//====================================================================================================================================================
221
222Bool_t AliMuonForwardTrackAnalysis::AnalyzeSingleMuon() {
223
224 if (fNTracksAnalyzedOfEvent>=fNTracksOfEvent) return kFALSE;
225
7e3dd1af 226 fMFTTrackWithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(fNTracksAnalyzedOfEvent);
227 fMFTTrackWithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(fNTracksAnalyzedOfEvent);
228
e806e863 229 fNTracksAnalyzedOfEvent++;
7e3dd1af 230
231 Bool_t passedCut = kFALSE;
232 if (fUseBransonForCut) passedCut = PassedCutSingleMuon(fMFTTrackWithBranson);
233 else passedCut = PassedCutSingleMuon(fMFTTrackWithoutBranson);
234 if (!passedCut) return kTRUE;
235
236 if (fUseBransonForKinematics) fMFTTrack = fMFTTrackWithBranson;
237 else fMFTTrack = fMFTTrackWithoutBranson;
e806e863 238 if (fMatchTrigger && !fMFTTrack->GetMatchTrigger()) return kTRUE;
239 fMCRefTrack = fMFTTrack->GetMCTrackRef();
240
d8c2cc3e 241 if (fOption!=kPionsKaons && !fMCRefTrack) return kTRUE;
242 if (fOption!=kPionsKaons && fMFTTrack->GetNWrongClustersMC()>fMaxNWrongClustersMC) return kTRUE;
243
244 if (fMCRefTrack) {
245 fHistZOriginSingleMuonsMC -> Fill(-1.*fMCRefTrack->Vz());
246 fHistZROriginSingleMuonsMC -> Fill(-1.*fMCRefTrack->Vz(), TMath::Sqrt(fMCRefTrack->Vx()*fMCRefTrack->Vx()+fMCRefTrack->Vy()*fMCRefTrack->Vy()));
247 }
e806e863 248
e806e863 249 AliMUONTrackParam *param = fMFTTrack->GetTrackParamAtMFTCluster(0);
7e3dd1af 250 AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
e806e863 251
252 TLorentzVector pMu;
253 Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
254 Double_t energy = TMath::Sqrt(param->P()*param->P() + mMu*mMu);
255 pMu.SetPxPyPzE(param->Px(), param->Py(), param->Pz(), energy);
256
e806e863 257 TMatrixD cov(5,5);
258 cov = param->GetCovariances();
259
260 fHistErrorSingleMuonsX -> Fill(1.e4*TMath::Sqrt(cov(0,0)));
261 fHistErrorSingleMuonsY -> Fill(1.e4*TMath::Sqrt(cov(2,2)));
262
7e3dd1af 263 Double_t dX = fMFTTrack->GetOffsetX(fPrimaryVtxX, fPrimaryVtxZ);
264 Double_t dY = fMFTTrack->GetOffsetY(fPrimaryVtxY, fPrimaryVtxZ);
e806e863 265
7e3dd1af 266 Double_t offset = fMFTTrack->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
267 Double_t weightedOffset = fMFTTrack->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
e806e863 268
269 // AliDebug(2, Form("pdg code = %d\n", fMCRefTrack->GetPdgCode()));
270
271 fHistOffsetSingleMuonsX -> Fill(1.e4*dX);
272 fHistOffsetSingleMuonsY -> Fill(1.e4*dY);
cd2f51d2 273
d8c2cc3e 274 if (fOption!=kPionsKaons) fHistSingleMuonsPtRapidityMC-> Fill(fMCRefTrack->Y(), fMCRefTrack->Pt());
b5ab1ac4 275 fHistOffsetSingleMuons -> Fill(1.e4*offset);
276 fHistWOffsetSingleMuons -> Fill(weightedOffset);
277 Double_t chi2OverNdf = fMFTTrack->GetGlobalChi2()/Double_t(fMFTTrack->GetNMFTClusters()+fMFTTrack->GetNMUONClusters());
278 fHistSingleMuonsOffsetChi2 -> Fill(1.e4*offset, chi2OverNdf);
e806e863 279
280 fNTracksAnalyzed++;
d8c2cc3e 281 fNTracksAnalyzedOfEventAfterCut++;
e806e863 282
283 return kTRUE;
284
285}
286
287//====================================================================================================================================================
288
289Bool_t AliMuonForwardTrackAnalysis::AnalyzeMuonPair() {
290
291 if (fNPairsAnalyzedOfEvent>=fNPairsOfEvent) return kFALSE;
292
7e3dd1af 293 fMFTTrackPairWithBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithBranson->At(fNPairsAnalyzedOfEvent);
294 fMFTTrackPairWithoutBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithoutBranson->At(fNPairsAnalyzedOfEvent);
e806e863 295
7e3dd1af 296 fNPairsAnalyzedOfEvent++;
297
298 Bool_t passedCut = kFALSE;
299 if (fUseBransonForCut) passedCut = PassedCutMuonPair(fMFTTrackPairWithBranson);
300 else passedCut = PassedCutMuonPair(fMFTTrackPairWithoutBranson);
301
302 if (!passedCut) return kTRUE;
303
304 if (fUseBransonForKinematics) fMFTTrackPair = fMFTTrackPairWithBranson;
305 else fMFTTrackPair = fMFTTrackPairWithoutBranson;
e806e863 306
7e3dd1af 307 if ( fMFTTrackPair->GetTrack(0)->GetNWrongClustersMC()>fMaxNWrongClustersMC ||
308 fMFTTrackPair->GetTrack(1)->GetNWrongClustersMC()>fMaxNWrongClustersMC ) return kTRUE;
e806e863 309
7e3dd1af 310 if (fOption==kResonanceOnly && !fMFTTrackPair->IsResonance()) return kTRUE;
311
789c2e43 312 if (fOption==kResonanceOnly && fMFTTrackPair->GetCharge() != 0) return kTRUE;
313
314 Double_t pca[3] = {0};
315 fMFTTrackPair -> GetPointOfClosestApproach(pca);
316 Double_t distancePrimaryVtxPCA = TMath::Sqrt(TMath::Power(fPrimaryVtxX-pca[0],2)+
317 TMath::Power(fPrimaryVtxY-pca[1],2)+
318 TMath::Power(fPrimaryVtxZ-pca[2],2));
319
d8c2cc3e 320 fMFTTrackPair -> SetKinem(fPrimaryVtxZ);
321
789c2e43 322 Int_t ptBin = fPtAxisDimuons->FindBin(fMFTTrackPair->GetPt());
cd2f51d2 323
324 if (1<=ptBin && ptBin<=fNPtDimuBins) {
d8c2cc3e 325 fHistMassMuonPairs[ptBin] -> Fill(fMFTTrackPair->GetMass());
cd2f51d2 326 fHistMassMuonPairsWithoutMFT[ptBin] -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
789c2e43 327 if (fOption!=kPionsKaons) fHistMassMuonPairsMC[ptBin] -> Fill(fMFTTrackPair->GetMassMC());
328 fHistWOffsetMuonPairsAtPrimaryVtx[ptBin] -> Fill(fMFTTrackPair->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
329 fHistWOffsetMuonPairsAtPCA[ptBin] -> Fill(fMFTTrackPair->GetWeightedOffset(pca[0], pca[1], pca[2]));
330 fHistDistancePrimaryVtxPCA[ptBin] -> Fill(distancePrimaryVtxPCA*1.e4);
331 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[ptBin] -> Fill(fMFTTrackPair->GetWeightedOffset(pca[0], pca[1], pca[2]), distancePrimaryVtxPCA*1.e4);
332 if (fEvalDimuonVtxResolution) {
333 fHistDimuonVtxResolutionX[ptBin]->Fill(pca[0]*1.e4);
334 fHistDimuonVtxResolutionY[ptBin]->Fill(pca[1]*1.e4);
335 fHistDimuonVtxResolutionZ[ptBin]->Fill(pca[2]*1.e4);
336 }
cd2f51d2 337 }
d8c2cc3e 338 fHistMassMuonPairs[0] -> Fill(fMFTTrackPair->GetMass());
cd2f51d2 339 fHistMassMuonPairsWithoutMFT[0] -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
789c2e43 340 if (fOption!=kPionsKaons) fHistMassMuonPairsMC[0] -> Fill(fMFTTrackPair->GetMassMC());
341 fHistWOffsetMuonPairsAtPrimaryVtx[0] -> Fill(fMFTTrackPair->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
342 fHistWOffsetMuonPairsAtPCA[0] -> Fill(fMFTTrackPair->GetWeightedOffset(pca[0], pca[1], pca[2]));
343 fHistDistancePrimaryVtxPCA[0] -> Fill(distancePrimaryVtxPCA*1.e4);
344 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[0] -> Fill(fMFTTrackPair->GetWeightedOffset(pca[0], pca[1], pca[2]), distancePrimaryVtxPCA*1.e4);
345 if (fEvalDimuonVtxResolution) {
346 fHistDimuonVtxResolutionX[0]->Fill(pca[0]*1.e4);
347 fHistDimuonVtxResolutionY[0]->Fill(pca[1]*1.e4);
348 fHistDimuonVtxResolutionZ[0]->Fill(pca[2]*1.e4);
349 }
cd2f51d2 350
789c2e43 351 if (fMFTTrackPair->GetCharge() == 0) {
352 fHistMassMuonPairsVsPt -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
353 fHistMassMuonPairsWithoutMFTVsPt -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
354 }
355 else if (fMFTTrackPair->GetCharge() == -2) {
356 fHistMassMuonPairsVsPtLSm -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
357 fHistMassMuonPairsWithoutMFTVsPtLSm -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
358 }
359 else if (fMFTTrackPair->GetCharge() == 2) {
360 fHistMassMuonPairsVsPtLSp -> Fill(fMFTTrackPair->GetMass(), fMFTTrackPair->GetPt());
361 fHistMassMuonPairsWithoutMFTVsPtLSp -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ), fMFTTrackPair->GetPt());
362 }
363
364 if (fOption!=kPionsKaons) fHistRapidityPtMuonPairs -> Fill(fMFTTrackPair->GetRapidityMC(), fMFTTrackPair->GetPt());
365 else fHistRapidityPtMuonPairs -> Fill(fMFTTrackPair->GetRapidity(), fMFTTrackPair->GetPt());
e806e863 366
d8c2cc3e 367 AliDebug(1, Form("mass = %f MC = %f", fMFTTrackPair->GetMass(), fMFTTrackPair->GetMassMC()));
368
369 fNPairsAnalyzedOfEventAfterCut++;
e806e863 370
371 return kTRUE;
372
373}
374
375//====================================================================================================================================================
376
377void AliMuonForwardTrackAnalysis::BuildMuonPairs() {
378
789c2e43 379 Int_t nMuonPairs = 0;
380
381 for (Int_t iTrack=0; iTrack<fMuonForwardTracksWithBranson->GetEntriesFast(); iTrack++) {
e806e863 382 for (Int_t jTrack=0; jTrack<iTrack; jTrack++) {
789c2e43 383
7e3dd1af 384 AliMuonForwardTrack *track0_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(iTrack);
385 AliMuonForwardTrack *track1_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(jTrack);
789c2e43 386
7e3dd1af 387 AliMuonForwardTrack *track0_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(iTrack);
388 AliMuonForwardTrack *track1_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(jTrack);
e806e863 389
7e3dd1af 390 if (fMatchTrigger) if (!track0_WithBranson->GetMatchTrigger() || !track1_WithBranson->GetMatchTrigger()) continue;
d8c2cc3e 391 if (fOption!=kPionsKaons && (!track0_WithBranson->GetMCTrackRef() || !track1_WithBranson->GetMCTrackRef())) continue;
e806e863 392
789c2e43 393 new ((*fMuonForwardTrackPairsWithBranson)[nMuonPairs]) AliMuonForwardTrackPair(track0_WithBranson, track1_WithBranson);
394 new ((*fMuonForwardTrackPairsWithoutBranson)[nMuonPairs]) AliMuonForwardTrackPair(track0_WithoutBranson, track1_WithoutBranson);
395
396 AliMuonForwardTrackPair *trackPairWithBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithBranson->At(nMuonPairs);
397 if (!(fOption==kResonanceOnly && !trackPairWithBranson->IsResonance())) nMuonPairs++;
398
e806e863 399 }
400 }
401
402}
403
404//====================================================================================================================================================
405
7e3dd1af 406Bool_t AliMuonForwardTrackAnalysis::PassedCutSingleMuon(AliMuonForwardTrack *track) {
407
408 AliMUONTrackParam *param = track->GetTrackParamAtMFTCluster(0);
409 AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
410
411 if (track->Pt()<fPtMinSingleMuons) return kFALSE;
412
413 if (fCutOnOffsetChi2) {
414 Double_t offset = 1.e4*track->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
415 Double_t chi2OverNdf = track->GetGlobalChi2() / Double_t(track->GetNMFTClusters()+track->GetNMUONClusters());
416 offset /= fScaleOffset;
417 chi2OverNdf /= fScaleChi2;
418 offset -= fCenterOffset;
419 chi2OverNdf -= fCenterChi2;
420 // printf("cut on offset and chi2: returning %d\n", TMath::Sqrt(offset*offset + chi2OverNdf*chi2OverNdf)>fRadiusCut);
421 if (TMath::Sqrt(offset*offset + chi2OverNdf*chi2OverNdf) > fRadiusCut) return kFALSE;
422 }
423
424 return kTRUE;
425
426}
427
428//====================================================================================================================================================
429
430Bool_t AliMuonForwardTrackAnalysis::PassedCutMuonPair(AliMuonForwardTrackPair *pair) {
431
432 return PassedCutSingleMuon(pair->GetTrack(0)) && PassedCutSingleMuon(pair->GetTrack(1));
433
434}
435
436//====================================================================================================================================================
437
e806e863 438void AliMuonForwardTrackAnalysis::Terminate(Char_t *outputFileName) {
439
e806e863 440 TFile *fileOut = new TFile(Form("%s/%s",fOutputDir.Data(),outputFileName), "recreate");
441
442 printf("Writing output objects to file %s\n", fileOut->GetName());
443
444 fHistOffsetSingleMuonsX -> Write();
445 fHistOffsetSingleMuonsY -> Write();
446 fHistOffsetSingleMuons -> Write();
447 fHistWOffsetSingleMuons -> Write();
448 fHistErrorSingleMuonsX -> Write();
449 fHistErrorSingleMuonsY -> Write();
450
d8c2cc3e 451 fHistSingleMuonsPtRapidityMC -> Write();
b5ab1ac4 452 fHistSingleMuonsOffsetChi2 -> Write();
e806e863 453
d8c2cc3e 454 fHistZOriginSingleMuonsMC -> Write();
455 fHistZROriginSingleMuonsMC -> Write();
456
cd2f51d2 457 for (Int_t iPtBin=0; iPtBin<fNPtDimuBins+1; iPtBin++) {
789c2e43 458 fHistMassMuonPairs[iPtBin] -> Write();
459 fHistMassMuonPairsWithoutMFT[iPtBin] -> Write();
460 fHistMassMuonPairsMC[iPtBin] -> Write();
461 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] -> Write();
462 fHistWOffsetMuonPairsAtPCA[iPtBin] -> Write();
463 fHistDistancePrimaryVtxPCA[iPtBin] -> Write();
464 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] -> Write();
465 if (fEvalDimuonVtxResolution) {
466 fHistDimuonVtxResolutionX[iPtBin] -> Write();
467 fHistDimuonVtxResolutionY[iPtBin] -> Write();
468 fHistDimuonVtxResolutionZ[iPtBin] -> Write();
469 }
cd2f51d2 470 }
471
789c2e43 472 fHistMassMuonPairsVsPt -> Write();
473 fHistMassMuonPairsWithoutMFTVsPt -> Write();
474 fHistMassMuonPairsVsPtLSp -> Write();
475 fHistMassMuonPairsWithoutMFTVsPtLSp -> Write();
476 fHistMassMuonPairsVsPtLSm -> Write();
477 fHistMassMuonPairsWithoutMFTVsPtLSm -> Write();
478
d8c2cc3e 479 fHistRapidityPtMuonPairs -> Write();
e806e863 480
481 fileOut -> Close();
482
483}
484
485//====================================================================================================================================================
486
487void AliMuonForwardTrackAnalysis::BookHistos() {
488
489 fHistOffsetSingleMuonsX = new TH1D("fHistOffsetSingleMuonsX", "Offset for single muons along X", 200, -1000, 1000);
490 fHistOffsetSingleMuonsY = new TH1D("fHistOffsetSingleMuonsY", "Offset for single muons along Y", 200, -1000, 1000);
491 fHistErrorSingleMuonsX = new TH1D("fHistErrorSingleMuonsX", "Coordinate Error for single muons along X", 200, 0, 1000);
492 fHistErrorSingleMuonsY = new TH1D("fHistErrorSingleMuonsY", "Coordinate Error for single muons along Y", 200, 0, 1000);
b5ab1ac4 493 fHistOffsetSingleMuons = new TH1D("fHistOffsetSingleMuons", "Offset for single muons", 200, 0, 2000);
e806e863 494 fHistWOffsetSingleMuons = new TH1D("fHistWOffsetSingleMuons", "Weighted Offset for single muons", 300, 0, 15);
495
d8c2cc3e 496 fHistSingleMuonsPtRapidityMC = new TH2D("fHistSingleMuonsPtRapidityMC", "Phase Space for single muons", 100, -4.5, -2., 100, 0., 10.);
789c2e43 497 fHistSingleMuonsOffsetChi2 = new TH2D("fHistSingleMuonsOffsetChi2", "Offset vs #chi^{2}/ndf for single muons", 400, 0, 4000, 200, 0, 10);
d8c2cc3e 498 fHistZOriginSingleMuonsMC = new TH1D("fHistZOriginSingleMuonsMC", "Z origin for single muons (from MC)", 1000, 0., 500.);
499 fHistZROriginSingleMuonsMC = new TH2D("fHistZROriginSingleMuonsMC", "Z-R origin for single muons (from MC)", 1000, 0., 500., 1000, 0., 100.);
e806e863 500
501 fHistOffsetSingleMuonsX -> SetXTitle("Offset(X) [#mum]");
502 fHistOffsetSingleMuonsY -> SetXTitle("Offset(Y) [#mum]");
503 fHistErrorSingleMuonsX -> SetXTitle("Err. on track position at z_{vtx} (X) [#mum]");
504 fHistErrorSingleMuonsY -> SetXTitle("Err. on track position at z_{vtx} (Y) [#mum]");
505 fHistOffsetSingleMuons -> SetXTitle("Offset [#mum]");
506 fHistWOffsetSingleMuons -> SetXTitle("Weighted Offset");
507
d8c2cc3e 508 fHistSingleMuonsPtRapidityMC -> SetXTitle("y^{#mu}");
509 fHistSingleMuonsPtRapidityMC -> SetYTitle("p_{T}^{#mu} [GeV/c]");
b5ab1ac4 510 fHistSingleMuonsOffsetChi2 -> SetXTitle("Offset [#mum]");
511 fHistSingleMuonsOffsetChi2 -> SetYTitle("#chi^{2}/ndf");
d8c2cc3e 512
513 fHistZOriginSingleMuonsMC -> SetXTitle("Z [cm]");
514 fHistZROriginSingleMuonsMC -> SetXTitle("Z [cm]");
515 fHistZROriginSingleMuonsMC -> SetXTitle("R [cm]");
e806e863 516
517 fHistOffsetSingleMuonsX -> Sumw2();
518 fHistOffsetSingleMuonsY -> Sumw2();
519 fHistErrorSingleMuonsX -> Sumw2();
520 fHistErrorSingleMuonsY -> Sumw2();
521 fHistOffsetSingleMuons -> Sumw2();
522 fHistWOffsetSingleMuons -> Sumw2();
523
d8c2cc3e 524 fHistZOriginSingleMuonsMC -> Sumw2();
525 fHistZROriginSingleMuonsMC -> Sumw2();
526
527 fHistSingleMuonsPtRapidityMC -> Sumw2();
cd2f51d2 528 fHistSingleMuonsOffsetChi2 -> Sumw2();
e806e863 529
530 //--------------------------------------------
531
cd2f51d2 532 for (Int_t iPtBin=0; iPtBin<=fNPtDimuBins+1; iPtBin++) {
533
534 if (!iPtBin) {
cd2f51d2 535 fHistMassMuonPairs[iPtBin] = new TH1D(Form("fHistMassMuonPairs_%d",iPtBin),
536 "Dimuon Mass (MUON+MFT) (All p_{T}^{#mu#mu})",
537 fNMassBins, fMassMin, fMassMax);
538 fHistMassMuonPairsWithoutMFT[iPtBin] = new TH1D(Form("fHistMassMuonPairsWithoutMFT_%d",iPtBin),
539 "Dimuon Mass (MUON only) (All p_{T}^{#mu#mu})",
540 fNMassBins, fMassMin, fMassMax);
541 fHistMassMuonPairsMC[iPtBin] = new TH1D(Form("fHistMassMuonPairsMC_%d",iPtBin),
542 "Dimuon Mass (MC) (All p_{T}^{#mu#mu})",
543 fNMassBins, fMassMin, fMassMax);
789c2e43 544 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] = new TH1D(Form("fHistWOffsetMuonPairsAtPrimaryVtx_%d",iPtBin),
545 "Weighted Offset for Muon Pairs at Primary Vertex (All p_{T}^{#mu#mu})",
546 300, 0, 60);
547 fHistWOffsetMuonPairsAtPCA[iPtBin] = new TH1D(Form("fHistWOffsetMuonPairsAtPCA_%d",iPtBin),
548 "Weighted Offset for Muon Pairs at PCA (All p_{T}^{#mu#mu})",
549 300, 0, 60);
550 fHistDistancePrimaryVtxPCA[iPtBin] = new TH1D(Form("fHistDistancePrimaryVtxPCA_%d",iPtBin),
551 "Distance between PCA and primary vertex (All p_{T}^{#mu#mu})",
552 1000, 0, 50000);
553 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] = new TH2D(Form("fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA_%d",iPtBin),
554 "DistancePrimaryVtxPCA vs WOffsetMuonPairsAtPCA (All p_{T}^{#mu#mu})",
555 300, 0, 60, 1000, 0, 50000);
556 if (fEvalDimuonVtxResolution) {
557 fHistDimuonVtxResolutionX[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionX_%d",iPtBin),
558 "Dimuon Vtx offset along X (All p_{T}^{#mu#mu})",
559 2000, -1000, 1000);
560 fHistDimuonVtxResolutionY[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionY_%d",iPtBin),
561 "Dimuon Vtx offset along Y (All p_{T}^{#mu#mu})",
562 2000, -1000, 1000);
563 fHistDimuonVtxResolutionZ[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionZ_%d",iPtBin),
564 "Dimuon Vtx offset along Z (All p_{T}^{#mu#mu})",
565 2000, -10000, 10000);
566 }
cd2f51d2 567 }
e806e863 568
cd2f51d2 569 else {
570 Double_t ptMin = fPtAxisDimuons->GetBinLowEdge(iPtBin);
571 Double_t ptMax = fPtAxisDimuons->GetBinUpEdge(iPtBin);
cd2f51d2 572 fHistMassMuonPairs[iPtBin] = new TH1D(Form("fHistMassMuonPairs_%d",iPtBin),
573 Form("Dimuon Mass (MUON+MFT) (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
574 fNMassBins, fMassMin, fMassMax);
575 fHistMassMuonPairsWithoutMFT[iPtBin] = new TH1D(Form("fHistMassMuonPairsWithoutMFT_%d",iPtBin),
576 Form("Dimuon Mass (MUON only) (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
577 fNMassBins, fMassMin, fMassMax);
578 fHistMassMuonPairsMC[iPtBin] = new TH1D(Form("fHistMassMuonPairsMC_%d",iPtBin),
579 Form("Dimuon Mass (MC) (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
789c2e43 580 fNMassBins, fMassMin, fMassMax);
581 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] = new TH1D(Form("fHistWOffsetMuonPairsAtPrimaryVtx_%d",iPtBin),
582 Form("Weighted Offset for Muon Pairs at Primary Vertex (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",
583 ptMin,ptMax), 300, 0, 60);
584 fHistWOffsetMuonPairsAtPCA[iPtBin] = new TH1D(Form("fHistWOffsetMuonPairsAtPCA_%d",iPtBin),
585 Form("Weighted Offset for Muon Pairs at PCA (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",
586 ptMin,ptMax), 300, 0, 60);
587 fHistDistancePrimaryVtxPCA[iPtBin] = new TH1D(Form("fHistDistancePrimaryVtxPCA_%d",iPtBin),
588 Form("Distance between PCA and primary vertex (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",
589 ptMin,ptMax), 1000, 0, 50000);
590 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] = new TH2D(Form("fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA_%d",iPtBin),
591 Form("DistancePrimaryVtxPCA vs WOffsetMuonPairsAtPCA (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",
592 ptMin,ptMax), 300, 0, 60, 1000, 0, 50000);
593 if (fEvalDimuonVtxResolution) {
594 fHistDimuonVtxResolutionX[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionX_%d",iPtBin),
595 Form("Dimuon Vtx offset along X (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
596 2000, -1000, 1000);
597 fHistDimuonVtxResolutionY[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionY_%d",iPtBin),
598 Form("Dimuon Vtx offset along Y (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
599 2000, -1000, 1000);
600 fHistDimuonVtxResolutionZ[iPtBin] = new TH1D(Form("fHistDimuonVtxResolutionZ_%d",iPtBin),
601 Form("Dimuon Vtx offset along Z (%3.1f < p_{T}^{#mu#mu} < %3.1f GeV/c)",ptMin,ptMax),
602 2000, -10000, 10000);
603 }
cd2f51d2 604 }
605
789c2e43 606 fHistMassMuonPairs[iPtBin] -> SetXTitle("Mass [GeV/c^{2}]");
607 fHistMassMuonPairsWithoutMFT[iPtBin] -> SetXTitle("Mass [GeV/c^{2}]");
608 fHistMassMuonPairsMC[iPtBin] -> SetXTitle("Mass [GeV/c^{2}]");
609 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] -> SetXTitle("Weighted Offset at Primary Vtx");
610 fHistWOffsetMuonPairsAtPCA[iPtBin] -> SetXTitle("Weighted Offset at PCA");
611 fHistDistancePrimaryVtxPCA[iPtBin] -> SetXTitle("PCA - Primary Vtx [#mum]");
612 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] -> SetXTitle("Weighted Offset at PCA");
613 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] -> SetYTitle("PCA - Primary Vtx [#mum]");
614 if (fEvalDimuonVtxResolution) {
615 fHistDimuonVtxResolutionX[iPtBin] -> SetXTitle("Offset [#mum]");
616 fHistDimuonVtxResolutionY[iPtBin] -> SetXTitle("Offset [#mum]");
617 fHistDimuonVtxResolutionZ[iPtBin] -> SetXTitle("Offset [#mum]");
618 }
cd2f51d2 619
789c2e43 620 fHistMassMuonPairs[iPtBin] -> Sumw2();
621 fHistMassMuonPairsWithoutMFT[iPtBin] -> Sumw2();
622 fHistMassMuonPairsMC[iPtBin] -> Sumw2();
623 fHistWOffsetMuonPairsAtPrimaryVtx[iPtBin] -> Sumw2();
624 fHistWOffsetMuonPairsAtPCA[iPtBin] -> Sumw2();
625 fHistDistancePrimaryVtxPCA[iPtBin] -> Sumw2();
626 fHistDistancePrimaryVtxPCAvsWOffsetMuonPairsAtPCA[iPtBin] -> Sumw2();
627 if (fEvalDimuonVtxResolution) {
628 fHistDimuonVtxResolutionX[iPtBin] -> Sumw2();
629 fHistDimuonVtxResolutionY[iPtBin] -> Sumw2();
630 fHistDimuonVtxResolutionZ[iPtBin] -> Sumw2();
631 }
cd2f51d2 632
633 }
634
789c2e43 635 fHistMassMuonPairsVsPt = new TH2D("fHistMassMuonPairsVsPt", "Dimuon Mass (MUON+MFT) vs p_{T}", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
636 fHistMassMuonPairsWithoutMFTVsPt = new TH2D("fHistMassMuonPairsWithoutMFTVsPt", "Dimuon Mass (MUON only) vs p_{T}", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
637 fHistMassMuonPairsVsPtLSp = new TH2D("fHistMassMuonPairsVsPtLSp", "Dimuon Mass (MUON+MFT) vs p_{T} Like sign ++", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
638 fHistMassMuonPairsWithoutMFTVsPtLSp = new TH2D("fHistMassMuonPairsWithoutMFTVsPtLSp", "Dimuon Mass (MUON only) vs p_{T} Like sign ++", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
639 fHistMassMuonPairsVsPtLSm = new TH2D("fHistMassMuonPairsVsPtLSm", "Dimuon Mass (MUON+MFT) vs p_{T} Like sign --", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
640 fHistMassMuonPairsWithoutMFTVsPtLSm = new TH2D("fHistMassMuonPairsWithoutMFTVsPtLSm", "Dimuon Mass (MUON only) vs p_{T} Like sign --", fNMassBins, fMassMin, fMassMax, 20, 0., 10.);
641
642 fHistMassMuonPairsVsPt -> SetXTitle("Mass [GeV/c^{2}]");
643 fHistMassMuonPairsWithoutMFTVsPt -> SetXTitle("Mass [GeV/c^{2}]");
644 fHistMassMuonPairsVsPt -> SetYTitle("p_{T} [GeV/c]");
645 fHistMassMuonPairsWithoutMFTVsPt -> SetYTitle("p_{T} [GeV/c]");
646 fHistMassMuonPairsVsPt -> Sumw2();
647 fHistMassMuonPairsWithoutMFTVsPt -> Sumw2();
648 fHistMassMuonPairsVsPtLSp -> SetXTitle("Mass [GeV/c^{2}]");
649 fHistMassMuonPairsWithoutMFTVsPtLSp -> SetXTitle("Mass [GeV/c^{2}]");
650 fHistMassMuonPairsVsPtLSp -> SetYTitle("p_{T} [GeV/c]");
651 fHistMassMuonPairsWithoutMFTVsPtLSp -> SetYTitle("p_{T} [GeV/c]");
652 fHistMassMuonPairsVsPtLSp -> Sumw2();
653 fHistMassMuonPairsWithoutMFTVsPtLSp -> Sumw2();
654 fHistMassMuonPairsVsPtLSm -> SetXTitle("Mass [GeV/c^{2}]");
655 fHistMassMuonPairsWithoutMFTVsPtLSm -> SetXTitle("Mass [GeV/c^{2}]");
656 fHistMassMuonPairsVsPtLSm -> SetYTitle("p_{T} [GeV/c]");
657 fHistMassMuonPairsWithoutMFTVsPtLSm -> SetYTitle("p_{T} [GeV/c]");
658 fHistMassMuonPairsVsPtLSm -> Sumw2();
659 fHistMassMuonPairsWithoutMFTVsPtLSm -> Sumw2();
660
d8c2cc3e 661 if (fOption==kPionsKaons) fHistRapidityPtMuonPairs = new TH2D("fHistRapidityPtMuonPairs", "Dimuon Phase Space (rec)", 20, -4.5, -2., 20, 0., 10.);
662 else fHistRapidityPtMuonPairs = new TH2D("fHistRapidityPtMuonPairs", "Dimuon Phase Space (MC)", 100, -4.5, -2., 100, 0., 10.);
663 fHistRapidityPtMuonPairs -> SetXTitle("Rapidity");
664 fHistRapidityPtMuonPairs -> SetYTitle("p_{T} [GeV/c]");
665 fHistRapidityPtMuonPairs -> Sumw2();
e806e863 666
667}
668
669//====================================================================================================================================================
670