FMD/analysis/AliFMDAnalysisTaskDndeta.cxx

   1
   2 #include <TROOT.h>
   3 #include <TSystem.h>
   4 #include <TInterpreter.h>
   5 #include <TChain.h>
   6 #include <TFile.h>
   7 #include <TList.h>
   8 #include <iostream>
   9 #include "TH1F.h"
  10 #include "TH2F.h"
  11 #include "AliFMDAnalysisTaskDndeta.h"
  12 #include "AliAnalysisManager.h"
  13 #include "AliESDFMD.h"
  14 #include "AliESDEvent.h"
  15 #include "AliAODEvent.h"
  16 #include "AliAODHandler.h"
  17 #include "AliMCEventHandler.h"
  18 #include "AliStack.h"
  19 #include "AliLog.h"
  20 #include "AliESDVertex.h"
  21 #include "TMath.h"
  22 #include "AliFMDAnaParameters.h"
  23 #include "AliFMDGeometry.h"
  24 #include "AliGenEventHeader.h"
  25 #include "TDatabasePDG.h"
  26 #include "TParticlePDG.h"
  27 ClassImp(AliFMDAnalysisTaskDndeta)
  28
  29
  30 AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta()
  31 : fDebug(0),
  32   fOutputList(0),
  33   fInputList(0),
  34   fArray(0),
  35   fInputArray(0),
  36   fVertexString(0x0),
  37   fNevents(),
  38   fNMCevents(),
  39   fStandalone(kTRUE),
  40   fMCevent(0)
  41 {
  42   // Default constructor
  43   DefineInput (0, TList::Class());
  44   DefineOutput(0, TList::Class());
  45 }
  46 //_____________________________________________________________________
  47 AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta(const char* name, Bool_t SE):
  48     AliAnalysisTask(name, "Density"),
  49     fDebug(0),
  50     fOutputList(0),
  51     fInputList(0),
  52     fArray(),
  53     fInputArray(0),
  54     fVertexString(0x0),
  55     fNevents(),
  56     fNMCevents(),
  57     fStandalone(kTRUE),
  58     fMCevent(0)
  59 {
  60   fStandalone = SE;
  61   if(fStandalone) {
  62     DefineInput (0, TList::Class());
  63     DefineInput(1, TObjString::Class());
  64     DefineOutput(0, TList::Class());
  65
  66   }
  67 }
  68 //_____________________________________________________________________
  69 void AliFMDAnalysisTaskDndeta::CreateOutputObjects()
  70 {
  71   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
  72
  73   fArray.SetName("FMD");
  74   fArray.SetOwner();
  75
  76   if(!fOutputList)
  77     fOutputList = new TList();
  78   fOutputList->SetName("BackgroundCorrected");
  79
  80
  81   TH2F* hMult = 0;
  82   TH1F* hPrimVertexBin = 0;
  83
  84
  85   TH2F* hBg = pars->GetBackgroundCorrection(1, 'I', 0);
  86   TH1F* hPrimary = new TH1F("hMultvsEta","hMultvsEta",
  87                             hBg->GetNbinsX(),
  88                             hBg->GetXaxis()->GetXmin(),
  89                             hBg->GetXaxis()->GetXmax());
  90   hPrimary->Sumw2();
  91   fOutputList->Add(hPrimary);
  92   Int_t nVtxbins = pars->GetNvtxBins();
  93
  94
  95   for(Int_t det =1; det<=3;det++)
  96     {
  97       TObjArray* detArray = new TObjArray();
  98       detArray->SetName(Form("FMD%d",det));
  99       fArray.AddAtAndExpand(detArray,det);
 100       Int_t nRings = (det==1 ? 1 : 2);
 101       for(Int_t ring = 0;ring<nRings;ring++)
 102         {
 103           Char_t ringChar = (ring == 0 ? 'I' : 'O');
 104           Int_t  nSec     = (ring == 0 ? 20 : 40);
 105
 106           TObjArray* vtxArray = new TObjArray();
 107           vtxArray->SetName(Form("FMD%d%c",det,ringChar));
 108           detArray->AddAtAndExpand(vtxArray,ring);
 109           for(Int_t i = 0; i< nVtxbins; i++) {
 110             TH2F* hBg = pars->GetBackgroundCorrection(det, ringChar, i);
 111             hMult  = new TH2F(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),
 112                               hBg->GetNbinsX(),
 113                               hBg->GetXaxis()->GetXmin(),
 114                               hBg->GetXaxis()->GetXmax(),
 115                               nSec, 0, 2*TMath::Pi());
 116             hMult->Sumw2();
 117             fOutputList->Add(hMult);
 118             vtxArray->AddAtAndExpand(hMult,i);
 119
 120           }
 121         }
 122     }
 123
 124   for(Int_t i = 0; i< nVtxbins; i++) {
 125
 126     hPrimVertexBin = new TH1F(Form("primmult_vtxbin%d",i),
 127                               Form("primmult_vtxbin%d",i),
 128                               hBg->GetNbinsX(),
 129                               hBg->GetXaxis()->GetXmin(),
 130                               hBg->GetXaxis()->GetXmax());
 131     hPrimVertexBin->Sumw2();
 132     fOutputList->Add(hPrimVertexBin);
 133
 134   }
 135
 136   fNevents.SetBins(nVtxbins,0,nVtxbins);
 137   fNevents.SetName("nEvents");
 138   fNMCevents.SetBins(nVtxbins,0,nVtxbins);
 139   fNMCevents.SetName("nMCEvents");
 140   fOutputList->Add(&fNevents);
 141   fOutputList->Add(&fNMCevents);
 142
 143 }
 144 //_____________________________________________________________________
 145 void AliFMDAnalysisTaskDndeta::ConnectInputData(Option_t */*option*/)
 146 {
 147   if(fStandalone) {
 148     fInputList   = (TList*)GetInputData(0);
 149     fVertexString = (TObjString*)GetInputData(1);
 150   }
 151 }
 152 //_____________________________________________________________________
 153 void AliFMDAnalysisTaskDndeta::Exec(Option_t */*option*/)
 154 {
 155   Int_t vtxbin   = fVertexString->GetString().Atoi();
 156   fNevents.Fill(vtxbin);
 157   for(UShort_t det=1;det<=3;det++) {
 158     //TObjArray* detInputArray = (TObjArray*)fInputArray->At(det);
 159     TObjArray* detArray = (TObjArray*)fArray.At(det);
 160     Int_t nRings = (det==1 ? 1 : 2);
 161     for (UShort_t ir = 0; ir < nRings; ir++) {
 162       Char_t ringChar = (ir == 0 ? 'I' : 'O');
 163       //TObjArray* vtxInputArray = (TObjArray*)detInputArray->At(ir);
 164       TObjArray* vtxArray = (TObjArray*)detArray->At(ir);
 165       TH2F* hMultTotal = (TH2F*)vtxArray->At(vtxbin);
 166
 167
 168       TH2F* hMultInput = (TH2F*)fInputList->FindObject(Form("mult_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
 169
 170       hMultTotal->Add(hMultInput);
 171
 172
 173     }
 174   }
 175
 176   if(fMCevent)
 177     ProcessPrimary();
 178
 179   if(fStandalone) {
 180     PostData(0, fOutputList);
 181   }
 182
 183 }
 184 //_____________________________________________________________________
 185 void AliFMDAnalysisTaskDndeta::Terminate(Option_t */*option*/) {
 186
 187
 188   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
 189
 190   Int_t nVtxbins = pars->GetNvtxBins();
 191
 192   for(UShort_t det=1;det<=3;det++) {
 193     TObjArray* detArray = (TObjArray*)fArray.At(det);
 194     Int_t nRings = (det==1 ? 1 : 2);
 195     for (UShort_t ir = 0; ir < nRings; ir++) {
 196       TObjArray* vtxArray = (TObjArray*)detArray->At(ir);
 197       Char_t ringChar = (ir == 0 ? 'I' : 'O');
 198       for(Int_t i =0; i<nVtxbins; i++) {
 199         TH2F* hMultTotal = (TH2F*)vtxArray->At(i);
 200         TH1D* hMultProj   = hMultTotal->ProjectionX(Form("dNdeta_FMD%d%c_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY());
 201         fOutputList->Add(hMultProj);
 202       }
 203     }
 204   }
 205
 206 }
 207 //_____________________________________________________________________
 208 void AliFMDAnalysisTaskDndeta::ProcessPrimary() {
 209
 210   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
 211
 212   AliMCParticle* particle = 0;
 213   AliStack* stack = fMCevent->Stack();
 214
 215   TH1F* hPrimary = (TH1F*)fOutputList->FindObject("hMultvsEta");
 216
 217   Bool_t firstTrack = kTRUE;
 218   Int_t nTracks = fMCevent->GetNumberOfTracks();
 219   for(Int_t i = 0 ;i<nTracks;i++) {
 220     particle = fMCevent->GetTrack(i);
 221     if(!particle)
 222       continue;
 223     if(TMath::Abs(particle->Zv()) > pars->GetVtxCutZ())
 224       continue;
 225
 226     if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
 227       hPrimary->Fill(particle->Eta());
 228       Double_t delta           = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
 229       Double_t vertexBinDouble = (particle->Zv() + pars->GetVtxCutZ()) / delta;
 230       Int_t    vertexBin       = (Int_t)vertexBinDouble;
 231
 232       TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin));
 233       hPrimVtxBin->Fill(particle->Eta());
 234       if(firstTrack) {
 235         fNMCevents.Fill(vertexBin);
 236         firstTrack = kFALSE;
 237       }
 238     }
 239
 240   }
 241
 242
 243
 244
 245
 246 }
 247 //_____________________________________________________________________
 248 //
 249 //
 250 // EOF