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.
10 // Author: J.Otwinowski 04/02/2008
11 //------------------------------------------------------------------------------
15 // after running comparison task, read the file, and get component
16 gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
18 TFile f("Output.root");
19 //AliComparisonDEdx * compObj = (AliComparisonDEdx*)f.Get("AliComparisonDEdx");
20 AliComparisonDEdx * compObj = (AliComparisonDEdx*)cOutput->FindObject("AliComparisonDEdx");
22 // Analyse comparison data
25 // the output histograms/graphs will be stored in the folder "folderDEdx"
26 compObj->GetAnalysisFolder()->ls("*");
28 // user can save whole comparison object (or only folder with anlysed histograms)
29 // in the seperate output file (e.g.)
30 TFile fout("Analysed_DEdx.root"."recreate");
31 compObj->Write(); // compObj->GetAnalysisFolder()->Write();
40 #include "AliComparisonDEdx.h"
41 #include "AliESDEvent.h"
42 #include "AliESDRecInfo.h"
44 #include "AliMCInfo.h"
45 #include "AliMCInfoCuts.h"
46 #include "AliMathBase.h"
47 #include "AliRecInfoCuts.h"
48 #include "AliTreeDraw.h"
52 ClassImp(AliComparisonDEdx)
54 //_____________________________________________________________________________
55 AliComparisonDEdx::AliComparisonDEdx():
56 // TNamed("AliComparisonDEdx","AliComparisonDEdx"),
57 AliComparisonObject("AliComparisonDEdx"),
69 // default constructor
72 //_____________________________________________________________________________
73 AliComparisonDEdx::AliComparisonDEdx(Char_t* name="AliComparisonDEdx", Char_t* title="AliComparisonDEdx",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
74 AliComparisonObject(name,title),
88 SetAnalysisMode(analysisMode);
89 SetHptGenerator(hptGenerator);
94 //_____________________________________________________________________________
95 AliComparisonDEdx::~AliComparisonDEdx()
98 if(fDeDxHisto) delete fDeDxHisto; fDeDxHisto=0;
100 if(fCutsRC) delete fCutsRC; fCutsRC=0;
101 if(fCutsMC) delete fCutsMC; fCutsMC=0;
103 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
106 //_____________________________________________________________________________
107 void AliComparisonDEdx::Init()
113 Double_t binsP[32] = {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.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
114 Double_t pMin = 0., pMax = 10.;
116 if(IsHptGenerator() == kTRUE) {
118 pMin = 0.; pMax = 100.;
121 //signal:alpha:y:z:snp:tgl:ncls:pid:p
122 Int_t binsQA[9] = {600,50, 50, 50, 50, 50, 80, 5, nPBins};
123 Double_t xminQA[9] = {0, -4,-20,-250, -1, -2, 0, 0., pMin};
124 Double_t xmaxQA[9] = {300, 4, 20, 250, 1, 2, 160, 5., pMax};
126 fDeDxHisto = new THnSparseF("fDeDxHisto","signal:alpha:y:z:snp:tgl:ncls:pid:momentum",9,binsQA,xminQA,xmaxQA);
127 if(!IsHptGenerator()) fDeDxHisto->SetBinEdges(8,binsP);
129 fDeDxHisto->GetAxis(0)->SetTitle("signal");
130 fDeDxHisto->GetAxis(1)->SetTitle("alpha (rad)");
131 fDeDxHisto->GetAxis(2)->SetTitle("y (cm)");
132 fDeDxHisto->GetAxis(3)->SetTitle("z (cm)");
133 fDeDxHisto->GetAxis(4)->SetTitle("snp");
134 fDeDxHisto->GetAxis(5)->SetTitle("tgl");
135 fDeDxHisto->GetAxis(6)->SetTitle("ncls");
136 fDeDxHisto->GetAxis(6)->SetTitle("pid");
137 fDeDxHisto->GetAxis(7)->SetTitle("p (GeV/c)");
142 AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
144 AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
147 fAnalysisFolder = CreateFolder("folderDEdx","Analysis de/dx Folder");
150 //_____________________________________________________________________________
151 void AliComparisonDEdx::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
153 // Fill dE/dx comparison information
155 // Check selection cuts
156 if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
159 if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0;
160 if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1;
161 if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2;
162 if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3;
163 if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4;
165 if (!infoRC->GetESDtrack()) return;
166 if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
168 Float_t dedx = infoRC->GetESDtrack()->GetTPCsignal();
169 Int_t ncls = infoRC->GetESDtrack()->GetTPCNcls();
171 const AliExternalTrackParam *innerParam = 0;
172 if ((innerParam = infoRC->GetESDtrack()->GetInnerParam()) == 0) return;
174 Double_t pt = innerParam->Pt();
175 Double_t lam = TMath::ATan2(innerParam->Pz(),innerParam->Pt());
176 Double_t p = pt/TMath::Cos(lam);
177 Double_t alpha = innerParam->GetAlpha();
178 Double_t y = innerParam->GetY();
179 Double_t z = innerParam->GetZ();
180 Double_t snp = innerParam->GetSnp();
181 Double_t tgl = innerParam->GetTgl();
183 Double_t vDeDxHisto[9] = {dedx,alpha,y,z,snp,tgl,ncls,pid,p};
184 fDeDxHisto->Fill(vDeDxHisto);
187 //_____________________________________________________________________________
188 void AliComparisonDEdx::ProcessTPCITS(AliMCInfo* const /*infoMC*/, AliESDRecInfo *const /*infoRC*/)
190 // Fill dE/dx comparison information
192 AliDebug(AliLog::kWarning, "Warning: Not implemented");
195 //_____________________________________________________________________________
196 void AliComparisonDEdx::ProcessConstrained(AliMCInfo* const /*infoMC*/, AliESDRecInfo *const /*infoRC*/)
198 // Fill dE/dx comparison information
200 AliDebug(AliLog::kWarning, "Warning: Not implemented");
203 //_____________________________________________________________________________
204 Long64_t AliComparisonDEdx::Merge(TCollection* const list)
206 // Merge list of objects (needed by PROOF)
214 TIterator* iter = list->MakeIterator();
217 // collection of generated histograms
219 while((obj = iter->Next()) != 0)
221 AliComparisonDEdx* entry = dynamic_cast<AliComparisonDEdx*>(obj);
222 if (entry == 0) continue;
224 fDeDxHisto->Add(entry->fDeDxHisto);
231 //_____________________________________________________________________________
232 void AliComparisonDEdx::Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
234 // Process comparison information
236 if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);
237 else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);
238 else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);
240 printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
245 //_____________________________________________________________________________
246 TH1F* AliComparisonDEdx::MakeResol(TH2F * his, Int_t integ, Bool_t type)
248 // Make resolution histograms
250 if (!gPad) new TCanvas;
251 hisr = AliTreeDraw::CreateResHistoI(his,&hism,integ);
252 if (type) return hism;
257 //_____________________________________________________________________________
258 void AliComparisonDEdx::Analyse()
260 // Analyze comparison information and store output histograms
261 // in the folder "folderDEdx"
263 TH1::AddDirectory(kFALSE);
264 TObjArray *aFolderObj = new TObjArray;
270 // write results in the folder
271 TCanvas * c = new TCanvas("can","TPC dedx");
274 hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,0);
275 hiss->SetXTitle("Tan(#theta)");
276 hiss->SetYTitle("#sigma_{dEdx}");
278 hiss->SetName("TPCdEdxResolTan");
280 aFolderObj->Add(hiss);
282 hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,1);
283 hiss->SetXTitle("Tan(#theta)");
284 hiss->SetYTitle("<dEdx>");
286 hiss->SetName("TPCdEdxMeanTan");
288 aFolderObj->Add(hiss);
290 gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,4);
291 gr->GetXaxis()->SetTitle("Tan(#theta)");
292 gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
293 gr->GetZaxis()->SetTitle("<dEdx>");
294 gr->SetName("TPCdEdxMeanTanPt_1");
295 gr->GetHistogram()->Draw("colz");
297 gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_1");
298 aFolderObj->Add(gr->GetHistogram());
300 gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,5);
301 gr->GetXaxis()->SetTitle("Tan(#theta)");
302 gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
303 gr->GetZaxis()->SetTitle("#sigma_{dEdx}");
304 gr->SetName("TPCdEdxMeanTanPt_2");
305 gr->GetHistogram()->Draw("colz");
307 gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_2");
308 aFolderObj->Add(gr->GetHistogram());
311 // export objects to analysis folder
312 fAnalysisFolder = ExportToFolder(aFolderObj);
314 // delete only TObjrArray
315 if(aFolderObj) delete aFolderObj;
318 //_____________________________________________________________________________
319 TFolder* AliComparisonDEdx::ExportToFolder(TObjArray * array)
321 // recreate folder avery time and export objects to new one
323 AliComparisonDEdx * comp=this;
324 TFolder *folder = comp->GetAnalysisFolder();
327 TFolder *newFolder = 0;
329 Int_t size = array->GetSize();
332 // get name and title from old folder
333 name = folder->GetName();
334 title = folder->GetTitle();
340 newFolder = CreateFolder(name.Data(),title.Data());
341 newFolder->SetOwner();
343 // add objects to folder
345 newFolder->Add(array->At(i));
354 //_____________________________________________________________________________
355 TFolder* AliComparisonDEdx::CreateFolder(TString name,TString title) {
356 // create folder for analysed histograms
358 folder = new TFolder(name.Data(),title.Data());