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