]>
Commit | Line | Data |
---|---|---|
3baa4bfd | 1 | //------------------------------------------------------------------------------ |
2 | // Implementation of AliComparisonDEdx class. It keeps information from | |
3 | // comparison of reconstructed and MC particle tracks. In addtion, | |
4 | // it keeps selection cuts used during comparison. The comparison | |
5 | // information is stored in the ROOT histograms. Analysis of these | |
6 | // histograms can be done by using Analyse() class function. The result of | |
7 | // the analysis (histograms/graphs) are stored in the folder which is | |
8 | // a data of AliComparisonDEdx. | |
9 | // | |
10 | // Author: J.Otwinowski 04/02/2008 | |
11 | //------------------------------------------------------------------------------ | |
12 | ||
13 | /* | |
14 | ||
15 | // after running comparison task, read the file, and get component | |
b4126c69 | 16 | gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C"); |
17 | LoadMyLibs(); | |
3baa4bfd | 18 | TFile f("Output.root"); |
19 | AliComparisonDEdx * compObj = (AliComparisonDEdx*)f.Get("AliComparisonDEdx"); | |
20 | ||
b4126c69 | 21 | // Analyse comparison data |
3baa4bfd | 22 | compObj->Analyse(); |
23 | ||
24 | // the output histograms/graphs will be stored in the folder "folderDEdx" | |
25 | compObj->GetAnalysisFolder()->ls("*"); | |
26 | ||
27 | // user can save whole comparison object (or only folder with anlysed histograms) | |
28 | // in the seperate output file (e.g.) | |
29 | TFile fout("Analysed_DEdx.root"."recreate"); | |
30 | compObj->Write(); // compObj->GetAnalysisFolder()->Write(); | |
31 | fout.Close(); | |
32 | ||
33 | */ | |
34 | ||
35 | #include <iostream> | |
36 | ||
37 | #include "TFile.h" | |
38 | #include "TCint.h" | |
39 | #include "TH3F.h" | |
40 | #include "TH2F.h" | |
41 | #include "TF1.h" | |
42 | #include "TProfile.h" | |
43 | #include "TProfile2D.h" | |
44 | #include "TGraph2D.h" | |
45 | #include "TCanvas.h" | |
46 | #include "TGraph.h" | |
47 | // | |
48 | #include "AliESDEvent.h" | |
49 | #include "AliESD.h" | |
50 | #include "AliESDfriend.h" | |
51 | #include "AliESDfriendTrack.h" | |
52 | #include "AliRecInfoCuts.h" | |
53 | #include "AliMCInfoCuts.h" | |
54 | #include "AliLog.h" | |
55 | // | |
56 | #include "AliMathBase.h" | |
57 | #include "AliTreeDraw.h" | |
58 | //#include "TStatToolkit.h" | |
59 | ||
60 | #include "AliMCInfo.h" | |
61 | #include "AliESDRecInfo.h" | |
62 | #include "AliComparisonDEdx.h" | |
63 | ||
64 | using namespace std; | |
65 | ||
66 | ClassImp(AliComparisonDEdx) | |
67 | ||
68 | //_____________________________________________________________________________ | |
69 | AliComparisonDEdx::AliComparisonDEdx(): | |
70 | // TNamed("AliComparisonDEdx","AliComparisonDEdx"), | |
71 | AliComparisonObject("AliComparisonDEdx"), | |
72 | ||
73 | // dEdx | |
74 | fTPCSignalNormTan(0), | |
75 | fTPCSignalNormSPhi(0), | |
76 | fTPCSignalNormTPhi(0), | |
77 | // | |
78 | fTPCSignalNormTanSPhi(0), | |
79 | fTPCSignalNormTanTPhi(0), | |
80 | fTPCSignalNormTanSPt(0), | |
81 | ||
82 | // Cuts | |
83 | fCutsRC(0), | |
84 | fCutsMC(0), | |
85 | fMCPtMin(0), | |
86 | fMCAbsTanThetaMax(0), | |
87 | fMCPdgCode(0), | |
88 | ||
89 | // histogram folder | |
90 | fAnalysisFolder(0) | |
91 | { | |
92 | Init(); | |
93 | } | |
94 | ||
95 | //_____________________________________________________________________________ | |
96 | AliComparisonDEdx::~AliComparisonDEdx(){ | |
97 | ||
98 | if(fTPCSignalNormTan) delete fTPCSignalNormTan; fTPCSignalNormTan=0; | |
99 | if(fTPCSignalNormSPhi) delete fTPCSignalNormSPhi; fTPCSignalNormSPhi=0; | |
100 | if(fTPCSignalNormTPhi) delete fTPCSignalNormTPhi; fTPCSignalNormTPhi=0; | |
101 | // | |
102 | if(fTPCSignalNormTanSPhi) delete fTPCSignalNormTanSPhi; fTPCSignalNormTanSPhi=0; | |
103 | if(fTPCSignalNormTanTPhi) delete fTPCSignalNormTanTPhi; fTPCSignalNormTanTPhi=0; | |
104 | if(fTPCSignalNormTanSPt) delete fTPCSignalNormTanSPt; fTPCSignalNormTanSPt=0; | |
105 | ||
106 | if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0; | |
107 | } | |
108 | ||
109 | //_____________________________________________________________________________ | |
110 | void AliComparisonDEdx::Init() | |
111 | { | |
112 | // Init histograms | |
113 | ||
114 | // TPC dEdx | |
115 | fTPCSignalNormTan = new TH2F("CdEdxTan","CdEdxTan",50, -2,2, 40,30,70); | |
116 | fTPCSignalNormTan->SetXTitle("tan(#theta)"); | |
117 | fTPCSignalNormTan->SetYTitle("rec. dE/dx / calc. dE/dx"); | |
118 | ||
119 | fTPCSignalNormSPhi = new TH2F("CdEdxSPhi","CdEdxSPhi",10,0.0,1,40,30,70); | |
120 | fTPCSignalNormSPhi->SetXTitle("sin(#phi)"); | |
121 | fTPCSignalNormSPhi->SetYTitle("rec. dE/dx / calc. dE/dx"); | |
122 | ||
123 | fTPCSignalNormTPhi = new TH2F("CdEdxTPhi","CdEdxTPhi",10,0.0,2,40,30,70); | |
124 | fTPCSignalNormTPhi->SetXTitle("tan(#phi)"); | |
125 | fTPCSignalNormTPhi->SetYTitle("rec. dE/dx / calc. dE/dx"); | |
126 | ||
127 | fTPCSignalNormTanSPhi= new TH3F("CdEdxTanSPhi","CdEdxTanSPhi",20, -2,2, 10,0.0 ,1, 40,30,70); | |
128 | fTPCSignalNormTanSPhi->SetXTitle("tan(#theta)"); | |
129 | fTPCSignalNormTanSPhi->SetYTitle("sin(#phi)"); | |
130 | fTPCSignalNormTanSPhi->SetZTitle("rec. dE/dx / calc. dE/dx"); | |
131 | ||
132 | fTPCSignalNormTanTPhi= new TH3F("CdEdxTanTPhi","CdEdxTanTPhi",20, -2,2, 10,0.0 ,1, 40,30,70); | |
133 | fTPCSignalNormTanTPhi->SetXTitle("tan(#theta)"); | |
134 | fTPCSignalNormTanTPhi->SetYTitle("tan(#phi)"); | |
135 | fTPCSignalNormTanTPhi->SetZTitle("rec. dE/dx / calc. dE/dx"); | |
136 | ||
137 | fTPCSignalNormTanSPt= new TH3F("CdEdxTanSPt","CdEdxTanSPt",20, -2,2, 10,0.3 ,3, 40,30,70); | |
138 | fTPCSignalNormTanSPt->SetXTitle("tan(#theta)"); | |
139 | fTPCSignalNormTanSPt->SetYTitle("#sqrt{p_{t}}"); | |
140 | fTPCSignalNormTanSPt->SetZTitle("rec. dE/dx / calc. dE/dx"); | |
141 | ||
142 | // Init cuts | |
143 | if(!fCutsMC) | |
144 | AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object"); | |
145 | if(!fCutsRC) | |
146 | AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object"); | |
147 | ||
148 | // init folder | |
149 | fAnalysisFolder = CreateFolder("folderDEdx","Analysis de/dx Folder"); | |
150 | } | |
151 | ||
152 | //_____________________________________________________________________________ | |
153 | void AliComparisonDEdx::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC){ | |
154 | ||
155 | // Fill dE/dx comparison information | |
156 | ||
157 | Float_t mcpt = infoMC->GetParticle().Pt(); | |
158 | Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5); | |
159 | Float_t mprim = infoMC->GetPrim(); | |
160 | ||
161 | // distance to Prim. vertex | |
162 | const Double_t* dv = infoMC->GetVDist(); | |
163 | ||
164 | Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz(); | |
165 | ||
166 | // Check selection cuts | |
167 | if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; | |
168 | if (!isPrim) return; | |
169 | if (infoRC->GetStatus(1)!=3) return; | |
170 | if (!infoRC->GetESDtrack()) return; | |
171 | if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return; | |
172 | if (!infoRC->GetESDtrack()->GetConstrainedParam()) return; | |
173 | //if (mprim>1.4) return; | |
174 | //if (mprim<0.5) return; | |
175 | if (mprim > fCutsMC->GetMaxTPCSignal()) return; | |
176 | if (mprim < fCutsMC->GetMinTPCSignal()) return; | |
177 | if (infoRC->GetESDtrack()->GetTPCsignalN()<fCutsRC->GetMinTPCsignalN()) return; | |
178 | // | |
179 | Float_t ratio = infoRC->GetESDtrack()->GetTPCsignal()/infoMC->GetPrim(); | |
180 | Float_t sphi = infoRC->GetESDtrack()->GetInnerParam()->GetSnp(); | |
181 | Float_t tphi = sphi/TMath::Sqrt(1-sphi*sphi); | |
182 | ||
183 | if (TMath::Abs(infoMC->GetParticle().GetPdgCode()) != GetMCPdgCode()) return; | |
184 | //if (mcpt>0.5) { | |
185 | if (mcpt > GetMCPtMin()) { | |
186 | fTPCSignalNormTan->Fill(tantheta,ratio); // only subset | |
187 | } | |
188 | ||
189 | //if (TMath::Abs(tantheta)<0.5){ | |
190 | if (TMath::Abs(tantheta) < GetMCAbsTanThetaMax()){ | |
191 | fTPCSignalNormSPhi->Fill(sphi,ratio); // only subset | |
192 | fTPCSignalNormTPhi->Fill(tphi,ratio); // only subset | |
193 | } | |
194 | fTPCSignalNormTanSPhi->Fill(tantheta,sphi,ratio); | |
195 | fTPCSignalNormTanTPhi->Fill(tantheta,tphi,ratio); | |
196 | fTPCSignalNormTanSPt->Fill(tantheta,TMath::Sqrt(mcpt),ratio); | |
197 | } | |
198 | ||
199 | //_____________________________________________________________________________ | |
200 | Long64_t AliComparisonDEdx::Merge(TCollection* list) | |
201 | { | |
202 | // Merge list of objects (needed by PROOF) | |
203 | ||
204 | if (!list) | |
205 | return 0; | |
206 | ||
207 | if (list->IsEmpty()) | |
208 | return 1; | |
209 | ||
210 | TIterator* iter = list->MakeIterator(); | |
211 | TObject* obj = 0; | |
212 | ||
213 | // collection of generated histograms | |
214 | Int_t count=0; | |
215 | while((obj = iter->Next()) != 0) | |
216 | { | |
217 | AliComparisonDEdx* entry = dynamic_cast<AliComparisonDEdx*>(obj); | |
218 | if (entry == 0) continue; | |
219 | ||
220 | fTPCSignalNormTan->Add(entry->fTPCSignalNormTan); | |
221 | fTPCSignalNormSPhi->Add(entry->fTPCSignalNormSPhi); | |
222 | fTPCSignalNormTPhi->Add(entry->fTPCSignalNormTPhi); | |
223 | // | |
224 | fTPCSignalNormTanSPhi->Add(entry->fTPCSignalNormTanSPhi); | |
225 | fTPCSignalNormTanTPhi->Add(entry->fTPCSignalNormTanTPhi); | |
226 | fTPCSignalNormTanSPt->Add(entry->fTPCSignalNormTanSPt); | |
227 | ||
228 | count++; | |
229 | } | |
230 | ||
231 | return count; | |
232 | } | |
233 | ||
234 | //_____________________________________________________________________________ | |
235 | void AliComparisonDEdx::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC) | |
236 | { | |
237 | // Process comparison information | |
238 | Process(infoMC,infoRC); | |
239 | } | |
240 | ||
241 | //_____________________________________________________________________________ | |
242 | TH1F* AliComparisonDEdx::MakeResol(TH2F * his, Int_t integ, Bool_t type) | |
243 | { | |
244 | // Make resolution histograms | |
245 | TH1F *hisr, *hism; | |
246 | if (!gPad) new TCanvas; | |
247 | hisr = AliTreeDraw::CreateResHistoI(his,&hism,integ); | |
248 | if (type) return hism; | |
249 | else | |
250 | return hisr; | |
251 | } | |
252 | ||
253 | //_____________________________________________________________________________ | |
254 | void AliComparisonDEdx::Analyse() | |
255 | { | |
b4126c69 | 256 | // Analyze comparison information and store output histograms |
3baa4bfd | 257 | // in the folder "folderDEdx" |
258 | // | |
259 | ||
260 | TH1::AddDirectory(kFALSE); | |
261 | ||
262 | AliComparisonDEdx * comp=this; | |
b4126c69 | 263 | TObjArray *aFolderObj = new TObjArray; |
3baa4bfd | 264 | |
265 | TH1F *hiss=0; | |
266 | TGraph2D * gr=0; | |
267 | ||
3baa4bfd | 268 | // write results in the folder |
269 | TCanvas * c = new TCanvas("can","TPC dedx"); | |
270 | c->cd(); | |
271 | ||
272 | hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,0); | |
273 | hiss->SetXTitle("Tan(#theta)"); | |
274 | hiss->SetYTitle("#sigma_{dEdx}"); | |
275 | hiss->Draw(); | |
276 | hiss->SetName("TPCdEdxResolTan"); | |
277 | ||
b4126c69 | 278 | aFolderObj->Add(hiss); |
3baa4bfd | 279 | // |
280 | hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,1); | |
281 | hiss->SetXTitle("Tan(#theta)"); | |
282 | hiss->SetYTitle("<dEdx>"); | |
283 | hiss->Draw(); | |
284 | hiss->SetName("TPCdEdxMeanTan"); | |
285 | ||
b4126c69 | 286 | aFolderObj->Add(hiss); |
3baa4bfd | 287 | // |
288 | gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,4); | |
289 | gr->GetXaxis()->SetTitle("Tan(#theta)"); | |
290 | gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}"); | |
291 | gr->GetZaxis()->SetTitle("<dEdx>"); | |
292 | gr->SetName("TPCdEdxMeanTanPt_1"); | |
293 | gr->GetHistogram()->Draw("colz"); | |
294 | ||
decf0997 | 295 | gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_1"); |
b4126c69 | 296 | aFolderObj->Add(gr->GetHistogram()); |
3baa4bfd | 297 | // |
298 | gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,5); | |
299 | gr->GetXaxis()->SetTitle("Tan(#theta)"); | |
300 | gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}"); | |
301 | gr->GetZaxis()->SetTitle("#sigma_{dEdx}"); | |
302 | gr->SetName("TPCdEdxMeanTanPt_2"); | |
303 | gr->GetHistogram()->Draw("colz"); | |
304 | ||
decf0997 | 305 | gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_2"); |
b4126c69 | 306 | aFolderObj->Add(gr->GetHistogram()); |
307 | ||
308 | // export objects to analysis folder | |
309 | fAnalysisFolder = ExportToFolder(aFolderObj); | |
310 | ||
311 | // delete only TObjrArray | |
312 | if(aFolderObj) delete aFolderObj; | |
313 | } | |
314 | ||
315 | //_____________________________________________________________________________ | |
316 | TFolder* AliComparisonDEdx::ExportToFolder(TObjArray * array) | |
317 | { | |
318 | // recreate folder avery time and export objects to new one | |
319 | // | |
320 | AliComparisonDEdx * comp=this; | |
321 | TFolder *folder = comp->GetAnalysisFolder(); | |
322 | ||
323 | TString name, title; | |
324 | TFolder *newFolder = 0; | |
325 | Int_t i = 0; | |
326 | Int_t size = array->GetSize(); | |
327 | ||
328 | if(folder) { | |
329 | // get name and title from old folder | |
330 | name = folder->GetName(); | |
331 | title = folder->GetTitle(); | |
332 | ||
333 | // delete old one | |
334 | delete folder; | |
3baa4bfd | 335 | |
b4126c69 | 336 | // create new one |
337 | newFolder = CreateFolder(name.Data(),title.Data()); | |
338 | newFolder->SetOwner(); | |
339 | ||
340 | // add objects to folder | |
341 | while(i < size) { | |
342 | newFolder->Add(array->At(i)); | |
343 | i++; | |
344 | } | |
345 | } | |
346 | ||
347 | return newFolder; | |
3baa4bfd | 348 | } |
349 | ||
b4126c69 | 350 | |
3baa4bfd | 351 | //_____________________________________________________________________________ |
352 | TFolder* AliComparisonDEdx::CreateFolder(TString name,TString title) { | |
353 | // create folder for analysed histograms | |
354 | TFolder *folder = 0; | |
355 | folder = new TFolder(name.Data(),title.Data()); | |
356 | ||
357 | return folder; | |
358 | } | |
359 |