]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FORWARD/analysis/AliFMDAnalysisTaskSharing.cxx
removing the setting of AOD track references for TPC only tracks
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis / AliFMDAnalysisTaskSharing.cxx
CommitLineData
3bb122c7 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 <TMath.h>
78f6f750 10//#include "AliFMDDebug.h"
3bb122c7 11#include "AliFMDAnalysisTaskSharing.h"
12#include "AliAnalysisManager.h"
13#include "AliESDFMD.h"
78f6f750 14//#include "AliFMDGeometry.h"
3bb122c7 15#include "AliMCEventHandler.h"
b64db9b1 16#include "AliMCEvent.h"
3bb122c7 17#include "AliStack.h"
18#include "AliESDVertex.h"
6289b3e8 19#include "AliMultiplicity.h"
3bb122c7 20#include "AliFMDAnaParameters.h"
0210c863 21#include "TH1F.h"
22#include "TObjString.h"
bb0a45c3 23//#include "/home/canute/ALICE/AliRoot/PWG0/AliPWG0Helper.h"
78f6f750 24//#include "AliFMDParameters.h"
b64db9b1 25#include "AliGenEventHeader.h"
da0805e2 26#include "AliGenPythiaEventHeader.h"
b64db9b1 27#include "AliHeader.h"
28#include "AliStack.h"
29#include "AliMCParticle.h"
0a2f2742 30#include "AliFMDStripIndex.h"
da0805e2 31#include "AliESDVZERO.h"
32#include "AliESDtrack.h"
507687cd 33#include "AliGenDPMjetEventHeader.h"
879c0f42 34#include "AliLog.h"
3bb122c7 35
0210c863 36// This is the task to do the FMD sharing or hit merging.
37// It reads the input ESDFMD data and posts an ESDFMD object to
38// the tasks that must be performed after this task ie.
39// Density, BackgroundCorrection and Dndeta.
40// Author: Hans Hjersing Dalsgaard, hans.dalsgaard@cern.ch
41
42
3bb122c7 43ClassImp(AliFMDAnalysisTaskSharing)
44
45//_____________________________________________________________________
46AliFMDAnalysisTaskSharing::AliFMDAnalysisTaskSharing()
47: fDebug(0),
48 fESD(0x0),
8dc7c4c2 49 foutputESDFMD(),
c78bc12b 50 fSharedThis(kFALSE),
7c3e5162 51 fSharedPrev(kFALSE),
a38d1871 52 fDiagList(0),
7c3e5162 53 fStandalone(kTRUE),
bb8a464f 54 fEsdVertex(0),
0a2f2742 55 fStatus(kTRUE),
aa303f0c 56 fLastTrackByStrip(0),
57 fLastOrbit(0)
3bb122c7 58{
59 // Default constructor
60 DefineInput (0, AliESDEvent::Class());
7c3e5162 61 DefineOutput(0, AliESDFMD::Class());
62 DefineOutput(1, AliESDVertex::Class());
63 DefineOutput(2, AliESDEvent::Class());
64 DefineOutput(3, TList::Class());
3bb122c7 65}
66//_____________________________________________________________________
1b418b63 67AliFMDAnalysisTaskSharing::AliFMDAnalysisTaskSharing(const char* name,
68 Bool_t SE):
3bb122c7 69 AliAnalysisTask(name, "AnalysisTaskFMD"),
70 fDebug(0),
71 fESD(0x0),
8dc7c4c2 72 foutputESDFMD(),
c78bc12b 73 fSharedThis(kFALSE),
7c3e5162 74 fSharedPrev(kFALSE),
a38d1871 75 fDiagList(0),
7c3e5162 76 fStandalone(kTRUE),
bb8a464f 77 fEsdVertex(0),
0a2f2742 78 fStatus(kTRUE),
da0805e2 79 fLastTrackByStrip(0),
aa303f0c 80 fLastOrbit(0)
3bb122c7 81{
0210c863 82 // named constructor
7c3e5162 83 fStandalone = SE;
84 if(fStandalone) {
85 DefineInput (0, AliESDEvent::Class());
86 DefineOutput(0, AliESDFMD::Class());
87 DefineOutput(1, AliESDVertex::Class());
88 DefineOutput(2, AliESDEvent::Class());
89 DefineOutput(3, TList::Class());
90 }
3bb122c7 91}
92//_____________________________________________________________________
daedf077 93Float_t AliFMDAnalysisTaskSharing::GetVtxEfficiencyFromData(){
94
95 TH1F* hEventsVtx = (TH1F*)fDiagList->FindObject("hEventsVtx");
96 TH1F* hEventsTr = (TH1F*)fDiagList->FindObject("hEventsTr");
97
1b418b63 98 if(hEventsTr->GetEntries() != 0 &&
99 hEventsVtx->GetEntries() != 0 &&
100 hEventsTr->GetEntries() != hEventsVtx->GetEntries())
daedf077 101 return hEventsVtx->GetEntries() / hEventsTr->GetEntries();
102 else return -1;
103
104}
059c7c6b 105//_____________________________________________________________________
1b418b63 106Float_t AliFMDAnalysisTaskSharing::GetNSDVtxEfficiencyFromData()
107{
059c7c6b 108
109 TH1F* hEventsNSDVtx = (TH1F*)fDiagList->FindObject("hEventsNSDVtx");
1b418b63 110 TH1F* hEventsNSD = (TH1F*)fDiagList->FindObject("hEventsNSD");
059c7c6b 111
1b418b63 112 if(hEventsNSD->GetEntries() != 0 &&
113 hEventsNSDVtx->GetEntries() != 0 &&
114 hEventsNSD->GetEntries() != hEventsNSDVtx->GetEntries())
059c7c6b 115 return hEventsNSDVtx->GetEntries() / hEventsNSD->GetEntries();
116 else return -1;
117
118}
119
120
daedf077 121//_____________________________________________________________________
3bb122c7 122void AliFMDAnalysisTaskSharing::CreateOutputObjects()
123{
0210c863 124 // Create the output objects
7c3e5162 125 if(!foutputESDFMD)
126 foutputESDFMD = new AliESDFMD();
127
128 if(!fEsdVertex)
129 fEsdVertex = new AliESDVertex();
130 //Diagnostics
a38d1871 131 if(!fDiagList)
132 fDiagList = new TList();
133
134 fDiagList->SetName("Sharing diagnostics");
b64db9b1 135
136 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
1b418b63 137 TH2F* hBg = pars->GetBackgroundCorrection(1, 'I', 0);
138 TH1F* hPrimary = new TH1F("hMultvsEtaNoCuts","hMultvsEtaNoCuts",
139 hBg->GetNbinsX(),
140 hBg->GetXaxis()->GetXmin(),
141 hBg->GetXaxis()->GetXmax());
b64db9b1 142 hPrimary->Sumw2();
143 fDiagList->Add(hPrimary);
059c7c6b 144
145 TH1F* hPrimaryNSD = new TH1F("hMultvsEtaNSDNoCuts","hMultvsEtaNSDNoCuts",
146 hBg->GetNbinsX(),
147 hBg->GetXaxis()->GetXmin(),
148 hBg->GetXaxis()->GetXmax());
149 hPrimaryNSD->Sumw2();
150 fDiagList->Add(hPrimaryNSD);
151
da0805e2 152 TH1F* hXvtx = new TH1F("hXvtx","x vertex distribution",100,-2,2);
153 TH1F* hYvtx = new TH1F("hYvtx","y vertex distribution",100,-2,2);
1b418b63 154 TH1F* hZvtx = new TH1F("hZvtx","z vertex distribution",
155 4*pars->GetNvtxBins(),
156 -4*pars->GetVtxCutZ(),
157 +4*pars->GetVtxCutZ());
f55d559b 158
da0805e2 159 fDiagList->Add(hXvtx);
160 fDiagList->Add(hYvtx);
f55d559b 161 fDiagList->Add(hZvtx);
162
b64db9b1 163 TH1F* hPrimVertexBin = 0;
0a2f2742 164 TH1F* hHits = 0;
b64db9b1 165 for(Int_t i = 0; i< pars->GetNvtxBins(); i++) {
166
167 hPrimVertexBin = new TH1F(Form("primmult_NoCuts_vtxbin%d",i),
168 Form("primmult_NoCuts_vtxbin%d",i),
169 hBg->GetNbinsX(),
170 hBg->GetXaxis()->GetXmin(),
171 hBg->GetXaxis()->GetXmax());
172 hPrimVertexBin->Sumw2();
173 fDiagList->Add(hPrimVertexBin);
174
175 }
176
7c3e5162 177 for(Int_t det = 1; det<=3; det++) {
178 Int_t nRings = (det==1 ? 1 : 2);
179
180 for(Int_t iring = 0;iring<nRings; iring++) {
181 Char_t ringChar = (iring == 0 ? 'I' : 'O');
182 TH1F* hEdist = new TH1F(Form("Edist_before_sharing_FMD%d%c", det, ringChar),
183 Form("Edist_before_sharing_FMD%d%c", det, ringChar),
bb8a464f 184 1000,0,25);
0210c863 185 TH1F* hEdistAfter = new TH1F(Form("Edist_after_sharing_FMD%d%c", det, ringChar),
7c3e5162 186 Form("Edist_after_sharing_FMD%d%c", det, ringChar),
bb8a464f 187 1000,0,25);
188
189
b64db9b1 190 //TH1F* hNstripsHit = new TH1F(Form("N_strips_hit_FMD%d%c",det,ringChar),
191 // Form("N_strips_hit_FMD%d%c",det,ringChar),
192 // 25,0,25);
a38d1871 193 fDiagList->Add(hEdist);
0210c863 194 fDiagList->Add(hEdistAfter);
b64db9b1 195 //fDiagList->Add(hNstripsHit);
0a2f2742 196
197 for(Int_t i = 0; i< pars->GetNvtxBins(); i++) {
198 hHits = new TH1F(Form("hMCHits_nocuts_FMD%d%c_vtxbin%d",det,ringChar,i),Form("hMCHits_FMD%d%c_vtxbin%d",det,ringChar,i),
199 hBg->GetNbinsX(),
200 hBg->GetXaxis()->GetXmin(),
201 hBg->GetXaxis()->GetXmax());
202 hHits->Sumw2();
203 fDiagList->Add(hHits);
7c3e5162 204
0a2f2742 205 }
206
7c3e5162 207 }
208 }
da0805e2 209 TH2F* hCorrelationFMDSPDhits = new TH2F("hCorrelationFMDSPDhits","hCorrelationFMDSPDhits;SPD;FMD ",100,0,200,100,0,500);
210 TH2F* hCorrelationFMDSPD = new TH2F("hCorrelationFMDSPD","hCorrelationFMDSPD;SPD;FMD ",100,0,200,100,0,500);
211 TH2F* hCorrelationFMD1SPD = new TH2F("hCorrelationFMD1SPD","hCorrelationFMD1SPD;SPD;FMD1 ",100,0,200,100,0,200);
212 TH2F* hCorrelationFMD2ISPD = new TH2F("hCorrelationFMD2ISPD","hCorrelationFMD2ISPD;SPD;FMD2I ",100,0,200,100,0,200);
213 TH2F* hCorrelationFMD2OSPD = new TH2F("hCorrelationFMD2OSPD","hCorrelationFMD2OSPD;SPD;FMD2O ",100,0,200,100,0,200);
214 TH2F* hCorrelationFMD3ISPD = new TH2F("hCorrelationFMD3ISPD","hCorrelationFMD3ISPD;SPD;FMD3I ",100,0,200,100,0,200);
215 TH2F* hCorrelationFMD3OSPD = new TH2F("hCorrelationFMD3OSPD","hCorrelationFMD3OSPD;SPD;FMD3O ",100,0,200,100,0,200);
216 TH2F* hCorrelationFMDGoodtracks = new TH2F("hCorrelationFMDGoodtracks","hCorrelationGoodtracks;good tracks;FMD ",100,0,200,100,0,200);
217 TH2F* hCorrelationFMDBadtracks = new TH2F("hCorrelationFMDBadtracks","hCorrelationBadtracks;bad tracks;FMD ",100,0,200,100,0,200);
218 TH2F* hCorrelationGoodbadtracks = new TH2F("hCorrelationGoodbadtracks","hCorrelationGoodbadtracks;good tracks;bad tracks ",100,0,200,100,0,200);
219 TH2F* hCorrelationSPDTracklets = new TH2F("hCorrelationSPDTracklets","hCorrelationSPDTracklets;hits ; tracklets ",100,0,500,100,0,200);
220 TH2F* hCorrelationClustersTracklets = new TH2F("hCorrelationClustersTracklets","hCorrelationClustersTracklets;clusters ; tracklets ",500,0,500,100,0,100);
221 TH2F* hCorrelationHitsRadius = new TH2F("hCorrelationHitsRadius","hCorrelationHitsRadius;hits ; radius ",100,0,500,100,0,10);
222 TH2F* hCorrelationHitsX = new TH2F("hCorrelationHitsX","hCorrelationHitsX;hits ; X ",100,0,500,100,-5,5);
223 TH2F* hCorrelationHitsY = new TH2F("hCorrelationHitsY","hCorrelationHitsY;hits ; Y ",100,0,500,100,-5,5);
224 fDiagList->Add(hCorrelationHitsRadius);
225 fDiagList->Add(hCorrelationHitsX);
226 fDiagList->Add(hCorrelationHitsY);
227 fDiagList->Add(hCorrelationFMDSPD);
228 fDiagList->Add(hCorrelationFMD1SPD);
229 fDiagList->Add(hCorrelationFMD2ISPD);
230 fDiagList->Add(hCorrelationFMD2OSPD);
231 fDiagList->Add(hCorrelationFMD3ISPD);
232 fDiagList->Add(hCorrelationFMD3OSPD);
233 fDiagList->Add(hCorrelationFMDGoodtracks);
234 fDiagList->Add(hCorrelationFMDBadtracks);
235 fDiagList->Add(hCorrelationGoodbadtracks);
236fDiagList->Add(hCorrelationFMDSPDhits);
237 fDiagList->Add(hCorrelationClustersTracklets);
238 fDiagList->Add(hCorrelationSPDTracklets);
239 TH2F* hCorrelationFMD23 = new TH2F("hCorrelationFMD23","hCorrelationFMD23;FMD2 ;FMD3 ",100,0,500,100,0,500);
240 TH2F* hCorrelationFMD2diff23 = new TH2F("hCorrelationFMD2diff23","hCorrelationFMD2diff23;FMD2 ;diff FMD23 ",100,0,100,100,0,100);
241 TH2F* hCorrelationFMD3diff23 = new TH2F("hCorrelationFMD3diff23","hCorrelationFMD3diff23;FMD3 ;diff FMD23 ",100,0,100,100,0,100);
242 TH2F* hCorrelationFMD1diff13 = new TH2F("hCorrelationFMD1diff13","hCorrelationFMD1diff13;FMD1 ;diff FMD13 ",100,0,100,100,0,100);
243 TH2F* hCorrelationFMD1diff12 = new TH2F("hCorrelationFMD1diff12","hCorrelationFMD1diff12;FMD1 ;diff FMD12 ",100,0,100,100,0,100);
244 TH2F* hCorrelationFMD12 = new TH2F("hCorrelationFMD12","hCorrelationFMD12;FMD1 ;FMD2 ",100,0,500,100,0,500);
245 TH2F* hCorrelationFMDBgCand = new TH2F("hCorrelationFMDBgCand","hCorrelationFMDBgCand;Bg Tr ;FMD ",100,0,100,500,0,500);
246
247 TH2F* hCorrelationFMDFlatTr = new TH2F("hCorrelationFMDFlatTr","hCorrelationFMDFlatTr;Bg Tr ;FMD ",100,0,100,500,0,500);
248 TH2F* hCorrelationFMDRatioFlatTr = new TH2F("hCorrelationFMDRatioFlatTr","hCorrelationFMDRatioFlatTr;Bg Tr ;FMD ",100,0,1,500,0,500);
249 fDiagList->Add(hCorrelationFMDBgCand);
250 fDiagList->Add(hCorrelationFMDFlatTr);
251 fDiagList->Add(hCorrelationFMDRatioFlatTr);
252 TH2F* hCorrelationFMDBgCandRelative = new TH2F("hCorrelationFMDBgCandRelative","hCorrelationFMDBgCandRelative;Bg Tr ;FMD ",100,0,2,500,0,500);
253 fDiagList->Add(hCorrelationFMDBgCandRelative);
254 fDiagList->Add(hCorrelationFMD2diff23);
255 fDiagList->Add(hCorrelationFMD3diff23);
256 fDiagList->Add(hCorrelationFMD1diff13);
257 fDiagList->Add(hCorrelationFMD1diff12);
258 fDiagList->Add(hCorrelationFMD23);
259 fDiagList->Add(hCorrelationFMD12);
260 TH2F* hTimeCorrelation = new TH2F("hCorrelationTime","hCorrelationTime ; time ; FMD hits",500,0,500,100,0,200);
261 fDiagList->Add(hTimeCorrelation);
262 TH1F* hHitDistribution = new TH1F("hFMDHitDistribution","hFMDHitDistribution ; FMD hits",500,0,500);
263
264 TH1F* hHitDistributionFMD1 = new TH1F("hFMDHitDistributionFMD1","hFMDHitDistributionFMD1 ; FMD1 hits",500,0,500);
265 TH1F* hHitDistributionFMD2I = new TH1F("hFMDHitDistributionFMD2I","hFMDHitDistributionFMD2I ; FMD2I hits",500,0,500);
266 TH1F* hHitDistributionFMD2O = new TH1F("hFMDHitDistributionFMD2O","hFMDHitDistributionFMD2O ; FMD2O hits",500,0,500);
267 TH1F* hHitDistributionFMD3I = new TH1F("hFMDHitDistributionFMD3I","hFMDHitDistributionFMD3I ; FMD3I hits",500,0,500);
268 TH1F* hHitDistributionFMD3O = new TH1F("hFMDHitDistributionFMD3O","hFMDHitDistributionFMD3O ; FMD3O hits",500,0,500);
269 TH1F* hTrVtxDistribution = new TH1F("hTrVtxDistribution","hTrVtxDistribution ; TrVtx",200,-500,500);
270 TH1F* hTrEtaDistribution = new TH1F("hTrEtaDistribution","hTrEtaDistribution ; TrEta",200,-9,9);
271 TH1F* hTrEtaDistribution2 = new TH1F("hTrEtaDistribution2","hTrEtaDistribution2 ; TrEta",200,-9,9);
272
273 TH1F* hFlatTracks = new TH1F("hFlatTracks","hFlatTracks ; Horizontal tracks",100,0,100);
aa303f0c 274
275 TH1F* hEnergyOfParticles = new TH1F("hEnergyOfParticles","hEnergyOfParticles",1000000,0,10);
276 fDiagList->Add(hEnergyOfParticles);
da0805e2 277 fDiagList->Add(hTrVtxDistribution);
278 fDiagList->Add(hTrEtaDistribution);
279 fDiagList->Add(hTrEtaDistribution2);
280 fDiagList->Add(hFlatTracks);
281 fDiagList->Add(hHitDistribution);
282 fDiagList->Add(hHitDistributionFMD1);
283 fDiagList->Add(hHitDistributionFMD2I);
284 fDiagList->Add(hHitDistributionFMD2O);
285 fDiagList->Add(hHitDistributionFMD3I);
286 fDiagList->Add(hHitDistributionFMD3O);
b64db9b1 287 TH1F* nMCevents = new TH1F("nMCEventsNoCuts","nMCEventsNoCuts",pars->GetNvtxBins(),0,pars->GetNvtxBins());
059c7c6b 288 TH1F* nMCeventsNSD= new TH1F("nMCEventsNSDNoCuts","nMCEventsNSDNoCuts",pars->GetNvtxBins(),0,pars->GetNvtxBins());
b64db9b1 289 fDiagList->Add(nMCevents);
059c7c6b 290 fDiagList->Add(nMCeventsNSD);
daedf077 291
059c7c6b 292 TH1F* hEventsVtx = new TH1F("hEventsVtx","hEventsVtx",pars->GetNvtxBins(),0,pars->GetNvtxBins());
293 TH1F* hEventsTr = new TH1F("hEventsTr","hEventsTr",pars->GetNvtxBins(),0,pars->GetNvtxBins());
294 TH1F* hEventsNSD = new TH1F("hEventsNSD","hEventsNSD",pars->GetNvtxBins(),0,pars->GetNvtxBins());
295 TH1F* hEventsNSDVtx = new TH1F("hEventsNSDVtx","hEventsNSDVtx",pars->GetNvtxBins(),0,pars->GetNvtxBins());
daedf077 296 fDiagList->Add(hEventsVtx);
297 fDiagList->Add(hEventsTr);
059c7c6b 298 fDiagList->Add(hEventsNSD);
299 fDiagList->Add(hEventsNSDVtx);
daedf077 300
3bb122c7 301}
302//_____________________________________________________________________
303void AliFMDAnalysisTaskSharing::ConnectInputData(Option_t */*option*/)
304{
0210c863 305 // connect the input data
7c3e5162 306 if(fStandalone)
307 fESD = (AliESDEvent*)GetInputData(0);
3bb122c7 308}
309//_____________________________________________________________________
310void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
311{
0210c863 312 //perform analysis on one event
3bb122c7 313 AliESD* old = fESD->GetAliESDOld();
314 if (old) {
315 fESD->CopyFromOldESD();
316 }
317
7c3e5162 318 foutputESDFMD->Clear();
541c19ed 319
320
da0805e2 321 Int_t delta = fESD->GetOrbitNumber() - fLastOrbit;
322 fLastOrbit = fESD->GetOrbitNumber();
323
324
b64db9b1 325 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
059c7c6b 326 Double_t vertex[3]={0,0,0};
f55d559b 327 Bool_t vtxStatus = pars->GetVertex(fESD,vertex);
bb8a464f 328 fEsdVertex->SetXYZ(vertex);
b64db9b1 329
059c7c6b 330 // std::cout<<vtxStatus<<" "<<vertex[0]<<" "<<vertex[1]<<" "<<vertex[2]<<std::endl;
331
b64db9b1 332 // Process primaries here to get true MC distribution
333 if(pars->GetProcessPrimary())
334 ProcessPrimary();
da0805e2 335 const AliMultiplicity* testmult = fESD->GetMultiplicity();
336 Int_t nTrackLets = testmult->GetNumberOfTracklets();
337 TH2F* hCorrelationClustersTracklets = (TH2F*)fDiagList->FindObject("hCorrelationClustersTracklets");
338 hCorrelationClustersTracklets->Fill(testmult->GetNumberOfSingleClusters(),nTrackLets);
b64db9b1 339
f55d559b 340
6289b3e8 341
059c7c6b 342 Bool_t isTriggered = pars->IsEventTriggered(AliFMDAnaParameters::kMB1);
343 Bool_t nsd = pars->IsEventTriggered(AliFMDAnaParameters::kNSD);
344
daedf077 345 Double_t delta2 = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
346 Double_t vertexBinDouble = (vertex[2] + pars->GetVtxCutZ()) / delta2;
347 Int_t vtxbin = (Int_t)vertexBinDouble;
059c7c6b 348 TH1F* hEventsVtx = (TH1F*)fDiagList->FindObject("hEventsVtx");
349 TH1F* hEventsTr = (TH1F*)fDiagList->FindObject("hEventsTr");
350 TH1F* hEventsNSD = (TH1F*)fDiagList->FindObject("hEventsNSD");
351 TH1F* hEventsNSDVtx = (TH1F*)fDiagList->FindObject("hEventsNSDVtx");
352
df2a9c32 353 if( TMath::Abs(vertex[2]) > pars->GetVtxCutZ()) {
354 fStatus = kFALSE;
355 return;
356 }
059c7c6b 357
daedf077 358 if(vtxStatus) hEventsVtx->Fill(vtxbin);
359 if(isTriggered) hEventsTr->Fill(vtxbin);
da0805e2 360
059c7c6b 361 if(vtxStatus && nsd) hEventsNSDVtx->Fill(vtxbin);
362 if(nsd) hEventsNSD->Fill(vtxbin);
363
da0805e2 364 if(!isTriggered || !vtxStatus ) {
bb8a464f 365 fStatus = kFALSE;
366 return;
367 }
368 else
369 fStatus = kTRUE;
b64db9b1 370
da0805e2 371 TH1F* hXvtx = (TH1F*)fDiagList->FindObject("hXvtx");
059c7c6b 372 if(vtxStatus ) hXvtx->Fill(vertex[0]);
da0805e2 373 TH1F* hYvtx = (TH1F*)fDiagList->FindObject("hYvtx");
059c7c6b 374 if(vtxStatus ) hYvtx->Fill(vertex[1]);
f55d559b 375 TH1F* hZvtx = (TH1F*)fDiagList->FindObject("hZvtx");
376 hZvtx->Fill(vertex[2]);
da0805e2 377 //const AliMultiplicity* testmult = fESD->GetMultiplicity();
378 //std::cout<<vertex[2]<<std::endl;
379 //Int_t nTrackLets = testmult->GetNumberOfTracklets();
6289b3e8 380
df2a9c32 381
da0805e2 382
6289b3e8 383 if(nTrackLets < 1000) foutputESDFMD->SetUniqueID(kTRUE);
384 else foutputESDFMD->SetUniqueID(kFALSE);
b64db9b1 385
3bb122c7 386 AliESDFMD* fmd = fESD->GetFMDData();
059c7c6b 387
3bb122c7 388 if (!fmd) return;
da0805e2 389 Int_t nHits[3][2] = {{0,0},{0,0},{0,0}};
032f0117 390 //Int_t nDead = 0;
3bb122c7 391 for(UShort_t det=1;det<=3;det++) {
392 Int_t nRings = (det==1 ? 1 : 2);
393 for (UShort_t ir = 0; ir < nRings; ir++) {
394 Char_t ring = (ir == 0 ? 'I' : 'O');
395 UShort_t nsec = (ir == 0 ? 20 : 40);
396 UShort_t nstr = (ir == 0 ? 512 : 256);
7c3e5162 397
a38d1871 398 TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_before_sharing_FMD%d%c",det,ring));
7c3e5162 399
8dc7c4c2 400 for(UShort_t sec =0; sec < nsec; sec++) {
401 fSharedThis = kFALSE;
402 fSharedPrev = kFALSE;
0210c863 403
3bb122c7 404 for(UShort_t strip = 0; strip < nstr; strip++) {
7c3e5162 405 foutputESDFMD->SetMultiplicity(det,ring,sec,strip,0.);
3bb122c7 406 Float_t mult = fmd->Multiplicity(det,ring,sec,strip);
032f0117 407 //if(mult == AliESDFMD::kInvalidMult)
408 // nDead++;
409 // if(mult > 0)
410 // std::cout<<mult<<std::endl;
bb8a464f 411
032f0117 412 if(mult == AliESDFMD::kInvalidMult || mult == 0 || mult > 20) continue;
da0805e2 413
78f6f750 414 //Double_t eta = fmd->Eta(det,ring,sec,strip);
415 Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex[2]);
da0805e2 416
7c3e5162 417 hEdist->Fill(mult);
032f0117 418
419
bb8a464f 420 if(fmd->IsAngleCorrected())
cc066cb9 421 mult = mult/TMath::Cos(Eta2Theta(eta));
0210c863 422 Float_t prevE = 0;
423 Float_t nextE = 0;
8dc823cc 424 if(strip != 0)
bb8a464f 425 if(fmd->Multiplicity(det,ring,sec,strip-1) != AliESDFMD::kInvalidMult) {
0210c863 426 prevE = fmd->Multiplicity(det,ring,sec,strip-1);
bb8a464f 427 if(fmd->IsAngleCorrected())
0210c863 428 prevE = prevE/TMath::Cos(Eta2Theta(fmd->Eta(det,ring,sec,strip-1)));
bb8a464f 429 }
8dc823cc 430 if(strip != nstr - 1)
bb8a464f 431 if(fmd->Multiplicity(det,ring,sec,strip+1) != AliESDFMD::kInvalidMult) {
0210c863 432 nextE = fmd->Multiplicity(det,ring,sec,strip+1);
bb8a464f 433 if(fmd->IsAngleCorrected())
0210c863 434 nextE = nextE/TMath::Cos(Eta2Theta(fmd->Eta(det,ring,sec,strip+1)));
bb8a464f 435 }
8dc823cc 436
0210c863 437 Float_t mergedEnergy = GetMultiplicityOfStrip(mult,eta,prevE,nextE,det,ring,sec,strip);
032f0117 438
da0805e2 439 if(mergedEnergy > 0.3 )
440 nHits[det-1][ir]++;
0210c863 441 foutputESDFMD->SetMultiplicity(det,ring,sec,strip,mergedEnergy);
bb8a464f 442 foutputESDFMD->SetEta(det,ring,sec,strip,eta);
7c3e5162 443
3bb122c7 444 }
445 }
446 }
447 }
032f0117 448
449 //std::cout<<nDead<<std::endl;
450
da0805e2 451 //cluster cut
452 //if(testmult->GetNumberOfSingleClusters() > 15 + nTrackLets)
453 // {fStatus = kFALSE; std::cout<<"FMD : "<<nHits[0][0]<<" "<<nHits[1][0]<<" "<<nHits[1][1]<<" "<<nHits[2][0]<<" "<<nHits[2][1]<<" tracks "<<testmult->GetNumberOfSingleClusters()<<" "<<nTrackLets<<std::endl; return;}
032f0117 454
da0805e2 455 TH2F* hCorrelationFMD23 = (TH2F*)fDiagList->FindObject("hCorrelationFMD23");
456 TH2F* hCorrelationFMD12 = (TH2F*)fDiagList->FindObject("hCorrelationFMD12");
457 TH2F* hCorrelationFMD2diff23 = (TH2F*)fDiagList->FindObject("hCorrelationFMD2diff23");
458 TH2F* hCorrelationFMD3diff23 = (TH2F*)fDiagList->FindObject("hCorrelationFMD3diff23");
459 TH2F* hCorrelationFMD1diff13 = (TH2F*)fDiagList->FindObject("hCorrelationFMD1diff13");
460 TH2F* hCorrelationFMD1diff12 = (TH2F*)fDiagList->FindObject("hCorrelationFMD1diff12");
461
462 TH2F* hCorrelationFMDSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMDSPD");
463 TH2F* hCorrelationFMD1SPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD1SPD");
464 TH2F* hCorrelationFMD2ISPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD2ISPD");
465 TH2F* hCorrelationFMD2OSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD2OSPD");
466 TH2F* hCorrelationFMD3ISPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD3ISPD");
467 TH2F* hCorrelationFMD3OSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD3OSPD");
468
469 TH2F* hCorrelationFMDSPDhits = (TH2F*)fDiagList->FindObject("hCorrelationFMDSPDhits");
470
471 TH2F* hCorrelationSPDTracklets = (TH2F*)fDiagList->FindObject("hCorrelationSPDTracklets");
472 TH2F* hTimeCorrelation = (TH2F*)fDiagList->FindObject("hCorrelationTime");
473 TH2F* hHitsRadius = (TH2F*)fDiagList->FindObject("hCorrelationHitsRadius");
474 TH2F* hHitsX = (TH2F*)fDiagList->FindObject("hCorrelationHitsX");
475 TH2F* hHitsY = (TH2F*)fDiagList->FindObject("hCorrelationHitsY");
476 TH1F* hFMDHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistribution");
477 TH1F* hFMD1HitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD1");
478 TH1F* hFMD2IHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2I");
479 TH1F* hFMD2OHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2O");
480 TH1F* hFMD3IHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3I");
481 TH1F* hFMD3OHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3O");
482 TH1F* hCorrelationFMDGoodtracks = (TH1F*)fDiagList->FindObject("hCorrelationFMDGoodtracks");
483 TH1F* hCorrelationFMDBadtracks = (TH1F*)fDiagList->FindObject("hCorrelationFMDBadtracks");
484 TH1F* hCorrelationGoodbadtracks = (TH1F*)fDiagList->FindObject("hCorrelationGoodbadtracks");
485 TH2F* hCorrelationFMDBgCand = (TH2F*)fDiagList->FindObject("hCorrelationFMDBgCand");
486 TH2F* hCorrelationFMDBgCandRelative = (TH2F*)fDiagList->FindObject("hCorrelationFMDBgCandRelative");
487
488 TH2F* hCorrelationFMDFlatTr = (TH2F*)fDiagList->FindObject("hCorrelationFMDFlatTr");
489 TH2F* hCorrelationFMDRatioFlatTr = (TH2F*)fDiagList->FindObject("hCorrelationFMDRatioFlatTr");
490 TH1F* hTrVtxDistribution = (TH1F*)fDiagList->FindObject("hTrVtxDistribution");
491 TH1F* hTrEtaDistribution = (TH1F*)fDiagList->FindObject("hTrEtaDistribution");
492 TH1F* hTrEtaDistribution2 = (TH1F*)fDiagList->FindObject("hTrEtaDistribution2");
493 hCorrelationFMDSPD->Fill(nTrackLets,nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
494 TH1F* hFlatTracks = (TH1F*)fDiagList->FindObject("hFlatTracks");
495 hCorrelationFMD1SPD->Fill(nTrackLets,nHits[0][0]);
496 hCorrelationFMD2ISPD->Fill(nTrackLets,nHits[1][0]);
497 hCorrelationFMD2OSPD->Fill(nTrackLets,nHits[1][1]);
498 hCorrelationFMD3ISPD->Fill(nTrackLets,nHits[2][0]);
499 hCorrelationFMD3OSPD->Fill(nTrackLets,nHits[2][1]);
500 hCorrelationFMDSPDhits->Fill(testmult->GetNumberOfFiredChips(0),nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
501 hCorrelationSPDTracklets->Fill(testmult->GetNumberOfFiredChips(0),nTrackLets);
502
503 hTimeCorrelation->Fill(delta,nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
504 hCorrelationFMD23->Fill(nHits[1][0]+nHits[1][1],nHits[2][0]+nHits[2][1]);
505 hCorrelationFMD12->Fill(nHits[0][0],nHits[1][0]+nHits[1][1]);
506
507 // if(TMath::Abs(nHits[1]-nHits[2]) > 15 && (nHits[1]+nHits[2]) > 35) {fStatus = kFALSE; std::cout<<"difference : "<<TMath::Abs(nHits[1]-nHits[2])<<std::endl; return;}
508
509 // if(testmult->GetNumberOfFiredChips(0))
510 // if(testmult->GetNumberOfFiredChips(0) > 15 && ((Float_t)nTrackLets / (Float_t)testmult->GetNumberOfFiredChips(0)) < 0.4)
511 // {fStatus = kFALSE; std::cout<<nTrackLets<<" "<<testmult->GetNumberOfFiredChips(0)<<" "<<nHits[0]<<" "<<nHits[1]<<" "<<nHits[2]<<std::endl; return;}
512
513
514 Float_t diff23 = (Float_t)TMath::Abs(nHits[2][1] + nHits[2][0] - nHits[1][1] - nHits[1][0]);
515
516 Float_t diff13 = TMath::Abs(nHits[2][1] + nHits[2][0] - nHits[1][1] - nHits[1][0] - nHits[0][0]);
517 Float_t diff12 = TMath::Abs(nHits[1][0] - nHits[0][0]);
518
519 hCorrelationFMD1diff12->Fill(nHits[0][0], diff12);
520 hCorrelationFMD1diff13->Fill(nHits[0][0], diff13);
521 hCorrelationFMD2diff23->Fill(nHits[1][1], diff23);
522 hCorrelationFMD3diff23->Fill(nHits[2][1], diff23);
523
524 //
525 Float_t nTotalFMDhits = nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1] ;
526 Float_t radius = TMath::Sqrt(TMath::Power(vertex[0] + 0.03715,2) + TMath::Power(vertex[1] - 0.1659,2));
527
528 if(vertex[1] !=0 || vertex[1] !=0) {
529 hHitsRadius->Fill(nTotalFMDhits,radius);
530 hHitsX->Fill(nTotalFMDhits,vertex[0]);
531 hHitsY->Fill(nTotalFMDhits,vertex[1]); }
3bb122c7 532
da0805e2 533 hFMDHitDistribution->Fill(nTotalFMDhits);
534 hFMD1HitDistribution->Fill(nHits[0][0]);
535 hFMD2IHitDistribution->Fill(nHits[1][0]);
536 hFMD2OHitDistribution->Fill(nHits[1][1]);
537 hFMD3IHitDistribution->Fill(nHits[2][0]);
538 hFMD3OHitDistribution->Fill(nHits[2][0]);
539
540 // if(radius > 0.5) {fStatus = kFALSE; std::cout<<"FMD : "<<nTotalFMDhits<<std::endl; foutputESDFMD->Clear(); return;}
541 //if(TMath::Abs(vertex[1] - 0.1659) > 0.1 ) {fStatus = kFALSE; std::cout<<"FMD : "<<nTotalFMDhits<<std::endl; foutputESDFMD->Clear(); return;}
542
543 if(nTrackLets < pars->GetLowSPDLimit() || nTrackLets > pars->GetHighSPDLimit())
544 {fStatus = kFALSE; std::cout<<nTrackLets<<" "<<" "<<nHits[0][0]<<" "<<nHits[1][0]<<" "<<nHits[1][1]<<" "<<nHits[2][0]<<" "<<nHits[2][1]<<std::endl; return;}
545
546
547 AliESDtrack* track = 0;
548 Int_t ntracks = fESD->GetNumberOfTracks();
549 Float_t ngood =0, nbad = 0;
550 //std::cout<<" Primary vtx : "<<vertex[0]<<" "<<vertex[1]<<" "<<vertex[2]<<" "<<nTotalFMDhits<<std::endl;
551 Int_t nBgCandidates = 0;
552 Float_t nFlat = 0;
553 for(Int_t i=0;i<ntracks;i++) {
554 track = fESD->GetTrack(i);
555 //std::cout<<track->GetX()-vertex[0]<<" "<<track->GetY()-vertex[1]<<" "<<track->GetZ()-vertex[2]<<std::endl;
556 //std::cout<<track->GetX()<<" "<<track->GetY()<<" "<<track->GetZ()<<std::endl;
557 hTrVtxDistribution->Fill(track->GetZ());
558
559 if(TMath::Abs( track->GetZ()) > 50 && TMath::Abs(track->GetZ()) < 300) { // && TMath::Abs( track->GetY()) < 1)
560 nBgCandidates++;
561 hTrEtaDistribution->Fill(track->Eta());
562 }
563
564
565
566 if(TMath::Abs(track->GetX()-vertex[0]) > 0.3 || TMath::Abs(track->GetY()-vertex[1]) > 0.3 || TMath::Abs(track->GetZ()-vertex[2]) > 0.3) {
567 nbad++;
568 hTrEtaDistribution2->Fill(track->Eta()); }
569 else ngood++;
570
571 if(TMath::Abs(track->Pt()) < 0.1)
572 nFlat++;
573 }
574
575 Float_t ratioFlat = 0;
576 if(fESD->GetNumberOfTracks())
577 ratioFlat = nFlat/(Float_t)fESD->GetNumberOfTracks();
578 hCorrelationFMDFlatTr->Fill(nFlat,nTotalFMDhits);
579 hCorrelationFMDRatioFlatTr->Fill(ratioFlat,nTotalFMDhits);
580 hFlatTracks->Fill(nFlat);
b64db9b1 581
da0805e2 582 // std::cout<<fESD->GetT0zVertex()<<" "<<vertex[2]<<std::endl;
583 Float_t ratioBg = 0;
584 //if(fESD->GetNumberOfTracks() > 0)
585
586 if(fESD->GetNumberOfTracks() > 0)
587 ratioBg = (Float_t)nBgCandidates/(Float_t)fESD->GetNumberOfTracks();
588 hCorrelationFMDBgCand->Fill(nBgCandidates,nTotalFMDhits);
589 hCorrelationFMDBgCandRelative->Fill(ratioBg,nTotalFMDhits);
590
591
aa303f0c 592 // Float_t ratio = (nbad > 0 ? ngood / nbad : 0);
da0805e2 593
594 hCorrelationFMDGoodtracks->Fill(ngood,nTotalFMDhits);
595 hCorrelationFMDBadtracks->Fill(nbad,nTotalFMDhits);
596 hCorrelationGoodbadtracks->Fill(ngood,nbad);
032f0117 597
7c3e5162 598 if(fStandalone) {
599 PostData(0, foutputESDFMD);
600 PostData(1, fEsdVertex);
601 PostData(2, fESD);
a38d1871 602 PostData(3, fDiagList);
7c3e5162 603 }
3bb122c7 604}
daedf077 605
8dc823cc 606//_____________________________________________________________________
607Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
bb8a464f 608 Float_t eta,
0210c863 609 Float_t prevE,
610 Float_t nextE,
bb8a464f 611 UShort_t det,
612 Char_t ring,
41bad769 613 UShort_t /*sec*/,
614 UShort_t /*strip*/) {
0210c863 615 //analyse and perform sharing on one strip
8dc823cc 616 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
6289b3e8 617
0210c863 618 Float_t mergedEnergy = 0;
cc066cb9 619 //Float_t nParticles = 0;
da0805e2 620 Float_t cutLow = 0.3;//0.15;
879c0f42 621
622 Double_t mpv = pars->GetMPV(det,ring,eta);
da0805e2 623
879c0f42 624 Double_t w = pars->GetSigma(det,ring,eta);
625 if (mpv > 100)
626 AliError(Form("FMD%d%c, eta=%f, MPV=%f w=%f", det, ring, eta, mpv, w));
627
628 Float_t cutHigh = mpv - 2 * w;
aa303f0c 629
36f7d459 630 // if(ring == 'I')
631 // cutLow = 0.1;
cc066cb9 632
f58a4769 633 //cutLow = 0;
634 //AliFMDParameters* recopars = AliFMDParameters::Instance();
635 //cutLow = (5*recopars->GetPedestalWidth(det,ring,sec,strip))/(recopars->GetPulseGain(det,ring,sec,strip)*recopars->GetDACPerMIP());
da0805e2 636 //if(foutputESDFMD->GetUniqueID() == kFALSE ) {
637
638 if(mult > 12 || mult < cutLow)
639 {
640 // std::cout<<"rejecting hit in FMD "<<det<<" "<<ring<<std::endl;
641 fSharedThis = kFALSE;
642 fSharedPrev = kFALSE;
643 return 0;
644 }
645
cc066cb9 646
647
f58a4769 648
6289b3e8 649 // Float_t cutPart = pars->GetMPV(det,ring,eta) - 5*pars->GetSigma(det,ring,eta);
0210c863 650 Float_t totalE = mult;
bb8a464f 651
cc066cb9 652
cc066cb9 653 //std::cout<<det<<ring<<" "<<sec<<" "<<strip<<" "<<cutLow<<std::endl;
c78bc12b 654 if(fSharedThis) {
655 fSharedThis = kFALSE;
656 fSharedPrev = kTRUE;
8dc823cc 657 return 0.;
658 }
659
cc066cb9 660 /* if(mult < 0.33*pars->GetMPV(det,ring,eta)) {
6289b3e8 661 fSharedThis = kFALSE;
662 fSharedPrev = kFALSE;
663 return 0;
cc066cb9 664 }*/
0210c863 665 if(mult<nextE && nextE>cutHigh && foutputESDFMD->GetUniqueID() == kTRUE)
cc066cb9 666 {
667 fSharedThis = kFALSE;
668 fSharedPrev = kFALSE;
669 return 0;
670 }
da0805e2 671
bb8a464f 672
0210c863 673 if(prevE > cutLow && prevE < cutHigh && !fSharedPrev ) {
674 totalE += prevE;
8dc823cc 675 }
676
0210c863 677 if(nextE > cutLow && nextE < cutHigh ) {
678 totalE += nextE;
c78bc12b 679 fSharedThis = kTRUE;
8dc823cc 680 }
da0805e2 681 totalE = totalE*TMath::Cos(Eta2Theta(eta));
a38d1871 682 TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
da0805e2 683 if(totalE > cutLow)
684 hEdist->Fill(totalE);
685
bb8a464f 686
0210c863 687 if(totalE > 0) {
cc066cb9 688
0210c863 689 mergedEnergy = totalE;
c78bc12b 690 fSharedPrev = kTRUE;
cc066cb9 691 // if(det == 1 && ring =='I')
c78bc12b 692 }
da0805e2 693 else{
694 fSharedThis = kFALSE;
695 fSharedPrev = kFALSE;
c78bc12b 696 }
da0805e2 697
0210c863 698 // mergedEnergy = mult;
8dc823cc 699
da0805e2 700
701 /* }
702else {
703 TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
704 if(mult > cutLow)
705 fEtotal+=mult;
706 if(mult < cutLow) {
707 mergedEnergy = fEtotal;
708 fEtotal = 0;
709 hEdist->Fill(mergedEnergy);
710
711 }
712
713 }*/
714
0210c863 715 return mergedEnergy;
cc066cb9 716 //}
8dc823cc 717}
da0805e2 718//_____________________________________________________________________
719void AliFMDAnalysisTaskSharing::Terminate(Option_t* /* option*/) {
720
721 TH1F* hFMDHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistribution");
722 TH1F* hFMD1HitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD1");
723 TH1F* hFMD2IHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2I");
724 TH1F* hFMD2OHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2O");
725 TH1F* hFMD3IHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3I");
726 TH1F* hFMD3OHitDistribution = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3O");
727
728 TH1F* hZvtx = (TH1F*)fDiagList->FindObject("hZvtx");
729 for(UShort_t det=1;det<=3;det++) {
730 Int_t nRings = (det==1 ? 1 : 2);
731 for (UShort_t ir = 0; ir < nRings; ir++) {
732 Char_t ring = (ir == 0 ? 'I' : 'O');
733
734 TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_before_sharing_FMD%d%c",det,ring));
735 TH1F* hEdistAfter = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
736 if(hZvtx->GetEntries()) {
737 hEdist->Scale(1./(Float_t)hZvtx->GetEntries());
738 hEdistAfter->Scale(1./(Float_t)hZvtx->GetEntries());
739 }
740
741 }
742
743 }
744 TH1F* hFlatTracks = (TH1F*)fDiagList->FindObject("hFlatTracks");
745 TH1F* hTrVtxDistribution = (TH1F*)fDiagList->FindObject("hTrVtxDistribution");
746 TH1F* hTrEtaDistribution = (TH1F*)fDiagList->FindObject("hTrEtaDistribution");
747 TH1F* hTrEtaDistribution2 = (TH1F*)fDiagList->FindObject("hTrEtaDistribution2");
748 if(hZvtx->GetEntries()) {
749 hFMDHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
750 hFMD1HitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
751 hFMD2IHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
752 hFMD2OHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
753 hFMD3IHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
754 hFMD3OHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
755 hFlatTracks->Scale(1./(Float_t)hZvtx->GetEntries());
756 hTrVtxDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
757 hTrEtaDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
758 hTrEtaDistribution2->Scale(1./(Float_t)hZvtx->GetEntries());
759 }
760
059c7c6b 761 //TH1F* hEventsVtx = (TH1F*)fDiagList->FindObject("hEventsVtx");
762 //TH1F* hEventsTr = (TH1F*)fDiagList->FindObject("hEventsTr");
daedf077 763 //hEventsVtx->Divide(hEventsTr);
da0805e2 764
765
766}
bb8a464f 767//_____________________________________________________________________
0210c863 768Float_t AliFMDAnalysisTaskSharing::Eta2Theta(Float_t eta) const{
769 //convert the eta of a strip to a theta
bb8a464f 770 Float_t theta = 2*TMath::ATan(TMath::Exp(-1*eta));
771
772 if(eta < 0)
773 theta = theta-TMath::Pi();
774
cc066cb9 775 // std::cout<<"From eta2Theta: "<<theta<<" "<<eta<<std::endl;
bb8a464f 776 return theta;
777
778
779
780}
b64db9b1 781
782
783
bb8a464f 784//_____________________________________________________________________
b64db9b1 785void AliFMDAnalysisTaskSharing::ProcessPrimary() {
aa303f0c 786 //Get the unspoiled MC dN/deta before event cuts
df2a9c32 787
b64db9b1 788 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
789 AliMCEvent* mcEvent = eventHandler->MCEvent();
790 if(!mcEvent)
791 return;
0a2f2742 792 fLastTrackByStrip.Reset(-1);
b64db9b1 793 AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
794
795 AliMCParticle* particle = 0;
796
797 AliStack* stack = mcEvent->Stack();
798
799 TH1F* hPrimary = (TH1F*)fDiagList->FindObject("hMultvsEtaNoCuts");
059c7c6b 800 TH1F* hPrimaryNSD = (TH1F*)fDiagList->FindObject("hMultvsEtaNSDNoCuts");
aa303f0c 801 TH1F* hEnergyOfParticles = (TH1F*)fDiagList->FindObject("hEnergyOfParticles");
b64db9b1 802 AliHeader* header = mcEvent->Header();
803 AliGenEventHeader* genHeader = header->GenEventHeader();
bb8a464f 804
850ad94d 805 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
507687cd 806 AliGenDPMjetEventHeader* dpmHeader = dynamic_cast<AliGenDPMjetEventHeader*>(header->GenEventHeader());
807
808
0f6643af 809 Bool_t nsd = kTRUE;
bb8a464f 810
0f6643af 811
507687cd 812 if (!pythiaGenHeader && !dpmHeader) {
813 std::cout<<" no pythia or dpm header!"<<std::endl;
0f6643af 814 }
507687cd 815 else {
816 if(pythiaGenHeader) {
817 Int_t pythiaType = pythiaGenHeader->ProcessType();
818
819 if(pythiaType==92||pythiaType==93)
820 nsd = kFALSE;
821
822 }
823 if(dpmHeader) {
824 Int_t processType = dpmHeader->ProcessType();
825 if(processType == 5 || processType == 6)
826 nsd = kFALSE;
827
828 }
850ad94d 829 }
507687cd 830
b64db9b1 831 TArrayF vertex;
832 genHeader->PrimaryVertex(vertex);
0a2f2742 833
b64db9b1 834 if(TMath::Abs(vertex.At(2)) > pars->GetVtxCutZ())
835 return;
836
837 Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
838 Double_t vertexBinDouble = (vertex.At(2) + pars->GetVtxCutZ()) / delta;
839 Int_t vertexBin = (Int_t)vertexBinDouble;
840
841 Bool_t firstTrack = kTRUE;
059c7c6b 842 Bool_t firstTrackNSD = kTRUE;
bb8a464f 843
b64db9b1 844 Int_t nTracks = stack->GetNprimary();
0a2f2742 845 if(pars->GetProcessHits())
846 nTracks = stack->GetNtrack();
b64db9b1 847 TH1F* nMCevents = (TH1F*)fDiagList->FindObject("nMCEventsNoCuts");
059c7c6b 848 TH1F* nMCeventsNSD = (TH1F*)fDiagList->FindObject("nMCEventsNSDNoCuts");
b64db9b1 849 for(Int_t i = 0 ;i<nTracks;i++) {
7aad0c47 850 particle = (AliMCParticle*) mcEvent->GetTrack(i);
b64db9b1 851 if(!particle)
852 continue;
853
854 if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
855 hPrimary->Fill(particle->Eta());
059c7c6b 856 if(nsd) {
857 hPrimaryNSD->Fill(particle->Eta());
858
859 if(firstTrackNSD) {
860 nMCeventsNSD->Fill(vertexBin);
861 firstTrackNSD = kFALSE;
862 }
863 }
b64db9b1 864
865 TH1F* hPrimVtxBin = (TH1F*)fDiagList->FindObject(Form("primmult_NoCuts_vtxbin%d",vertexBin));
866 hPrimVtxBin->Fill(particle->Eta());
867
868 if(firstTrack) {
869 nMCevents->Fill(vertexBin);
870 firstTrack = kFALSE;
871 }
872
873 }
0a2f2742 874 if(pars->GetProcessHits()) {
875
876 for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) {
877
878 AliTrackReference* ref = particle->GetTrackReference(j);
879 UShort_t det,sec,strip;
880 Char_t ring;
881 if(ref->DetectorId() != AliTrackReference::kFMD)
882 continue;
aa303f0c 883 if(particle->Charge() != 0)
884 hEnergyOfParticles->Fill(particle->E());
885
0a2f2742 886 AliFMDStripIndex::Unpack(ref->UserId(),det,ring,sec,strip);
887 Float_t thisStripTrack = fLastTrackByStrip(det,ring,sec,strip);
888 if(particle->Charge() != 0 && i != thisStripTrack ) {
889 //Double_t x,y,z;
890
891 Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));//-1*TMath::Log(TMath::Tan(0.5*theta));
892 TH1F* hHits = (TH1F*)fDiagList->FindObject(Form("hMCHits_nocuts_FMD%d%c_vtxbin%d",det,ring,vertexBin));
893
894
895 hHits->Fill(eta);
896
897 Float_t nstrips = (ring =='O' ? 256 : 512);
898
899 fLastTrackByStrip(det,ring,sec,strip) = (Float_t)i;
900
901 if(strip >0)
902 fLastTrackByStrip(det,ring,sec,strip-1) = (Float_t)i;
903 if(strip < (nstrips - 1))
904 fLastTrackByStrip(det,ring,sec,strip+1) = (Float_t)i;
905
906 }
907
908
909 }
910
911
912 }
913
b64db9b1 914 }
915
916}
917
3bb122c7 918//_____________________________________________________________________
919//
920// EOF
921//
5be76c19 922// EOF