]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FORWARD/analysis/AliFMDAnalysisTaskGenerateCorrection.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding correction for FMD dead channels
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis / AliFMDAnalysisTaskGenerateCorrection.cxx
CommitLineData
d31dce33 1#include "AliFMDAnalysisTaskGenerateCorrection.h"
2#include "AliESDEvent.h"
3#include "iostream"
4#include "AliESDFMD.h"
5#include "TH2F.h"
6#include "AliTrackReference.h"
7#include "AliStack.h"
8#include "AliFMDAnaParameters.h"
9#include "AliFMDStripIndex.h"
10#include "AliStack.h"
11#include "AliMCParticle.h"
12#include "AliMCEvent.h"
13//#include "AliFMDGeometry.h"
14#include "TArray.h"
15#include "AliGenEventHeader.h"
04f1ff3d 16#include "AliMultiplicity.h"
d31dce33 17#include "AliHeader.h"
18#include "AliFMDAnaCalibBackgroundCorrection.h"
19#include "AliFMDAnaCalibEventSelectionEfficiency.h"
059c7c6b 20#include "AliGenPythiaEventHeader.h"
d31dce33 21//#include "AliCDBManager.h"
22//#include "AliCDBId.h"
23//#include "AliCDBMetaData.h"
24#include "TSystem.h"
25#include "TROOT.h"
26#include "TAxis.h"
507687cd 27#include "AliGenDPMjetEventHeader.h"
28
d31dce33 29ClassImp(AliFMDAnalysisTaskGenerateCorrection)
30
31//_____________________________________________________________________
32AliFMDAnalysisTaskGenerateCorrection::AliFMDAnalysisTaskGenerateCorrection():
33AliAnalysisTaskSE(),
34 fListOfHits(),
35 fListOfPrimaries(),
36 fListOfCorrection(),
37 fVertexBins(),
38 fLastTrackByStrip(0),
39 fHitsByStrip(0),
40 fZvtxCut(10),
41 fNvtxBins(10),
42 fNbinsEta(200),
9f55be54 43 fBackground(0),
44 fEventSelectionEff(0)
d31dce33 45{
46 // Default constructor
47}
48//_____________________________________________________________________
49AliFMDAnalysisTaskGenerateCorrection::AliFMDAnalysisTaskGenerateCorrection(const char* name):
50 AliAnalysisTaskSE(name),
51 fListOfHits(),
52 fListOfPrimaries(),
53 fListOfCorrection(),
54 fVertexBins(),
55 fLastTrackByStrip(0),
56 fHitsByStrip(0),
57 fZvtxCut(10),
58 fNvtxBins(10),
59 fNbinsEta(200),
9f55be54 60 fBackground(0),
61 fEventSelectionEff(0)
d31dce33 62{
63
64 DefineOutput(1, TList::Class());
65 DefineOutput(2, TList::Class());
66 DefineOutput(3, TH1F::Class());
67 DefineOutput(4, TList::Class());
68}
69//_____________________________________________________________________
70void AliFMDAnalysisTaskGenerateCorrection::UserCreateOutputObjects()
71{
72// Create the output containers
73//
74
75 std::cout<<"Creating output objects"<<std::endl;
507687cd 76
df2a9c32 77 // SelectCollisionCandidates(AliVEvent::kAny);
507687cd 78
04f1ff3d 79 for(Int_t v=0; v<fNvtxBins;v++) {
80
81 TH2F* hSPDhits = new TH2F(Form("hSPDhits_vtx%d",v),
82 Form("hSPDhits_vtx%d",v),
820f4fcc 83 fNbinsEta, -4,6, 20, 0,2*TMath::Pi());
04f1ff3d 84 hSPDhits->Sumw2();
85 fListOfHits.Add(hSPDhits);
86
f7356393 87 TH1F* hReadChannels = new TH1F(Form("hFMDReadChannels_vtx%d",v),
88 Form("hFMDReadChannels_vtx%d",v),
89 fNbinsEta,-4,6);
90 hReadChannels->Sumw2();
91 fListOfHits.Add(hReadChannels);
92 TH1F* hAllChannels = new TH1F(Form("hFMDAllChannels_vtx%d",v),
93 Form("hFMDAllChannels_vtx%d",v),
94 fNbinsEta,-4,6);
95 hAllChannels->Sumw2();
96 fListOfHits.Add(hAllChannels);
97
98
99
04f1ff3d 100 for(Int_t iring = 0; iring<2;iring++) {
101 Char_t ringChar = (iring == 0 ? 'I' : 'O');
102 Int_t nSec = (iring == 1 ? 40 : 20);
103
d31dce33 104 TH2F* hPrimary = new TH2F(Form("hPrimary_FMD_%c_vtx%d",ringChar,v),
105 Form("hPrimary_FMD_%c_vtx%d",ringChar,v),
820f4fcc 106 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
d31dce33 107 hPrimary->Sumw2();
108 fListOfPrimaries.Add(hPrimary);
507687cd 109 TH2F* hPrimaryNSD = new TH2F(Form("hPrimaryNSD_FMD_%c_vtx%d",ringChar,v),
110 Form("hPrimaryNSD_FMD_%c_vtx%d",ringChar,v),
111 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
112 hPrimaryNSD->Sumw2();
113 fListOfPrimaries.Add(hPrimaryNSD);
114
04f1ff3d 115
116
d31dce33 117 }
118 }
119
04f1ff3d 120
d31dce33 121
122
123 for(Int_t det =1; det<=3;det++) {
124 Int_t nRings = (det==1 ? 1 : 2);
125 for(Int_t ring = 0;ring<nRings;ring++) {
126 Int_t nSec = (ring == 1 ? 40 : 20);
127 Char_t ringChar = (ring == 0 ? 'I' : 'O');
128 TH1F* doubleHits = new TH1F(Form("DoubleHits_FMD%d%c",det,ringChar),
820f4fcc 129 Form("DoubleHits_FMD%d%c",det,ringChar),fNbinsEta, -4,6);
d31dce33 130 TH1F* allHits = new TH1F(Form("allHits_FMD%d%c",det,ringChar),
820f4fcc 131 Form("allHits_FMD%d%c",det,ringChar), fNbinsEta, -4,6);
d31dce33 132
133 doubleHits->Sumw2();
134 allHits->Sumw2();
135 fListOfHits.Add(allHits);
136 fListOfHits.Add(doubleHits);
137
138 for(Int_t v=0; v<fNvtxBins;v++) {
139 TH2F* hHits = new TH2F(Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
140 Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
820f4fcc 141 fNbinsEta, -4,6, nSec, 0, 2*TMath::Pi());
d31dce33 142 hHits->Sumw2();
143 fListOfHits.Add(hHits);
507687cd 144 TH2F* hHitsNSD = new TH2F(Form("hHitsNSD_FMD%d%c_vtx%d", det,ringChar,v),
145 Form("hHitsNSD_FMD%d%c_vtx%d", det,ringChar,v),
146 fNbinsEta, -4,6, nSec, 0, 2*TMath::Pi());
147 hHitsNSD->Sumw2();
148 fListOfHits.Add(hHitsNSD);
149
d31dce33 150
151 }
152 }
153 }
154
155 TH1F* hEventsSelected = new TH1F("EventsSelected","EventsSelected",fNvtxBins,0,fNvtxBins);
156 TH1F* hEventsAll = new TH1F("EventsAll","EventsAll",fNvtxBins,0,fNvtxBins);
059c7c6b 157 TH1F* hEventsAllNSD = new TH1F("EventsAllNSD","EventsAllNSD",fNvtxBins,0,fNvtxBins);
f55d559b 158 TH1F* hEventsSelectedVtx = new TH1F("EventsSelectedVtx","EventsSelectedVtx",fNvtxBins,0,fNvtxBins);
159 TH1F* hEventsSelectedTrigger = new TH1F("EventsSelectedTrigger","EventsSelectedTrigger",fNvtxBins,0,fNvtxBins);
f55d559b 160
059c7c6b 161 TH1F* hEventsSelectedNSDVtx = new TH1F("EventsSelectedNSDVtx","EventsSelectedNSDVtx",fNvtxBins,0,fNvtxBins);
162 TH1F* hEventsSelectedNSD = new TH1F("EventsSelectedNSD","EventsSelectedNSD",fNvtxBins,0,fNvtxBins);
163
850ad94d 164 TH1F* hXvtx = new TH1F("hXvtx","x vertex distribution",100,-2,2);
165 TH1F* hYvtx = new TH1F("hYvtx","y vertex distribution",100,-2,2);
166 TH1F* hZvtx = new TH1F("hZvtx","z vertex distribution",4*fNvtxBins,-4*fZvtxCut,4*fZvtxCut);
167
168 fListOfPrimaries.Add(hXvtx);
169 fListOfPrimaries.Add(hYvtx);
170 fListOfPrimaries.Add(hZvtx);
f55d559b 171
0b0a4ae5 172 hEventsSelected->Sumw2();
173 hEventsAll->Sumw2();
059c7c6b 174 hEventsAllNSD->Sumw2();
175
d31dce33 176 fListOfHits.Add(hEventsSelected);
f55d559b 177 fListOfHits.Add(hEventsSelectedVtx);
178 fListOfHits.Add(hEventsSelectedTrigger);
059c7c6b 179 fListOfHits.Add(hEventsSelectedNSDVtx);
180 fListOfHits.Add(hEventsSelectedNSD);
181
182
d31dce33 183 fListOfPrimaries.Add(hEventsAll);
059c7c6b 184 fListOfPrimaries.Add(hEventsAllNSD);
d31dce33 185
ab3e0abc 186 for(Int_t v=0; v<fNvtxBins;v++) {
187
188 for(Int_t iring = 0; iring<2;iring++) {
189 Char_t ringChar = (iring == 0 ? 'I' : 'O');
190 Int_t nSec = (iring == 1 ? 40 : 20);
191 TH2F* hAnalysed = new TH2F(Form("Analysed_FMD%c_vtx%d",ringChar,v),
192 Form("Analysed_FMD%c_vtx%d",ringChar,v),
820f4fcc 193 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
ab3e0abc 194
195 hAnalysed->Sumw2();
196 fListOfPrimaries.Add(hAnalysed);
197
059c7c6b 198 TH2F* hAnalysedNSD = new TH2F(Form("AnalysedNSD_FMD%c_vtx%d",ringChar,v),
199 Form("AnalysedNSD_FMD%c_vtx%d",ringChar,v),
200 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
201
202 hAnalysedNSD->Sumw2();
203 fListOfPrimaries.Add(hAnalysedNSD);
204
ab3e0abc 205 TH2F* hInel = new TH2F(Form("Inel_FMD%c_vtx%d",ringChar,v),
206 Form("Inel_FMD%c_vtx%d",ringChar,v),
820f4fcc 207 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
ab3e0abc 208
209 hInel->Sumw2();
210 fListOfPrimaries.Add(hInel);
211
059c7c6b 212 TH2F* hNSD = new TH2F(Form("NSD_FMD%c_vtx%d",ringChar,v),
213 Form("NSD_FMD%c_vtx%d",ringChar,v),
214 fNbinsEta, -4,6, nSec, 0,2*TMath::Pi());
215
216 hNSD->Sumw2();
217 fListOfPrimaries.Add(hNSD);
218
ab3e0abc 219 }
220 }
d31dce33 221
222 fVertexBins.SetName("VertexBins");
223 fVertexBins.GetXaxis()->Set(fNvtxBins,-1*fZvtxCut,fZvtxCut);
224
507687cd 225
226
d31dce33 227}
228//_____________________________________________________________________
229void AliFMDAnalysisTaskGenerateCorrection::Init()
230{
231 fLastTrackByStrip.Reset(-1);
232
233
234}
235//_____________________________________________________________________
236void AliFMDAnalysisTaskGenerateCorrection::UserExec(Option_t */*option*/)
237{
238
239 fLastTrackByStrip.Reset(-1);
240 fHitsByStrip.Reset(0);
241 AliMCEvent* mcevent = MCEvent();
242
243 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
244
245
246 AliESDEvent* esdevent = (AliESDEvent*)InputEvent();
059c7c6b 247
248 pars->SetTriggerStatus(esdevent);
249
d31dce33 250 Double_t esdvertex[3];
f55d559b 251 Bool_t vtxStatus = pars->GetVertex(esdevent,esdvertex);
d31dce33 252
f7356393 253 /* Double_t deltaEsd = 2*fZvtxCut/fNvtxBins;
254 Double_t vertexBinDoubleEsd = (esdvertex[2] + fZvtxCut) / deltaEsd;
255 Int_t vertexBinEsd = (Int_t)vertexBinDoubleEsd;*/
256
d31dce33 257 AliMCParticle* particle = 0;
258 AliStack* stack = mcevent->Stack();
259
260 UShort_t det,sec,strip;
261 Char_t ring;
262
263 Int_t nTracks = mcevent->GetNumberOfTracks();
264 AliHeader* header = mcevent->Header();
265 AliGenEventHeader* genHeader = header->GenEventHeader();
266
059c7c6b 267 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
507687cd 268 AliGenDPMjetEventHeader* dpmHeader = dynamic_cast<AliGenDPMjetEventHeader*>(header->GenEventHeader());
269
270 if (!pythiaGenHeader && !dpmHeader) {
271 std::cout<<" no pythia or dpm header! - NSD selection unusable"<<std::endl;
04d7f248 272 //return;
059c7c6b 273 }
059c7c6b 274 Bool_t nsd = kTRUE;
04d7f248 275 if(pythiaGenHeader) {
276 Int_t pythiaType = pythiaGenHeader->ProcessType();
277
278 if(pythiaType==92 || pythiaType==93)
279 nsd = kFALSE;
280 }
507687cd 281 if(dpmHeader) {
282 Int_t processType = dpmHeader->ProcessType();
283 if(processType == 5 || processType == 6)
284 nsd = kFALSE;
285
286 }
287
059c7c6b 288
d31dce33 289 TArrayF vertex;
290 genHeader->PrimaryVertex(vertex);
291
850ad94d 292 TH1F* hXvtx = (TH1F*)fListOfPrimaries.FindObject("hXvtx");
293 hXvtx->Fill(vertex.At(0));
294 TH1F* hYvtx = (TH1F*)fListOfPrimaries.FindObject("hYvtx");
295 hYvtx->Fill(vertex.At(1));
296 TH1F* hZvtx = (TH1F*)fListOfPrimaries.FindObject("hZvtx");
297 hZvtx->Fill(vertex.At(2));
298
299
059c7c6b 300
d31dce33 301
302 Double_t delta = 2*fZvtxCut/fNvtxBins;
303 Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta;
304 Int_t vertexBin = (Int_t)vertexBinDouble;
305
306 // Vertex determination correction
f55d559b 307 TH1F* hEventsSelected = (TH1F*)fListOfHits.FindObject("EventsSelected");
308 TH1F* hEventsSelectedVtx = (TH1F*)fListOfHits.FindObject("EventsSelectedVtx");
309 TH1F* hEventsSelectedTrigger = (TH1F*)fListOfHits.FindObject("EventsSelectedTrigger");
059c7c6b 310 TH1F* hEventsSelectedNSDVtx = (TH1F*)fListOfHits.FindObject("EventsSelectedNSDVtx");
311 TH1F* hEventsSelectedNSD = (TH1F*)fListOfHits.FindObject("EventsSelectedNSD");
f55d559b 312 TH1F* hEventsAll = (TH1F*)fListOfPrimaries.FindObject("EventsAll");
059c7c6b 313 TH1F* hEventsAllNSD = (TH1F*)fListOfPrimaries.FindObject("EventsAllNSD");
f55d559b 314
d31dce33 315 // TH1F* hTriggered = (TH1F*)fListOfHits.FindObject("Triggered");
316 // TH1F* hTriggeredAll = (TH1F*)fListOfPrimaries.FindObject("TriggeredAll");
317
318 Bool_t vtxFound = kTRUE;
f55d559b 319 if(!vtxStatus)
d31dce33 320 vtxFound = kFALSE;
d31dce33 321
4b8bdb60 322 //if(TMath::Abs(vertex.At(2)) > fZvtxCut) {
323 // vtxFound = kFALSE;
059c7c6b 324
4b8bdb60 325 //}
df2a9c32 326
327 if(TMath::Abs(vertex.At(2)) > fZvtxCut) {
328 return;
329 }
330
331
059c7c6b 332 Bool_t isTriggered = pars->IsEventTriggered(AliFMDAnaParameters::kMB1);
333 Bool_t isTriggeredNSD = pars->IsEventTriggered(AliFMDAnaParameters::kNSD);
f55d559b 334 if(vtxFound && isTriggered) hEventsSelected->Fill(vertexBin);
d31dce33 335
f55d559b 336 if(vtxFound) hEventsSelectedVtx->Fill(vertexBin);
337 if(isTriggered) hEventsSelectedTrigger->Fill(vertexBin);
059c7c6b 338
339 if(vtxFound && isTriggeredNSD) hEventsSelectedNSDVtx->Fill(vertexBin);
340 if(isTriggeredNSD) hEventsSelectedNSD->Fill(vertexBin);
341
342
d31dce33 343 hEventsAll->Fill(vertexBin);
059c7c6b 344 if(nsd) hEventsAllNSD->Fill(vertexBin);
70d74659 345
f7356393 346 // FMD dead channels
347
348 TH1F* hFMDReadChannels = (TH1F*)fListOfHits.FindObject(Form("hFMDReadChannels_vtx%d",vertexBin));
349 TH1F* hFMDAllChannels = (TH1F*)fListOfHits.FindObject(Form("hFMDAllChannels_vtx%d",vertexBin));
350
351 AliESDFMD* fmd = esdevent->GetFMDData();
352
353 for(UShort_t d=1;d<=3;d++) {
354 Int_t nRings = (d==1 ? 1 : 2);
355 for (UShort_t ir = 0; ir < nRings; ir++) {
356 Char_t r = (ir == 0 ? 'I' : 'O');
357 UShort_t nsec = (ir == 0 ? 20 : 40);
358 UShort_t nstr = (ir == 0 ? 512 : 256);
359 for(UShort_t s =0; s < nsec; s++) {
360 for(UShort_t str = 0; str < nstr; str++) {
361 Float_t mult = fmd->Multiplicity(d,r,s,str);
362
363 hFMDAllChannels->Fill(pars->GetEtaFromStrip(d,r,s,str,esdvertex[2]));
364
365 if(mult != AliESDFMD::kInvalidMult)
366 hFMDReadChannels->Fill(pars->GetEtaFromStrip(d,r,s,str,esdvertex[2]));
367 }
368 }
369 }
370 }
371 // End, FMD dead channels
df2a9c32 372
d31dce33 373 for(Int_t i = 0 ;i<nTracks;i++) {
7aad0c47 374 particle = (AliMCParticle*) mcevent->GetTrack(i);
d31dce33 375
376 if(!particle)
377 continue;
378
379 if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
380
df2a9c32 381 if(vtxFound && isTriggered) {
382 TH2F* hPrimaryInner = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",'I',vertexBin));
383 TH2F* hPrimaryOuter = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",'O',vertexBin));
384 hPrimaryInner->Fill(particle->Eta(),particle->Phi());
385 hPrimaryOuter->Fill(particle->Eta(),particle->Phi());
386 if( isTriggeredNSD) {
387 TH2F* hPrimaryInnerNSD = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimaryNSD_FMD_%c_vtx%d",'I',vertexBin));
388 TH2F* hPrimaryOuterNSD = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimaryNSD_FMD_%c_vtx%d",'O',vertexBin));
389 hPrimaryInnerNSD->Fill(particle->Eta(),particle->Phi());
390 hPrimaryOuterNSD->Fill(particle->Eta(),particle->Phi());
391 }
507687cd 392 }
ab3e0abc 393 TH2F* hAnalysedInner = (TH2F*)fListOfPrimaries.FindObject( Form("Analysed_FMD%c_vtx%d",'I',vertexBin));
394 TH2F* hAnalysedOuter = (TH2F*)fListOfPrimaries.FindObject( Form("Analysed_FMD%c_vtx%d",'O',vertexBin));
059c7c6b 395 TH2F* hAnalysedNSDInner = (TH2F*)fListOfPrimaries.FindObject( Form("AnalysedNSD_FMD%c_vtx%d",'I',vertexBin));
396 TH2F* hAnalysedNSDOuter = (TH2F*)fListOfPrimaries.FindObject( Form("AnalysedNSD_FMD%c_vtx%d",'O',vertexBin));
ab3e0abc 397 TH2F* hInelInner = (TH2F*)fListOfPrimaries.FindObject( Form("Inel_FMD%c_vtx%d",'I',vertexBin));
398 TH2F* hInelOuter = (TH2F*)fListOfPrimaries.FindObject( Form("Inel_FMD%c_vtx%d",'O',vertexBin));
059c7c6b 399 TH2F* hNSDInner = (TH2F*)fListOfPrimaries.FindObject( Form("NSD_FMD%c_vtx%d",'I',vertexBin));
400 TH2F* hNSDOuter = (TH2F*)fListOfPrimaries.FindObject( Form("NSD_FMD%c_vtx%d",'O',vertexBin));
ab3e0abc 401
820f4fcc 402 //if(vtxFound && isTriggered) {
059c7c6b 403 if(isTriggeredNSD) {
404 hAnalysedNSDInner->Fill(particle->Eta(),particle->Phi());
405 hAnalysedNSDOuter->Fill(particle->Eta(),particle->Phi());
406 }
df2a9c32 407 if(isTriggered && vtxFound) {
ab3e0abc 408 hAnalysedInner->Fill(particle->Eta(),particle->Phi());
409 hAnalysedOuter->Fill(particle->Eta(),particle->Phi());
410 }
411 hInelInner->Fill(particle->Eta(),particle->Phi());
412 hInelOuter->Fill(particle->Eta(),particle->Phi());
413
059c7c6b 414 if(nsd) {
415 hNSDInner->Fill(particle->Eta(),particle->Phi());
416 hNSDOuter->Fill(particle->Eta(),particle->Phi());
417 }
d31dce33 418 }
df2a9c32 419 if(!vtxFound || !isTriggered) continue;
d31dce33 420
421 for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) {
422
423 AliTrackReference* ref = particle->GetTrackReference(j);
424
425 if(ref->DetectorId() != AliTrackReference::kFMD)
426 continue;
ab3e0abc 427
d31dce33 428 AliFMDStripIndex::Unpack(ref->UserId(),det,ring,sec,strip);
429 Float_t thisStripTrack = fLastTrackByStrip(det,ring,sec,strip);
430 if(particle->Charge() != 0 && i != thisStripTrack ) {
431
432 Float_t phi = pars->GetPhiFromSector(det,ring,sec);
433 Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));
434
435 TH2F* hHits = (TH2F*)fListOfHits.FindObject(Form("hHits_FMD%d%c_vtx%d", det,ring,vertexBin));
436 hHits->Fill(eta,phi);
507687cd 437 if( isTriggeredNSD) {
438 TH2F* hHitsNSD = (TH2F*)fListOfHits.FindObject(Form("hHitsNSD_FMD%d%c_vtx%d", det,ring,vertexBin));
439 hHitsNSD->Fill(eta,phi);
440 }
441
d31dce33 442 Float_t nstrips = (ring =='O' ? 256 : 512);
443 fHitsByStrip(det,ring,sec,strip) +=1;
444 TH1F* allHits = (TH1F*)fListOfHits.FindObject(Form("allHits_FMD%d%c",det,ring));
445 TH1F* doubleHits = (TH1F*)fListOfHits.FindObject(Form("DoubleHits_FMD%d%c",det,ring));
446
447 if(fHitsByStrip(det,ring,sec,strip) == 1)
448 allHits->Fill(eta);
449
450 doubleHits->Fill(eta);
451
452 fLastTrackByStrip(det,ring,sec,strip) = (Float_t)i;
453 if(strip >0)
454 fLastTrackByStrip(det,ring,sec,strip-1) = (Float_t)i;
455 if(strip < (nstrips - 1))
456 fLastTrackByStrip(det,ring,sec,strip+1) = (Float_t)i;
457 }
458 }
459
460 }
04f1ff3d 461
04f1ff3d 462
df2a9c32 463 if(vtxFound && isTriggered) {
464 //SPD part HHD
465 TH2F* hSPDMult = (TH2F*)fListOfHits.FindObject(Form("hSPDhits_vtx%d", vertexBin));
466
467 const AliMultiplicity* spdmult = esdevent->GetMultiplicity();
468 for(Int_t j = 0; j< spdmult->GetNumberOfTracklets();j++)
469 hSPDMult->Fill(spdmult->GetEta(j),spdmult->GetPhi(j));
d31dce33 470
df2a9c32 471 for(Int_t j = 0; j< spdmult->GetNumberOfSingleClusters();j++)
472 hSPDMult->Fill(-TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.)),spdmult->GetPhiSingle(j));
473
474 }
475
d31dce33 476 PostData(1, &fListOfHits);
477 PostData(2, &fListOfPrimaries);
478 PostData(3, &fVertexBins);
479}
480//_____________________________________________________________________
481void AliFMDAnalysisTaskGenerateCorrection::Terminate(Option_t */*option*/)
482{
483 /* TH1F* allHits = (TH1F*)fListOfHits.FindObject("allHits");
484 TH1F* doubleHits = (TH1F*)fListOfHits.FindObject("DoubleHits");
485
486 doubleHits->Divide(allHits);
487 GenerateCorrection();
488 PostData(1, &fListOfHits);
489 PostData(4, &fListOfCorrection);*/
490
491}
492//_____________________________________________________________________
493void AliFMDAnalysisTaskGenerateCorrection::GenerateCorrection() {
494
495 fBackground = new AliFMDAnaCalibBackgroundCorrection();
496 fEventSelectionEff = new AliFMDAnaCalibEventSelectionEfficiency();
497
ab3e0abc 498 //Event selection
f55d559b 499 TH1F* hEventsSelected = (TH1F*)fListOfHits.FindObject("EventsSelected");
500 TH1F* hEventsSelectedVtx = (TH1F*)fListOfHits.FindObject("EventsSelectedVtx");
501 TH1F* hEventsSelectedTrigger = (TH1F*)fListOfHits.FindObject("EventsSelectedTrigger");
059c7c6b 502 TH1F* hEventsSelectedNSDVtx = (TH1F*)fListOfHits.FindObject("EventsSelectedNSDVtx");
503 TH1F* hEventsSelectedNSD = (TH1F*)fListOfHits.FindObject("EventsSelectedNSD");
504
f55d559b 505 TH1F* hEventsAll = (TH1F*)fListOfPrimaries.FindObject("EventsAll");
059c7c6b 506 TH1F* hEventsAllNSD = (TH1F*)fListOfPrimaries.FindObject("EventsAllNSD");
820f4fcc 507 TH1F* hEventsSelectedVtxDivByTr = (TH1F*)hEventsSelectedVtx->Clone("hEventsSelectedVtxDivByTr");
059c7c6b 508 TH1F* hEventsSelectedNSDVtxDivByNSD = (TH1F*)hEventsSelectedNSDVtx->Clone("hEventsSelectedNSDVtxDivByNSD");
509
820f4fcc 510 fListOfHits.Add(hEventsSelectedVtxDivByTr);
059c7c6b 511 fListOfHits.Add(hEventsSelectedNSDVtxDivByNSD);
512 hEventsSelectedNSDVtxDivByNSD->Divide(hEventsSelectedNSD);
820f4fcc 513 hEventsSelectedVtxDivByTr->Divide(hEventsSelectedTrigger);
059c7c6b 514
515 for(Int_t v=0;v<fNvtxBins ; v++) {
ab3e0abc 516 TH2F* hAnalysedInner = (TH2F*)fListOfPrimaries.FindObject(Form("Analysed_FMD%c_vtx%d",'I',v));
517 TH2F* hAnalysedOuter = (TH2F*)fListOfPrimaries.FindObject(Form("Analysed_FMD%c_vtx%d",'O',v));
059c7c6b 518 TH2F* hAnalysedNSDInner = (TH2F*)fListOfPrimaries.FindObject( Form("AnalysedNSD_FMD%c_vtx%d",'I',v));
519 TH2F* hAnalysedNSDOuter = (TH2F*)fListOfPrimaries.FindObject( Form("AnalysedNSD_FMD%c_vtx%d",'O',v));
520
ab3e0abc 521 TH2F* hInelInner = (TH2F*)fListOfPrimaries.FindObject(Form("Inel_FMD%c_vtx%d",'I',v));
522 TH2F* hInelOuter = (TH2F*)fListOfPrimaries.FindObject(Form("Inel_FMD%c_vtx%d",'O',v));
059c7c6b 523 TH2F* hNSDInner = (TH2F*)fListOfPrimaries.FindObject( Form("NSD_FMD%c_vtx%d",'I',v));
524 TH2F* hNSDOuter = (TH2F*)fListOfPrimaries.FindObject( Form("NSD_FMD%c_vtx%d",'O',v));
820f4fcc 525 // hAnalysedInner->Scale((hEventsSelected->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelected->GetBinContent(v+1)));
526 //hAnalysedOuter->Scale((hEventsSelected->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelected->GetBinContent(v+1)));
527
528 hAnalysedInner->Scale((hEventsSelectedTrigger->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelectedTrigger->GetBinContent(v+1)));
529
530 hAnalysedOuter->Scale((hEventsSelectedTrigger->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelectedTrigger->GetBinContent(v+1)));
059c7c6b 531 hAnalysedNSDInner->Scale((hEventsSelectedNSD->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelectedNSD->GetBinContent(v+1)));
532
533 hAnalysedNSDOuter->Scale((hEventsSelectedNSD->GetBinContent(v+1) == 0 ? 0 : 1/hEventsSelectedNSD->GetBinContent(v+1)));
534
820f4fcc 535
059c7c6b 536
ab3e0abc 537 hInelInner->Scale((hEventsAll->GetBinContent(v+1) == 0 ? 0 : 1/hEventsAll->GetBinContent(v+1)));
538 hInelOuter->Scale((hEventsAll->GetBinContent(v+1) == 0 ? 0 : 1/hEventsAll->GetBinContent(v+1)));
820f4fcc 539
059c7c6b 540 hNSDInner->Scale((hEventsAllNSD->GetBinContent(v+1) == 0 ? 0 : 1/hEventsAllNSD->GetBinContent(v+1)));
541 hNSDOuter->Scale((hEventsAllNSD->GetBinContent(v+1) == 0 ? 0 : 1/hEventsAllNSD->GetBinContent(v+1)));
542
820f4fcc 543
ab3e0abc 544 hAnalysedInner->Divide(hInelInner);
545 hAnalysedOuter->Divide(hInelOuter);
820f4fcc 546
059c7c6b 547 hAnalysedNSDInner->Divide(hNSDInner);
548 hAnalysedNSDOuter->Divide(hNSDOuter);
549
550 fEventSelectionEff->SetCorrection("INEL",v,'I',hAnalysedInner);
551 fEventSelectionEff->SetCorrection("INEL",v,'O',hAnalysedOuter);
552 //NSD
553 fEventSelectionEff->SetCorrection("NSD",v,'I',hAnalysedNSDInner);
554 fEventSelectionEff->SetCorrection("NSD",v,'O',hAnalysedNSDOuter);
555
f7356393 556 //FMD dead channels
557 TH1F* hFMDReadChannels = (TH1F*)fListOfHits.FindObject(Form("hFMDReadChannels_vtx%d",v));
558
559 TH1F* hFMDAllChannels = (TH1F*)fListOfHits.FindObject(Form("hFMDAllChannels_vtx%d",v));
560 if(hFMDReadChannels) {
561 TH1F* hFMDDeadCorrection = (TH1F*)hFMDReadChannels->Clone("hFMDDeadCorrection");
562 hFMDDeadCorrection->Divide(hFMDAllChannels);
563 fBackground->SetFMDDeadCorrection(v,hFMDDeadCorrection);
564 fListOfCorrection.Add(hFMDDeadCorrection);
565 }
566 else AliWarning("No Dead Channel Correction generated");
059c7c6b 567
541c19ed 568 }
ab3e0abc 569
820f4fcc 570 Float_t vtxEff = 1;
571 if(hEventsSelectedTrigger->GetEntries())
572 vtxEff = hEventsSelectedVtx->GetEntries() / hEventsSelectedTrigger->GetEntries();
573 fEventSelectionEff->SetVtxToTriggerRatio(vtxEff);
ab3e0abc 574
b85ea106 575 // hEventsAll->Divide(hEventsAll,hEventsSelected,1,1,"B");
f55d559b 576 hEventsSelectedVtx->Divide(hEventsAll);
577 hEventsSelectedTrigger->Divide(hEventsAll);
0b0a4ae5 578
059c7c6b 579 hEventsSelectedNSDVtx->Divide(hEventsAllNSD);
580 hEventsSelectedNSD->Divide(hEventsAllNSD);
581
0b0a4ae5 582 for(Int_t i = 1; i<=hEventsSelected->GetNbinsX(); i++) {
583 if(hEventsSelected->GetBinContent(i) == 0 )
584 continue;
059c7c6b 585 Float_t b = hEventsSelected->GetBinContent(i);
586 // Float_t b = hEventsSelected->GetBinContent(i)*hEventsSelectedTrigger->GetBinContent(i);
b85ea106 587 Float_t db = hEventsSelected->GetBinError(i);
0b0a4ae5 588 Float_t sum = hEventsAll->GetBinContent(i);
589 Float_t dsum = hEventsAll->GetBinError(i);
b85ea106 590 Float_t a = sum-b;
591 Float_t da = TMath::Sqrt(TMath::Power(db,2) + TMath::Power(dsum,2));
0b0a4ae5 592
b85ea106 593 Float_t cor = sum / b;
594 Float_t ecor = TMath::Sqrt(TMath::Power(da,2) + TMath::Power(a/(b*db),2)) / b;
0b0a4ae5 595
b85ea106 596 hEventsAll->SetBinContent(i,cor);
597 hEventsAll->SetBinError(i,ecor);
0b0a4ae5 598
599 }
d31dce33 600
ab3e0abc 601 //TH1F* hEventTest = (TH1F*)hEventsAll->Clone("hEventTest");
602
603
604
b85ea106 605 fEventSelectionEff->SetCorrection(hEventsAll);
d31dce33 606
607 for(Int_t det= 1; det <=3; det++) {
608 Int_t nRings = (det==1 ? 1 : 2);
609
610 for(Int_t iring = 0; iring<nRings; iring++) {
611 Char_t ring = (iring == 0 ? 'I' : 'O');
612 TH1F* allHits = (TH1F*)fListOfHits.FindObject(Form("allHits_FMD%d%c",det,ring));
613 TH1F* doubleHits = (TH1F*)fListOfHits.FindObject(Form("DoubleHits_FMD%d%c",det,ring));
b85ea106 614 allHits->Divide(doubleHits);
615
616 fBackground->SetDoubleHitCorrection(det,ring,allHits);
617
d31dce33 618 for(Int_t vertexBin=0;vertexBin<fNvtxBins ;vertexBin++) {
619 TH2F* hHits = (TH2F*)fListOfHits.FindObject(Form("hHits_FMD%d%c_vtx%d", det,ring,vertexBin));
507687cd 620 TH2F* hPrimary = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",ring,vertexBin));
621 TH2F* hCorrection = (TH2F*)hHits->Clone(Form("FMD%d%c_vtxbin_%d_correction",det,ring,vertexBin));
d31dce33 622 hCorrection->Divide(hPrimary);
507687cd 623
624 TH2F* hHitsNSD = (TH2F*)fListOfHits.FindObject(Form("hHitsNSD_FMD%d%c_vtx%d", det,ring,vertexBin));
625 TH2F* hPrimaryNSD = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimaryNSD_FMD_%c_vtx%d",ring,vertexBin));
df2a9c32 626 TH2F* hCorrectionNSD = 0;
627 if(hHitsNSD && hPrimaryNSD) {
628 hCorrectionNSD = (TH2F*)hHitsNSD->Clone(Form("FMDNSD%d%c_vtxbin_%d_correction",det,ring,vertexBin));
629 hCorrectionNSD->Divide(hPrimaryNSD);
630 }
0b0a4ae5 631
d31dce33 632 hCorrection->SetTitle(hCorrection->GetName());
633 fListOfCorrection.Add(hCorrection);
634 fBackground->SetBgCorrection(det,ring,vertexBin,hCorrection);
df2a9c32 635
636 if(hCorrectionNSD) {
637 fBackground->SetNSDBgCorrection(det,ring,vertexBin,hCorrectionNSD);
638 hCorrectionNSD->SetTitle(hCorrectionNSD->GetName());
639 fListOfCorrection.Add(hCorrectionNSD);
640 }
d31dce33 641
642 }
643
644 }
645 }
04f1ff3d 646 for(Int_t vertexBin=0;vertexBin<fNvtxBins ;vertexBin++) {
647 TH2F* hPrimary = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",'I',vertexBin));
648 TH2F* hSPDMult = (TH2F*)fListOfHits.FindObject(Form("hSPDhits_vtx%d", vertexBin));
11e702b3 649 if(!hSPDMult) continue;
650
04f1ff3d 651 TH2F* hCorrection = (TH2F*)hSPDMult->Clone(Form("SPD_vtxbin_%d_correction",vertexBin));
652 hCorrection->SetTitle(hCorrection->GetName());
653 fListOfCorrection.Add(hCorrection);
654 hCorrection->Divide(hPrimary);
f7356393 655
04f1ff3d 656
657 TH1F* hAlive = new TH1F(Form("hAliveSPD_vtxbin%d",vertexBin),Form("hAliveSPD_vtxbin%d",vertexBin),hSPDMult->GetNbinsX(),hSPDMult->GetXaxis()->GetXmin(), hSPDMult->GetXaxis()->GetXmax());
658 TH1F* hPresent = new TH1F(Form("hPresentSPD_vtxbin%d",vertexBin),Form("hPresentSPD_vtxbin%d",vertexBin),hSPDMult->GetNbinsX(),hSPDMult->GetXaxis()->GetXmin(), hSPDMult->GetXaxis()->GetXmax());
659 for(Int_t xx = 1; xx <=hSPDMult->GetNbinsX(); xx++) {
660
df2a9c32 661
662
04f1ff3d 663 for(Int_t yy = 1; yy <=hSPDMult->GetNbinsY(); yy++) {
df2a9c32 664 if(TMath::Abs(hCorrection->GetXaxis()->GetBinCenter(xx)) > 1.9) {
665 hCorrection->SetBinContent(xx,yy,0.);
666 hCorrection->SetBinError(xx,yy,0.);
667 }
668
669 if(hCorrection->GetBinContent(xx,yy) > 0.9) {
670 hAlive->Fill(hCorrection->GetXaxis()->GetBinCenter(xx) );
671 }
672 else {
673 hCorrection->SetBinContent(xx,yy,0.);
674 hCorrection->SetBinError(xx,yy,0.);
675 }
676
04f1ff3d 677 hPresent->Fill(hCorrection->GetXaxis()->GetBinCenter(xx));
678
679 }
680 }
681 TH1F* hDeadCorrection = (TH1F*)hAlive->Clone(Form("hSPDDeadCorrection_vtxbin%d",vertexBin));
682 hDeadCorrection->Divide(hPresent);
683 fBackground->SetSPDDeadCorrection(vertexBin,hDeadCorrection);
684 fListOfCorrection.Add(hDeadCorrection);
ab3e0abc 685
f7356393 686 fBackground->SetBgCorrection(0,'Q',vertexBin,hCorrection);
687 }
ab3e0abc 688
689
d31dce33 690 TAxis refAxis(fNvtxBins,-1*fZvtxCut,fZvtxCut);
691 fBackground->SetRefAxis(&refAxis);
692
693}
694//_____________________________________________________________________
695void AliFMDAnalysisTaskGenerateCorrection::ReadFromFile(const Char_t* filename, Bool_t storeInOCDB, Int_t /*runNo*/) {
0b0a4ae5 696
d31dce33 697 TFile infile(filename);
698 TH1F* hVertex = (TH1F*)infile.Get("VertexBins");
699 fZvtxCut = hVertex->GetXaxis()->GetXmax();
700 fNvtxBins = hVertex->GetXaxis()->GetNbins();
701 fVertexBins.SetName("VertexBins");
702 fVertexBins.GetXaxis()->Set(fNvtxBins,-1*fZvtxCut,fZvtxCut);
703
704 TList* listOfHits = (TList*)infile.Get("Hits");
705 TList* listOfPrim = (TList*)infile.Get("Primaries");
706
f55d559b 707 TH1F* hEventsSelected = (TH1F*)listOfHits->FindObject("EventsSelected");
708 TH1F* hEventsSelectedVtx = (TH1F*)listOfHits->FindObject("EventsSelectedVtx");
709 TH1F* hEventsSelectedTrigger = (TH1F*)listOfHits->FindObject("EventsSelectedTrigger");
059c7c6b 710 TH1F* hEventsSelectedNSDVtx = (TH1F*)listOfHits->FindObject("EventsSelectedNSDVtx");
711 TH1F* hEventsSelectedNSD = (TH1F*)listOfHits->FindObject("EventsSelectedNSD");
f55d559b 712 TH1F* hEventsAll = (TH1F*)listOfPrim->FindObject("EventsAll");
059c7c6b 713 TH1F* hEventsAllNSD = (TH1F*)listOfPrim->FindObject("EventsAllNSD");
f55d559b 714
d31dce33 715 fListOfHits.Add(hEventsSelected);
f55d559b 716 fListOfHits.Add(hEventsSelectedVtx);
059c7c6b 717 fListOfHits.Add(hEventsSelectedNSD);
718 fListOfHits.Add(hEventsSelectedNSDVtx);
f55d559b 719 fListOfHits.Add(hEventsSelectedTrigger);
d31dce33 720 fListOfPrimaries.Add(hEventsAll);
059c7c6b 721 fListOfPrimaries.Add(hEventsAllNSD);
d31dce33 722
850ad94d 723 TH1F* hXvtx = (TH1F*)listOfPrim->FindObject("hXvtx");
724 TH1F* hYvtx = (TH1F*)listOfPrim->FindObject("hYvtx");
725 TH1F* hZvtx = (TH1F*)listOfPrim->FindObject("hZvtx");
726 fListOfPrimaries.Add(hXvtx);
727 fListOfPrimaries.Add(hYvtx);
728 fListOfPrimaries.Add(hZvtx);
729
d31dce33 730 for(Int_t det =1; det<=3;det++)
059c7c6b 731 {
732 Int_t nRings = (det==1 ? 1 : 2);
733 for(Int_t ring = 0;ring<nRings;ring++)
734 {
735 Char_t ringChar = (ring == 0 ? 'I' : 'O');
736 TH1F* allHits = (TH1F*)listOfHits->FindObject(Form("allHits_FMD%d%c",det,ringChar));
737 TH1F* doubleHits = (TH1F*)listOfHits->FindObject(Form("DoubleHits_FMD%d%c",det,ringChar));
738 fListOfHits.Add(allHits);
739 fListOfHits.Add(doubleHits);
740 for(Int_t v=0; v<fNvtxBins;v++)
741 {
742
743 TH2F* hHits = (TH2F*)listOfHits->FindObject(Form("hHits_FMD%d%c_vtx%d", det,ringChar,v));
744 fListOfHits.Add(hHits);
507687cd 745 TH2F* hHitsNSD = (TH2F*)listOfHits->FindObject(Form("hHitsNSD_FMD%d%c_vtx%d", det,ringChar,v));
df2a9c32 746 if(hHitsNSD)
747 fListOfHits.Add(hHitsNSD);
059c7c6b 748 }
749 }
750 }
04f1ff3d 751 for(Int_t v=0; v<fNvtxBins;v++) {
752 TH2F* hSPDHits = (TH2F*)listOfHits->FindObject(Form("hSPDhits_vtx%d", v));
753 fListOfHits.Add(hSPDHits);
f7356393 754 TH1F* hFMDReadChannels = (TH1F*)listOfHits->FindObject(Form("hFMDReadChannels_vtx%d",v));
755 TH1F* hFMDAllChannels = (TH1F*)listOfHits->FindObject(Form("hFMDAllChannels_vtx%d",v));
04f1ff3d 756
f7356393 757 if(hFMDReadChannels && hFMDAllChannels) {
758 fListOfHits.Add(hFMDReadChannels);
759 fListOfHits.Add(hFMDAllChannels);
760 }
04f1ff3d 761 for(Int_t iring = 0; iring<2;iring++) {
762 Char_t ringChar = (iring == 0 ? 'I' : 'O');
763
d31dce33 764
765 TH2F* hPrimary = (TH2F*)listOfPrim->FindObject( Form("hPrimary_FMD_%c_vtx%d",ringChar,v));
507687cd 766 TH2F* hPrimaryNSD = (TH2F*)listOfPrim->FindObject( Form("hPrimaryNSD_FMD_%c_vtx%d",ringChar,v));
767
ab3e0abc 768 TH2F* hAnalysed = (TH2F*)listOfPrim->FindObject(Form("Analysed_FMD%c_vtx%d",ringChar,v));
059c7c6b 769 TH2F* hAnalysedNSD = (TH2F*)listOfPrim->FindObject(Form("AnalysedNSD_FMD%c_vtx%d",ringChar,v));
ab3e0abc 770 TH2F* hInel = (TH2F*)listOfPrim->FindObject(Form("Inel_FMD%c_vtx%d",ringChar,v));
059c7c6b 771 TH2F* hNSD = (TH2F*)listOfPrim->FindObject(Form("NSD_FMD%c_vtx%d",ringChar,v));
d31dce33 772
ab3e0abc 773 fListOfPrimaries.Add(hPrimary);
df2a9c32 774 if(hPrimaryNSD)
775 fListOfPrimaries.Add(hPrimaryNSD);
ab3e0abc 776 fListOfPrimaries.Add(hAnalysed);
777 fListOfPrimaries.Add(hInel);
059c7c6b 778 fListOfPrimaries.Add(hAnalysedNSD);
779 fListOfPrimaries.Add(hNSD);
d31dce33 780 }
781 }
782 GenerateCorrection();
783
784 TFile fout("backgroundFromFile.root","recreate");
785 fListOfHits.Write();
786 fListOfPrimaries.Write();
787 fListOfCorrection.Write();
788 fVertexBins.Write();
789
790 fout.Close();
0b0a4ae5 791 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
d31dce33 792 if(storeInOCDB) {
0b0a4ae5 793 TFile fbg(pars->GetPath(pars->GetBackgroundID()),"RECREATE");
d31dce33 794 fBackground->Write(AliFMDAnaParameters::GetBackgroundID());
795 fbg.Close();
0b0a4ae5 796 TFile feselect(pars->GetPath(pars->GetEventSelectionEffID()),"RECREATE");
d31dce33 797 fEventSelectionEff->Write(AliFMDAnaParameters::GetEventSelectionEffID());
798 feselect.Close();
799
800 }
801
802
803
804
805}
806//_____________________________________________________________________
807//
808// EOF
809//