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