]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/ChargedHadrons/dNdPt/macros/GenerateCorrMatr_pp.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / dNdPt / macros / GenerateCorrMatr_pp.C
CommitLineData
62e3b4b6 1Double_t GenerateCorrMatr_pp(const char* mcfile, const char* mctask, const char* idstring ,const char* outfile, const char* gifdir = 0)
2//Double_t GenerateCorrMatr_pp()
3{
4//tmp setting
5// const char* mcfile = "/d/alice09/mknichel/train/V007.MC_pp/2011-05-05_2347.7147/mergedPeriods/MC_pp/7TeV/LHC10f6a/mknichel_dNdPtpp_TPCITS.root";
6// const char* mctask = "mknichel_dNdPtpp_TPCITS";
7// const char* idstring = "TPCITS";
8// const char* outfile = "/u/mknichel/alice/dNdPt_pp/temp/corrMatr_LHC10f6a_TPCITS.root";
9// const char* gifdir = "/u/mknichel/alice/dNdPt_pp/temp";
10
11TString fname (mcfile);
12
13TString id (idstring);
14TString taskname = "mknichel_dNdPtpp_" + id;
15TString objname = "dNdPtAnalysis_" + id;
16
17 // settings vor zVertex cut (event and track level)
18 Double_t zVert = 10.0;
19
20 // setting on eta cut (track level)
21 Double_t eta = 0.8;
22
23 // strangeness scaling factor (for secondaries from strange decays)
24 Double_t sscale = 2.0; //ignored
25
26
27 // Set pt binning
28 const Int_t ptNbins = 31;
29 Double_t bins[32] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0,2.4,2.8, 3.0, 4.0, 50.0, 100.0 };
30 Double_t* binsPt = new Double_t[32];
31 for (int ibin=0; ibin<32; ibin++) {binsPt[ibin] = bins[ibin];}
32
33
34
35 Int_t ptNbinsTrackEventCorr = 36;
36 Int_t etaNbins = 30;
37 Int_t zvNbins = 12;
38
39 Double_t binsMultDefault[28] = {-0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5,
40 9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5,
41 19.5,20.5, 21.5, 22.5, 23.5, 24.5, 29.5, 149.5};
42
43 Double_t binsPtTrackEventCorrDefault[37] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,3.0,4.0,50.0};
44
45 Double_t binsPtDefault[69] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.5,5.0,5.5,6.0,6.5,7.0,8.0,9.0,10.0,11.0,12.0,13.0,14.0,15.0,16.0,18.0,20.0,22.0,24.0,26.0,28.0,30.0,32.0,34.0,36.0,40.0,45.0,50.0};
46
47 Double_t binsEtaDefault[31] = {-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5};
48
49 Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
50
51 // binning for corrections (note difference for efficiency and secondaries)
52 if (fname.Contains("2.76TeV") || fname.Contains("900GeV")) {
53
54 const Int_t ptNbinsSecCorr = 28;
55 const Int_t ptNbinsEffCorr = 31;
56 const Int_t etaNbinsCorr = 4;
57 const Int_t zvNbinsCorr = 4;
58 const Int_t multNbinsCorr = 2;
59 Double_t ptBinsSecCorr[29] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.2, 1.4, 1.6, 2.0, 2.4, 3.0, 50.0, 100.0 };
60 Double_t ptBinsEffCorr[32] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.2, 1.4, 1.6, 2.0, 2.4, 2.8, 3.0, 3.4, 4.0, 50.0, 100.0 };
61 Double_t etaBinsCorr[5] = {-0.8,-0.4,0.,0.4,0.8};
62
63 Double_t zvBinsCorr[5] = {-10.,-5.,0.,5.,10.};
64 Double_t multBinsCorr[3] = {-0.5,0.5,500.};
65
66
67 } else {
68
69 const Int_t ptNbinsSecCorr = 38;
70 const Int_t ptNbinsEffCorr = 50;
71 const Int_t etaNbinsCorr = 8;
72 const Int_t zvNbinsCorr = 4;
73 const Int_t multNbinsCorr = 2;
74 Double_t ptBinsSecCorr[39] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,3.0,3.4,4.0,50.0,100.0};
75 Double_t ptBinsEffCorr[51] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6, 3.8, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 10.0, 50.0, 100.0 };
76 Double_t etaBinsCorr[9] = {-0.8,-0.6,-0.4,-0.2,0.,0.2,0.4,0.6,0.8};
77
78 Double_t zvBinsCorr[5] = {-10.,-5.,0.,5.,10.};
79 Double_t multBinsCorr[3] = {-0.5,0.5,500.};
80
81 }
82
83
84 Int_t binsTrackEffCorr[3]={zvNbinsCorr,ptNbinsEffCorr,etaNbinsCorr};
85 Int_t binsTrackSecCorr[3]={zvNbinsCorr,ptNbinsSecCorr,etaNbinsCorr};
86
87 THnSparse* hSTrackEffCorr = new THnSparseF("hSTrackEffCorr","Zv:pT:eta",3,binsTrackEffCorr);
88 hSTrackEffCorr->SetBinEdges(0,zvBinsCorr);
89 hSTrackEffCorr->SetBinEdges(1,ptBinsEffCorr);
90 hSTrackEffCorr->SetBinEdges(2,etaBinsCorr);
91// hSTrackEffCorr->SetBinEdges(3,multBinsCorr);
92 hSTrackEffCorr->GetAxis(0)->SetTitle("Zv (cm)");
93 hSTrackEffCorr->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
94 hSTrackEffCorr->GetAxis(2)->SetTitle("#eta");
95// hSTrackEffCorr->GetAxis(3)->SetTitle("mult MB");
96 hSTrackEffCorr->Sumw2();
97
98 THnSparse* hSTrackSecCorr = new THnSparseF("hSTrackSecCorr","Zv:pT:eta",3,binsTrackSecCorr);
99 hSTrackSecCorr->SetBinEdges(0,zvBinsCorr);
100 hSTrackSecCorr->SetBinEdges(1,ptBinsSecCorr);
101 hSTrackSecCorr->SetBinEdges(2,etaBinsCorr);
102// hSTrackSecCorr->SetBinEdges(3,multBinsCorr);
103 hSTrackSecCorr->GetAxis(0)->SetTitle("Zv (cm)");
104 hSTrackSecCorr->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
105 hSTrackSecCorr->GetAxis(2)->SetTitle("#eta");
106// hSTrackSecCorr->GetAxis(3)->SetTitle("mult MB");
107 hSTrackSecCorr->Sumw2();
108
109
110 // make plots nicer
111 gROOT->SetStyle("Plain");
112 gStyle->SetPalette(1);
113
114 // array for all correction matrices
115 TObjArray* CorrMatr = new TObjArray();
116
117 // array for all control histograms
118 TObjArray* ContHist = new TObjArray();
119
120
121 // load mc information
122 TFile* fmc = TFile::Open(mcfile,"READ");
123 if (!fmc) return -1;
124 TList* lmc = dynamic_cast<TList*>(fmc->Get(taskname.Data()));
125 if (!lmc) return -1;
126 AlidNdPtAnalysis *obj = dynamic_cast<AlidNdPtAnalysis*>(lmc->FindObject(objname.Data()));
127 if (!obj) return -1;
128
129 //Event statistics
130 THnSparse *fRecEventMatrix = obj->GetRecEventMatrix(); //all reconstructed events
131 TH2D* h2RecEventAll = fRecEventMatrix->Projection(0,1)->Clone("h2RecEventAll");
132 fRecEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
133 TH2D* h2RecEvent = fRecEventMatrix->Projection(0,1)->Clone("h2RecEvent");
134 Double_t MCReconstructedEvents = h2RecEvent->Integral();
135 Double_t MCReconstructedEventsAll = h2RecEventAll->Integral();
136 ContHist->Add(h2RecEvent);
137 ContHist->Add(h2RecEventAll);
138
139 THnSparse *fTriggerEventMatrix = obj->GetTriggerEventMatrix(); //all triggered events
140 Double_t MCTriggeredEventsAll0mult = obj->GetTriggerEventMatrix()->Projection(1)->Integral(1,1); // 0 contributers to vertex
141 TH2D* h2TriggerEventAll = fTriggerEventMatrix->Projection(0,1)->Clone("h2TriggerEventAll");
142 fTriggerEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
143 Double_t MCTriggeredEvents0mult = obj->GetTriggerEventMatrix()->Projection(1)->Integral(1,1); // 0 contributers to vertex
144 TH2D* h2TriggerEvent = fTriggerEventMatrix->Projection(0,1)->Clone("h2TriggerEvent");
145 Double_t MCTriggeredEvents = h2TriggerEvent->Integral();
146 Double_t MCTriggeredEventsAll = h2TriggerEventAll->Integral();
147 ContHist->Add(h2TriggerEvent);
148 ContHist->Add(h2TriggerEventAll);
149
150 THnSparse *fGenEventMatrix = obj->GetGenEventMatrix(); //all generated events
151 TH2D* h2GenEventAll = fGenEventMatrix->Projection(0,1)->Clone("h2GenEventAll");
152 fGenEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
153 TH2D* h2GenEvent = fGenEventMatrix->Projection(0,1)->Clone("h2GenEvent");
154 Double_t MCGeneratedEvents = h2GenEvent->Integral();
155 Double_t MCGeneratedEventsAll = h2GenEventAll->Integral();
156 ContHist->Add(h2RecEvent);
157 ContHist->Add(h2RecEventAll);
158
159 THnSparse *fRecEventHist = obj->GetRecEventHist();
160 THnSparse *fRecEventHist2 = obj->GetRecEventHist2();
161
162 printf("=== generated MC events for correction matrices %lf ===\n",MCGeneratedEvents);
163 printf("=== triggered MC events for correction matrices %lf ===\n",MCTriggeredEvents);
164 printf("=== recontructed MC events for correction matrices %lf ===\n",MCReconstructedEvents);
165 printf("\n");
166 printf("=== cut on the zVertex +- %lf ===\n",zVert);
167 printf("=== generated MC events (before zVertex cut) %lf ===\n",MCGeneratedEventsAll);
168 printf("=== triggered MC events (before zVertex cut) %lf ===\n",MCTriggeredEventsAll);
169 printf("=== recontructed MC events (before zVertex cut) %lf ===\n",MCReconstructedEventsAll);
170
171 TH1D* h1MCGeneratedEvents = new TH1D("h1MCGeneratedEvents","h1MCGeneratedEvents",1,0,1);
172 TH1D* h1MCTriggeredEvents = new TH1D("h1MCTriggeredEvents","h1MCTriggeredEvents",1,0,1);
173 TH1D* h1MCReconstructedEvents = new TH1D("h1MCReconstructedEvents","h1MCReconstructedEvents",1,0,1);
174 TH1D* h1MCGeneratedEventsAll = new TH1D("h1MCGeneratedEventsAll","h1MCGeneratedEventsAll",1,0,1);
175 TH1D* h1MCTriggeredEventsAll = new TH1D("h1MCTriggeredEventsAll","h1MCTriggeredEventsAll",1,0,1);
176 TH1D* h1MCReconstructedEventsAll = new TH1D("h1MCReconstructedEventsAll","h1MCReconstructedEventsAll",1,0,1);
177 TH1D* h1MCTriggeredEventsAll0mult = new TH1D("h1MCTriggeredEventsAll0mult","h1MCTriggeredEventsAll0mult",1,0,1);
178 TH1D* h1MCTriggeredEvents0mult = new TH1D("h1MCTriggeredEvents0mult","h1MCTriggeredEvents0mult",1,0,1);
179
180 h1MCGeneratedEvents->Fill(0,MCGeneratedEvents);
181 h1MCGeneratedEvents->SetEntries(MCGeneratedEvents);
182 h1MCTriggeredEvents->Fill(0,MCTriggeredEvents);
183 h1MCTriggeredEvents->SetEntries(MCTriggeredEvents);
184 h1MCReconstructedEvents->Fill(0,MCReconstructedEvents);
185 h1MCReconstructedEvents->SetEntries(MCReconstructedEvents);
186 h1MCGeneratedEventsAll->Fill(0,MCGeneratedEventsAll);
187 h1MCGeneratedEventsAll->SetEntries(MCGeneratedEventsAll);
188 h1MCTriggeredEventsAll->Fill(0,MCTriggeredEventsAll);
189 h1MCTriggeredEventsAll->SetEntries(MCTriggeredEventsAll);
190 h1MCReconstructedEventsAll->Fill(0,MCReconstructedEventsAll);
191 h1MCReconstructedEventsAll->SetEntries(MCReconstructedEventsAll);
192 h1MCTriggeredEventsAll0mult->Fill(0,MCTriggeredEventsAll0mult);
193 h1MCTriggeredEventsAll0mult->SetEntries(MCTriggeredEventsAll0mult);
194 h1MCTriggeredEvents0mult->Fill(0,MCTriggeredEvents0mult);
195 h1MCTriggeredEvents0mult->SetEntries(MCTriggeredEvents0mult);
196
197
198 ContHist->Add(h1MCGeneratedEvents);
199 ContHist->Add(h1MCTriggeredEvents);
200 ContHist->Add(h1MCReconstructedEvents);
201 ContHist->Add(h1MCGeneratedEventsAll);
202 ContHist->Add(h1MCTriggeredEventsAll);
203 ContHist->Add(h1MCReconstructedEventsAll);
204 ContHist->Add(h1MCTriggeredEvents0mult);
205 ContHist->Add(h1MCTriggeredEventsAll0mult);
206
207
208 // efficienfy and correction matrices for tigger and vertex efficiency
209 TH2D* h2EventTriggerEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEventAll,h2GenEventAll,"h2EventTriggerEffAll");
210 TH2D* h2EventTriggerCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2GenEventAll,h2TriggerEventAll,"h2EventTriggerCorrAll");
211 TH2D* h2EventTriggerEff = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEvent,h2GenEvent,"h2EventTriggerEff");
212 TH2D* h2EventTriggerCorr = AlidNdPtHelper::GenerateCorrMatrix(h2GenEvent,h2TriggerEvent,"h2EventTriggerCorr");
213
214 TH2D* h2EventRecEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2RecEventAll,h2TriggerEventAll,"h2EventRecEffAll");
215 TH2D* h2EventRecCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEventAll,h2RecEventAll,"h2EventRecCorrAll");
216 TH2D* h2EventRecEff = AlidNdPtHelper::GenerateCorrMatrix(h2RecEvent,h2TriggerEvent,"h2EventRecEff");
217 TH2D* h2EventRecCorr = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEvent,h2RecEvent,"h2EventRecCorr");
218
219 TH2D* h2EventEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2RecEventAll,h2GenEventAll,"h2EventEffAll");
220 TH2D* h2EventCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2GenEventAll,h2RecEventAll,"h2EventCorrAll");
221 TH2D* h2EventEff = AlidNdPtHelper::GenerateCorrMatrix(h2RecEvent,h2GenEvent,"h2EventEff");
222 TH2D* h2EventCorr = AlidNdPtHelper::GenerateCorrMatrix(h2GenEvent,h2RecEvent,"h2EventCorr");
223
224 CorrMatr->Add(h2EventTriggerEffAll);
225 CorrMatr->Add(h2EventTriggerCorrAll);
226 CorrMatr->Add(h2EventTriggerEff);
227 CorrMatr->Add(h2EventTriggerCorr);
228 CorrMatr->Add(h2EventRecEffAll);
229 CorrMatr->Add(h2EventRecCorrAll);
230 CorrMatr->Add(h2EventRecEff);
231 CorrMatr->Add(h2EventRecCorr);
232 CorrMatr->Add(h2EventEffAll);
233 CorrMatr->Add(h2EventCorrAll);
234 CorrMatr->Add(h2EventEff);
235 CorrMatr->Add(h2EventCorr);
236
237 THnSparse* hSMultCorrelation = obj->GetEventMultCorrelationMatrix()->Clone("hSMultCorrelation");
238 CorrMatr->Add(hSMultCorrelation);
239 TH2D* h2MultCorrelation = hSMultCorrelation->Projection(0,1)->Clone("h2MultCorrelation");
240 CorrMatr->Add(h2MultCorrelation);
241
242 // all recontructed
243 THnSparse *fRecTrackMatrix = obj->GetRecTrackMatrix();
244 fRecTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
245 fRecTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
246 TH3D* h3RecTrack = fRecTrackMatrix->Projection(0,1,2)->Clone("h3RecTrack");
247 TH2D* h2RecTrack_zv_pt = fRecTrackMatrix->Projection(0,1)->Clone("h2RecTrack_zv_pt");
248 TH2D* h2RecTrack_zv_eta = fRecTrackMatrix->Projection(0,2)->Clone("h2RecTrack_zv_eta");
249 TH2D* h2RecTrack_pt_eta = fRecTrackMatrix->Projection(1,2)->Clone("h2RecTrack_pt_eta");
250 TH1D* h1RecTrack_zv = fRecTrackMatrix->Projection(0)->Clone("h1RecTrack_zv");
251 TH1D* h1RecTrack_pt = fRecTrackMatrix->Projection(1)->Clone("h1RecTrack_pt");
252 TH1D* h1RecTrack_eta = fRecTrackMatrix->Projection(2)->Clone("h1RecTrack_eta");
253 Double_t MCReconstructedTracks = h3RecTrack->Integral();
254
255 ContHist->Add(h3RecTrack);
256 ContHist->Add(h2RecTrack_zv_pt);
257 ContHist->Add(h2RecTrack_zv_eta);
258 ContHist->Add(h2RecTrack_pt_eta);
259 ContHist->Add(h1RecTrack_zv);
260 ContHist->Add(h1RecTrack_pt);
261 ContHist->Add(h1RecTrack_eta);
262
263 // recontructed primary tracks
264 THnSparse *fRecPrimTrackMatrix = obj->GetRecPrimTrackMatrix();
265 THnSparse *fRecTrackMatrixScaled = fRecPrimTrackMatrix->Clone("fRecTrackMatrixScaled"); //used later for secondaries scaling
266 fRecPrimTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
267 fRecPrimTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
268 TH3D* h3RecPrimTrack = fRecPrimTrackMatrix->Projection(0,1,2)->Clone("h3RecPrimTrack");
269 TH2D* h2RecPrimTrack_zv_pt = fRecPrimTrackMatrix->Projection(0,1)->Clone("h2RecPrimTrack_zv_pt");
270 TH2D* h2RecPrimTrack_zv_eta = fRecPrimTrackMatrix->Projection(0,2)->Clone("h2RecPrimTrack_zv_eta");
271 TH2D* h2RecPrimTrack_pt_eta = fRecPrimTrackMatrix->Projection(1,2)->Clone("h2RecPrimTrack_pt_eta");
272 TH1D* h1RecPrimTrack_zv = fRecPrimTrackMatrix->Projection(0)->Clone("h1RecPrimTrack_zv");
273 TH1D* h1RecPrimTrack_pt = fRecPrimTrackMatrix->Projection(1)->Clone("h1RecPrimTrack_pt");
274 TH1D* h1RecPrimTrack_eta = fRecPrimTrackMatrix->Projection(2)->Clone("h1RecPrimTrack_eta");
275 Double_t MCReconstructedPrimTracks = h3RecPrimTrack->Integral();
276
277 ContHist->Add(h3RecPrimTrack);
278 ContHist->Add(h2RecPrimTrack_zv_pt);
279 ContHist->Add(h2RecPrimTrack_zv_eta);
280 ContHist->Add(h2RecPrimTrack_pt_eta);
281 ContHist->Add(h1RecPrimTrack_zv);
282 ContHist->Add(h1RecPrimTrack_pt);
283 ContHist->Add(h1RecPrimTrack_eta);
284
285 // recontructed secondary tracks
286 THnSparse *fRecSecTrackMatrix = obj->GetRecSecTrackMatrix();
287 THnSparse *fRecSecTrackMatrixScaled = fRecSecTrackMatrix->Clone("fRecSecTrackMatrixScaled"); //used later for secondaries scaling
288 fRecSecTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
289 fRecSecTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
290 TH3D* h3RecSecTrack = fRecSecTrackMatrix->Projection(0,1,2)->Clone("h3RecSecTrack");
291 TH2D* h2RecSecTrack_zv_pt = fRecSecTrackMatrix->Projection(0,1)->Clone("h2RecSecTrack_zv_pt");
292 TH2D* h2RecSecTrack_zv_eta = fRecSecTrackMatrix->Projection(0,2)->Clone("h2RecSecTrack_zv_eta");
293 TH2D* h2RecSecTrack_pt_eta = fRecSecTrackMatrix->Projection(1,2)->Clone("h2RecSecTrack_pt_eta");
294 TH1D* h1RecSecTrack_zv = fRecSecTrackMatrix->Projection(0)->Clone("h1RecSecTrack_zv");
295 TH1D* h1RecSecTrack_pt = fRecSecTrackMatrix->Projection(1)->Clone("h1RecSecTrack_pt");
296 TH1D* h1RecSecTrack_eta = fRecSecTrackMatrix->Projection(2)->Clone("h1RecSecTrack_eta");
297 Double_t MCReconstructedSecTracks = h3RecSecTrack->Integral();
298
299 ContHist->Add(h3RecSecTrack);
300 ContHist->Add(h2RecSecTrack_zv_pt);
301 ContHist->Add(h2RecSecTrack_zv_eta);
302 ContHist->Add(h2RecSecTrack_pt_eta);
303 ContHist->Add(h1RecSecTrack_zv);
304 ContHist->Add(h1RecSecTrack_pt);
305 ContHist->Add(h1RecSecTrack_eta);
306
307 // generated primary tracks
308 THnSparse *fGenPrimTrackMatrix = obj->GetGenPrimTrackMatrix();
309 fGenPrimTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
310 fGenPrimTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
311 TH3D* h3GenPrimTrack = fGenPrimTrackMatrix->Projection(0,1,2)->Clone("h3GenPrimTrack");
312 TH2D* h2GenPrimTrack_zv_Pt = fGenPrimTrackMatrix->Projection(0,1)->Clone("h2GenPrimTrack_zv_pt");
313 TH2D* h2GenPrimTrack_zv_eta = fGenPrimTrackMatrix->Projection(0,2)->Clone("h2GenPrimTrack_zv_eta");
314 TH2D* h2GenPrimTrack_pt_eta = fGenPrimTrackMatrix->Projection(1,2)->Clone("h2GenPrimTrack_pt_eta");
315 TH1D* h1GenPrimTrack_zv = fGenPrimTrackMatrix->Projection(0)->Clone("h1GenPrimTrack_zv");
316 TH1D* h1GenPrimTrack_pt = fGenPrimTrackMatrix->Projection(1)->Clone("h1GenPrimTrack_pt");
317 TH1D* h1GenPrimTrack_eta = fGenPrimTrackMatrix->Projection(2)->Clone("h1GenPrimTrack_eta");
318 Double_t MCGeneratedPrimTracks = h3GenPrimTrack->Integral();
319
320 ContHist->Add(h3GenPrimTrack);
321 ContHist->Add(h2GenPrimTrack_zv_pt);
322 ContHist->Add(h2GenPrimTrack_zv_eta);
323 ContHist->Add(h2GenPrimTrack_pt_eta);
324 ContHist->Add(h1GenPrimTrack_zv);
325 ContHist->Add(h1GenPrimTrack_pt);
326 ContHist->Add(h1GenPrimTrack_eta);
327 printf("\n");
328 printf("==============================================================\n");
329 printf("=== recontructed MC tracks %lf ===\n",MCReconstructedTracks);
330 printf("=== recontructed MC secondary tracks %lf ===\n",MCReconstructedSecTracks);
331 printf("=== recontructed MC primary tracks %lf ===\n",MCReconstructedPrimTracks);
332 printf("=== generated MC primary track %lf ===\n",MCGeneratedPrimTracks);
333 printf("==============================================================\n");
334 printf("\n");
335
336 // mc truth histogram (for self-consistency check)
337 TH1D* dNdPt_MC = (TH1D*) h1GenPrimTrack_pt->Clone("dNdPt_MC");
338 TH1D* dNdEta_MC =(TH1D*) h1GenPrimTrack_eta->Clone("dNdEta_MC");
339 // normalization and finalization
340 // 1/N_evt 1/(2 pi pt) 1/width 1/etarange
341 for (int ii=1; ii <= dNdPt_MC->GetNbinsX() ;ii++) {
342 Double_t pt = dNdPt_MC->GetBinCenter(ii);
343 Double_t width = dNdPt_MC->GetBinWidth(ii);
344 Double_t val = dNdPt_MC->GetBinContent(ii);
345 Double_t err = dNdPt_MC->GetBinError(ii);
346 Double_t cval = (val)/(width * 2.0 * TMath::Pi() * 1.6 * MCGeneratedEvents * pt);
347 Double_t cerr = (err)/(width * 2.0 * TMath::Pi() * 1.6 * MCGeneratedEvents * pt);
348 dNdPt_MC->SetBinContent(ii,cval);
349 dNdPt_MC->SetBinError(ii,cerr);
350 }
351 dNdPt_MC->SetMarkerStyle(21);
352 dNdPt_MC->SetTitle("; p_{T} (GeV/c) ; 1/N_{evt} 1/(2#pi p_{T}) (d^{2}N_{ch})/(d#eta dp_{T})^{-2}");
353 dNdEta_MC->Scale(1,"width");
354 dNdEta_MC->Scale(1./MCGeneratedEvents);
355 ContHist->Add(dNdPt_MC);
356 ContHist->Add(dNdEta_MC);
357
358
359 // Rebin for corrections (pt)
360 cout << "rebinning...." << endl;
361
362 TH1D* h1RecPrimTrack_pt_Rebin = h1RecPrimTrack_pt->Rebin(ptNbins,"h1RecPrimTrack_pt_Rebin",binsPt);
363 TH1D* h1GenPrimTrack_pt_Rebin = h1GenPrimTrack_pt->Rebin(ptNbins,"h1GenPrimTrack_pt_Rebin",binsPt);
364 ContHist->Add(h1RecPrimTrack_pt_Rebin);
365 ContHist->Add(h1GenPrimTrack_pt_Rebin);
366
367
368// THnSparse *fSparseTriggerTrackEvent = obj->GetTriggerTrackEventMatrix();//Tracks from triggered events
369// THnSparse *fSparseVtxTrackEvent = obj->GetRecTrackEventMatrix();//Tracks from events with rec. vtx
370// THnSparse *fSparseGenTrackEvent = obj->GetGenTrackEventMatrix();//generated TrackEvent matrix
371
372 // tracking efficiencies + corrections
373 TH2D* h2TrackEff_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_zv_pt,h2GenPrimTrack_zv_pt,"h2TrackEff_zv_pt");
374 TH2D* h2TrackCorr_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_zv_pt,h2RecPrimTrack_zv_pt,"h2TrackCorr_zv_pt");
375 TH2D* h2TrackEff_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_zv_eta,h2GenPrimTrack_zv_eta,"h2TrackEff_zv_eta");
376 TH2D* h2TrackCorr_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_zv_eta,h2RecPrimTrack_zv_eta,"h2TrackCorr_zv_eta");
377 TH2D* h2TrackEff_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_pt_eta,h2GenPrimTrack_pt_eta,"h2TrackEff_pt_eta");
378 TH2D* h2TrackCorr_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_pt_eta,h2RecPrimTrack_pt_eta,"h2TrackCorr_pt_eta");
379
380
381 TH1D* h1TrackEff_zv = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_zv,h1GenPrimTrack_zv,"h1TrackEff_zv");
382 TH1D* h1TrackCorr_zv = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_zv,h1RecPrimTrack_zv,"h1TrackCorr_zv");
383 TH1D* h1TrackEff_pt = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_pt_Rebin,h1GenPrimTrack_pt_Rebin,"h1TrackEff_pt");
384 TH1D* h1TrackCorr_pt = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_pt_Rebin,h1RecPrimTrack_pt_Rebin,"h1TrackCorr_pt");
385 TH1D* h1TrackEff_eta = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_eta,h1GenPrimTrack_eta,"h1TrackEff_eta");
386 TH1D* h1TrackCorr_eta = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_eta,h1RecPrimTrack_eta,"h1TrackCorr_eta");
387 CorrMatr->Add(h2TrackEff_zv_pt);
388 CorrMatr->Add(h2TrackCorr_zv_pt);
389 CorrMatr->Add(h2TrackEff_zv_eta);
390 CorrMatr->Add(h2TrackCorr_zv_eta);
391 CorrMatr->Add(h2TrackEff_pt_eta);
392 CorrMatr->Add(h2TrackCorr_pt_eta);
393 CorrMatr->Add(h1TrackEff_zv);
394 CorrMatr->Add(h1TrackCorr_zv);
395 CorrMatr->Add(h1TrackEff_pt);
396 CorrMatr->Add(h1TrackCorr_pt);
397 CorrMatr->Add(h1TrackEff_eta);
398 CorrMatr->Add(h1TrackCorr_eta);
399
400 // scale the secondaries before calculating correction matrices
401 for (Long64_t i = 0; i < fRecSecTrackMatrixScaled->GetNbins(); i++) {
402 Int_t c[3];
403 Double_t val = fRecSecTrackMatrixScaled->GetBinContent(i,c);
404 Double_t err = fRecSecTrackMatrixScaled->GetBinError(i);
405 Double_t pt = fRecSecTrackMatrixScaled->GetAxis(1)->GetBinCenter(c[1]);
406 Double_t scale = GetStrangenessCorrFactor(pt,sscale);
407// Double_t scale = AlidNdPtHelper::GetStrangenessCorrFactor(pt);
408 fRecSecTrackMatrixScaled->SetBinContent(c,val*scale);
409 fRecSecTrackMatrixScaled->SetBinError(c,err*scale);
410 }
411
412 // for correct determination of secondaries contamination, also the total total tracks have to be scaled
413 // this is done by taking primaries and adding the scaled secondaries
414 fRecTrackMatrixScaled->Add(fRecSecTrackMatrixScaled);
415
416 fRecSecTrackMatrixScaled->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
417 fRecSecTrackMatrixScaled->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
418
419 TH3D* h3RecSecTrackScaled = fRecSecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecSecTrackScaled");
420 TH2D* h2RecSecTrackScaled_zv_pt = fRecSecTrackMatrixScaled->Projection(0,1)->Clone("h2RecSecTrackScaled_zv_pt");
421 TH2D* h2RecSecTrackScaled_zv_eta = fRecSecTrackMatrixScaled->Projection(0,2)->Clone("h2RecSecTrackScaled_zv_eta");
422 TH2D* h2RecSecTrackScaled_pt_eta = fRecSecTrackMatrixScaled->Projection(1,2)->Clone("h2RecSecTrackScaled_pt_eta");
423 TH1D* h1RecSecTrackScaled_zv = fRecSecTrackMatrixScaled->Projection(0)->Clone("h1RecSecTrackScaled_zv");
424 TH1D* h1RecSecTrackScaled_pt = fRecSecTrackMatrixScaled->Projection(1)->Clone("h1RecSecTrackScaled_pt");
425 TH1D* h1RecSecTrackScaled_eta = fRecSecTrackMatrixScaled->Projection(2)->Clone("h1RecSecTrackScaled_eta");
426
427 ContHist->Add(h3RecSecTrackScaled);
428 ContHist->Add(h2RecSecTrackScaled_zv_pt);
429 ContHist->Add(h2RecSecTrackScaled_zv_eta);
430 ContHist->Add(h2RecSecTrackScaled_pt_eta);
431 ContHist->Add(h1RecSecTrackScaled_zv);
432 ContHist->Add(h1RecSecTrackScaled_pt);
433 ContHist->Add(h1RecSecTrackScaled_eta);
434
435 fRecTrackMatrixScaled->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
436 fRecTrackMatrixScaled->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
437
438 TH3D* h3RecTrackScaled = fRecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecTrackScaled");
439 TH2D* h2RecTrackScaled_zv_pt = fRecTrackMatrixScaled->Projection(0,1)->Clone("h2RecTrackScaled_zv_pt");
440 TH2D* h2RecTrackScaled_zv_eta = fRecTrackMatrixScaled->Projection(0,2)->Clone("h2RecTrackScaled_zv_eta");
441 TH2D* h2RecTrackScaled_pt_eta = fRecTrackMatrixScaled->Projection(1,2)->Clone("h2RecTrackScaled_pt_eta");
442 TH1D* h1RecTrackScaled_zv = fRecTrackMatrixScaled->Projection(0)->Clone("h1RecTrackScaled_zv");
443 TH1D* h1RecTrackScaled_pt = fRecTrackMatrixScaled->Projection(1)->Clone("h1RecTrackScaled_pt");
444 TH1D* h1RecTrackScaled_eta = fRecTrackMatrixScaled->Projection(2)->Clone("h1RecTrackScaled_eta");
445
446 ContHist->Add(h3RecTrackScaled);
447 ContHist->Add(h2RecTrackScaled_zv_pt);
448 ContHist->Add(h2RecTrackScaled_zv_eta);
449 ContHist->Add(h2RecTrackScaled_pt_eta);
450 ContHist->Add(h1RecTrackScaled_zv);
451 ContHist->Add(h1RecTrackScaled_pt);
452 ContHist->Add(h1RecTrackScaled_eta);
453
454 // Rebin for corrections (pt)
455 TH1D* h1RecTrackScaled_pt_Rebin = h1RecTrackScaled_pt->Rebin(ptNbins,"h1RecTrackScaled_pt_Rebin",binsPt);
456 TH1D* h1RecSecTrackScaled_pt_Rebin = h1RecSecTrackScaled_pt->Rebin(ptNbins,"h1RecSecTrackScaled_pt_Rebin",binsPt);
457 ContHist->Add(h1RecTrackScaled_pt_Rebin);
458 ContHist->Add(h1RecSecTrackScaled_pt_Rebin);
459
460 // create histograms for secondaries contamination and correction
461
462 TH2D* h2SecCont_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_zv_pt,h2RecTrackScaled_zv_pt,"h2SecCont_zv_pt");
463 TH2D* h2SecCorr_zv_pt = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_zv_pt,h2RecTrackScaled_zv_pt,"h2SecCorr_zv_pt");
464 TH2D* h2SecCont_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_zv_eta,h2RecTrackScaled_zv_eta,"h2SecCont_zv_eta");
465 TH2D* h2SecCorr_zv_eta = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_zv_eta,h2RecTrackScaled_zv_eta,"h2SecCorr_zv_eta");
466 TH2D* h2SecCont_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_pt_eta,h2RecTrackScaled_pt_eta,"h2SecCont_pt_eta");
467 TH2D* h2SecCorr_pt_eta = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_pt_eta,h2RecTrackScaled_pt_eta,"h2SecCorr_pt_eta");
468 TH1D* h1SecCont_zv = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_zv,h1RecTrackScaled_zv,"h1SecCont_zv");
469 TH1D* h1SecCorr_zv = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_zv,h1RecTrackScaled_zv,"h1SecCorr_zv");
470 TH1D* h1SecCont_pt = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_pt_Rebin,h1RecTrackScaled_pt_Rebin,"h1SecCont_pt");
471 TH1D* h1SecCorr_pt = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_pt_Rebin,h1RecTrackScaled_pt_Rebin,"h1SecCorr_pt");
472 TH1D* h1SecCont_eta = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_eta,h1RecTrackScaled_eta,"h1SecCont_eta");
473 TH1D* h1SecCorr_eta = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_eta,h1RecTrackScaled_eta,"h1SecCorr_eta");
474
475 CorrMatr->Add(h2SecCont_zv_pt);
476 CorrMatr->Add(h2SecCorr_zv_pt);
477 CorrMatr->Add(h2SecCont_zv_eta);
478 CorrMatr->Add(h2SecCorr_zv_eta);
479 CorrMatr->Add(h2SecCont_pt_eta);
480 CorrMatr->Add(h2SecCorr_pt_eta);
481 CorrMatr->Add(h1SecCont_zv);
482 CorrMatr->Add(h1SecCorr_zv);
483 CorrMatr->Add(h1SecCont_pt);
484 CorrMatr->Add(h1SecCorr_pt);
485 CorrMatr->Add(h1SecCont_eta);
486 CorrMatr->Add(h1SecCorr_eta);
487
488 //vertex distribution for events with no vertex
489 fGenEventMatrix->GetAxis(0)->SetRange();
490 fGenEventMatrix->GetAxis(1)->SetRange(1,1);
491 TH1D* h1GenEventMatrix_bin0_zv = (TH1D*) fGenEventMatrix->Projection(0)->Clone("h1GenEventMatrix_bin0_zv");
492 ContHist->Add(h1GenEventMatrix_bin0_zv);
493
494 fTriggerEventMatrix->GetAxis(0)->SetRange();
495 fTriggerEventMatrix->GetAxis(1)->SetRange(1,1);
496 TH1D* h1TriggerEventMatrix_bin0_zv = (TH1D*) fTriggerEventMatrix->Projection(0)->Clone("h1TriggerEventMatrix_bin0_zv");
497 ContHist->Add(h1TriggerEventMatrix_bin0_zv);
498
499
500 fRecEventMatrix->GetAxis(0)->SetRange();
501 fRecEventMatrix->GetAxis(1)->SetRange();
502 TH1D* h1RecEventMatrix_zv = (TH1D*) fRecEventMatrix->Projection(0)->Clone("h1RecEventMatrix_zv");
503 h1RecEventMatrix_zv->Scale(1./h1RecEventMatrix_zv->Integral());
504 //ContHist->Add(h1TriggerEventMatrix_zv);
505
506 TH1D* h1TriggerEff_bin0_zv = AlidNdPtHelper::GenerateCorrMatrix(h1TriggerEventMatrix_bin0_zv,h1GenEventMatrix_bin0_zv,"h1TriggerEff_bin0_zv");
507 TH1D* h1TriggerCorr_bin0_zv = AlidNdPtHelper::GenerateCorrMatrix(h1GenEventMatrix_bin0_zv,h1TriggerEventMatrix_bin0_zv,"h1TriggerCorr_bin0_zv");
508 CorrMatr->Add(h1TriggerEff_bin0_zv);
509 CorrMatr->Add(h1TriggerCorr_bin0_zv);
510
511 TH1D* h1Ratio_zv = (TH1D*) h1GenEventMatrix_bin0_zv->Clone("h1Ratio_zv");
512 h1Ratio_zv->Scale(1./h1Ratio_zv->Integral());
513 h1Ratio_zv->Divide(h1RecEventMatrix_zv);
514 CorrMatr->Add(h1Ratio_zv);
515
516
517 // aded to correct sperately for zvert shape of triggered and untriggered bin0 events
518 fRecEventMatrix->GetAxis(0)->SetRange();
519 fRecEventMatrix->GetAxis(1)->SetRange();
520 fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0,0);
521 fGenEventMatrix->GetAxis(1)->SetRangeUser(0,0);
522
523 TH1D* rec_zv = fRecEventMatrix->Projection(0)->Clone("rec_zv");
524 TH1D* trig0_zv = fTriggerEventMatrix->Projection(0)->Clone("trig0_zv");
525 TH1D* corr_shape_trig0_zv = AlidNdPtHelper::GenerateCorrMatrix(trig0_zv,rec_zv,"corr_shape_trig0_zv");
526 CorrMatr->Add(corr_shape_trig0_zv);
527
528 TH1D* gen0_zv = fGenEventMatrix->Projection(0)->Clone("gen0_zv");
529 TH1D* gen0notrig_zv = gen0_zv->Clone("gen0notrig_zv");
530 gen0notrig_zv->Add(trig0_zv,-1.);
531 TH1D* corr_shape_notrig0_zv = AlidNdPtHelper::GenerateCorrMatrix(gen0notrig_zv,rec_zv,"corr_shape_notrig0_zv");
532 CorrMatr->Add(corr_shape_notrig0_zv);
533
534
535 // efficienfy and correction for THnSparse
536 THnSparse* hSRecSecTrackMatrixScaled = AlidNdPtHelper::RebinTHnSparse(fRecSecTrackMatrixScaled,hSTrackSecCorr,"hSRecSecTrackMatrixScaled");
537 THnSparse* hSRecTrackMatrixScaled = AlidNdPtHelper::RebinTHnSparse(fRecTrackMatrixScaled,hSTrackSecCorr,"hSRecTrackMatrixScaled");
538 THnSparse* hSRecPrimTrackMatrix = AlidNdPtHelper::RebinTHnSparse(fRecPrimTrackMatrix,hSTrackEffCorr,"hSRecPrimTrackMatrix");
539 THnSparse* hSGenPrimTrackMatrix = AlidNdPtHelper::RebinTHnSparse(fGenPrimTrackMatrix,hSTrackEffCorr,"hSGenPrimTrackMatrix");
540 ContHist->Add(hSRecSecTrackMatrixScaled);
541 ContHist->Add(hSRecTrackMatrixScaled);
542 ContHist->Add(hSRecPrimTrackMatrix);
543 ContHist->Add(hSGenPrimTrackMatrix);
544
545 THnSparse* hSTrackEff = AlidNdPtHelper::GenerateCorrMatrix(hSRecPrimTrackMatrix,hSGenPrimTrackMatrix,"hSTrackEff");
546 THnSparse* hSTrackCorr = AlidNdPtHelper::GenerateCorrMatrix(hSGenPrimTrackMatrix,hSRecPrimTrackMatrix,"hSTrackCorr");
547 THnSparse* hSSecCont = AlidNdPtHelper::GenerateCorrMatrix(hSRecSecTrackMatrixScaled,hSRecTrackMatrixScaled,"hSSecCont");
548 THnSparse* hSSecCorr = AlidNdPtHelper::GenerateContCorrMatrix(hSRecSecTrackMatrixScaled,hSRecTrackMatrixScaled,"hSSecCorr");
549 CorrMatr->Add(hSTrackEff);
550 CorrMatr->Add(hSTrackCorr);
551 CorrMatr->Add(hSSecCont);
552 CorrMatr->Add(hSSecCorr);
553
554 // create th3 from thnsparse, used for corrections
555 TH3D* h3RecPrimTrack_Rebin = hSRecPrimTrackMatrix->Projection(0,1,2)->Clone("h3RecPrimTrack_Rebin");
556 TH3D* h3GenPrimTrack_Rebin = hSGenPrimTrackMatrix->Projection(0,1,2)->Clone("h3GenPrimTrack_Rebin");
557 TH3D* h3RecSecTrackScaled_Rebin = hSRecSecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecSecTrackScaled_Rebin");
558 TH3D* h3RecTrackScaled_Rebin = hSRecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecTrackScaled_Rebin");
559
560 TH3D* h3TrackEff = AlidNdPtHelper::GenerateCorrMatrix(h3RecPrimTrack_Rebin,h3GenPrimTrack_Rebin,"h3TrackEff");
561 TH3D* h3TrackCorr = AlidNdPtHelper::GenerateCorrMatrix(h3GenPrimTrack_Rebin,h3RecPrimTrack_Rebin,"h3TrackCorr");
562 TH3D* h3SecCont = AlidNdPtHelper::GenerateCorrMatrix(h3RecSecTrackScaled_Rebin,h3RecTrackScaled_Rebin,"h3SecCont");
563 TH3D* h3SecCorr = AlidNdPtHelper::GenerateContCorrMatrix(h3RecSecTrackScaled_Rebin,h3RecTrackScaled_Rebin,"h3SecCorr");
564 CorrMatr->Add(h3TrackEff);
565 CorrMatr->Add(h3TrackCorr);
566 CorrMatr->Add(h3SecCont);
567 CorrMatr->Add(h3SecCorr);
568
569 // check for trigger/vertex bias on track level
570 THnSparse* hSRecTrackEventMatrix = (THnSparse*) obj->GetRecTrackEventMatrix()->Clone("hSRecTrackEventMatrix");
571 THnSparse* hSGenTrackEventMatrix = (THnSparse*) obj->GetGenTrackEventMatrix()->Clone("hSGenTrackEventMatrix");
572 THnSparse* hSTriggerTrackEventMatrix = (THnSparse*) obj->GetTriggerTrackEventMatrix()->Clone("hSTriggerTrackEventMatrix");
573
574 hSRecTrackEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
575 hSRecTrackEventMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
576 hSGenTrackEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
577 hSGenTrackEventMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
578 hSTriggerTrackEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer
579 hSTriggerTrackEventMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta
580
581
582 TH1D* h1TriggerBiasCorr_zv = AlidNdPtHelper::GenerateCorrMatrix(hSTriggerTrackEventMatrix->Projection(0),hSGenTrackEventMatrix->Projection(0),"h1TriggerBiasCorr_zv");
583 TH1D* h1TriggerBiasCorr_pt = AlidNdPtHelper::GenerateCorrMatrix(hSTriggerTrackEventMatrix->Projection(1),hSGenTrackEventMatrix->Projection(1),"h1TriggerBiasCorr_pt");
584 TH1D* h1TriggerBiasCorr_eta = AlidNdPtHelper::GenerateCorrMatrix(hSTriggerTrackEventMatrix->Projection(2),hSGenTrackEventMatrix->Projection(2),"h1TriggerBiasCorr_eta");
585
586 TH1D* h1VertexBiasCorr_zv = AlidNdPtHelper::GenerateCorrMatrix(hSRecTrackEventMatrix->Projection(0),hSTriggerTrackEventMatrix->Projection(0),"h1VertexBiasCorr_zv");
587 TH1D* h1VertexBiasCorr_pt = AlidNdPtHelper::GenerateCorrMatrix(hSRecTrackEventMatrix->Projection(1),hSTriggerTrackEventMatrix->Projection(1),"h1VertexBiasCorr_pt");
588 TH1D* h1VertexBiasCorr_eta = AlidNdPtHelper::GenerateCorrMatrix(hSRecTrackEventMatrix->Projection(2),hSTriggerTrackEventMatrix->Projection(2),"h1VertexBiasCorr_eta");
589
590 ContHist->Add(h1TriggerBiasCorr_zv);
591 ContHist->Add(h1TriggerBiasCorr_pt);
592 ContHist->Add(h1TriggerBiasCorr_eta);
593 ContHist->Add(h1VertexBiasCorr_zv);
594 ContHist->Add(h1VertexBiasCorr_pt);
595 ContHist->Add(h1VertexBiasCorr_eta);
596
597
598
599 // plot pictures and save to gifdir
600 for (i=0; i < CorrMatr->LastIndex(); i++) {
601 TCanvas* ctmp = PlotHist(CorrMatr->At(i),idstring);
602 if (gifdir && ctmp) {
603 TString gif(gifdir);
604 gif += '/';
605 gif += ctmp->GetName();
606 gif += ".gif";
607 ctmp->SaveAs(gif.Data(),"gif");
608 delete ctmp;
609 }
610 }
611 for (i=0; i < ContHist->LastIndex(); i++) {
612 TCanvas* ctmp = PlotHist(ContHist->At(i),idstring);
613 if (gifdir && ctmp) {
614 TString gif(gifdir);
615 gif += '/';
616 gif += ctmp->GetName();
617 gif += ".gif";
618 ctmp->SaveAs(gif.Data(),"gif");
619 delete ctmp;
620 }
621 }
622
623 // save all correction matrices and control histograms to file
624 if (!outfile) { return; }
625 TFile *out = TFile::Open(outfile,"RECREATE");
626 CorrMatr->Write();
627 ContHist->Write();
628 out->Close();
629
630 return MCReconstructedEvents;
631
632}
633
634
635//_____________________________________________________________________________
636Double_t GetStrangenessCorrFactor(Double_t pt, Double_t s)
637{
638 // data driven correction factor for secondaries (PbPb)
639
640 if (pt <= 0.17) return 1.0;
641 if (pt <= 0.4) return GetLinearInterpolationValue(0.17,1.0,0.4,1.07, pt);
642 if (pt <= 0.6) return GetLinearInterpolationValue(0.4,1.07,0.6,1.25, pt);
643 if (pt <= 1.2) return GetLinearInterpolationValue(0.6,1.25,1.2,1.5, pt);
644 return 1.5;
645}
646
647
648//___________________________________________________________________________
649Double_t GetLinearInterpolationValue(const Double_t x1,const Double_t y1,const Double_t x2,const Double_t y2, const Double_t pt)
650{
651 //
652 // linear interpolation
653 //
654 return ((y2-y1)/(x2-x1))*pt+(y2-(((y2-y1)/(x2-x1))*x2));
655}
656
657//___________________________________________________________________________
658TCanvas* PlotHist(TObject* hobj, const char* label=0)
659{
660 TH1* h = dynamic_cast<TH1*>(hobj);
661 if (!h) return 0;
662 if (h->GetDimension() > 2) return 0;
663 h->SetStats(0);
664 if ( TString(h->GetName()).Contains("Events")) { h->SetStats(1); }
665 TString t(label);
666 if (label) t += "_";
667 t += h->GetName();
668 h->SetTitle(t.Data());
669 TCanvas* c = new TCanvas(t.Data(),t.Data());
670 if (h->GetDimension() >= 1) {
671 TString xlabel(h->GetXaxis()->GetTitle());
672 if (xlabel.Contains("Pt")) { c->SetLogx(); h->GetXaxis()->SetRangeUser(0.1 , 100.); }
673 }
674 if (h->GetDimension() == 1) {
675 if (xlabel.Contains("p_{T}")) { c->SetLogx(); c->SetLogy(); h->GetXaxis()->SetRangeUser(0.1 , 100.); }
676 }
677 if (h->GetDimension() == 2) {
678 TString ylabel(h->GetYaxis()->GetTitle());
679 if (ylabel.Contains("Pt")) { c->SetLogy(); h->GetYaxis()->SetRangeUser(0.1 , 100.); }
680 if (ylabel.Contains("p_{T}")) { c->SetLogy(); h->GetYaxis()->SetRangeUser(0.1 , 100.); }
681 h->Draw("COLZ");
682 }
683 if (h->GetDimension() == 1) {
684 h->Draw();
685 }
686 return c;
687
688}
689
690Int_t CheckLoadLibrary(const char* library)
691{
692 // checks if a library is already loaded, if not loads the library
693
9a3036c4 694 if (strlen(gSystem->GetLibraries(library, "", kFALSE)) > 0)
62e3b4b6 695 return 1;
696
697 return gSystem->Load(library);
698}