]>
Commit | Line | Data |
---|---|---|
1 | Double_t GenerateCorrMatr(const char* mcfile, const char* mctask, const char* idstring ,const char* outfile, const char* gifdir = 0) | |
2 | //Double_t GenerateCorrMatr() | |
3 | { | |
4 | ||
5 | //tmp setting | |
6 | // const char* mcfile = "/lustre/alice/train/V005.PbPb_MC/2010-11-29_0108.4251/mergedPeriods/MC/PbPb/LHC10h8/mknichel_dNdPtPbPb_gt_v0_c0_syst4.root"; | |
7 | // const char* mctask = "mknichel_dNdPtPbPb_gt_v0_c0_syst4"; | |
8 | // const char* idstring = "gt_v0_c0_syst4"; | |
9 | // const char* outfile = "/u/mknichel/alice/dNdPt/2010-11-29_0318/corrMatr_gt_v0_c0_syst4.root"; | |
10 | // const char* gifdir = "/u/mknichel/alice/dNdPt/2010-11-29"; | |
11 | ||
12 | // settings vor zVertex cut (event and track level) | |
13 | Double_t zVert = 10.0; | |
14 | ||
15 | // setting on eta cut (track level) | |
16 | Double_t eta = 0.8; | |
17 | ||
18 | // strangeness scaling factor (for secondaries from strange decays) | |
19 | Double_t sscale = 1.4; | |
20 | ||
21 | //load required libraries | |
22 | //load required libraries | |
23 | gSystem->AddIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/ -I$ALICE_ROOT/include -I$ALICE_ROOT/STEER -I$ALICE_ROOT/ANALYSIS -I$ALICE_ROOT/PWG0 -I$ALICE_ROOT/PWGPP -I$ALICE_ROOT/PWG2 -I$ALICE_ROOT/PWG3 -I$ALICE_ROOT/PWG3/vertexingHF -I$ALICE_ROOT/PWG4 -I$ALICE_ROOT/CORRFW -I$ALICE_ROOT/TPC -I$ALICE_ROOT/TRD -I$ALICE_ROOT/PWG3/muon -I$ALICE_ROOT/JETAN -I$ALICE_ROOT/ANALYSIS/Tender"); | |
24 | ||
25 | gSystem->Load("libCore"); | |
26 | gSystem->Load("libPhysics"); | |
27 | gSystem->Load("libMinuit"); | |
28 | gSystem->Load("libGui"); | |
29 | gSystem->Load("libXMLParser"); | |
30 | ||
31 | gSystem->Load("libGeom"); | |
32 | gSystem->Load("libVMC"); | |
33 | ||
34 | gSystem->Load("libNet"); | |
35 | ||
36 | ||
37 | gSystem->Load("libSTEERBase"); | |
38 | gSystem->Load("libESD"); | |
39 | gSystem->Load("libCDB"); | |
40 | gSystem->Load("libRAWDatabase"); | |
41 | gSystem->Load("libRAWDatarec"); | |
42 | gSystem->Load("libANALYSIS"); | |
43 | ||
44 | ||
45 | ||
46 | gSystem->Load("libANALYSIS.so"); | |
47 | gSystem->Load("libANALYSISalice.so"); | |
48 | gSystem->Load("libTENDER.so"); | |
49 | gSystem->Load("libCORRFW.so"); | |
50 | gSystem->Load("libPWG0base.so"); | |
51 | gSystem->Load("libPWG0dep"); | |
52 | gSystem->Load("libPWG0selectors.so"); | |
53 | ||
54 | ||
55 | ||
56 | ||
57 | ||
58 | // make plots nicer | |
59 | gROOT->SetStyle("Plain"); | |
60 | gStyle->SetPalette(1); | |
61 | ||
62 | // array for all correction matrices | |
63 | TObjArray* CorrMatr = new TObjArray(); | |
64 | ||
65 | // array for all control histograms | |
66 | TObjArray* ContHist = new TObjArray(); | |
67 | ||
68 | ||
69 | // load mc information | |
70 | TFile* fmc = TFile::Open(mcfile,"READ"); | |
71 | if (!fmc) return -1; | |
72 | TList* lmc = dynamic_cast<TList*>(fmc->Get(mctask)); | |
73 | if (!lmc) return -1; | |
74 | AlidNdPtAnalysisPbPb *obj = dynamic_cast<AlidNdPtAnalysisPbPb*>(lmc->FindObject("dNdPtAnalysisPbPb")); | |
75 | if (!obj) return -1; | |
76 | ||
77 | //Event statistics | |
78 | THnSparse *fRecEventMatrix = obj->GetRecEventMatrix(); //all reconstructed events | |
79 | TH2D* h2RecEventAll = fRecEventMatrix->Projection(0,1)->Clone("h2RecEventAll"); | |
80 | fRecEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
81 | TH2D* h2RecEvent = fRecEventMatrix->Projection(0,1)->Clone("h2RecEvent"); | |
82 | Double_t MCReconstructedEvents = h2RecEvent->Integral(); | |
83 | Double_t MCReconstructedEventsAll = h2RecEventAll->Integral(); | |
84 | ContHist->Add(h2RecEvent); | |
85 | ContHist->Add(h2RecEventAll); | |
86 | ||
87 | THnSparse *fTriggerEventMatrix = obj->GetTriggerEventMatrix(); //all triggered events | |
88 | TH2D* h2TriggerEventAll = fTriggerEventMatrix->Projection(0,1)->Clone("h2TriggerEventAll"); | |
89 | fTriggerEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
90 | TH2D* h2TriggerEvent = fTriggerEventMatrix->Projection(0,1)->Clone("h2TriggerEvent"); | |
91 | Double_t MCTriggeredEvents = h2TriggerEvent->Integral(); | |
92 | Double_t MCTriggeredEventsAll = h2TriggerEventAll->Integral(); | |
93 | ContHist->Add(h2TriggerEvent); | |
94 | ContHist->Add(h2TriggerEventAll); | |
95 | ||
96 | THnSparse *fGenEventMatrix = obj->GetGenEventMatrix(); //all generated events | |
97 | TH2D* h2GenEventAll = fGenEventMatrix->Projection(0,1)->Clone("h2GenEventAll"); | |
98 | fGenEventMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
99 | TH2D* h2GenEvent = fGenEventMatrix->Projection(0,1)->Clone("h2GenEvent"); | |
100 | Double_t MCGeneratedEvents = h2GenEvent->Integral(); | |
101 | Double_t MCGeneratedEventsAll = h2GenEventAll->Integral(); | |
102 | ContHist->Add(h2RecEvent); | |
103 | ContHist->Add(h2RecEventAll); | |
104 | ||
105 | printf("=== generated MC events for correction matrices %lf ===\n",MCGeneratedEvents); | |
106 | printf("=== triggered MC events for correction matrices %lf ===\n",MCTriggeredEvents); | |
107 | printf("=== recontructed MC events for correction matrices %lf ===\n",MCReconstructedEvents); | |
108 | printf("\n"); | |
109 | printf("=== cut on the zVertex +- %lf ===\n",zVert); | |
110 | printf("=== generated MC events (before zVertex cut) %lf ===\n",MCGeneratedEventsAll); | |
111 | printf("=== triggered MC events (before zVertex cut) %lf ===\n",MCTriggeredEventsAll); | |
112 | printf("=== recontructed MC events (before zVertex cut) %lf ===\n",MCReconstructedEventsAll); | |
113 | ||
114 | TH1D* h1MCGeneratedEvents = new TH1D("h1MCGeneratedEvents","h1MCGeneratedEvents",1,0,1); | |
115 | TH1D* h1MCTriggeredEvents = new TH1D("h1MCTriggeredEvents","h1MCTriggeredEvents",1,0,1); | |
116 | TH1D* h1MCReconstructedEvents = new TH1D("h1MCReconstructedEvents","h1MCReconstructedEvents",1,0,1); | |
117 | TH1D* h1MCGeneratedEventsAll = new TH1D("h1MCGeneratedEventsAll","h1MCGeneratedEventsAll",1,0,1); | |
118 | TH1D* h1MCTriggeredEventsAll = new TH1D("h1MCTriggeredEventsAll","h1MCTriggeredEventsAll",1,0,1); | |
119 | TH1D* h1MCReconstructedEventsAll = new TH1D("h1MCReconstructedEventsAll","h1MCReconstructedEventsAll",1,0,1); | |
120 | ||
121 | h1MCGeneratedEvents->Fill(0,MCGeneratedEvents); | |
122 | h1MCGeneratedEvents->SetEntries(MCGeneratedEvents); | |
123 | h1MCTriggeredEvents->Fill(0,MCTriggeredEvents); | |
124 | h1MCTriggeredEvents->SetEntries(MCTriggeredEvents); | |
125 | h1MCReconstructedEvents->Fill(0,MCReconstructedEvents); | |
126 | h1MCReconstructedEvents->SetEntries(MCReconstructedEvents); | |
127 | h1MCGeneratedEventsAll->Fill(0,MCGeneratedEventsAll); | |
128 | h1MCGeneratedEventsAll->SetEntries(MCGeneratedEventsAll); | |
129 | h1MCTriggeredEventsAll->Fill(0,MCTriggeredEventsAll); | |
130 | h1MCTriggeredEventsAll->SetEntries(MCTriggeredEventsAll); | |
131 | h1MCReconstructedEventsAll->Fill(0,MCReconstructedEventsAll); | |
132 | h1MCReconstructedEventsAll->SetEntries(MCReconstructedEventsAll); | |
133 | ||
134 | ContHist->Add(h1MCGeneratedEvents); | |
135 | ContHist->Add(h1MCTriggeredEvents); | |
136 | ContHist->Add(h1MCReconstructedEvents); | |
137 | ContHist->Add(h1MCGeneratedEventsAll); | |
138 | ContHist->Add(h1MCTriggeredEventsAll); | |
139 | ContHist->Add(h1MCReconstructedEventsAll); | |
140 | ||
141 | // efficienfy and correction matrices for tigger and vertex efficiency | |
142 | TH2D* h2EventTriggerEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEventAll,h2GenEventAll,"h2EventTriggerEffAll"); | |
143 | TH2D* h2EventTriggerCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2GenEventAll,h2TriggerEventAll,"h2EventTriggerCorrAll"); | |
144 | TH2D* h2EventTriggerEff = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEvent,h2GenEvent,"h2EventTriggerEff"); | |
145 | TH2D* h2EventTriggerCorr = AlidNdPtHelper::GenerateCorrMatrix(h2GenEvent,h2TriggerEvent,"h2EventTriggerCorr"); | |
146 | ||
147 | TH2D* h2EventRecEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2RecEventAll,h2TriggerEventAll,"h2EventRecEffAll"); | |
148 | TH2D* h2EventRecCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEventAll,h2RecEventAll,"h2EventRecCorrAll"); | |
149 | TH2D* h2EventRecEff = AlidNdPtHelper::GenerateCorrMatrix(h2RecEvent,h2TriggerEvent,"h2EventRecEff"); | |
150 | TH2D* h2EventRecCorr = AlidNdPtHelper::GenerateCorrMatrix(h2TriggerEvent,h2RecEvent,"h2EventRecCorr"); | |
151 | ||
152 | TH2D* h2EventEffAll = AlidNdPtHelper::GenerateCorrMatrix(h2RecEventAll,h2GenEventAll,"h2EventEffAll"); | |
153 | TH2D* h2EventCorrAll = AlidNdPtHelper::GenerateCorrMatrix(h2GenEventAll,h2RecEventAll,"h2EventCorrAll"); | |
154 | TH2D* h2EventEff = AlidNdPtHelper::GenerateCorrMatrix(h2RecEvent,h2GenEvent,"h2EventEff"); | |
155 | TH2D* h2EventCorr = AlidNdPtHelper::GenerateCorrMatrix(h2GenEvent,h2RecEvent,"h2EventCorr"); | |
156 | ||
157 | CorrMatr->Add(h2EventTriggerEffAll); | |
158 | CorrMatr->Add(h2EventTriggerCorrAll); | |
159 | CorrMatr->Add(h2EventTriggerEff); | |
160 | CorrMatr->Add(h2EventTriggerCorr); | |
161 | CorrMatr->Add(h2EventRecEffAll); | |
162 | CorrMatr->Add(h2EventRecCorrAll); | |
163 | CorrMatr->Add(h2EventRecEff); | |
164 | CorrMatr->Add(h2EventRecCorr); | |
165 | CorrMatr->Add(h2EventEffAll); | |
166 | CorrMatr->Add(h2EventCorrAll); | |
167 | CorrMatr->Add(h2EventEff); | |
168 | CorrMatr->Add(h2EventCorr); | |
169 | ||
170 | ||
171 | ||
172 | // all recontructed | |
173 | THnSparse *fRecTrackMatrix = obj->GetRecTrackMatrix(); | |
174 | fRecTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
175 | fRecTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
176 | TH3D* h3RecTrack = fRecTrackMatrix->Projection(0,1,2)->Clone("h3RecTrack"); | |
177 | TH2D* h2RecTrack_zv_pt = fRecTrackMatrix->Projection(0,1)->Clone("h2RecTrack_zv_pt"); | |
178 | TH2D* h2RecTrack_zv_eta = fRecTrackMatrix->Projection(0,2)->Clone("h2RecTrack_zv_eta"); | |
179 | TH2D* h2RecTrack_pt_eta = fRecTrackMatrix->Projection(1,2)->Clone("h2RecTrack_pt_eta"); | |
180 | TH1D* h1RecTrack_zv = fRecTrackMatrix->Projection(0)->Clone("h1RecTrack_zv"); | |
181 | TH1D* h1RecTrack_pt = fRecTrackMatrix->Projection(1)->Clone("h1RecTrack_pt"); | |
182 | TH1D* h1RecTrack_eta = fRecTrackMatrix->Projection(2)->Clone("h1RecTrack_eta"); | |
183 | Double_t MCReconstructedTracks = h3RecTrack->Integral(); | |
184 | ||
185 | ContHist->Add(h3RecTrack); | |
186 | ContHist->Add(h2RecTrack_zv_pt); | |
187 | ContHist->Add(h2RecTrack_zv_eta); | |
188 | ContHist->Add(h2RecTrack_pt_eta); | |
189 | ContHist->Add(h1RecTrack_zv); | |
190 | ContHist->Add(h1RecTrack_pt); | |
191 | ContHist->Add(h1RecTrack_eta); | |
192 | ||
193 | // recontructed primary tracks | |
194 | THnSparse *fRecPrimTrackMatrix = obj->GetRecPrimTrackMatrix(); | |
195 | THnSparse *fRecTrackMatrixScaled = fRecPrimTrackMatrix->Clone("fRecTrackMatrixScaled"); //used later for secondaries scaling | |
196 | fRecPrimTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
197 | fRecPrimTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
198 | TH3D* h3RecPrimTrack = fRecPrimTrackMatrix->Projection(0,1,2)->Clone("h3RecPrimTrack"); | |
199 | TH2D* h2RecPrimTrack_zv_pt = fRecPrimTrackMatrix->Projection(0,1)->Clone("h2RecPrimTrack_zv_pt"); | |
200 | TH2D* h2RecPrimTrack_zv_eta = fRecPrimTrackMatrix->Projection(0,2)->Clone("h2RecPrimTrack_zv_eta"); | |
201 | TH2D* h2RecPrimTrack_pt_eta = fRecPrimTrackMatrix->Projection(1,2)->Clone("h2RecPrimTrack_pt_eta"); | |
202 | TH1D* h1RecPrimTrack_zv = fRecPrimTrackMatrix->Projection(0)->Clone("h1RecPrimTrack_zv"); | |
203 | TH1D* h1RecPrimTrack_pt = fRecPrimTrackMatrix->Projection(1)->Clone("h1RecPrimTrack_pt"); | |
204 | TH1D* h1RecPrimTrack_eta = fRecPrimTrackMatrix->Projection(2)->Clone("h1RecPrimTrack_eta"); | |
205 | Double_t MCReconstructedPrimTracks = h3RecPrimTrack->Integral(); | |
206 | ||
207 | ContHist->Add(h3RecPrimTrack); | |
208 | ContHist->Add(h2RecPrimTrack_zv_pt); | |
209 | ContHist->Add(h2RecPrimTrack_zv_eta); | |
210 | ContHist->Add(h2RecPrimTrack_pt_eta); | |
211 | ContHist->Add(h1RecPrimTrack_zv); | |
212 | ContHist->Add(h1RecPrimTrack_pt); | |
213 | ContHist->Add(h1RecPrimTrack_eta); | |
214 | ||
215 | // recontructed secondary tracks | |
216 | THnSparse *fRecSecTrackMatrix = obj->GetRecSecTrackMatrix(); | |
217 | THnSparse *fRecSecTrackMatrixScaled = fRecSecTrackMatrix->Clone("fRecSecTrackMatrixScaled"); //used later for secondaries scaling | |
218 | fRecSecTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
219 | fRecSecTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
220 | TH3D* h3RecSecTrack = fRecSecTrackMatrix->Projection(0,1,2)->Clone("h3RecSecTrack"); | |
221 | TH2D* h2RecSecTrack_zv_pt = fRecSecTrackMatrix->Projection(0,1)->Clone("h2RecSecTrack_zv_pt"); | |
222 | TH2D* h2RecSecTrack_zv_eta = fRecSecTrackMatrix->Projection(0,2)->Clone("h2RecSecTrack_zv_eta"); | |
223 | TH2D* h2RecSecTrack_pt_eta = fRecSecTrackMatrix->Projection(1,2)->Clone("h2RecSecTrack_pt_eta"); | |
224 | TH1D* h1RecSecTrack_zv = fRecSecTrackMatrix->Projection(0)->Clone("h1RecSecTrack_zv"); | |
225 | TH1D* h1RecSecTrack_pt = fRecSecTrackMatrix->Projection(1)->Clone("h1RecSecTrack_pt"); | |
226 | TH1D* h1RecSecTrack_eta = fRecSecTrackMatrix->Projection(2)->Clone("h1RecSecTrack_eta"); | |
227 | Double_t MCReconstructedSecTracks = h3RecSecTrack->Integral(); | |
228 | ||
229 | ContHist->Add(h3RecSecTrack); | |
230 | ContHist->Add(h2RecSecTrack_zv_pt); | |
231 | ContHist->Add(h2RecSecTrack_zv_eta); | |
232 | ContHist->Add(h2RecSecTrack_pt_eta); | |
233 | ContHist->Add(h1RecSecTrack_zv); | |
234 | ContHist->Add(h1RecSecTrack_pt); | |
235 | ContHist->Add(h1RecSecTrack_eta); | |
236 | ||
237 | // generated primary tracks | |
238 | THnSparse *fGenPrimTrackMatrix = obj->GetGenPrimTrackMatrix(); | |
239 | fGenPrimTrackMatrix->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
240 | fGenPrimTrackMatrix->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
241 | TH3D* h3GenPrimTrack = fGenPrimTrackMatrix->Projection(0,1,2)->Clone("h3GenPrimTrack"); | |
242 | TH2D* h2GenPrimTrack_zv_Pt = fGenPrimTrackMatrix->Projection(0,1)->Clone("h2GenPrimTrack_zv_pt"); | |
243 | TH2D* h2GenPrimTrack_zv_eta = fGenPrimTrackMatrix->Projection(0,2)->Clone("h2GenPrimTrack_zv_eta"); | |
244 | TH2D* h2GenPrimTrack_pt_eta = fGenPrimTrackMatrix->Projection(1,2)->Clone("h2GenPrimTrack_pt_eta"); | |
245 | TH1D* h1GenPrimTrack_zv = fGenPrimTrackMatrix->Projection(0)->Clone("h1GenPrimTrack_zv"); | |
246 | TH1D* h1GenPrimTrack_pt = fGenPrimTrackMatrix->Projection(1)->Clone("h1GenPrimTrack_pt"); | |
247 | TH1D* h1GenPrimTrack_eta = fGenPrimTrackMatrix->Projection(2)->Clone("h1GenPrimTrack_eta"); | |
248 | Double_t MCGeneratedPrimTracks = h3GenPrimTrack->Integral(); | |
249 | ||
250 | ContHist->Add(h3GenPrimTrack); | |
251 | ContHist->Add(h2GenPrimTrack_zv_pt); | |
252 | ContHist->Add(h2GenPrimTrack_zv_eta); | |
253 | ContHist->Add(h2GenPrimTrack_pt_eta); | |
254 | ContHist->Add(h1GenPrimTrack_zv); | |
255 | ContHist->Add(h1GenPrimTrack_pt); | |
256 | ContHist->Add(h1GenPrimTrack_eta); | |
257 | printf("\n"); | |
258 | printf("==============================================================\n"); | |
259 | printf("=== recontructed MC tracks %lf ===\n",MCReconstructedTracks); | |
260 | printf("=== recontructed MC secondary tracks %lf ===\n",MCReconstructedSecTracks); | |
261 | printf("=== recontructed MC primary tracks %lf ===\n",MCReconstructedPrimTracks); | |
262 | printf("=== generated MC primary track %lf ===\n",MCGeneratedPrimTracks); | |
263 | printf("==============================================================\n"); | |
264 | printf("\n"); | |
265 | ||
266 | // THnSparse *fSparseTriggerTrackEvent = obj->GetTriggerTrackEventMatrix();//Tracks from triggered events | |
267 | // THnSparse *fSparseVtxTrackEvent = obj->GetRecTrackEventMatrix();//Tracks from events with rec. vtx | |
268 | // THnSparse *fSparseGenTrackEvent = obj->GetGenTrackEventMatrix();//generated TrackEvent matrix | |
269 | ||
270 | // tracking efficiencies + corrections | |
271 | TH3D* h3TrackEff = AlidNdPtHelper::GenerateCorrMatrix(h3RecPrimTrack,h3GenPrimTrack,"h3TrackEff"); | |
272 | TH3D* h3TrackCorr = AlidNdPtHelper::GenerateCorrMatrix(h3GenPrimTrack,h3RecPrimTrack,"h3TrackCorr"); | |
273 | ||
274 | TH2D* h2TrackEff_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_zv_pt,h2GenPrimTrack_zv_pt,"h2TrackEff_zv_pt"); | |
275 | TH2D* h2TrackCorr_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_zv_pt,h2RecPrimTrack_zv_pt,"h2TrackCorr_zv_pt"); | |
276 | TH2D* h2TrackEff_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_zv_eta,h2GenPrimTrack_zv_eta,"h2TrackEff_zv_eta"); | |
277 | TH2D* h2TrackCorr_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_zv_eta,h2RecPrimTrack_zv_eta,"h2TrackCorr_zv_eta"); | |
278 | TH2D* h2TrackEff_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecPrimTrack_pt_eta,h2GenPrimTrack_pt_eta,"h2TrackEff_pt_eta"); | |
279 | TH2D* h2TrackCorr_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2GenPrimTrack_pt_eta,h2RecPrimTrack_pt_eta,"h2TrackCorr_pt_eta"); | |
280 | ||
281 | TH1D* h1TrackEff_zv = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_zv,h1GenPrimTrack_zv,"h1TrackEff_zv"); | |
282 | TH1D* h1TrackCorr_zv = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_zv,h1RecPrimTrack_zv,"h1TrackCorr_zv"); | |
283 | TH1D* h1TrackEff_pt = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_pt,h1GenPrimTrack_pt,"h1TrackEff_pt"); | |
284 | TH1D* h1TrackCorr_pt = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_pt,h1RecPrimTrack_pt,"h1TrackCorr_pt"); | |
285 | TH1D* h1TrackEff_eta = AlidNdPtHelper::GenerateCorrMatrix(h1RecPrimTrack_eta,h1GenPrimTrack_eta,"h1TrackEff_eta"); | |
286 | TH1D* h1TrackCorr_eta = AlidNdPtHelper::GenerateCorrMatrix(h1GenPrimTrack_eta,h1RecPrimTrack_eta,"h1TrackCorr_eta"); | |
287 | ||
288 | CorrMatr->Add(h3TrackEff); | |
289 | CorrMatr->Add(h3TrackCorr); | |
290 | CorrMatr->Add(h2TrackEff_zv_pt); | |
291 | CorrMatr->Add(h2TrackCorr_zv_pt); | |
292 | CorrMatr->Add(h2TrackEff_zv_eta); | |
293 | CorrMatr->Add(h2TrackCorr_zv_eta); | |
294 | CorrMatr->Add(h2TrackEff_pt_eta); | |
295 | CorrMatr->Add(h2TrackCorr_pt_eta); | |
296 | CorrMatr->Add(h1TrackEff_zv); | |
297 | CorrMatr->Add(h1TrackCorr_zv); | |
298 | CorrMatr->Add(h1TrackEff_pt); | |
299 | CorrMatr->Add(h1TrackCorr_pt); | |
300 | CorrMatr->Add(h1TrackEff_eta); | |
301 | CorrMatr->Add(h1TrackCorr_eta); | |
302 | ||
303 | // scale the secondaries before calculating correction matrices | |
304 | for (Long64_t i = 0; i < fRecSecTrackMatrixScaled->GetNbins(); i++) { | |
305 | Int_t c[3]; | |
306 | Double_t val = fRecSecTrackMatrixScaled->GetBinContent(i,c); | |
307 | Double_t err = fRecSecTrackMatrixScaled->GetBinError(i); | |
308 | Double_t pt = fRecSecTrackMatrixScaled->GetAxis(1)->GetBinCenter(c[1]); | |
309 | Double_t scale = GetStrangenessCorrFactorPbPb(pt, sscale); | |
310 | fRecSecTrackMatrixScaled->SetBinContent(c,val*scale); | |
311 | fRecSecTrackMatrixScaled->SetBinError(c,err*scale); | |
312 | } | |
313 | ||
314 | // for correct determination of secondaries contamination, also the total total tracks have to be scaled | |
315 | // this is done by taking primaries and adding the scaled secondaries | |
316 | fRecTrackMatrixScaled->Add(fRecSecTrackMatrixScaled); | |
317 | ||
318 | fRecSecTrackMatrixScaled->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
319 | fRecSecTrackMatrixScaled->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
320 | ||
321 | TH3D* h3RecSecTrackScaled = fRecSecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecSecTrackScaled"); | |
322 | TH2D* h2RecSecTrackScaled_zv_pt = fRecSecTrackMatrixScaled->Projection(0,1)->Clone("h2RecSecTrackScaled_zv_pt"); | |
323 | TH2D* h2RecSecTrackScaled_zv_eta = fRecSecTrackMatrixScaled->Projection(0,2)->Clone("h2RecSecTrackScaled_zv_eta"); | |
324 | TH2D* h2RecSecTrackScaled_pt_eta = fRecSecTrackMatrixScaled->Projection(1,2)->Clone("h2RecSecTrackScaled_pt_eta"); | |
325 | TH1D* h1RecSecTrackScaled_zv = fRecSecTrackMatrixScaled->Projection(0)->Clone("h1RecSecTrackScaled_zv"); | |
326 | TH1D* h1RecSecTrackScaled_pt = fRecSecTrackMatrixScaled->Projection(1)->Clone("h1RecSecTrackScaled_pt"); | |
327 | TH1D* h1RecSecTrackScaled_eta = fRecSecTrackMatrixScaled->Projection(2)->Clone("h1RecSecTrackScaled_eta"); | |
328 | ||
329 | ContHist->Add(h3RecSecTrackScaled); | |
330 | ContHist->Add(h2RecSecTrackScaled_zv_pt); | |
331 | ContHist->Add(h2RecSecTrackScaled_zv_eta); | |
332 | ContHist->Add(h2RecSecTrackScaled_pt_eta); | |
333 | ContHist->Add(h1RecSecTrackScaled_zv); | |
334 | ContHist->Add(h1RecSecTrackScaled_pt); | |
335 | ContHist->Add(h1RecSecTrackScaled_eta); | |
336 | ||
337 | fRecTrackMatrixScaled->GetAxis(0)->SetRangeUser(-zVert, zVert-0.01);//zVer | |
338 | fRecTrackMatrixScaled->GetAxis(2)->SetRangeUser(-eta, eta-0.01);//eta | |
339 | ||
340 | TH3D* h3RecTrackScaled = fRecTrackMatrixScaled->Projection(0,1,2)->Clone("h3RecTrackScaled"); | |
341 | TH2D* h2RecTrackScaled_zv_pt = fRecTrackMatrixScaled->Projection(0,1)->Clone("h2RecTrackScaled_zv_pt"); | |
342 | TH2D* h2RecTrackScaled_zv_eta = fRecTrackMatrixScaled->Projection(0,2)->Clone("h2RecTrackScaled_zv_eta"); | |
343 | TH2D* h2RecTrackScaled_pt_eta = fRecTrackMatrixScaled->Projection(1,2)->Clone("h2RecTrackScaled_pt_eta"); | |
344 | TH1D* h1RecTrackScaled_zv = fRecTrackMatrixScaled->Projection(0)->Clone("h1RecTrackScaled_zv"); | |
345 | TH1D* h1RecTrackScaled_pt = fRecTrackMatrixScaled->Projection(1)->Clone("h1RecTrackScaled_pt"); | |
346 | TH1D* h1RecTrackScaled_eta = fRecTrackMatrixScaled->Projection(2)->Clone("h1RecTrackScaled_eta"); | |
347 | ||
348 | ContHist->Add(h3RecTrackScaled); | |
349 | ContHist->Add(h2RecTrackScaled_zv_pt); | |
350 | ContHist->Add(h2RecTrackScaled_zv_eta); | |
351 | ContHist->Add(h2RecTrackScaled_pt_eta); | |
352 | ContHist->Add(h1RecTrackScaled_zv); | |
353 | ContHist->Add(h1RecTrackScaled_pt); | |
354 | ContHist->Add(h1RecTrackScaled_eta); | |
355 | ||
356 | // create histograms for secondaries contamination and correction | |
357 | ||
358 | TH3D* h3SecCont = AlidNdPtHelper::GenerateCorrMatrix(h3RecSecTrackScaled,h3RecTrackScaled,"h3SecCont"); | |
359 | TH3D* h3SecCorr = AlidNdPtHelper::GenerateContCorrMatrix(h3RecSecTrackScaled,h3RecTrackScaled,"h3SecCorr"); | |
360 | TH2D* h2SecCont_zv_pt = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_zv_pt,h2RecTrackScaled_zv_pt,"h2SecCont_zv_pt"); | |
361 | TH2D* h2SecCorr_zv_pt = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_zv_pt,h2RecTrackScaled_zv_pt,"h2SecCorr_zv_pt"); | |
362 | TH2D* h2SecCont_zv_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_zv_eta,h2RecTrackScaled_zv_eta,"h2SecCont_zv_eta"); | |
363 | TH2D* h2SecCorr_zv_eta = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_zv_eta,h2RecTrackScaled_zv_eta,"h2SecCorr_zv_eta"); | |
364 | TH2D* h2SecCont_pt_eta = AlidNdPtHelper::GenerateCorrMatrix(h2RecSecTrackScaled_pt_eta,h2RecTrackScaled_pt_eta,"h2SecCont_pt_eta"); | |
365 | TH2D* h2SecCorr_pt_eta = AlidNdPtHelper::GenerateContCorrMatrix(h2RecSecTrackScaled_pt_eta,h2RecTrackScaled_pt_eta,"h2SecCorr_pt_eta"); | |
366 | TH1D* h1SecCont_zv = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_zv,h1RecTrackScaled_zv,"h1SecCont_zv"); | |
367 | TH1D* h1SecCorr_zv = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_zv,h1RecTrackScaled_zv,"h1SecCorr_zv"); | |
368 | TH1D* h1SecCont_pt = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_pt,h1RecTrackScaled_pt,"h1SecCont_pt"); | |
369 | TH1D* h1SecCorr_pt = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_pt,h1RecTrackScaled_pt,"h1SecCorr_pt"); | |
370 | TH1D* h1SecCont_eta = AlidNdPtHelper::GenerateCorrMatrix(h1RecSecTrackScaled_eta,h1RecTrackScaled_eta,"h1SecCont_eta"); | |
371 | TH1D* h1SecCorr_eta = AlidNdPtHelper::GenerateContCorrMatrix(h1RecSecTrackScaled_eta,h1RecTrackScaled_eta,"h1SecCorr_eta"); | |
372 | ||
373 | CorrMatr->Add(h3SecCont); | |
374 | CorrMatr->Add(h3SecCorr); | |
375 | CorrMatr->Add(h2SecCont_zv_pt); | |
376 | CorrMatr->Add(h2SecCorr_zv_pt); | |
377 | CorrMatr->Add(h2SecCont_zv_eta); | |
378 | CorrMatr->Add(h2SecCorr_zv_eta); | |
379 | CorrMatr->Add(h2SecCont_pt_eta); | |
380 | CorrMatr->Add(h2SecCorr_pt_eta); | |
381 | CorrMatr->Add(h1SecCont_zv); | |
382 | CorrMatr->Add(h1SecCorr_zv); | |
383 | CorrMatr->Add(h1SecCont_pt); | |
384 | CorrMatr->Add(h1SecCorr_pt); | |
385 | CorrMatr->Add(h1SecCont_eta); | |
386 | CorrMatr->Add(h1SecCorr_eta); | |
387 | ||
388 | // plot pictures and save to gifdir | |
389 | for (i=0; i < CorrMatr->LastIndex(); i++) { | |
390 | TCanvas* ctmp = PlotHist(CorrMatr->At(i),idstring); | |
391 | if (gifdir && ctmp) { | |
392 | TString gif(gifdir); | |
393 | gif += '/'; | |
394 | gif += ctmp->GetName(); | |
395 | gif += ".gif"; | |
396 | ctmp->SaveAs(gif.Data(),"gif"); | |
397 | delete ctmp; | |
398 | } | |
399 | } | |
400 | for (i=0; i < ContHist->LastIndex(); i++) { | |
401 | TCanvas* ctmp = PlotHist(ContHist->At(i),idstring); | |
402 | if (gifdir && ctmp) { | |
403 | TString gif(gifdir); | |
404 | gif += '/'; | |
405 | gif += ctmp->GetName(); | |
406 | gif += ".gif"; | |
407 | ctmp->SaveAs(gif.Data(),"gif"); | |
408 | delete ctmp; | |
409 | } | |
410 | } | |
411 | ||
412 | // save all correction matrices and control histograms to file | |
413 | if (!outfile) { return; } | |
414 | TFile *out = TFile::Open(outfile,"RECREATE"); | |
415 | CorrMatr->Write(); | |
416 | ContHist->Write(); | |
417 | out->Close(); | |
418 | ||
419 | return MCReconstructedEvents; | |
420 | ||
421 | } | |
422 | ||
423 | ||
424 | //_____________________________________________________________________________ | |
425 | Double_t GetStrangenessCorrFactorPbPb(Double_t pt, Double_t s) | |
426 | { | |
427 | // data driven correction factor for secondaries (PbPb) | |
428 | ||
429 | if (pt <= 0.25) return 1.0; | |
430 | if (pt <= 0.5) return GetLinearInterpolationValue(0.25,1.0,0.5,0.60+0.40*s, pt); | |
431 | if (pt <= 1.0) return GetLinearInterpolationValue(0.5,0.60+0.40*s,1.0,0.53+0.47*s, pt); | |
432 | if (pt <= 2.0) return GetLinearInterpolationValue(1.0,0.53+0.47*s,2.0,0.44+0.56*s, pt); | |
433 | if (pt <= 5.0) return GetLinearInterpolationValue(2.0,0.44+0.56*s,5.0,0.33+0.67*s, pt); | |
434 | return 0.33+0.67*s; | |
435 | } | |
436 | ||
437 | ||
438 | //___________________________________________________________________________ | |
439 | Double_t GetLinearInterpolationValue(const Double_t x1,const Double_t y1,const Double_t x2,const Double_t y2, const Double_t pt) | |
440 | { | |
441 | // | |
442 | // linear interpolation | |
443 | // | |
444 | return ((y2-y1)/(x2-x1))*pt+(y2-(((y2-y1)/(x2-x1))*x2)); | |
445 | } | |
446 | ||
447 | //___________________________________________________________________________ | |
448 | TCanvas* PlotHist(TObject* hobj, const char* label=0) | |
449 | { | |
450 | TH1* h = dynamic_cast<TH1*>(hobj); | |
451 | if (!h) return 0; | |
452 | if (h->GetDimension() > 2) return 0; | |
453 | h->SetStats(0); | |
454 | if ( TString(h->GetName()).Contains("Events")) { h->SetStats(1); } | |
455 | TString t(label); | |
456 | if (label) t += "_"; | |
457 | t += h->GetName(); | |
458 | h->SetTitle(t.Data()); | |
459 | TCanvas* c = new TCanvas(t.Data(),t.Data()); | |
460 | if (h->GetDimension() >= 1) { | |
461 | TString xlabel(h->GetXaxis()->GetTitle()); | |
462 | if (xlabel.Contains("Pt")) { c->SetLogx(); h->GetXaxis()->SetRangeUser(0.1 , 50.); } | |
463 | } | |
464 | if (h->GetDimension() == 1) { | |
465 | if (xlabel.Contains("p_{T}")) { c->SetLogx(); c->SetLogy(); h->GetXaxis()->SetRangeUser(0.1 , 50.); } | |
466 | } | |
467 | if (h->GetDimension() == 2) { | |
468 | TString ylabel(h->GetYaxis()->GetTitle()); | |
469 | if (ylabel.Contains("Pt")) { c->SetLogy(); h->GetYaxis()->SetRangeUser(0.1 , 50.); } | |
470 | if (ylabel.Contains("p_{T}")) { c->SetLogy(); h->GetYaxis()->SetRangeUser(0.1 , 50.); } | |
471 | h->Draw("COLZ"); | |
472 | } | |
473 | if (h->GetDimension() == 1) { | |
474 | h->Draw(); | |
475 | } | |
476 | return c; | |
477 | ||
478 | } | |
479 | ||
480 | Int_t CheckLoadLibrary(const char* library) | |
481 | { | |
482 | // checks if a library is already loaded, if not loads the library | |
483 | ||
484 | if (strlen(gSystem->GetLibraries(Form("%s.so", library), "", kFALSE)) > 0) | |
485 | return 1; | |
486 | ||
487 | return gSystem->Load(library); | |
488 | } |