All FMD corrections are now divided into the analysis + adding new corrections
[u/mrichter/AliRoot.git] / FMD / analysis / AliFMDAnalysisTaskDndeta.cxx
CommitLineData
7c3e5162 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>
bb8a464f 9#include "TH1F.h"
7c3e5162 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"
78f6f750 23//#include "AliFMDGeometry.h"
bb8a464f 24#include "AliGenEventHeader.h"
cc066cb9 25#include "AliHeader.h"
78f6f750 26//#include "TDatabasePDG.h"
27//#include "TParticlePDG.h"
cc066cb9 28#include "AliFMDStripIndex.h"
7c3e5162 29ClassImp(AliFMDAnalysisTaskDndeta)
30
31
32AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta()
33: fDebug(0),
34 fOutputList(0),
35 fInputList(0),
7c3e5162 36 fVertexString(0x0),
37 fNevents(),
bb8a464f 38 fNMCevents(),
39 fStandalone(kTRUE),
b64db9b1 40 fLastTrackByStrip(0)
7c3e5162 41{
42 // Default constructor
43 DefineInput (0, TList::Class());
44 DefineOutput(0, TList::Class());
45}
46//_____________________________________________________________________
47AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta(const char* name, Bool_t SE):
48 AliAnalysisTask(name, "Density"),
49 fDebug(0),
50 fOutputList(0),
51 fInputList(0),
7c3e5162 52 fVertexString(0x0),
53 fNevents(),
bb8a464f 54 fNMCevents(),
55 fStandalone(kTRUE),
b64db9b1 56 fLastTrackByStrip(0)
7c3e5162 57{
58 fStandalone = SE;
59 if(fStandalone) {
60 DefineInput (0, TList::Class());
61 DefineInput(1, TObjString::Class());
62 DefineOutput(0, TList::Class());
63
64 }
65}
66//_____________________________________________________________________
67void AliFMDAnalysisTaskDndeta::CreateOutputObjects()
68{
69 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
70
7c3e5162 71 if(!fOutputList)
72 fOutputList = new TList();
73 fOutputList->SetName("BackgroundCorrected");
74
75
76 TH2F* hMult = 0;
cc066cb9 77 TH1F* hHits = 0;
bb8a464f 78 TH1F* hPrimVertexBin = 0;
7c3e5162 79
bb8a464f 80
5ea1e0a9 81 TH2F* hBgTmp = pars->GetBackgroundCorrection(1, 'I', 0);
bb8a464f 82 TH1F* hPrimary = new TH1F("hMultvsEta","hMultvsEta",
5ea1e0a9 83 hBgTmp->GetNbinsX(),
84 hBgTmp->GetXaxis()->GetXmin(),
85 hBgTmp->GetXaxis()->GetXmax());
bb8a464f 86 hPrimary->Sumw2();
87 fOutputList->Add(hPrimary);
7c3e5162 88 Int_t nVtxbins = pars->GetNvtxBins();
b85ea106 89 TH2F* hBg = 0;
90 for(Int_t i = 0; i< nVtxbins; i++) {
91
92 for(Int_t det =1; det<=3;det++)
93 {
94 Int_t nRings = (det==1 ? 1 : 2);
95 for(Int_t ring = 0;ring<nRings;ring++)
96 {
97 Char_t ringChar = (ring == 0 ? 'I' : 'O');
98 Int_t nSec = (ring == 0 ? 20 : 40);
99
100
101
102 hBg = pars->GetBackgroundCorrection(det, ringChar, i);
7c3e5162 103 hMult = new TH2F(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),
104 hBg->GetNbinsX(),
105 hBg->GetXaxis()->GetXmin(),
106 hBg->GetXaxis()->GetXmax(),
107 nSec, 0, 2*TMath::Pi());
cc066cb9 108
18d4b9ae 109 hHits = new TH1F(Form("hMCHits_FMD%d%c_vtxbin%d",det,ringChar,i),Form("hMCHits_FMD%d%c_vtxbin%d",det,ringChar,i),
cc066cb9 110 hBg->GetNbinsX(),
111 hBg->GetXaxis()->GetXmin(),
112 hBg->GetXaxis()->GetXmax());
b85ea106 113 hHits->Sumw2();
114 fOutputList->Add(hHits);
cc066cb9 115
7c3e5162 116 hMult->Sumw2();
117 fOutputList->Add(hMult);
b85ea106 118
7c3e5162 119 }
b85ea106 120 }
121 }
7c3e5162 122
bb8a464f 123 for(Int_t i = 0; i< nVtxbins; i++) {
124
125 hPrimVertexBin = new TH1F(Form("primmult_vtxbin%d",i),
126 Form("primmult_vtxbin%d",i),
5ea1e0a9 127 hBgTmp->GetNbinsX(),
128 hBgTmp->GetXaxis()->GetXmin(),
129 hBgTmp->GetXaxis()->GetXmax());
bb8a464f 130 hPrimVertexBin->Sumw2();
131 fOutputList->Add(hPrimVertexBin);
132
133 }
134
7c3e5162 135 fNevents.SetBins(nVtxbins,0,nVtxbins);
136 fNevents.SetName("nEvents");
bb8a464f 137 fNMCevents.SetBins(nVtxbins,0,nVtxbins);
138 fNMCevents.SetName("nMCEvents");
7c3e5162 139 fOutputList->Add(&fNevents);
bb8a464f 140 fOutputList->Add(&fNMCevents);
7c3e5162 141
142}
143//_____________________________________________________________________
144void AliFMDAnalysisTaskDndeta::ConnectInputData(Option_t */*option*/)
145{
146 if(fStandalone) {
147 fInputList = (TList*)GetInputData(0);
7c3e5162 148 }
149}
150//_____________________________________________________________________
151void AliFMDAnalysisTaskDndeta::Exec(Option_t */*option*/)
152{
b64db9b1 153 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
9f55be54 154 fVertexString = (TObjString*)fInputList->At(0);
155
7c3e5162 156 Int_t vtxbin = fVertexString->GetString().Atoi();
9f55be54 157
7c3e5162 158 fNevents.Fill(vtxbin);
9f55be54 159
7c3e5162 160 for(UShort_t det=1;det<=3;det++) {
7c3e5162 161 Int_t nRings = (det==1 ? 1 : 2);
162 for (UShort_t ir = 0; ir < nRings; ir++) {
163 Char_t ringChar = (ir == 0 ? 'I' : 'O');
9f55be54 164
b64db9b1 165 TH2F* hMultTotal = (TH2F*)fOutputList->FindObject(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
166
7c3e5162 167
168 TH2F* hMultInput = (TH2F*)fInputList->FindObject(Form("mult_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
1282ce49 169
7c3e5162 170 hMultTotal->Add(hMultInput);
9f55be54 171
7c3e5162 172 }
173 }
bb8a464f 174
b64db9b1 175 if(pars->GetProcessPrimary())
bb8a464f 176 ProcessPrimary();
177
7c3e5162 178 if(fStandalone) {
179 PostData(0, fOutputList);
180 }
181
182}
183//_____________________________________________________________________
184void AliFMDAnalysisTaskDndeta::Terminate(Option_t */*option*/) {
185
1282ce49 186
7c3e5162 187 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
188
189 Int_t nVtxbins = pars->GetNvtxBins();
190
191 for(UShort_t det=1;det<=3;det++) {
7c3e5162 192 Int_t nRings = (det==1 ? 1 : 2);
193 for (UShort_t ir = 0; ir < nRings; ir++) {
1282ce49 194 Char_t ringChar = (ir == 0 ? 'I' : 'O');
7c3e5162 195 for(Int_t i =0; i<nVtxbins; i++) {
b64db9b1 196
b64db9b1 197 TH2F* hMultTotal = (TH2F*)fOutputList->FindObject(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i));
f1eff0e8 198 TH2F* hMultTrVtx = (TH2F*)hMultTotal->Clone(Form("dNdeta_FMD%d%c_TrVtx_vtxbin%d",det,ringChar,i));
b64db9b1 199
b85ea106 200 hMultTotal->Scale(1/pars->GetEventSelectionEfficiency(i));
201
1282ce49 202 TH1D* hMultProj = hMultTotal->ProjectionX(Form("dNdeta_FMD%d%c_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY());
f1eff0e8 203 TH1D* hMultProjTrVtx = hMultTrVtx->ProjectionX(Form("dNdeta_FMD%d%c_TrVtx_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY());
b64db9b1 204 fOutputList->Add(hMultTrVtx);
1282ce49 205 fOutputList->Add(hMultProj);
b64db9b1 206 fOutputList->Add(hMultProjTrVtx);
7c3e5162 207 }
208 }
209 }
1282ce49 210
bb8a464f 211}
212//_____________________________________________________________________
213void AliFMDAnalysisTaskDndeta::ProcessPrimary() {
cc066cb9 214
b64db9b1 215 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
216 AliMCEvent* mcEvent = eventHandler->MCEvent();
217 if(!mcEvent)
218 return;
219
cc066cb9 220 fLastTrackByStrip.Reset(-1);
221
bb8a464f 222 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
223
224 AliMCParticle* particle = 0;
b64db9b1 225 AliStack* stack = mcEvent->Stack();
bb8a464f 226
227 TH1F* hPrimary = (TH1F*)fOutputList->FindObject("hMultvsEta");
b64db9b1 228 AliHeader* header = mcEvent->Header();
cc066cb9 229 AliGenEventHeader* genHeader = header->GenEventHeader();
230
231 TArrayF vertex;
232 genHeader->PrimaryVertex(vertex);
233 if(TMath::Abs(vertex.At(2)) > pars->GetVtxCutZ())
234 return;
235 Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
236 Double_t vertexBinDouble = (vertex.At(2) + pars->GetVtxCutZ()) / delta;
237 Int_t vertexBin = (Int_t)vertexBinDouble;
bb8a464f 238
239 Bool_t firstTrack = kTRUE;
bb0a45c3 240
241 // we loop over the primaries only unless we need the hits (diagnostics running slowly)
242 Int_t nTracks = stack->GetNprimary();
b64db9b1 243 if(pars->GetProcessHits())
bb0a45c3 244 nTracks = stack->GetNtrack();
245
bb8a464f 246 for(Int_t i = 0 ;i<nTracks;i++) {
7aad0c47 247 particle = (AliMCParticle*) mcEvent->GetTrack(i);
bb8a464f 248 if(!particle)
249 continue;
bb0a45c3 250
bb8a464f 251 if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
252 hPrimary->Fill(particle->Eta());
cc066cb9 253
6289b3e8 254
bb8a464f 255 TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin));
256 hPrimVtxBin->Fill(particle->Eta());
257 if(firstTrack) {
258 fNMCevents.Fill(vertexBin);
259 firstTrack = kFALSE;
260 }
b3546e91 261
bb0a45c3 262 }
f1eff0e8 263 if(pars->GetProcessHits()) {
b85ea106 264
bb0a45c3 265 for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) {
cc066cb9 266
55cadbf9 267 AliTrackReference* ref = particle->GetTrackReference(j);
268 UShort_t det,sec,strip;
269 Char_t ring;
270 if(ref->DetectorId() != AliTrackReference::kFMD)
271 continue;
272 AliFMDStripIndex::Unpack(ref->UserId(),det,ring,sec,strip);
273 Float_t thisStripTrack = fLastTrackByStrip(det,ring,sec,strip);
274 if(particle->Charge() != 0 && i != thisStripTrack ) {
275 //Double_t x,y,z;
276
277 Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));//-1*TMath::Log(TMath::Tan(0.5*theta));
18d4b9ae 278 TH1F* hHits = (TH1F*)fOutputList->FindObject(Form("hMCHits_FMD%d%c_vtxbin%d",det,ring,vertexBin));
b85ea106 279
280
55cadbf9 281 hHits->Fill(eta);
b85ea106 282
55cadbf9 283 Float_t nstrips = (ring =='O' ? 256 : 512);
284
285 fLastTrackByStrip(det,ring,sec,strip) = (Float_t)i;
cc066cb9 286
55cadbf9 287 if(strip >0)
288 fLastTrackByStrip(det,ring,sec,strip-1) = (Float_t)i;
289 if(strip < (nstrips - 1))
290 fLastTrackByStrip(det,ring,sec,strip+1) = (Float_t)i;
291
292 }
b85ea106 293
294
b3546e91 295 }
b85ea106 296
297
bb0a45c3 298 }
b85ea106 299
300
301
bb0a45c3 302 }
7c3e5162 303}
304//_____________________________________________________________________
305//
306//
307// EOF