]>
Commit | Line | Data |
---|---|---|
c04c80e6 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
27de2dfb | 15 | |
16 | /* $Id$ */ | |
17 | ||
c04c80e6 | 18 | // |
19 | // Class for spectrum correction | |
20 | // Subtraction of hadronic background, Unfolding of the data and | |
21 | // Renormalization done here | |
22 | // The following containers have to be set: | |
23 | // - Correction framework container for real data | |
24 | // - Correction framework container for MC (Efficiency Map) | |
25 | // - Correction framework container for background coming from data | |
26 | // - Correction framework container for background coming from MC | |
27 | // | |
245877d0 | 28 | // Author: |
29 | // Raphaelle Bailhache <R.Bailhache@gsi.de> | |
30 | // Markus Fasel <M.Fasel@gsi.de> | |
c04c80e6 | 31 | // |
32 | ||
33 | #include <TArrayD.h> | |
34 | #include <TH1.h> | |
35 | #include <TList.h> | |
36 | #include <TObjArray.h> | |
37 | #include <TROOT.h> | |
38 | #include <TCanvas.h> | |
39 | #include <TLegend.h> | |
40 | #include <TStyle.h> | |
41 | #include <TMath.h> | |
42 | #include <TAxis.h> | |
43 | #include <TGraphErrors.h> | |
44 | #include <TFile.h> | |
45 | #include <TPad.h> | |
67fe7bd0 | 46 | #include <TH2D.h> |
c04c80e6 | 47 | |
3a72645a | 48 | #include "AliPID.h" |
c04c80e6 | 49 | #include "AliCFContainer.h" |
50 | #include "AliCFDataGrid.h" | |
51 | #include "AliCFEffGrid.h" | |
52 | #include "AliCFGridSparse.h" | |
53 | #include "AliCFUnfolding.h" | |
54 | #include "AliLog.h" | |
55 | ||
56 | #include "AliHFEspectrum.h" | |
3a72645a | 57 | #include "AliHFEcuts.h" |
58 | #include "AliHFEcontainer.h" | |
c04c80e6 | 59 | |
60 | ClassImp(AliHFEspectrum) | |
61 | ||
62 | //____________________________________________________________ | |
63 | AliHFEspectrum::AliHFEspectrum(const char *name): | |
64 | TNamed(name, ""), | |
65 | fCFContainers(NULL), | |
66 | fTemporaryObjects(NULL), | |
67 | fCorrelation(NULL), | |
68 | fBackground(NULL), | |
3a72645a | 69 | fInclusiveSpectrum(kTRUE), |
c04c80e6 | 70 | fDumpToFile(kFALSE), |
3a72645a | 71 | fNbDimensions(1), |
c04c80e6 | 72 | fNEvents(0), |
73 | fStepMC(-1), | |
74 | fStepTrue(-1), | |
75 | fStepData(-1), | |
3a72645a | 76 | fStepBeforeCutsV0(-1), |
77 | fStepAfterCutsV0(-1), | |
c04c80e6 | 78 | fStepGuessedUnfolding(-1), |
3a72645a | 79 | fNumberOfIterations(5), |
80 | fDebugLevel(0) | |
c04c80e6 | 81 | { |
82 | // | |
83 | // Default constructor | |
84 | // | |
85 | } | |
86 | ||
87 | //____________________________________________________________ | |
88 | AliHFEspectrum::~AliHFEspectrum(){ | |
89 | // | |
90 | // Destructor | |
91 | // | |
92 | if(fCFContainers) delete fCFContainers; | |
93 | if(fTemporaryObjects){ | |
94 | fTemporaryObjects->Clear(); | |
95 | delete fTemporaryObjects; | |
96 | } | |
97 | } | |
3a72645a | 98 | //____________________________________________________________ |
6555e2ad | 99 | Bool_t AliHFEspectrum::Init(const AliHFEcontainer *datahfecontainer, const AliHFEcontainer *mchfecontainer, const AliHFEcontainer *v0hfecontainer){ |
3a72645a | 100 | // |
101 | // Init what we need for the correction: | |
102 | // | |
103 | // Raw spectrum, hadron contamination | |
104 | // MC efficiency maps, correlation matrix | |
105 | // V0 efficiency if wanted | |
106 | // | |
107 | // This for a given dimension. | |
108 | // If no inclusive spectrum, then take only efficiency map for beauty electron | |
109 | // and the appropriate correlation matrix | |
110 | // | |
111 | ||
112 | // Get the requested format | |
113 | Int_t dims[3]; | |
114 | switch(fNbDimensions){ | |
115 | case 1: dims[0] = 0; | |
116 | break; | |
117 | case 2: for(Int_t i = 0; i < 2; i++) dims[i] = i; | |
118 | break; | |
119 | case 3: for(Int_t i = 0; i < 3; i++) dims[i] = i; | |
120 | break; | |
121 | default: | |
122 | AliError("Container with this number of dimensions not foreseen (yet)"); | |
123 | return kFALSE; | |
124 | }; | |
125 | ||
126 | // Data container: raw spectrum + hadron contamination | |
127 | AliCFContainer *datacontainer = datahfecontainer->GetCFContainer("recTrackContReco"); | |
128 | AliCFContainer *contaminationcontainer = datahfecontainer->GetCFContainer("hadronicBackground"); | |
129 | if((!datacontainer) || (!contaminationcontainer)) return kFALSE; | |
130 | ||
131 | AliCFContainer *datacontainerD = GetSlicedContainer(datacontainer, fNbDimensions, dims); | |
132 | AliCFContainer *contaminationcontainerD = GetSlicedContainer(contaminationcontainer, fNbDimensions, dims); | |
133 | if((!datacontainerD) || (!contaminationcontainerD)) return kFALSE; | |
134 | SetContainer(datacontainerD,AliHFEspectrum::kDataContainer); | |
135 | SetContainer(contaminationcontainerD,AliHFEspectrum::kBackgroundData); | |
136 | ||
137 | // MC container: ESD/MC efficiency maps + MC/MC efficiency maps | |
138 | AliCFContainer *mccontaineresd = 0x0; | |
139 | AliCFContainer *mccontainermc = 0x0; | |
140 | if(fInclusiveSpectrum) { | |
141 | mccontaineresd = mchfecontainer->MakeMergedCFContainer("sumesd","sumesd","MCTrackCont:recTrackContReco"); | |
142 | mccontainermc = mchfecontainer->MakeMergedCFContainer("summc","summc","MCTrackCont:recTrackContMC"); | |
143 | } | |
144 | else { | |
145 | mccontaineresd = mchfecontainer->MakeMergedCFContainer("sumesd","sumesd","MCTrackCont:recTrackContReco:recTrackContDEReco"); | |
146 | mccontainermc = mchfecontainer->MakeMergedCFContainer("summc","summc","MCTrackCont:recTrackContMC:recTrackContDEMC"); | |
147 | } | |
148 | if((!mccontaineresd) || (!mccontainermc)) return kFALSE; | |
149 | ||
150 | Int_t source = -1; | |
151 | if(!fInclusiveSpectrum) source = 1; | |
152 | AliCFContainer *mccontaineresdD = GetSlicedContainer(mccontaineresd, fNbDimensions, dims, source); | |
153 | AliCFContainer *mccontainermcD = GetSlicedContainer(mccontainermc, fNbDimensions, dims, source); | |
154 | if((!mccontaineresdD) || (!mccontainermcD)) return kFALSE; | |
155 | SetContainer(mccontainermcD,AliHFEspectrum::kMCContainerMC); | |
156 | SetContainer(mccontaineresdD,AliHFEspectrum::kMCContainerESD); | |
157 | ||
158 | // MC container: correlation matrix | |
159 | THnSparseF *mccorrelation = 0x0; | |
bf892a6a | 160 | if(fInclusiveSpectrum) { |
161 | if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 2))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterPID"); | |
162 | if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 1))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterTOF"); | |
163 | if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepbeforePID"); | |
164 | } | |
3a72645a | 165 | else mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterDE"); |
166 | if(!mccorrelation) return kFALSE; | |
167 | THnSparseF *mccorrelationD = GetSlicedCorrelation(mccorrelation, fNbDimensions, dims); | |
168 | if(!mccorrelationD) { | |
169 | printf("No correlation\n"); | |
170 | return kFALSE; | |
171 | } | |
172 | SetCorrelation(mccorrelationD); | |
173 | ||
174 | // V0 container Electron, pt eta phi | |
175 | if(v0hfecontainer) { | |
176 | AliCFContainer *containerV0 = v0hfecontainer->GetCFContainer("taggedTrackContainerReco"); | |
177 | if(!containerV0) return kFALSE; | |
178 | AliCFContainer *containerV0Electron = GetSlicedContainer(containerV0, fNbDimensions, dims, AliPID::kElectron+1); | |
179 | if(!containerV0Electron) return kFALSE; | |
180 | SetContainer(containerV0Electron,AliHFEspectrum::kDataContainerV0); | |
181 | } | |
182 | ||
183 | if(fDebugLevel>0){ | |
184 | ||
185 | AliCFDataGrid *contaminationspectrum = (AliCFDataGrid *) ((AliCFDataGrid *)GetSpectrum(contaminationcontainer,1))->Clone(); | |
186 | contaminationspectrum->SetName("contaminationspectrum"); | |
187 | TCanvas * ccontaminationspectrum = new TCanvas("contaminationspectrum","contaminationspectrum",1000,700); | |
188 | ccontaminationspectrum->Divide(3,1); | |
189 | ccontaminationspectrum->cd(1); | |
6555e2ad | 190 | TH2D * contaminationspectrum2dpteta = (TH2D *) contaminationspectrum->Project(1,0); |
191 | TH2D * contaminationspectrum2dptphi = (TH2D *) contaminationspectrum->Project(2,0); | |
192 | TH2D * contaminationspectrum2detaphi = (TH2D *) contaminationspectrum->Project(1,2); | |
3a72645a | 193 | contaminationspectrum2dpteta->SetStats(0); |
194 | contaminationspectrum2dpteta->SetTitle(""); | |
195 | contaminationspectrum2dpteta->GetXaxis()->SetTitle("#eta"); | |
196 | contaminationspectrum2dpteta->GetYaxis()->SetTitle("p_{T} [GeV/c]"); | |
197 | contaminationspectrum2dptphi->SetStats(0); | |
198 | contaminationspectrum2dptphi->SetTitle(""); | |
199 | contaminationspectrum2dptphi->GetXaxis()->SetTitle("#phi [rad]"); | |
200 | contaminationspectrum2dptphi->GetYaxis()->SetTitle("p_{T} [GeV/c]"); | |
201 | contaminationspectrum2detaphi->SetStats(0); | |
202 | contaminationspectrum2detaphi->SetTitle(""); | |
203 | contaminationspectrum2detaphi->GetXaxis()->SetTitle("#eta"); | |
204 | contaminationspectrum2detaphi->GetYaxis()->SetTitle("#phi [rad]"); | |
205 | contaminationspectrum2dptphi->Draw("colz"); | |
206 | ccontaminationspectrum->cd(2); | |
207 | contaminationspectrum2dpteta->Draw("colz"); | |
208 | ccontaminationspectrum->cd(3); | |
209 | contaminationspectrum2detaphi->Draw("colz"); | |
210 | ||
211 | } | |
212 | ||
213 | ||
214 | return kTRUE; | |
215 | ||
216 | ||
217 | } | |
c04c80e6 | 218 | |
219 | //____________________________________________________________ | |
3a72645a | 220 | Bool_t AliHFEspectrum::Correct(Bool_t subtractcontamination){ |
c04c80e6 | 221 | // |
222 | // Correct the spectrum for efficiency and unfolding | |
223 | // with both method and compare | |
224 | // | |
225 | ||
226 | gStyle->SetPalette(1); | |
227 | gStyle->SetOptStat(1111); | |
228 | gStyle->SetPadBorderMode(0); | |
229 | gStyle->SetCanvasColor(10); | |
230 | gStyle->SetPadLeftMargin(0.13); | |
231 | gStyle->SetPadRightMargin(0.13); | |
67fe7bd0 | 232 | |
3a72645a | 233 | /////////////////////////// |
234 | // Check initialization | |
235 | /////////////////////////// | |
c04c80e6 | 236 | |
3a72645a | 237 | if((!GetContainer(kDataContainer)) || (!GetContainer(kMCContainerMC)) || (!GetContainer(kMCContainerESD))){ |
238 | AliInfo("You have to init before"); | |
239 | return kFALSE; | |
240 | } | |
241 | ||
242 | if((fStepTrue == 0) && (fStepMC == 0) && (fStepData == 0)) { | |
243 | AliInfo("You have to set the steps before: SetMCTruthStep, SetMCEffStep, SetStepToCorrect"); | |
244 | return kFALSE; | |
245 | } | |
246 | ||
247 | SetNumberOfIteration(50); | |
248 | SetStepGuessedUnfolding(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack); | |
249 | ||
250 | AliCFDataGrid *dataGridAfterFirstSteps = 0x0; | |
251 | ////////////////////////////////// | |
252 | // Subtract hadron background | |
253 | ///////////////////////////////// | |
67fe7bd0 | 254 | AliCFDataGrid *dataspectrumaftersubstraction = 0x0; |
3a72645a | 255 | if(subtractcontamination) { |
256 | dataspectrumaftersubstraction = SubtractBackground(kTRUE); | |
257 | dataGridAfterFirstSteps = dataspectrumaftersubstraction; | |
258 | } | |
259 | ||
260 | //////////////////////////////////////////////// | |
261 | // Correct for TPC efficiency from V0 | |
262 | /////////////////////////////////////////////// | |
263 | AliCFDataGrid *dataspectrumafterV0efficiencycorrection = 0x0; | |
264 | AliCFContainer *dataContainerV0 = GetContainer(kDataContainerV0); | |
265 | if(dataContainerV0){ | |
266 | dataspectrumafterV0efficiencycorrection = CorrectV0Efficiency(dataspectrumaftersubstraction); | |
267 | dataGridAfterFirstSteps = dataspectrumafterV0efficiencycorrection; | |
67fe7bd0 | 268 | } |
c04c80e6 | 269 | |
3a72645a | 270 | /////////////// |
c04c80e6 | 271 | // Unfold |
3a72645a | 272 | ////////////// |
273 | TList *listunfolded = Unfold(dataGridAfterFirstSteps); | |
c04c80e6 | 274 | if(!listunfolded){ |
275 | printf("Unfolded failed\n"); | |
3a72645a | 276 | return kFALSE; |
c04c80e6 | 277 | } |
278 | THnSparse *correctedspectrum = (THnSparse *) listunfolded->At(0); | |
279 | THnSparse *residualspectrum = (THnSparse *) listunfolded->At(1); | |
280 | if(!correctedspectrum){ | |
281 | AliError("No corrected spectrum\n"); | |
3a72645a | 282 | return kFALSE; |
c04c80e6 | 283 | } |
67fe7bd0 | 284 | if(!residualspectrum){ |
c04c80e6 | 285 | AliError("No residul spectrum\n"); |
3a72645a | 286 | return kFALSE; |
c04c80e6 | 287 | } |
67fe7bd0 | 288 | |
3a72645a | 289 | ///////////////////// |
c04c80e6 | 290 | // Simply correct |
3a72645a | 291 | //////////////////// |
292 | AliCFDataGrid *alltogetherCorrection = CorrectForEfficiency(dataGridAfterFirstSteps); | |
67fe7bd0 | 293 | |
3a72645a | 294 | |
67fe7bd0 | 295 | ////////// |
c04c80e6 | 296 | // Plot |
297 | ////////// | |
3a72645a | 298 | if(fDebugLevel > 0.0) { |
299 | ||
300 | TCanvas * ccorrected = new TCanvas("corrected","corrected",1000,700); | |
301 | ccorrected->Divide(2,1); | |
302 | ccorrected->cd(1); | |
303 | gPad->SetLogy(); | |
304 | TGraphErrors* correctedspectrumD = Normalize(correctedspectrum); | |
305 | correctedspectrumD->SetTitle(""); | |
306 | correctedspectrumD->GetYaxis()->SetTitleOffset(1.5); | |
307 | correctedspectrumD->GetYaxis()->SetRangeUser(0.000000001,1.0); | |
308 | correctedspectrumD->SetMarkerStyle(26); | |
309 | correctedspectrumD->SetMarkerColor(kBlue); | |
310 | correctedspectrumD->SetLineColor(kBlue); | |
311 | correctedspectrumD->Draw("AP"); | |
312 | TGraphErrors* alltogetherspectrumD = Normalize(alltogetherCorrection); | |
313 | alltogetherspectrumD->SetTitle(""); | |
314 | alltogetherspectrumD->GetYaxis()->SetTitleOffset(1.5); | |
315 | alltogetherspectrumD->GetYaxis()->SetRangeUser(0.000000001,1.0); | |
316 | alltogetherspectrumD->SetMarkerStyle(25); | |
317 | alltogetherspectrumD->SetMarkerColor(kBlack); | |
318 | alltogetherspectrumD->SetLineColor(kBlack); | |
319 | alltogetherspectrumD->Draw("P"); | |
320 | TLegend *legcorrected = new TLegend(0.4,0.6,0.89,0.89); | |
321 | legcorrected->AddEntry(correctedspectrumD,"Corrected","p"); | |
322 | legcorrected->AddEntry(alltogetherspectrumD,"Alltogether","p"); | |
323 | legcorrected->Draw("same"); | |
324 | ccorrected->cd(2); | |
325 | TH1D *correctedTH1D = correctedspectrum->Projection(0); | |
6555e2ad | 326 | TH1D *alltogetherTH1D = (TH1D *) alltogetherCorrection->Project(0); |
3a72645a | 327 | TH1D* ratiocorrected = (TH1D*)correctedTH1D->Clone(); |
328 | ratiocorrected->SetName("ratiocorrected"); | |
329 | ratiocorrected->SetTitle(""); | |
330 | ratiocorrected->GetYaxis()->SetTitle("Unfolded/DirectCorrected"); | |
331 | ratiocorrected->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
332 | ratiocorrected->Divide(correctedTH1D,alltogetherTH1D,1,1); | |
333 | ratiocorrected->SetStats(0); | |
334 | ratiocorrected->Draw(); | |
335 | ||
336 | ||
337 | // Dump to file if needed | |
338 | ||
339 | if(fDumpToFile) { | |
340 | ||
341 | TFile *out = new TFile("finalSpectrum.root","recreate"); | |
342 | out->cd(); | |
343 | // | |
344 | correctedspectrumD->SetName("UnfoldingCorrectedSpectrum"); | |
345 | correctedspectrumD->Write(); | |
346 | alltogetherspectrumD->SetName("AlltogetherSpectrum"); | |
347 | alltogetherspectrumD->Write(); | |
348 | ratiocorrected->SetName("RatioUnfoldingAlltogetherSpectrum"); | |
349 | ratiocorrected->Write(); | |
350 | // | |
351 | correctedspectrum->SetName("UnfoldingCorrectedNotNormalizedSpectrum"); | |
352 | correctedspectrum->Write(); | |
353 | alltogetherCorrection->SetName("AlltogetherCorrectedNotNormalizedSpectrum"); | |
354 | alltogetherCorrection->Write(); | |
355 | // | |
356 | out->Close(); delete out; | |
357 | } | |
358 | ||
359 | ||
360 | } | |
361 | ||
362 | ||
363 | ||
364 | ||
365 | return kTRUE; | |
366 | } | |
367 | //____________________________________________________________ | |
368 | AliCFDataGrid* AliHFEspectrum::SubtractBackground(Bool_t setBackground){ | |
369 | // | |
370 | // Apply background subtraction | |
371 | // | |
372 | ||
373 | // Raw spectrum | |
374 | AliCFContainer *dataContainer = GetContainer(kDataContainer); | |
375 | if(!dataContainer){ | |
376 | AliError("Data Container not available"); | |
377 | return NULL; | |
378 | } | |
379 | AliCFDataGrid *spectrumSubtracted = new AliCFDataGrid("spectrumSubtracted", "Data Grid for spectrum after Background subtraction", *dataContainer,fStepData); | |
380 | ||
381 | AliCFDataGrid *dataspectrumbeforesubstraction = (AliCFDataGrid *) ((AliCFDataGrid *)GetSpectrum(GetContainer(kDataContainer),fStepData))->Clone(); | |
382 | dataspectrumbeforesubstraction->SetName("dataspectrumbeforesubstraction"); | |
383 | ||
384 | // Background Estimate | |
385 | AliCFContainer *backgroundContainer = GetContainer(kBackgroundData); | |
386 | if(!backgroundContainer){ | |
387 | AliError("MC background container not found"); | |
388 | return NULL; | |
389 | } | |
390 | ||
391 | Int_t stepbackground = 1; | |
392 | if(!fInclusiveSpectrum) stepbackground = 2; | |
393 | AliCFDataGrid *backgroundGrid = new AliCFDataGrid("ContaminationGrid","ContaminationGrid",*backgroundContainer,stepbackground); | |
394 | ||
395 | // Subtract | |
396 | spectrumSubtracted->Add(backgroundGrid,-1.0); | |
397 | if(setBackground){ | |
398 | if(fBackground) delete fBackground; | |
399 | fBackground = backgroundGrid; | |
400 | } else delete backgroundGrid; | |
401 | ||
402 | ||
403 | if(fDebugLevel > 0) { | |
67fe7bd0 | 404 | |
405 | TCanvas * cbackgroundsubtraction = new TCanvas("backgroundsubtraction","backgroundsubtraction",1000,700); | |
3a72645a | 406 | cbackgroundsubtraction->Divide(3,1); |
67fe7bd0 | 407 | cbackgroundsubtraction->cd(1); |
3a72645a | 408 | gPad->SetLogy(); |
6555e2ad | 409 | TH1D *measuredTH1Daftersubstraction = (TH1D *) spectrumSubtracted->Project(0); |
410 | TH1D *measuredTH1Dbeforesubstraction = (TH1D *) dataspectrumbeforesubstraction->Project(0); | |
67fe7bd0 | 411 | CorrectFromTheWidth(measuredTH1Daftersubstraction); |
412 | CorrectFromTheWidth(measuredTH1Dbeforesubstraction); | |
413 | measuredTH1Daftersubstraction->SetStats(0); | |
414 | measuredTH1Daftersubstraction->SetTitle(""); | |
415 | measuredTH1Daftersubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]"); | |
416 | measuredTH1Daftersubstraction->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
417 | measuredTH1Daftersubstraction->SetMarkerStyle(25); | |
418 | measuredTH1Daftersubstraction->SetMarkerColor(kBlack); | |
419 | measuredTH1Daftersubstraction->SetLineColor(kBlack); | |
420 | measuredTH1Dbeforesubstraction->SetStats(0); | |
421 | measuredTH1Dbeforesubstraction->SetTitle(""); | |
422 | measuredTH1Dbeforesubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]"); | |
423 | measuredTH1Dbeforesubstraction->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
424 | measuredTH1Dbeforesubstraction->SetMarkerStyle(24); | |
425 | measuredTH1Dbeforesubstraction->SetMarkerColor(kBlue); | |
426 | measuredTH1Dbeforesubstraction->SetLineColor(kBlue); | |
427 | measuredTH1Daftersubstraction->Draw(); | |
428 | measuredTH1Dbeforesubstraction->Draw("same"); | |
429 | TLegend *legsubstraction = new TLegend(0.4,0.6,0.89,0.89); | |
3a72645a | 430 | legsubstraction->AddEntry(measuredTH1Dbeforesubstraction,"With hadron contamination","p"); |
431 | legsubstraction->AddEntry(measuredTH1Daftersubstraction,"Without hadron contamination ","p"); | |
67fe7bd0 | 432 | legsubstraction->Draw("same"); |
433 | cbackgroundsubtraction->cd(2); | |
3a72645a | 434 | gPad->SetLogy(); |
67fe7bd0 | 435 | TH1D* ratiomeasuredcontamination = (TH1D*)measuredTH1Dbeforesubstraction->Clone(); |
436 | ratiomeasuredcontamination->SetName("ratiomeasuredcontamination"); | |
437 | ratiomeasuredcontamination->SetTitle(""); | |
438 | ratiomeasuredcontamination->GetYaxis()->SetTitle("(with contamination - without contamination) / with contamination"); | |
439 | ratiomeasuredcontamination->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
440 | ratiomeasuredcontamination->Add(measuredTH1Daftersubstraction,-1.0); | |
441 | ratiomeasuredcontamination->Divide(measuredTH1Dbeforesubstraction); | |
442 | ratiomeasuredcontamination->SetStats(0); | |
443 | ratiomeasuredcontamination->SetMarkerStyle(26); | |
444 | ratiomeasuredcontamination->SetMarkerColor(kBlack); | |
445 | ratiomeasuredcontamination->SetLineColor(kBlack); | |
446 | ratiomeasuredcontamination->Draw(); | |
3a72645a | 447 | cbackgroundsubtraction->cd(3); |
6555e2ad | 448 | TH1D *measuredTH1background = (TH1D *) backgroundGrid->Project(0); |
3a72645a | 449 | CorrectFromTheWidth(measuredTH1background); |
450 | measuredTH1background->SetStats(0); | |
451 | measuredTH1background->SetTitle(""); | |
452 | measuredTH1background->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]"); | |
453 | measuredTH1background->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
454 | measuredTH1background->SetMarkerStyle(26); | |
455 | measuredTH1background->SetMarkerColor(kBlack); | |
456 | measuredTH1background->SetLineColor(kBlack); | |
457 | measuredTH1background->Draw(); | |
458 | ||
67fe7bd0 | 459 | } |
460 | ||
c04c80e6 | 461 | |
3a72645a | 462 | return spectrumSubtracted; |
c04c80e6 | 463 | } |
c04c80e6 | 464 | //____________________________________________________________ |
3a72645a | 465 | AliCFDataGrid *AliHFEspectrum::CorrectV0Efficiency(AliCFDataGrid* const bgsubpectrum){ |
466 | ||
c04c80e6 | 467 | // |
3a72645a | 468 | // Apply TPC pid efficiency correction from V0 |
c04c80e6 | 469 | // |
470 | ||
3a72645a | 471 | AliCFContainer *v0Container = GetContainer(kDataContainerV0); |
472 | if(!v0Container){ | |
473 | AliError("V0 Container not available"); | |
c04c80e6 | 474 | return NULL; |
475 | } | |
3a72645a | 476 | |
477 | // Efficiency | |
478 | AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*v0Container); | |
479 | efficiencyD->CalculateEfficiency(fStepAfterCutsV0,fStepBeforeCutsV0); | |
480 | ||
481 | // Data in the right format | |
482 | AliCFDataGrid *dataGrid = 0x0; | |
483 | if(bgsubpectrum) { | |
484 | dataGrid = bgsubpectrum; | |
c04c80e6 | 485 | } |
3a72645a | 486 | else { |
487 | ||
488 | AliCFContainer *dataContainer = GetContainer(kDataContainer); | |
489 | if(!dataContainer){ | |
490 | AliError("Data Container not available"); | |
491 | return NULL; | |
492 | } | |
c04c80e6 | 493 | |
3a72645a | 494 | dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData); |
495 | } | |
c04c80e6 | 496 | |
3a72645a | 497 | // Correct |
498 | AliCFDataGrid *result = (AliCFDataGrid *) dataGrid->Clone(); | |
499 | result->ApplyEffCorrection(*efficiencyD); | |
c04c80e6 | 500 | |
3a72645a | 501 | if(fDebugLevel > 0) { |
502 | ||
503 | TCanvas * cV0Efficiency = new TCanvas("V0Efficiency","V0Efficiency",1000,700); | |
504 | cV0Efficiency->Divide(2,1); | |
505 | cV0Efficiency->cd(1); | |
6555e2ad | 506 | TH1D *afterE = (TH1D *) result->Project(0); |
507 | TH1D *beforeE = (TH1D *) dataGrid->Project(0); | |
3a72645a | 508 | afterE->SetStats(0); |
509 | afterE->SetTitle(""); | |
510 | afterE->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]"); | |
511 | afterE->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
512 | afterE->SetMarkerStyle(25); | |
513 | afterE->SetMarkerColor(kBlack); | |
514 | afterE->SetLineColor(kBlack); | |
515 | beforeE->SetStats(0); | |
516 | beforeE->SetTitle(""); | |
517 | beforeE->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]"); | |
518 | beforeE->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
519 | beforeE->SetMarkerStyle(24); | |
520 | beforeE->SetMarkerColor(kBlue); | |
521 | beforeE->SetLineColor(kBlue); | |
522 | afterE->Draw(); | |
523 | beforeE->Draw("same"); | |
524 | TLegend *legV0efficiency = new TLegend(0.4,0.6,0.89,0.89); | |
525 | legV0efficiency->AddEntry(beforeE,"Before Efficiency correction","p"); | |
526 | legV0efficiency->AddEntry(afterE,"After Efficiency correction","p"); | |
527 | legV0efficiency->Draw("same"); | |
528 | cV0Efficiency->cd(2); | |
6555e2ad | 529 | TH1D* efficiencyDproj = (TH1D *) efficiencyD->Project(0); |
3a72645a | 530 | efficiencyDproj->SetTitle(""); |
531 | efficiencyDproj->SetStats(0); | |
532 | efficiencyDproj->SetMarkerStyle(25); | |
533 | efficiencyDproj->Draw(); | |
c04c80e6 | 534 | |
3a72645a | 535 | |
536 | } | |
537 | ||
538 | ||
539 | return result; | |
540 | ||
541 | } | |
c04c80e6 | 542 | //____________________________________________________________ |
3a72645a | 543 | TList *AliHFEspectrum::Unfold(AliCFDataGrid* const bgsubpectrum){ |
c04c80e6 | 544 | |
545 | // | |
546 | // Unfold and eventually correct for efficiency the bgsubspectrum | |
547 | // | |
548 | ||
3a72645a | 549 | AliCFContainer *mcContainer = GetContainer(kMCContainerMC); |
c04c80e6 | 550 | if(!mcContainer){ |
551 | AliError("MC Container not available"); | |
552 | return NULL; | |
553 | } | |
554 | ||
555 | if(!fCorrelation){ | |
556 | AliError("No Correlation map available"); | |
557 | return NULL; | |
558 | } | |
559 | ||
3a72645a | 560 | // Data |
c04c80e6 | 561 | AliCFDataGrid *dataGrid = 0x0; |
562 | if(bgsubpectrum) { | |
c04c80e6 | 563 | dataGrid = bgsubpectrum; |
c04c80e6 | 564 | } |
565 | else { | |
566 | ||
567 | AliCFContainer *dataContainer = GetContainer(kDataContainer); | |
568 | if(!dataContainer){ | |
569 | AliError("Data Container not available"); | |
570 | return NULL; | |
571 | } | |
572 | ||
3a72645a | 573 | dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData); |
c04c80e6 | 574 | } |
575 | ||
c04c80e6 | 576 | // Guessed |
3a72645a | 577 | AliCFDataGrid* guessedGrid = new AliCFDataGrid("guessed","",*mcContainer, fStepGuessedUnfolding); |
c04c80e6 | 578 | THnSparse* guessedTHnSparse = ((AliCFGridSparse*)guessedGrid->GetData())->GetGrid(); |
579 | ||
580 | // Efficiency | |
3a72645a | 581 | AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*mcContainer); |
c04c80e6 | 582 | efficiencyD->CalculateEfficiency(fStepMC,fStepTrue); |
583 | ||
584 | // Unfold | |
585 | ||
3a72645a | 586 | AliCFUnfolding unfolding("unfolding","",fNbDimensions,fCorrelation,efficiencyD->GetGrid(),dataGrid->GetGrid(),guessedTHnSparse); |
c04c80e6 | 587 | unfolding.SetMaxNumberOfIterations(fNumberOfIterations); |
588 | unfolding.UseSmoothing(); | |
589 | unfolding.Unfold(); | |
590 | ||
591 | // Results | |
592 | THnSparse* result = unfolding.GetUnfolded(); | |
593 | THnSparse* residual = unfolding.GetEstMeasured(); | |
594 | TList *listofresults = new TList; | |
595 | listofresults->SetOwner(); | |
596 | listofresults->AddAt((THnSparse*)result->Clone(),0); | |
597 | listofresults->AddAt((THnSparse*)residual->Clone(),1); | |
598 | ||
3a72645a | 599 | if(fDebugLevel > 0) { |
600 | ||
601 | TCanvas * cresidual = new TCanvas("residual","residual",1000,700); | |
602 | cresidual->Divide(2,1); | |
603 | cresidual->cd(1); | |
604 | gPad->SetLogy(); | |
605 | TGraphErrors* residualspectrumD = Normalize(residual); | |
606 | if(!residualspectrumD) { | |
607 | AliError("Number of Events not set for the normalization"); | |
608 | return kFALSE; | |
609 | } | |
610 | residualspectrumD->SetTitle(""); | |
611 | residualspectrumD->GetYaxis()->SetTitleOffset(1.5); | |
612 | residualspectrumD->GetYaxis()->SetRangeUser(0.000000001,1.0); | |
613 | residualspectrumD->SetMarkerStyle(26); | |
614 | residualspectrumD->SetMarkerColor(kBlue); | |
615 | residualspectrumD->SetLineColor(kBlue); | |
616 | residualspectrumD->Draw("AP"); | |
617 | AliCFDataGrid *dataGridBis = (AliCFDataGrid *) (dataGrid->Clone()); | |
618 | dataGridBis->SetName("dataGridBis"); | |
619 | TGraphErrors* measuredspectrumD = Normalize(dataGridBis); | |
620 | measuredspectrumD->SetTitle(""); | |
621 | measuredspectrumD->GetYaxis()->SetTitleOffset(1.5); | |
622 | measuredspectrumD->GetYaxis()->SetRangeUser(0.000000001,1.0); | |
623 | measuredspectrumD->SetMarkerStyle(25); | |
624 | measuredspectrumD->SetMarkerColor(kBlack); | |
625 | measuredspectrumD->SetLineColor(kBlack); | |
626 | measuredspectrumD->Draw("P"); | |
627 | TLegend *legres = new TLegend(0.4,0.6,0.89,0.89); | |
628 | legres->AddEntry(residualspectrumD,"Residual","p"); | |
629 | legres->AddEntry(measuredspectrumD,"Measured","p"); | |
630 | legres->Draw("same"); | |
631 | cresidual->cd(2); | |
632 | TH1D *residualTH1D = residual->Projection(0); | |
6555e2ad | 633 | TH1D *measuredTH1D = (TH1D *) dataGridBis->Project(0); |
3a72645a | 634 | TH1D* ratioresidual = (TH1D*)residualTH1D->Clone(); |
635 | ratioresidual->SetName("ratioresidual"); | |
636 | ratioresidual->SetTitle(""); | |
637 | ratioresidual->GetYaxis()->SetTitle("Folded/Measured"); | |
638 | ratioresidual->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
639 | ratioresidual->Divide(residualTH1D,measuredTH1D,1,1); | |
640 | ratioresidual->SetStats(0); | |
641 | ratioresidual->Draw(); | |
642 | ||
643 | } | |
644 | ||
c04c80e6 | 645 | return listofresults; |
646 | ||
647 | } | |
648 | ||
649 | //____________________________________________________________ | |
3a72645a | 650 | AliCFDataGrid *AliHFEspectrum::CorrectForEfficiency(AliCFDataGrid* const bgsubpectrum){ |
c04c80e6 | 651 | |
652 | // | |
653 | // Apply unfolding and efficiency correction together to bgsubspectrum | |
654 | // | |
655 | ||
3a72645a | 656 | AliCFContainer *mcContainer = GetContainer(kMCContainerESD); |
c04c80e6 | 657 | if(!mcContainer){ |
658 | AliError("MC Container not available"); | |
659 | return NULL; | |
660 | } | |
661 | ||
c04c80e6 | 662 | // Efficiency |
3a72645a | 663 | AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*mcContainer); |
c04c80e6 | 664 | efficiencyD->CalculateEfficiency(fStepMC,fStepTrue); |
665 | ||
666 | // Data in the right format | |
667 | AliCFDataGrid *dataGrid = 0x0; | |
668 | if(bgsubpectrum) { | |
c04c80e6 | 669 | dataGrid = bgsubpectrum; |
c04c80e6 | 670 | } |
671 | else { | |
3a72645a | 672 | |
c04c80e6 | 673 | AliCFContainer *dataContainer = GetContainer(kDataContainer); |
674 | if(!dataContainer){ | |
675 | AliError("Data Container not available"); | |
676 | return NULL; | |
677 | } | |
678 | ||
3a72645a | 679 | dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData); |
c04c80e6 | 680 | } |
681 | ||
682 | // Correct | |
683 | AliCFDataGrid *result = (AliCFDataGrid *) dataGrid->Clone(); | |
684 | result->ApplyEffCorrection(*efficiencyD); | |
c04c80e6 | 685 | |
3a72645a | 686 | if(fDebugLevel > 0) { |
687 | ||
bf892a6a | 688 | printf("Step MC: %d\n",fStepMC); |
689 | printf("Step tracking: %d\n",AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack); | |
690 | printf("Step MC true: %d\n",fStepTrue); | |
3a72645a | 691 | AliCFEffGrid *efficiencymcPID = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),fStepMC,AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack); |
692 | AliCFEffGrid *efficiencymctrackinggeo = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack,fStepTrue); | |
693 | AliCFEffGrid *efficiencymcall = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),fStepMC,fStepTrue); | |
694 | ||
695 | AliCFEffGrid *efficiencyesdall = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerESD),fStepMC,fStepTrue); | |
696 | ||
697 | TCanvas * cefficiency = new TCanvas("efficiency","efficiency",1000,700); | |
698 | cefficiency->cd(1); | |
6555e2ad | 699 | TH1D* efficiencymcPIDD = (TH1D *) efficiencymcPID->Project(0); |
3a72645a | 700 | efficiencymcPIDD->SetTitle(""); |
701 | efficiencymcPIDD->SetStats(0); | |
702 | efficiencymcPIDD->SetMarkerStyle(25); | |
703 | efficiencymcPIDD->Draw(); | |
6555e2ad | 704 | TH1D* efficiencymctrackinggeoD = (TH1D *) efficiencymctrackinggeo->Project(0); |
3a72645a | 705 | efficiencymctrackinggeoD->SetTitle(""); |
706 | efficiencymctrackinggeoD->SetStats(0); | |
707 | efficiencymctrackinggeoD->SetMarkerStyle(26); | |
708 | efficiencymctrackinggeoD->Draw("same"); | |
6555e2ad | 709 | TH1D* efficiencymcallD = (TH1D *) efficiencymcall->Project(0); |
3a72645a | 710 | efficiencymcallD->SetTitle(""); |
711 | efficiencymcallD->SetStats(0); | |
712 | efficiencymcallD->SetMarkerStyle(27); | |
713 | efficiencymcallD->Draw("same"); | |
6555e2ad | 714 | TH1D* efficiencyesdallD = (TH1D *) efficiencyesdall->Project(0); |
3a72645a | 715 | efficiencyesdallD->SetTitle(""); |
716 | efficiencyesdallD->SetStats(0); | |
717 | efficiencyesdallD->SetMarkerStyle(24); | |
718 | efficiencyesdallD->Draw("same"); | |
719 | TLegend *legeff = new TLegend(0.4,0.6,0.89,0.89); | |
720 | legeff->AddEntry(efficiencymcPIDD,"PID efficiency","p"); | |
721 | legeff->AddEntry(efficiencymctrackinggeoD,"Tracking geometry efficiency","p"); | |
722 | legeff->AddEntry(efficiencymcallD,"Overall efficiency","p"); | |
723 | legeff->AddEntry(efficiencyesdallD,"Overall efficiency ESD","p"); | |
724 | legeff->Draw("same"); | |
725 | } | |
726 | ||
c04c80e6 | 727 | return result; |
728 | ||
729 | } | |
3a72645a | 730 | |
c04c80e6 | 731 | //__________________________________________________________________________________ |
732 | TGraphErrors *AliHFEspectrum::Normalize(THnSparse * const spectrum) const { | |
733 | // | |
734 | // Normalize the spectrum to 1/(2*Pi*p_{T})*dN/dp_{T} (GeV/c)^{-2} | |
735 | // Give the final pt spectrum to be compared | |
736 | // | |
737 | ||
738 | if(fNEvents > 0) { | |
739 | ||
740 | TH1D* projection = spectrum->Projection(0); | |
741 | CorrectFromTheWidth(projection); | |
742 | TGraphErrors *graphError = NormalizeTH1(projection); | |
743 | return graphError; | |
744 | ||
745 | } | |
746 | ||
747 | return 0x0; | |
748 | ||
749 | ||
750 | } | |
751 | //__________________________________________________________________________________ | |
752 | TGraphErrors *AliHFEspectrum::Normalize(AliCFDataGrid * const spectrum) const { | |
753 | // | |
754 | // Normalize the spectrum to 1/(2*Pi*p_{T})*dN/dp_{T} (GeV/c)^{-2} | |
755 | // Give the final pt spectrum to be compared | |
756 | // | |
757 | if(fNEvents > 0) { | |
758 | ||
6555e2ad | 759 | TH1D* projection = (TH1D *) spectrum->Project(0); |
c04c80e6 | 760 | CorrectFromTheWidth(projection); |
761 | TGraphErrors *graphError = NormalizeTH1(projection); | |
762 | return graphError; | |
763 | ||
764 | } | |
765 | ||
766 | return 0x0; | |
767 | ||
768 | ||
769 | } | |
770 | //__________________________________________________________________________________ | |
771 | TGraphErrors *AliHFEspectrum::NormalizeTH1(TH1 *input) const { | |
772 | // | |
773 | // Normalize the spectrum to 1/(2*Pi*p_{T})*dN/dp_{T} (GeV/c)^{-2} | |
774 | // Give the final pt spectrum to be compared | |
775 | // | |
776 | if(fNEvents > 0) { | |
777 | ||
778 | TGraphErrors *spectrumNormalized = new TGraphErrors(input->GetNbinsX()); | |
779 | Double_t p = 0, dp = 0; Int_t point = 1; | |
780 | Double_t n = 0, dN = 0; | |
781 | Double_t nCorr = 0, dNcorr = 0; | |
782 | Double_t errdN = 0, errdp = 0; | |
783 | for(Int_t ibin = input->GetXaxis()->GetFirst(); ibin <= input->GetXaxis()->GetLast(); ibin++){ | |
784 | point = ibin - input->GetXaxis()->GetFirst(); | |
785 | p = input->GetXaxis()->GetBinCenter(ibin); | |
786 | dp = input->GetXaxis()->GetBinWidth(ibin)/2.; | |
787 | n = input->GetBinContent(ibin); | |
788 | dN = input->GetBinError(ibin); | |
789 | ||
790 | // New point | |
791 | nCorr = 0.5 * 1/1.6 * 1./(Double_t)(fNEvents) * 1/(2. * TMath::Pi() * p) * n; | |
792 | errdN = 1./(2. * TMath::Pi() * p); | |
793 | errdp = 1./(2. * TMath::Pi() * p*p) * n; | |
794 | dNcorr = 0.5 * 1/1.6 * 1./(Double_t)(fNEvents) * TMath::Sqrt(errdN * errdN * dN *dN + errdp *errdp * dp *dp); | |
795 | ||
796 | spectrumNormalized->SetPoint(point, p, nCorr); | |
797 | spectrumNormalized->SetPointError(point, dp, dNcorr); | |
798 | } | |
799 | spectrumNormalized->GetXaxis()->SetTitle("p_{T} [GeV/c]"); | |
800 | spectrumNormalized->GetYaxis()->SetTitle("#frac{1}{2 #pi p_{T}} #frac{dN}{dp_{T}} / [GeV/c]^{-2}"); | |
801 | spectrumNormalized->SetMarkerStyle(22); | |
802 | spectrumNormalized->SetMarkerColor(kBlue); | |
803 | spectrumNormalized->SetLineColor(kBlue); | |
804 | ||
805 | return spectrumNormalized; | |
806 | ||
807 | } | |
808 | ||
809 | return 0x0; | |
810 | ||
811 | ||
812 | } | |
813 | //____________________________________________________________ | |
814 | void AliHFEspectrum::SetContainer(AliCFContainer *cont, AliHFEspectrum::CFContainer_t type){ | |
815 | // | |
816 | // Set the container for a given step to the | |
817 | // | |
818 | if(!fCFContainers) fCFContainers = new TList; | |
819 | fCFContainers->AddAt(cont, type); | |
820 | } | |
821 | ||
822 | //____________________________________________________________ | |
823 | AliCFContainer *AliHFEspectrum::GetContainer(AliHFEspectrum::CFContainer_t type){ | |
824 | // | |
825 | // Get Correction Framework Container for given type | |
826 | // | |
827 | if(!fCFContainers) return NULL; | |
828 | return dynamic_cast<AliCFContainer *>(fCFContainers->At(type)); | |
829 | } | |
c04c80e6 | 830 | //____________________________________________________________ |
3a72645a | 831 | AliCFContainer *AliHFEspectrum::GetSlicedContainer(AliCFContainer *container, Int_t nDim, Int_t *dimensions,Int_t source) { |
c04c80e6 | 832 | // |
3a72645a | 833 | // Slice bin for a given source of electron |
c04c80e6 | 834 | // |
835 | ||
836 | Double_t *varMin = new Double_t[container->GetNVar()], | |
837 | *varMax = new Double_t[container->GetNVar()]; | |
838 | ||
839 | Double_t *binLimits; | |
840 | for(Int_t ivar = 0; ivar < container->GetNVar(); ivar++){ | |
841 | ||
842 | binLimits = new Double_t[container->GetNBins(ivar)+1]; | |
843 | container->GetBinLimits(ivar,binLimits); | |
3a72645a | 844 | if((ivar == 4) && ((source>= 0) && (source<container->GetNBins(ivar)))) { |
845 | varMin[ivar] = binLimits[source]; | |
846 | varMax[ivar] = binLimits[source]; | |
847 | } | |
848 | else { | |
849 | varMin[ivar] = binLimits[0]; | |
850 | varMax[ivar] = binLimits[container->GetNBins(ivar)]; | |
851 | } | |
852 | ||
c04c80e6 | 853 | delete[] binLimits; |
854 | } | |
855 | ||
856 | AliCFContainer *k = container->MakeSlice(nDim, dimensions, varMin, varMax); | |
857 | AddTemporaryObject(k); | |
858 | delete[] varMin; delete[] varMax; | |
859 | ||
860 | return k; | |
861 | ||
862 | } | |
863 | ||
864 | //_________________________________________________________________________ | |
3a72645a | 865 | THnSparseF *AliHFEspectrum::GetSlicedCorrelation(THnSparseF *correlationmatrix, Int_t nDim, Int_t *dimensions) const { |
c04c80e6 | 866 | // |
3a72645a | 867 | // Slice correlation |
c04c80e6 | 868 | // |
869 | ||
3a72645a | 870 | Int_t ndimensions = correlationmatrix->GetNdimensions(); |
c04c80e6 | 871 | printf("Number of dimension %d correlation map\n",ndimensions); |
872 | if(ndimensions < (2*nDim)) { | |
873 | AliError("Problem in the dimensions"); | |
874 | return NULL; | |
875 | } | |
876 | Int_t ndimensionsContainer = (Int_t) ndimensions/2; | |
877 | printf("Number of dimension %d container\n",ndimensionsContainer); | |
878 | ||
879 | Int_t *dim = new Int_t[nDim*2]; | |
880 | for(Int_t iter=0; iter < nDim; iter++){ | |
881 | dim[iter] = dimensions[iter]; | |
882 | dim[iter+nDim] = ndimensionsContainer + dimensions[iter]; | |
883 | printf("For iter %d: %d and iter+nDim %d: %d\n",iter,dimensions[iter],iter+nDim,ndimensionsContainer + dimensions[iter]); | |
884 | } | |
885 | ||
3a72645a | 886 | THnSparseF *k = (THnSparseF *) correlationmatrix->Projection(nDim*2,dim); |
c04c80e6 | 887 | |
888 | delete[] dim; | |
889 | return k; | |
890 | ||
891 | } | |
892 | //___________________________________________________________________________ | |
893 | void AliHFEspectrum::CorrectFromTheWidth(TH1D *h1) const { | |
894 | // | |
895 | // Correct from the width of the bins --> dN/dp_{T} (GeV/c)^{-1} | |
896 | // | |
897 | ||
898 | TAxis *axis = h1->GetXaxis(); | |
899 | Int_t nbinX = h1->GetNbinsX(); | |
900 | ||
901 | for(Int_t i = 1; i <= nbinX; i++) { | |
902 | ||
903 | Double_t width = axis->GetBinWidth(i); | |
904 | Double_t content = h1->GetBinContent(i); | |
905 | Double_t error = h1->GetBinError(i); | |
906 | h1->SetBinContent(i,content/width); | |
907 | h1->SetBinError(i,error/width); | |
908 | } | |
909 | ||
910 | } | |
911 | //____________________________________________________________ | |
912 | void AliHFEspectrum::AddTemporaryObject(TObject *o){ | |
913 | // | |
914 | // Emulate garbage collection on class level | |
915 | // | |
916 | if(!fTemporaryObjects) { | |
917 | fTemporaryObjects= new TList; | |
918 | fTemporaryObjects->SetOwner(); | |
919 | } | |
920 | fTemporaryObjects->Add(o); | |
921 | } | |
922 | ||
923 | //____________________________________________________________ | |
924 | void AliHFEspectrum::ClearObject(TObject *o){ | |
925 | // | |
926 | // Do a safe deletion | |
927 | // | |
928 | if(fTemporaryObjects){ | |
929 | if(fTemporaryObjects->FindObject(o)) fTemporaryObjects->Remove(o); | |
930 | delete o; | |
931 | } else{ | |
932 | // just delete | |
933 | delete o; | |
934 | } | |
935 | } | |
936 | //_________________________________________________________________________ | |
937 | TObject* AliHFEspectrum::GetSpectrum(AliCFContainer * const c, Int_t step) { | |
245877d0 | 938 | AliCFDataGrid* data = new AliCFDataGrid("data","",*c, step); |
c04c80e6 | 939 | return data; |
940 | } | |
941 | //_________________________________________________________________________ | |
3a72645a | 942 | TObject* AliHFEspectrum::GetEfficiency(AliCFContainer * const c, Int_t step, Int_t step0){ |
c04c80e6 | 943 | // |
944 | // Create efficiency grid and calculate efficiency | |
945 | // of step to step0 | |
946 | // | |
947 | TString name("eff"); | |
948 | name += step; | |
949 | name+= step0; | |
950 | AliCFEffGrid* eff = new AliCFEffGrid((const char*)name,"",*c); | |
951 | eff->CalculateEfficiency(step,step0); | |
952 | return eff; | |
953 | } |