PWG3/AliAnalysisTaskSingleMu.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //----------------------------------------------------------------------------
  17 //    Implementation of the single muon analysis class
  18 // An example of usage can be found in the macro runSingleMuAnalysis.C.
  19 //----------------------------------------------------------------------------
  20
  21 #define AliAnalysisTaskSingleMu_cxx
  22
  23 // ROOT includes
  24 #include "Riostream.h"
  25 #include "TChain.h"
  26 #include "TH1.h"
  27 #include "TCanvas.h"
  28 #include "TSystem.h"
  29 #include "TROOT.h"
  30 #include "TParticle.h"
  31 #include "TLorentzVector.h"
  32
  33 // STEER includes
  34 #include "AliLog.h"
  35
  36 #include "AliAODEvent.h"
  37 #include "AliAODTrack.h"
  38 #include "AliAODVertex.h"
  39 #include "AliAODInputHandler.h"
  40
  41 #include "AliAnalysisTask.h"
  42 #include "AliAnalysisDataSlot.h"
  43 #include "AliAnalysisManager.h"
  44 #include "AliAnalysisTaskSingleMu.h"
  45
  46 ClassImp(AliAnalysisTaskSingleMu)
  47
  48 //________________________________________________________________________
  49 AliAnalysisTaskSingleMu::AliAnalysisTaskSingleMu(const char *name) :
  50   AliAnalysisTask(name,""),
  51   fAOD(0),
  52   fOutputContainer(0)
  53 {
  54   //
  55   /// Constructor.
  56   //
  57   ResetHistos();
  58   // Input slot #0 works with an Ntuple
  59   DefineInput(0, TChain::Class());
  60   // Output slot #0 writes into a TObjArray container
  61   DefineOutput(0,  TObjArray::Class());
  62 }
  63
  64 //___________________________________________________________________________
  65 void AliAnalysisTaskSingleMu::ConnectInputData(Option_t *)
  66 {
  67   //
  68   /// Connect AOD here
  69   /// Called once
  70   //
  71
  72   TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
  73   if (!tree) {
  74     Printf("ERROR: Could not read chain from input slot 0");
  75   } else {
  76     /*
  77     // Disable all branches and enable only the needed ones
  78     // The next two lines are different when data produced as AliAODEvent is read
  79     tree->SetBranchStatus("*", kFALSE);
  80     tree->SetBranchStatus("fTracks.*", kTRUE);
  81     */
  82
  83     AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
  84
  85     if (!aodH) {
  86       Printf("ERROR: Could not get AODInputHandler");
  87     } else
  88       printf("   ConnectInputData of task %s\n", GetName());
  89       fAOD = aodH->GetEvent();
  90   }
  91 }
  92
  93 //___________________________________________________________________________
  94 void AliAnalysisTaskSingleMu::CreateOutputObjects()
  95 {
  96   //
  97   /// Create output histograms
  98   //
  99   printf("   CreateOutputObjects of task %s\n", GetName());
 100
 101   Int_t ptBins = 60;
 102   Float_t ptLow = 0., ptHigh = 30.;
 103   Char_t *ptName = "P_{t} (GeV/c)";
 104
 105   Int_t vzBins = 40;
 106   Float_t vzLow = -20., vzHigh = 20.;
 107   Char_t *vzName = "Vz (cm)";
 108
 109   TString baseName, histoName;
 110   fOutputContainer = new TObjArray(fgkNhistos*fgkNTrigCuts);
 111   fOutputContainer->SetName("SingleMuAnalysisContainer");
 112   Int_t iHisto = 0;
 113
 114   for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){
 115
 116     // 2D histos
 117     if(!fVzVsPt[iTrig]){
 118       baseName = "fVzVsPt";
 119       histoName = baseName + trigName[iTrig];
 120       fVzVsPt[iTrig] = new TH2F(histoName, histoName,
 121                                 ptBins, ptLow, ptHigh,
 122                                 vzBins, vzLow, vzHigh);
 123       fVzVsPt[iTrig]->GetXaxis()->SetTitle(ptName);
 124       fVzVsPt[iTrig]->GetYaxis()->SetTitle(vzName);
 125
 126       fOutputContainer->AddAt(fVzVsPt[iTrig], iHisto);
 127       iHisto++;
 128     }
 129   }
 130 }
 131
 132 //________________________________________________________________________
 133 void AliAnalysisTaskSingleMu::Exec(Option_t *)
 134 {
 135   //
 136   /// Main loop
 137   /// Called for each event
 138   //
 139
 140   TTree *tinput = (TTree*)GetInputData(0);
 141   tinput->GetReadEntry();
 142
 143   if (!fAOD) {
 144     Printf("ERROR: fAOD not available");
 145     return;
 146   }
 147
 148
 149   // Object declaration
 150   AliAODTrack *muonTrack = 0x0;
 151   TLorentzVector lorVec;
 152   Int_t trigMatch = -1;
 153
 154   Int_t nTracks = fAOD->GetNumberOfTracks();
 155   for (Int_t itrack = 0; itrack < nTracks; itrack++) {
 156     muonTrack = fAOD->GetTrack(itrack);
 157
 158     // Apply cuts
 159     if(!MuonPassesCuts(*muonTrack, lorVec, trigMatch)) continue;
 160
 161     for(Int_t iTrig=0; iTrig<=trigMatch; iTrig++){
 162       fVzVsPt[iTrig]->Fill(lorVec.Pt(), fAOD->GetPrimaryVertex()->GetZ());
 163     }
 164   }
 165
 166   // Post final data. It will be written to a file with option "RECREATE"
 167   PostData(0, fOutputContainer);
 168 }
 169
 170 //________________________________________________________________________
 171 void AliAnalysisTaskSingleMu::Terminate(Option_t *) {
 172   //
 173   /// Draw some histogram at the end.
 174   //
 175   if (!gROOT->IsBatch()) {
 176     TCanvas *c1 = new TCanvas("c1","Vz vs Pt",10,10,310,310);
 177     c1->SetFillColor(10); c1->SetHighLightColor(10);
 178     c1->SetLeftMargin(0.15); c1->SetBottomMargin(0.15);
 179     c1->Divide(2,2);
 180     for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){
 181       c1->cd(iTrig+1);
 182       fVzVsPt[iTrig]->DrawCopy("COLZ");
 183     }
 184   }
 185 }
 186
 187 //________________________________________________________________________
 188 void AliAnalysisTaskSingleMu::ResetHistos()
 189 {
 190   //
 191   /// Reset histograms
 192   //
 193   for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){
 194     fVzVsPt[iTrig] = 0x0;
 195   }
 196   trigName[kNoMatchTrig] = "NoMatchTrig";
 197   trigName[kAllPtTrig]   = "AllPtTrig";
 198   trigName[kLowPtTrig]   = "LowPtTrig";
 199   trigName[kHighPtTrig]  = "HighPtTrig";
 200 }
 201
 202
 203 //________________________________________________________________________
 204 Bool_t AliAnalysisTaskSingleMu::MuonPassesCuts(AliAODTrack &muonTrack,
 205                                                TLorentzVector &lorVec,
 206                                                Int_t &trigMatch)
 207 {
 208   //
 209   /// Fill lorentz vector and check cuts
 210   //
 211
 212   // Check if track is a muon
 213   if(muonTrack.GetMostProbablePID()!=AliAODTrack::kMuon) return kFALSE;
 214
 215   // Check if track is triggered
 216   trigMatch = kNoMatchTrig;
 217   if (muonTrack.MatchTriggerHighPt()) {
 218     trigMatch = kHighPtTrig;
 219   } else if (muonTrack.MatchTriggerLowPt()) {
 220     trigMatch = kLowPtTrig;
 221   } else if (muonTrack.MatchTriggerAnyPt()){
 222     trigMatch = kAllPtTrig;
 223   }
 224
 225   // Fill track parameters
 226   Double_t px = muonTrack.Px();
 227   Double_t py = muonTrack.Py();
 228   Double_t pz = muonTrack.Pz();
 229   Double_t p  = muonTrack.P();
 230
 231   const Double_t kMuonMass = 0.105658369;
 232
 233   Double_t energy = TMath::Sqrt(p*p + kMuonMass*kMuonMass);
 234   lorVec.SetPxPyPzE(px,py,pz,energy);
 235
 236   return kTRUE;
 237 }