PWGCF/Correlations/DPhi/AliAnalysisTaskContMC.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2009, 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 //         AliAnalysisTaskContMC class
  18 //-----------------------------------------------------------------
  19
  20 #include "TChain.h"
  21 #include "TTree.h"
  22 #include "TLegend.h"
  23 #include "TH1F.h"
  24 #include "TH2F.h"
  25 #include "TH3F.h"
  26 #include "TCanvas.h"
  27 #include "AliAnalysisTask.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliVTrack.h"
  30 #include "AliAODMCParticle.h"
  31 #include "AliVParticle.h"
  32 #include "AliAODEvent.h"
  33 #include "AliAODInputHandler.h"
  34 #include "AliAnalysisTaskContMC.h"
  35 #include "AliAnalysisTaskESDfilter.h"
  36 #include "AliAnalysisDataContainer.h"
  37 #include "AliHelperPID.h"
  38 #include "AliCentrality.h"
  39 #include "TProof.h"
  40 #include "AliPID.h"
  41 #include "AliVEvent.h"
  42 #include "AliPIDResponse.h"
  43 #include "AliStack.h"
  44 #include <TMCProcess.h>
  45
  46 #include <iostream>
  47
  48 using namespace std;
  49
  50 ClassImp(AliAnalysisTaskContMC)
  51
  52 //________________________________________________________________________
  53 AliAnalysisTaskContMC::AliAnalysisTaskContMC(const char *name) : AliAnalysisTaskSE(name), fAOD(0), fNSigmaPID(0), fIsMC(0), fOutput(0), fHistID(0)
  54 {
  55   // Default constructor
  56
  57
  58   DefineInput(0, TChain::Class());
  59   //DefineOutput(1, AliHelperPID::Class());
  60   DefineOutput(1, TList::Class());
  61
  62 }
  63 //________________________________________________________________________
  64 //________________________________________________________________________
  65 void AliAnalysisTaskContMC::UserCreateOutputObjects()
  66 {
  67   Printf("\n\n\n\n\n\n In CreateOutput Object:");
  68
  69   //create output objects
  70   Printf("NSigma: %.1f",fNSigmaPID->GetNSigmaCut());
  71   Printf("IsMC: %d",fNSigmaPID->GetisMC());
  72   fOutput = new TList();
  73   fOutput->SetOwner();
  74   fOutput->SetName("list");
  75
  76   fHistID = new TH3F("fHistID","fHistID",6,-1.5,4.5,4,-1.5,2.5,38,0.2,4);
  77   fOutput->Add(fHistID);
  78
  79   if(!fNSigmaPID)AliFatal("PID object should be set in the steering macro");
  80   fOutput->Add(fNSigmaPID);
  81
  82   //PostData(1, fNSigmaPID  );
  83   PostData(1, fOutput  );
  84
  85 }
  86
  87 //________________________________________________________________________
  88 void AliAnalysisTaskContMC::UserExec(Option_t *)
  89 {
  90   const Int_t npart=4;
  91   const Int_t pdgcode[npart+2]={-1,211,321,2212,11,13};
  92   const Int_t partid[npart+2]={-1,0,1,2,3,4};
  93   // main event loop
  94   fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);
  95   if (!fAOD) {
  96     AliWarning("ERROR: AliAODEvent not available \n");
  97     return;
  98   }
  99
 100   if (strcmp(fAOD->ClassName(), "AliAODEvent"))
 101     {
 102       AliFatal("Not processing AODs");
 103     }
 104   //MC Loop
 105   TClonesArray *arrayMC = 0;
 106   Printf("fIsMC: %d",fIsMC);
 107   if (fIsMC)
 108     {
 109       arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 110       if (!arrayMC) {
 111         AliFatal("Error: MC particles branch not found!\n");
 112       }
 113     }
 114   Double_t centrality = 0;
 115
 116   AliCentrality *centralityObj = ((AliVAODHeader*)fAOD->GetHeader())->GetCentralityP();
 117   if (centralityObj)
 118     {
 119       centrality = centralityObj->GetCentralityPercentileUnchecked("V0M");
 120       AliInfo(Form("Centrality is %f", centrality));
 121     }
 122   else
 123     {
 124       Printf("WARNING: Centrality object is 0");
 125       centrality = -1;
 126     }
 127   if (centrality < 0)
 128     return;
 129
 130
 131   //vertex selection
 132   Int_t nVertex = ((AliAODEvent*)fAOD)->GetNumberOfVertices();
 133   if( nVertex > 0 ) {
 134     AliAODVertex* vertex = (AliAODVertex*)((AliAODEvent*)fAOD)->GetPrimaryVertex();
 135     //Int_t nTracksPrim = vertex->GetNContributors();
 136     Double_t zVertex = vertex->GetZ();
 137     //10 cm cut
 138     if(TMath::Abs(zVertex)>10) return;
 139     //AliInfo(Form(" Vertex in = %f with %d particles by  %s data ...",zVertex,nTracksPrim,vertex->GetName()));
 140     // Reject TPC only vertex
 141     TString name(vertex->GetName());
 142     if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex"))return;
 143   }
 144
 145   //Int_t count=0;
 146   //track loop
 147   for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
 148     AliAODTrack* track = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iTracks));
 149     if(!track) AliFatal("Not a standard AOD");
 150     if(!(track->TestFilterBit(32)))continue;
 151     if (TMath::Abs(track->Eta()) > .8 || track->Pt() < .2) continue;
 152
 153     Int_t pdg=-999;
 154     Int_t isph=-999;
 155     if (fIsMC && arrayMC) {
 156       AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(TMath::Abs(track->GetLabel()));
 157       if (!partMC) {
 158         AliError("Cannot get MC particle");
 159         continue;
 160       }
 161       pdg=TMath::Abs(partMC->GetPdgCode());
 162       isph=partMC->IsPhysicalPrimary();
 163     }
 164     if(!isph)continue;
 165     //step 1, TOF Matching
 166     UInt_t status;
 167     status=track->GetStatus();
 168     if((status&AliVTrack::kTOFout)==0 || (status&AliVTrack::kTIME)==0)continue;
 169
 170     //step 2, combined PID
 171     Int_t IDTPCTOF=fNSigmaPID->GetParticleSpecies(track,1);
 172     if(IDTPCTOF==999)IDTPCTOF=-1;
 173     Int_t IDMC=-1;
 174     for(Int_t ipart=0;ipart<npart+2;ipart++)if(TMath::Abs(pdg)==pdgcode[ipart])IDMC=partid[ipart];
 175     fHistID->Fill(IDMC,IDTPCTOF,track->Pt());
 176   } // end loop on tracks
 177
 178   PostData(1,fOutput);
 179
 180 }
 181
 182 //_________________________________________________________________
 183 void   AliAnalysisTaskContMC::Terminate(Option_t *)
 184 {
 185   // Terminate analysis
 186   //
 187   fOutput = dynamic_cast<TList*>(GetOutputData(1));
 188   if (!fOutput) {
 189     printf("ERROR: fOutput not available\n");
 190     return;
 191   }
 192   fHistID = dynamic_cast<TH3F*>(fOutput->FindObject("fHistID"));
 193   fNSigmaPID = dynamic_cast<AliHelperPID*>(fOutput->FindObject("fNSigmaPID"));
 194 }