MFT/AliMuonForwardTrackAnalysis.cxx

   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"
  25 #include "TGraph.h"
  26 #include "AliMuonForwardTrackAnalysis.h"
  27
  28 ClassImp(AliMuonForwardTrackAnalysis)
  29
  30 //====================================================================================================================================================
  31
  32 AliMuonForwardTrackAnalysis::AliMuonForwardTrackAnalysis():
  33   TObject(),
  34   fInputDir(0),
  35   fOutputDir(0),
  36   fInputTreeWithBranson(0x0),
  37   fInputTreeWithoutBranson(0x0),
  38   fMuonForwardTracksWithBranson(0),
  39   fMuonForwardTrackPairsWithBranson(0),
  40   fMuonForwardTracksWithoutBranson(0),
  41   fMuonForwardTrackPairsWithoutBranson(0),
  42   fMFTTrackWithBranson(0),
  43   fMFTTrackWithoutBranson(0),
  44   fMFTTrack(0),
  45   fMFTTrackPairWithBranson(0),
  46   fMFTTrackPairWithoutBranson(0),
  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),
  57   fHistOffsetSingleMuonsX(0x0),
  58   fHistOffsetSingleMuonsY(0x0),
  59   fHistOffsetSingleMuons(0x0),
  60   fHistWOffsetSingleMuons(0x0),
  61   fHistErrorSingleMuonsX(0x0),
  62   fHistErrorSingleMuonsY(0x0),
  63   fHistOffsetSingleMuonsX_vsPtRapidity(0x0),
  64   fHistOffsetSingleMuonsY_vsPtRapidity(0x0),
  65   fHistSingleMuonsPtRapidity(0x0),
  66   fHistSingleMuonsOffsetChi2(0x0),
  67   fHistWOffsetMuonPairs(0x0),
  68   fHistMassMuonPairs(0x0),
  69   fHistMassMuonPairsWithoutMFT(0x0),
  70   fHistMassMuonPairsMC(0x0),
  71   fHistRapidityPtMuonPairsMC(0x0),
  72   fGraphSingleMuonsOffsetChi2(0x0),
  73   fNMassBins(1000),
  74   fMassMin(0),
  75   fMassMax(10),
  76   fSingleMuonAnalysis(1),
  77   fMuonPairAnalysis(1),
  78   fMatchTrigger(0),
  79   fOption(0),
  80   fXVertResMC(50.e-4),
  81   fYVertResMC(50.e-4),
  82   fZVertResMC(50.e-4),
  83   fPrimaryVtxX(0.),
  84   fPrimaryVtxY(0.),
  85   fPrimaryVtxZ(0.),
  86   fMaxNWrongClustersMC(999),
  87   fPtMinSingleMuons(0),
  88   fUseBransonForCut(kFALSE),
  89   fUseBransonForKinematics(kFALSE),
  90   fCutOnOffsetChi2(kFALSE),
  91   fCenterOffset(0.),
  92   fCenterChi2(1.),
  93   fScaleOffset(250.),
  94   fScaleChi2(4.5),
  95   fRadiusCut(1.)
  96 {
  97
  98   // default constructor
  99
 100   for (Int_t rapBin=0; rapBin<fNRapBinsOffsetSingleMuons; rapBin++) {
 101     for (Int_t ptBin=0; ptBin<fNPtBinsOffsetSingleMuons; ptBin++) {
 102       fHistOffsetSingleMuonsX_tmp[rapBin][ptBin] = 0x0;
 103       fHistOffsetSingleMuonsY_tmp[rapBin][ptBin] = 0x0;
 104     }
 105   }
 106
 107 }
 108
 109 //====================================================================================================================================================
 110
 111 Bool_t AliMuonForwardTrackAnalysis::Init(Char_t *inputFileName) {
 112
 113   BookHistos();
 114
 115   TFile *inputFile = new TFile(Form("%s/%s",fInputDir.Data(),inputFileName));
 116   if (!inputFile || !inputFile->IsOpen()) {
 117     AliError(Form("Error opening file %s", inputFileName));
 118     return kFALSE;
 119   }
 120   fInputTreeWithBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithBranson");
 121   if (!fInputTreeWithBranson) {
 122     AliError("Error reading input tree");
 123     return kFALSE;
 124   }
 125   fInputTreeWithoutBranson = (TTree*) inputFile->Get("AliMuonForwardTracksWithoutBranson");
 126   if (!fInputTreeWithoutBranson) {
 127     AliError("Error reading input tree");
 128     return kFALSE;
 129   }
 130
 131   if (fFirstEvent<0 || fLastEvent<0 || fFirstEvent>fLastEvent || fFirstEvent>=fInputTreeWithBranson->GetEntries()) {
 132     fFirstEvent = 0;
 133     fLastEvent  = fInputTreeWithBranson->GetEntries()-1;
 134   }
 135   else {
 136     fLastEvent = TMath::Min(fLastEvent, Int_t(fInputTreeWithBranson->GetEntries()-1));
 137   }
 138
 139   AliInfo(Form("Analysing events %d to %d", fFirstEvent, fLastEvent));
 140
 141   fMuonForwardTracksWithBranson = new TClonesArray("AliMuonForwardTrack");
 142   fInputTreeWithBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithBranson);
 143
 144   fMuonForwardTracksWithoutBranson = new TClonesArray("AliMuonForwardTrack");
 145   fInputTreeWithoutBranson->SetBranchAddress("tracks", &fMuonForwardTracksWithoutBranson);
 146
 147   TGeoManager::Import(Form("%s/geometry.root",fInputDir.Data()));
 148
 149   AliMUONTrackExtrap::SetField();
 150
 151   fMuonForwardTrackPairsWithBranson    = new TClonesArray("AliMuonForwardTrackPair");
 152   fMuonForwardTrackPairsWithoutBranson = new TClonesArray("AliMuonForwardTrackPair");
 153
 154   return kTRUE;
 155
 156 }
 157
 158 //====================================================================================================================================================
 159
 160 Bool_t AliMuonForwardTrackAnalysis::LoadNextEvent() {
 161
 162   if (fEv>fLastEvent) return kFALSE;
 163   if (fEv<fFirstEvent) { fEv++; return kTRUE; }
 164   fMuonForwardTracksWithBranson -> Delete();
 165   fMuonForwardTracksWithoutBranson -> Delete();
 166   fInputTreeWithBranson->GetEvent(fEv);
 167   fInputTreeWithoutBranson->GetEvent(fEv);
 168   AliInfo(Form("**** analyzing event # %4d (%3d tracks) ****", fEv, fMuonForwardTracksWithBranson->GetEntries()));
 169
 170   fPrimaryVtxX = gRandom->Gaus(0., fXVertResMC);
 171   fPrimaryVtxY = gRandom->Gaus(0., fYVertResMC);
 172   fPrimaryVtxZ = gRandom->Gaus(0., fZVertResMC);
 173
 174   if (fSingleMuonAnalysis) {
 175     fNTracksAnalyzedOfEvent = 0;
 176     fNTracksOfEvent = fMuonForwardTracksWithBranson->GetEntries();
 177     while (AnalyzeSingleMuon()) continue;
 178   }
 179
 180   if (fMuonPairAnalysis) {
 181     if (fMuonForwardTrackPairsWithBranson) {
 182       fMuonForwardTrackPairsWithBranson->Delete();
 183       fMuonForwardTrackPairsWithoutBranson->Delete();
 184     }
 185     BuildMuonPairs();
 186     fNPairsAnalyzedOfEvent = 0;
 187     fNPairsOfEvent = fMuonForwardTrackPairsWithBranson->GetEntries();
 188     while (AnalyzeMuonPair()) continue;
 189   }
 190
 191   fEv++;
 192
 193   return kTRUE;
 194
 195 }
 196
 197 //====================================================================================================================================================
 198
 199 Bool_t AliMuonForwardTrackAnalysis::AnalyzeSingleMuon() {
 200
 201   if (fNTracksAnalyzedOfEvent>=fNTracksOfEvent) return kFALSE;
 202
 203   fMFTTrackWithBranson    = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(fNTracksAnalyzedOfEvent);
 204   fMFTTrackWithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(fNTracksAnalyzedOfEvent);
 205
 206   fNTracksAnalyzedOfEvent++;
 207
 208   Bool_t passedCut = kFALSE;
 209   if (fUseBransonForCut) passedCut = PassedCutSingleMuon(fMFTTrackWithBranson);
 210   else passedCut = PassedCutSingleMuon(fMFTTrackWithoutBranson);
 211   if (!passedCut) return kTRUE;
 212
 213   if (fUseBransonForKinematics) fMFTTrack = fMFTTrackWithBranson;
 214   else fMFTTrack = fMFTTrackWithoutBranson;
 215   if (fMatchTrigger && !fMFTTrack->GetMatchTrigger()) return kTRUE;
 216   fMCRefTrack = fMFTTrack->GetMCTrackRef();
 217
 218   if (!fMCRefTrack) return kTRUE;
 219   if (fMFTTrack->GetNWrongClustersMC()>fMaxNWrongClustersMC) return kTRUE;
 220
 221   AliMUONTrackParam *param = fMFTTrack->GetTrackParamAtMFTCluster(0);
 222   AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
 223
 224   TLorentzVector pMu;
 225   Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
 226   Double_t energy = TMath::Sqrt(param->P()*param->P() + mMu*mMu);
 227   pMu.SetPxPyPzE(param->Px(), param->Py(), param->Pz(), energy);
 228
 229   TMatrixD cov(5,5);
 230   cov = param->GetCovariances();
 231
 232   fHistErrorSingleMuonsX -> Fill(1.e4*TMath::Sqrt(cov(0,0)));
 233   fHistErrorSingleMuonsY -> Fill(1.e4*TMath::Sqrt(cov(2,2)));
 234
 235   Double_t dX = fMFTTrack->GetOffsetX(fPrimaryVtxX, fPrimaryVtxZ);
 236   Double_t dY = fMFTTrack->GetOffsetY(fPrimaryVtxY, fPrimaryVtxZ);
 237
 238   Double_t offset = fMFTTrack->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
 239   Double_t weightedOffset = fMFTTrack->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
 240
 241   //  AliDebug(2, Form("pdg code = %d\n", fMCRefTrack->GetPdgCode()));
 242
 243   fHistOffsetSingleMuonsX -> Fill(1.e4*dX);
 244   fHistOffsetSingleMuonsY -> Fill(1.e4*dY);
 245   Int_t rapBin = fHistOffsetSingleMuonsX_vsPtRapidity->GetXaxis()->FindBin(pMu.Rapidity());
 246   Int_t ptBin  = fHistOffsetSingleMuonsX_vsPtRapidity->GetYaxis()->FindBin(pMu.Pt());
 247   //  AliDebug(1, Form("pt = %f (%d), rap = %f (%d)\n", pMu.Pt(), pMu.Rapidity(), ptBin, rapBin));
 248   if (0<rapBin && rapBin<=fNRapBinsOffsetSingleMuons && 0<ptBin && ptBin<=fNPtBinsOffsetSingleMuons) {
 249     fHistOffsetSingleMuonsX_tmp[rapBin-1][ptBin-1]->Fill(1.e4*dX);
 250     fHistOffsetSingleMuonsY_tmp[rapBin-1][ptBin-1]->Fill(1.e4*dY);
 251   }
 252   fHistSingleMuonsPtRapidity  -> Fill(pMu.Rapidity(), pMu.Pt());
 253   fHistOffsetSingleMuons      -> Fill(1.e4*offset);
 254   fHistWOffsetSingleMuons     -> Fill(weightedOffset);
 255   Double_t chi2OverNdf = fMFTTrack->GetGlobalChi2()/Double_t(fMFTTrack->GetNMFTClusters()+fMFTTrack->GetNMUONClusters());
 256   fHistSingleMuonsOffsetChi2  -> Fill(1.e4*offset, chi2OverNdf);
 257   fGraphSingleMuonsOffsetChi2 -> SetPoint(fGraphSingleMuonsOffsetChi2->GetN(),1.e4*offset, chi2OverNdf);
 258
 259   fNTracksAnalyzed++;
 260
 261   return kTRUE;
 262
 263 }
 264
 265 //====================================================================================================================================================
 266
 267 Bool_t AliMuonForwardTrackAnalysis::AnalyzeMuonPair() {
 268
 269   if (fNPairsAnalyzedOfEvent>=fNPairsOfEvent) return kFALSE;
 270
 271   fMFTTrackPairWithBranson    = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithBranson->At(fNPairsAnalyzedOfEvent);
 272   fMFTTrackPairWithoutBranson = (AliMuonForwardTrackPair*) fMuonForwardTrackPairsWithoutBranson->At(fNPairsAnalyzedOfEvent);
 273
 274   fNPairsAnalyzedOfEvent++;
 275
 276   Bool_t passedCut = kFALSE;
 277   if (fUseBransonForCut) passedCut = PassedCutMuonPair(fMFTTrackPairWithBranson);
 278   else passedCut = PassedCutMuonPair(fMFTTrackPairWithoutBranson);
 279
 280   if (!passedCut) return kTRUE;
 281
 282   if (fUseBransonForKinematics) fMFTTrackPair = fMFTTrackPairWithBranson;
 283   else fMFTTrackPair = fMFTTrackPairWithoutBranson;
 284
 285   if ( fMFTTrackPair->GetTrack(0)->GetNWrongClustersMC()>fMaxNWrongClustersMC ||
 286        fMFTTrackPair->GetTrack(1)->GetNWrongClustersMC()>fMaxNWrongClustersMC ) return kTRUE;
 287
 288   if (fOption==kResonanceOnly && !fMFTTrackPair->IsResonance()) return kTRUE;
 289
 290   fHistMassMuonPairs           -> Fill(fMFTTrackPair->GetMass(fPrimaryVtxZ));
 291   fHistWOffsetMuonPairs        -> Fill(fMFTTrackPair->GetWeightedOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
 292   fHistMassMuonPairsWithoutMFT -> Fill(fMFTTrackPair->GetMassWithoutMFT(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ));
 293   fHistMassMuonPairsMC         -> Fill(fMFTTrackPair->GetMassMC());
 294   fHistRapidityPtMuonPairsMC   -> Fill(fMFTTrackPair->GetRapidityMC(), fMFTTrackPair->GetPtMC());
 295
 296   AliDebug(1, Form("mass = %f   MC = %f", fMFTTrackPair->GetMass(fPrimaryVtxZ), fMFTTrackPair->GetMassMC()));
 297
 298   return kTRUE;
 299
 300 }
 301
 302 //====================================================================================================================================================
 303
 304 void AliMuonForwardTrackAnalysis::BuildMuonPairs() {
 305
 306   for (Int_t iTrack=0; iTrack<fMuonForwardTracksWithBranson->GetEntries(); iTrack++) {
 307     for (Int_t jTrack=0; jTrack<iTrack; jTrack++) {
 308
 309       AliMuonForwardTrack *track0_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(iTrack);
 310       AliMuonForwardTrack *track1_WithBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithBranson->At(jTrack);
 311
 312       AliMuonForwardTrack *track0_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(iTrack);
 313       AliMuonForwardTrack *track1_WithoutBranson = (AliMuonForwardTrack*) fMuonForwardTracksWithoutBranson->At(jTrack);
 314
 315       if (fMatchTrigger) if (!track0_WithBranson->GetMatchTrigger() || !track1_WithBranson->GetMatchTrigger()) continue;
 316       if (!track0_WithBranson->GetMCTrackRef() || !track1_WithBranson->GetMCTrackRef()) continue;
 317
 318       AliMuonForwardTrackPair *trackPairWithBranson    = new AliMuonForwardTrackPair(track0_WithBranson, track1_WithBranson);
 319       AliMuonForwardTrackPair *trackPairWithoutBranson = new AliMuonForwardTrackPair(track0_WithoutBranson, track1_WithoutBranson);
 320       if (fOption==kResonanceOnly && !trackPairWithBranson->IsResonance()) {
 321         delete trackPairWithBranson;
 322         delete trackPairWithoutBranson;
 323         continue;
 324       }
 325       new ((*fMuonForwardTrackPairsWithBranson)[fMuonForwardTrackPairsWithBranson->GetEntries()]) AliMuonForwardTrackPair(*trackPairWithBranson);
 326       new ((*fMuonForwardTrackPairsWithoutBranson)[fMuonForwardTrackPairsWithoutBranson->GetEntries()]) AliMuonForwardTrackPair(*trackPairWithoutBranson);
 327     }
 328   }
 329
 330 }
 331
 332 //====================================================================================================================================================
 333
 334 Bool_t AliMuonForwardTrackAnalysis::PassedCutSingleMuon(AliMuonForwardTrack *track) {
 335
 336   AliMUONTrackParam *param = track->GetTrackParamAtMFTCluster(0);
 337   AliMUONTrackExtrap::ExtrapToZCov(param, fPrimaryVtxZ);
 338
 339   if (track->Pt()<fPtMinSingleMuons) return kFALSE;
 340
 341   if (fCutOnOffsetChi2) {
 342     Double_t offset = 1.e4*track->GetOffset(fPrimaryVtxX, fPrimaryVtxY, fPrimaryVtxZ);
 343     Double_t chi2OverNdf = track->GetGlobalChi2() / Double_t(track->GetNMFTClusters()+track->GetNMUONClusters());
 344     offset /= fScaleOffset;
 345     chi2OverNdf /= fScaleChi2;
 346     offset -= fCenterOffset;
 347     chi2OverNdf -= fCenterChi2;
 348     //    printf("cut on offset and chi2: returning %d\n", TMath::Sqrt(offset*offset + chi2OverNdf*chi2OverNdf)>fRadiusCut);
 349     if (TMath::Sqrt(offset*offset + chi2OverNdf*chi2OverNdf) > fRadiusCut) return kFALSE;
 350   }
 351
 352   return kTRUE;
 353
 354 }
 355
 356 //====================================================================================================================================================
 357
 358 Bool_t AliMuonForwardTrackAnalysis::PassedCutMuonPair(AliMuonForwardTrackPair *pair) {
 359
 360   return PassedCutSingleMuon(pair->GetTrack(0)) && PassedCutSingleMuon(pair->GetTrack(1));
 361
 362 }
 363
 364 //====================================================================================================================================================
 365
 366 void AliMuonForwardTrackAnalysis::Terminate(Char_t *outputFileName) {
 367
 368   for (Int_t rapBin=0; rapBin<fNRapBinsOffsetSingleMuons; rapBin++) {
 369     for (Int_t ptBin=0; ptBin<fNPtBinsOffsetSingleMuons; ptBin++) {
 370       Int_t binMin_x = fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->FindBin(-3*fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->GetRMS());
 371       Int_t binMax_x = fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->FindBin(+3*fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->GetRMS());
 372       for (Int_t bin=1; bin<=fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->GetNbinsX(); bin++) {
 373         if (bin<binMin_x || bin>binMax_x) fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->SetBinContent(bin, 0);
 374       }
 375       Int_t binMin_y = fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->FindBin(-3*fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->GetRMS());
 376       Int_t binMax_y = fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->FindBin(+3*fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->GetRMS());
 377       for (Int_t bin=1; bin<=fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->GetNbinsX(); bin++) {
 378         if (bin<binMin_y || bin>binMax_y) fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->SetBinContent(bin, 0);
 379       }
 380       fHistOffsetSingleMuonsX_vsPtRapidity->SetBinContent(rapBin+1, ptBin+1, fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->GetRMS());
 381       fHistOffsetSingleMuonsX_vsPtRapidity->SetBinError(rapBin+1, ptBin+1, fHistOffsetSingleMuonsX_tmp[rapBin][ptBin]->GetRMSError());
 382       fHistOffsetSingleMuonsY_vsPtRapidity->SetBinContent(rapBin+1, ptBin+1, fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->GetRMS());
 383       fHistOffsetSingleMuonsY_vsPtRapidity->SetBinError(rapBin+1, ptBin+1, fHistOffsetSingleMuonsY_tmp[rapBin][ptBin]->GetRMSError());
 384     }
 385   }
 386
 387   TFile *fileOut = new TFile(Form("%s/%s",fOutputDir.Data(),outputFileName), "recreate");
 388
 389   printf("Writing output objects to file %s\n", fileOut->GetName());
 390
 391   fHistOffsetSingleMuonsX -> Write();
 392   fHistOffsetSingleMuonsY -> Write();
 393   fHistOffsetSingleMuons  -> Write();
 394   fHistWOffsetSingleMuons -> Write();
 395   fHistErrorSingleMuonsX  -> Write();
 396   fHistErrorSingleMuonsY  -> Write();
 397
 398   fHistOffsetSingleMuonsX_vsPtRapidity -> Write();
 399   fHistOffsetSingleMuonsY_vsPtRapidity -> Write();
 400
 401 //   for (Int_t rapBin=0; rapBin<fNRapBinsOffsetSingleMuons; rapBin++) {
 402 //     for (Int_t ptBin=0; ptBin<fNPtBinsOffsetSingleMuons; ptBin++) {
 403 //       fHistOffsetSingleMuonsX_tmp[rapBin][ptBin] -> Write();
 404 //       fHistOffsetSingleMuonsY_tmp[rapBin][ptBin] -> Write();
 405 //     }
 406 //   }
 407
 408   fHistSingleMuonsPtRapidity -> Write();
 409   fHistSingleMuonsOffsetChi2 -> Write();
 410
 411   fGraphSingleMuonsOffsetChi2 -> Write();
 412
 413   fHistWOffsetMuonPairs        -> Write();
 414   fHistMassMuonPairs           -> Write();
 415   fHistMassMuonPairsWithoutMFT -> Write();
 416   fHistMassMuonPairsMC         -> Write();
 417   fHistRapidityPtMuonPairsMC   -> Write();
 418
 419   fileOut -> Close();
 420
 421 }
 422
 423 //====================================================================================================================================================
 424
 425 void AliMuonForwardTrackAnalysis::BookHistos() {
 426
 427   fHistOffsetSingleMuonsX = new TH1D("fHistOffsetSingleMuonsX", "Offset for single muons along X",  200, -1000, 1000);
 428   fHistOffsetSingleMuonsY = new TH1D("fHistOffsetSingleMuonsY", "Offset for single muons along Y",  200, -1000, 1000);
 429   fHistErrorSingleMuonsX  = new TH1D("fHistErrorSingleMuonsX",  "Coordinate Error for single muons along X",  200, 0, 1000);
 430   fHistErrorSingleMuonsY  = new TH1D("fHistErrorSingleMuonsY",  "Coordinate Error for single muons along Y",  200, 0, 1000);
 431   fHistOffsetSingleMuons  = new TH1D("fHistOffsetSingleMuons",  "Offset for single muons",          200, 0, 2000);
 432   fHistWOffsetSingleMuons = new TH1D("fHistWOffsetSingleMuons", "Weighted Offset for single muons", 300, 0, 15);
 433
 434   fHistOffsetSingleMuonsX_vsPtRapidity = new TH2D("fHistOffsetSingleMuonsX_vsPtRapidity", "Offset for single muons along X",
 435                                                   10, -4, -2.5, 10, 0.5, 5.5);
 436   fHistOffsetSingleMuonsY_vsPtRapidity = new TH2D("fHistOffsetSingleMuonsY_vsPtRapidity", "Offset for single muons along Y",
 437                                                   10, -4, -2.5, 10, 0.5, 5.5);
 438
 439   for (Int_t rapBin=0; rapBin<fNRapBinsOffsetSingleMuons; rapBin++) {
 440     for (Int_t ptBin=0; ptBin<fNPtBinsOffsetSingleMuons; ptBin++) {
 441       fHistOffsetSingleMuonsX_tmp[rapBin][ptBin] = new TH1D(Form("fHistOffsetSingleMuonsX_tmp_%02d_%02d",rapBin,ptBin), "",  1000, -1000, 1000);
 442       fHistOffsetSingleMuonsY_tmp[rapBin][ptBin] = new TH1D(Form("fHistOffsetSingleMuonsY_tmp_%02d_%02d",rapBin,ptBin), "",  1000, -1000, 1000);
 443     }
 444   }
 445
 446   fHistSingleMuonsPtRapidity = new TH2D("fHistSingleMuonsPtRapidity", "Phase Space for single muons", 10, -4, -2.5, 10, 0.5, 5.5);
 447   fHistSingleMuonsOffsetChi2 = new TH2D("fHistSingleMuonsOffsetChi2", "Offset vs #chi^{2}/ndf for single muons", 400, 0, 4000, 100, 0, 20);
 448
 449   fHistOffsetSingleMuonsX -> SetXTitle("Offset(X)  [#mum]");
 450   fHistOffsetSingleMuonsY -> SetXTitle("Offset(Y)  [#mum]");
 451   fHistErrorSingleMuonsX  -> SetXTitle("Err. on track position at z_{vtx} (X)  [#mum]");
 452   fHistErrorSingleMuonsY  -> SetXTitle("Err. on track position at z_{vtx} (Y)  [#mum]");
 453   fHistOffsetSingleMuons  -> SetXTitle("Offset  [#mum]");
 454   fHistWOffsetSingleMuons -> SetXTitle("Weighted Offset");
 455
 456   fHistOffsetSingleMuonsX_vsPtRapidity -> SetXTitle("y^{#mu}");
 457   fHistOffsetSingleMuonsY_vsPtRapidity -> SetXTitle("y^{#mu}");
 458   fHistOffsetSingleMuonsX_vsPtRapidity -> SetYTitle("p_{T}^{#mu}  [GeV/c]");
 459   fHistOffsetSingleMuonsY_vsPtRapidity -> SetYTitle("p_{T}^{#mu}  [GeV/c]");
 460
 461   fHistSingleMuonsPtRapidity -> SetXTitle("y^{#mu}");
 462   fHistSingleMuonsPtRapidity -> SetYTitle("p_{T}^{#mu}  [GeV/c]");
 463   fHistSingleMuonsOffsetChi2 -> SetXTitle("Offset  [#mum]");
 464   fHistSingleMuonsOffsetChi2 -> SetYTitle("#chi^{2}/ndf");
 465
 466   fHistOffsetSingleMuonsX -> Sumw2();
 467   fHistOffsetSingleMuonsY -> Sumw2();
 468   fHistErrorSingleMuonsX  -> Sumw2();
 469   fHistErrorSingleMuonsY  -> Sumw2();
 470   fHistOffsetSingleMuons  -> Sumw2();
 471   fHistWOffsetSingleMuons -> Sumw2();
 472
 473   fHistOffsetSingleMuonsX_vsPtRapidity -> Sumw2();
 474   fHistOffsetSingleMuonsY_vsPtRapidity -> Sumw2();
 475   fHistSingleMuonsPtRapidity           -> Sumw2();
 476   fHistSingleMuonsOffsetChi2           -> Sumw2();
 477
 478   //--------------------------------------------
 479
 480   fGraphSingleMuonsOffsetChi2 = new TGraph();
 481   fGraphSingleMuonsOffsetChi2 -> SetName("fGraphSingleMuonsOffsetChi2");
 482   fGraphSingleMuonsOffsetChi2 -> SetTitle("fGraphSingleMuonsOffsetChi2");
 483
 484   //--------------------------------------------
 485
 486   fHistWOffsetMuonPairs        = new TH1D("fHistWOffsetMuonPairs",        "Weighted Offset for Muon Pairs", 300, 0, 60);
 487   fHistMassMuonPairs           = new TH1D("fHistMassMuonPairs",           "Dimuon Mass (MUON+MFT)",  fNMassBins, fMassMin, fMassMax);
 488   fHistMassMuonPairsWithoutMFT = new TH1D("fHistMassMuonPairsWithoutMFT", "Dimuon Mass (MUON only)", fNMassBins, fMassMin, fMassMax);
 489   fHistMassMuonPairsMC         = new TH1D("fHistMassMuonPairsMC",         "Dimuon Mass (MC)",        fNMassBins, fMassMin, fMassMax);
 490   fHistRapidityPtMuonPairsMC   = new TH2D("fHistRapidityPtMuonPairsMC",   "Dimuon Phase Space (MC)", 100, -4.5, -2., 100, 0., 10.);
 491
 492   fHistWOffsetMuonPairs        -> SetXTitle("Weighted Offset");
 493   fHistMassMuonPairs           -> SetXTitle("Mass  [GeV/c^{2}]");
 494   fHistMassMuonPairsWithoutMFT -> SetXTitle("Mass  [GeV/c^{2}]");
 495   fHistMassMuonPairsMC         -> SetXTitle("Mass  [GeV/c^{2}]");
 496   fHistRapidityPtMuonPairsMC   -> SetXTitle("Rapidity");
 497   fHistRapidityPtMuonPairsMC   -> SetYTitle("p_{T}  [GeV/c]");
 498
 499   fHistWOffsetMuonPairs        -> Sumw2();
 500   fHistMassMuonPairs           -> Sumw2();
 501   fHistMassMuonPairsWithoutMFT -> Sumw2();
 502   fHistMassMuonPairsMC         -> Sumw2();
 503   fHistRapidityPtMuonPairsMC   -> Sumw2();
 504
 505 }
 506
 507 //====================================================================================================================================================
 508