2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
17 // Class for spectrum correction
18 // Subtraction of hadronic background, Unfolding of the data and
19 // Renormalization done here
20 // The following containers have to be set:
21 // - Correction framework container for real data
22 // - Correction framework container for MC (Efficiency Map)
23 // - Correction framework container for background coming from data
24 // - Correction framework container for background coming from MC
27 // Raphaelle Bailhache <R.Bailhache@gsi.de>
28 // Markus Fasel <M.Fasel@gsi.de>
34 #include <TObjArray.h>
41 #include <TGraphErrors.h>
48 #include "AliCFContainer.h"
49 #include "AliCFDataGrid.h"
50 #include "AliCFEffGrid.h"
51 #include "AliCFGridSparse.h"
52 #include "AliCFUnfolding.h"
55 #include "AliHFEInclusiveSpectrum.h"
56 #include "AliHFEInclusiveSpectrumQA.h"
57 #include "AliHFEcuts.h"
58 #include "AliHFEcontainer.h"
59 #include "AliHFEtools.h"
61 ClassImp(AliHFEInclusiveSpectrum)
63 //____________________________________________________________
64 AliHFEInclusiveSpectrum::AliHFEInclusiveSpectrum(const char *name):
65 AliHFECorrectSpectrumBase(name),
67 fNoCentralitySelectionMC(kFALSE)
70 // Default constructor
73 fQA = new AliHFEInclusiveSpectrumQA("AliHFEInclusiveSpectrumQA");
76 //____________________________________________________________
77 AliHFEInclusiveSpectrum::AliHFEInclusiveSpectrum(const AliHFEInclusiveSpectrum &ref):
78 AliHFECorrectSpectrumBase(ref),
80 fNoCentralitySelectionMC(ref.fNoCentralitySelectionMC)
88 //____________________________________________________________
89 AliHFEInclusiveSpectrum &AliHFEInclusiveSpectrum::operator=(const AliHFEInclusiveSpectrum &ref){
91 // Assignment operator
97 //____________________________________________________________
98 void AliHFEInclusiveSpectrum::Copy(TObject &o) const {
100 // Copy into object o
102 AliHFEInclusiveSpectrum &target = dynamic_cast<AliHFEInclusiveSpectrum &>(o);
104 target.fNoCentralitySelectionMC = fNoCentralitySelectionMC;
108 //____________________________________________________________
109 AliHFEInclusiveSpectrum::~AliHFEInclusiveSpectrum(){
116 //____________________________________________________________
117 Bool_t AliHFEInclusiveSpectrum::Init(const AliHFEcontainer *datahfecontainer, const AliHFEcontainer *mchfecontainer, const AliHFEcontainer */*bghfecontainer*/, const AliHFEcontainer *v0hfecontainer,AliCFContainer *photoniccontainerD){
119 // Init what we need for the correction:
121 // Raw spectrum, hadron contamination
122 // MC efficiency maps, correlation matrix
123 // V0 efficiency if wanted
125 // This for a given dimension.
130 Bool_t centralitySelectionData = kTRUE, centralitySelectionMC = !fNoCentralitySelectionMC;
133 // Data container: raw spectrum + hadron contamination
135 AliCFContainer *datacontainer = datahfecontainer->GetCFContainer("recTrackContReco");
136 AliCFContainer *contaminationcontainer = datahfecontainer->GetCFContainer("hadronicBackground");
137 if((!datacontainer) || (!contaminationcontainer)) return kFALSE;
138 AliCFContainer *datacontainerD = GetSlicedContainer(datacontainer, fNbDimensions, fDims, -1, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionData);
139 AliCFContainer *contaminationcontainerD = GetSlicedContainer(contaminationcontainer, fNbDimensions, fDims, -1, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionData);
140 if((!datacontainerD) || (!contaminationcontainerD)) return kFALSE;
141 SetContainer(datacontainerD,AliHFECorrectSpectrumBase::kDataContainer);
142 SetContainer(contaminationcontainerD,AliHFECorrectSpectrumBase::kBackgroundData);
144 // Photonic Background
145 SetContainer(photoniccontainerD,AliHFECorrectSpectrumBase::kPhotonicBackground);
148 Int_t dimqa = datacontainer->GetNVar();
150 for(Int_t i = 0; i < dimqa; i++) dimsqa[i] = i;
151 AliCFContainer *datacontainerDQA = GetSlicedContainer(datacontainer, dimqa, dimsqa, -1, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionData);
152 fQA->AddResultAt(datacontainerDQA,AliHFEInclusiveSpectrumQA::kDataProjection);
155 // MC container: ESD/MC efficiency maps + MC/MC efficiency maps
157 AliCFContainer *mccontaineresd = 0x0;
158 AliCFContainer *mccontainermc = 0x0;
159 mccontaineresd = mchfecontainer->MakeMergedCFContainer("sumesd","sumesd","MCTrackCont:recTrackContReco");
160 mccontainermc = mchfecontainer->MakeMergedCFContainer("summc","summc","MCTrackCont:recTrackContMC");
161 if((!mccontaineresd) || (!mccontainermc)) return kFALSE;
163 AliCFContainer *mccontaineresdD = GetSlicedContainer(mccontaineresd, fNbDimensions, fDims, source, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionMC);
164 AliCFContainer *mccontainermcD = GetSlicedContainer(mccontainermc, fNbDimensions, fDims, source, fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionMC);
165 if((!mccontaineresdD) || (!mccontainermcD)) return kFALSE;
166 SetContainer(mccontainermcD,AliHFECorrectSpectrumBase::kMCContainerMC);
167 SetContainer(mccontaineresdD,AliHFECorrectSpectrumBase::kMCContainerESD);
170 // Correlation matrix
172 THnSparseF *mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepbeforePID");
173 if(!mccorrelation) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterPID");
174 if(!mccorrelation) return kFALSE;
175 THnSparseF *mccorrelationD = GetSlicedCorrelation(mccorrelation, fNbDimensions, fDims,fChargeChoosen,fTestCentralityLow,fTestCentralityHigh, centralitySelectionMC);
176 if(!mccorrelationD) {
177 printf("No correlation\n");
180 SetCorrelation(mccorrelationD);
182 fQA->AddResultAt(mccorrelation,AliHFEInclusiveSpectrumQA::kCMProjection);
185 // V0 container Electron, pt eta phi
188 AliCFContainer *containerV0 = v0hfecontainer->GetCFContainer("taggedTrackContainerReco");
189 if(!containerV0) return kFALSE;
190 AliCFContainer *containerV0Electron = GetSlicedContainer(containerV0, fNbDimensions, fDims, AliPID::kElectron,fChargeChoosen,fTestCentralityLow,fTestCentralityHigh);
191 if(!containerV0Electron) return kFALSE;
192 SetContainer(containerV0Electron,AliHFECorrectSpectrumBase::kDataContainerV0);
196 fQA->DrawProjections();
201 //____________________________________________________________
202 Bool_t AliHFEInclusiveSpectrum::Correct(Bool_t subtractcontamination, Bool_t subtractphotonic){
204 // Correct the spectrum for efficiency and unfolding
205 // with both method and compare
208 gStyle->SetPalette(1);
209 gStyle->SetOptStat(1111);
210 gStyle->SetPadBorderMode(0);
211 gStyle->SetCanvasColor(10);
212 gStyle->SetPadLeftMargin(0.13);
213 gStyle->SetPadRightMargin(0.13);
215 printf("Steps are: stepdata %d, stepMC %d, steptrue %d, stepV0after %d, stepV0before %d\n",fStepData,fStepMC,fStepTrue,fStepAfterCutsV0,fStepBeforeCutsV0);
217 ///////////////////////////
218 // Check initialization
219 ///////////////////////////
221 if((!GetContainer(kDataContainer)) || (!GetContainer(kMCContainerMC)) || (!GetContainer(kMCContainerESD))){
222 AliInfo("You have to init before");
226 if((fStepTrue < 0) && (fStepMC < 0) && (fStepData < 0)) {
227 AliInfo("You have to set the steps before: SetMCTruthStep, SetMCEffStep, SetStepToCorrect");
231 SetStepGuessedUnfolding(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack);
233 AliCFDataGrid *dataGridAfterFirstSteps = 0x0;
235 //////////////////////////////////
236 // Subtract hadron background
237 /////////////////////////////////
238 AliCFDataGrid *dataspectrumaftersubstraction = 0x0;
239 if(subtractcontamination) {
240 dataspectrumaftersubstraction = SubtractBackground();
241 dataGridAfterFirstSteps = dataspectrumaftersubstraction;
244 //////////////////////////////////
245 // Subtract Photonic background
246 /////////////////////////////////
247 AliCFDataGrid *dataspectrumafterphotonicsubstraction = 0x0;
248 if(subtractphotonic) {
249 dataspectrumafterphotonicsubstraction = SubtractPhotonicBackground();
250 dataGridAfterFirstSteps = dataspectrumafterphotonicsubstraction;
253 ////////////////////////////////////////////////
254 // Correct for TPC efficiency from V0 if any
255 ///////////////////////////////////////////////
256 AliCFDataGrid *dataspectrumafterV0efficiencycorrection = 0x0;
257 AliCFContainer *dataContainerV0 = GetContainer(kDataContainerV0);
259 dataspectrumafterV0efficiencycorrection = CorrectV0Efficiency(dataspectrumaftersubstraction);
260 dataGridAfterFirstSteps = dataspectrumafterV0efficiencycorrection;
263 //////////////////////////////////////////////////////////////////////////////
264 // Correct for efficiency parametrized (if TPC efficiency is parametrized)
265 /////////////////////////////////////////////////////////////////////////////
266 AliCFDataGrid *dataspectrumafterefficiencyparametrizedcorrection = 0x0;
267 if(fEfficiencyFunction){
268 dataspectrumafterefficiencyparametrizedcorrection = CorrectParametrizedEfficiency(dataGridAfterFirstSteps);
269 dataGridAfterFirstSteps = dataspectrumafterefficiencyparametrizedcorrection;
272 TGraphErrors* correctedspectrumD = 0x0;
273 TGraphErrors* alltogetherspectrumD = 0x0;
274 if(fStepMC>fStepTrue) {
278 THnSparse *correctedspectrum = Unfold(dataGridAfterFirstSteps);
279 if(!correctedspectrum){
280 AliError("No corrected spectrum\n");
284 /////////////////////
287 AliCFDataGrid *alltogetherCorrection = CorrectForEfficiency(dataGridAfterFirstSteps);
292 correctedspectrumD = Normalize(correctedspectrum);
293 alltogetherspectrumD = Normalize(alltogetherCorrection);
300 if(dataGridAfterFirstSteps) {
301 correctedspectrumD = Normalize(dataGridAfterFirstSteps);
302 alltogetherspectrumD = Normalize(dataGridAfterFirstSteps);
308 fQA->AddResultAt(correctedspectrumD,AliHFEInclusiveSpectrumQA::kFinalResultUnfolded);
309 fQA->AddResultAt(alltogetherspectrumD,AliHFEInclusiveSpectrumQA::kFinalResultDirectEfficiency);
315 //____________________________________________________________
316 AliCFDataGrid* AliHFEInclusiveSpectrum::SubtractBackground(){
318 // Apply background subtraction
322 AliCFContainer *dataContainer = GetContainer(kDataContainer);
324 AliError("Data Container not available");
327 printf("Step data: %d\n",fStepData);
328 AliCFDataGrid *spectrumSubtracted = new AliCFDataGrid("spectrumSubtracted", "Data Grid for spectrum after Background subtraction", *dataContainer,fStepData);
330 AliCFDataGrid *dataspectrumbeforesubstraction = (AliCFDataGrid *) ((AliCFDataGrid *)GetSpectrum(GetContainer(kDataContainer),fStepData))->Clone();
331 dataspectrumbeforesubstraction->SetName("dataspectrumbeforesubstraction");
334 // Background Estimate
335 AliCFContainer *backgroundContainer = GetContainer(kBackgroundData);
336 if(!backgroundContainer){
337 AliError("MC background container not found");
341 Int_t stepbackground = 1;
342 AliCFDataGrid *backgroundGrid = new AliCFDataGrid("ContaminationGrid","ContaminationGrid",*backgroundContainer,stepbackground);
345 spectrumSubtracted->Add(backgroundGrid,-1.0);
348 TH1D *subtractedspectrum = (TH1D *) spectrumSubtracted->Project(0);
349 AliHFEtools::NormaliseBinWidth(subtractedspectrum);
350 TH1D *rawspectrum = (TH1D *) dataspectrumbeforesubstraction->Project(0);
351 AliHFEtools::NormaliseBinWidth(rawspectrum);
352 fQA->AddResultAt(subtractedspectrum,AliHFEInclusiveSpectrumQA::kAfterSC);
353 fQA->AddResultAt(rawspectrum,AliHFEInclusiveSpectrumQA::kBeforeSC);
354 fQA->DrawSubtractContamination();
356 if(fNbDimensions>=2) {
357 fQA->AddResultAt((TObject *) spectrumSubtracted,AliHFEInclusiveSpectrumQA::kAfterSCND);
358 fQA->AddResultAt((TObject *) dataspectrumbeforesubstraction,AliHFEInclusiveSpectrumQA::kBeforeSCND);
359 fQA->AddResultAt((TObject *) backgroundGrid,AliHFEInclusiveSpectrumQA::kHadronContaminationND);
360 fQA->DrawSubtractContaminationND();
364 return spectrumSubtracted;
367 //____________________________________________________________
368 AliCFDataGrid* AliHFEInclusiveSpectrum::SubtractPhotonicBackground(){
370 // Apply Photonic background subtraction
373 printf("Photonic Background Subtraction \n");
376 AliCFContainer *dataContainer = GetContainer(kDataContainer);
378 AliError("Data Container not available");
381 printf("Step data: %d\n",fStepData);
382 AliCFDataGrid *spectrumPhotonicSubtracted = new AliCFDataGrid("spectrumPhotonicSubtracted", "Data Grid for spectrum after Photonic Background subtraction", *dataContainer,fStepData);
384 AliCFDataGrid *dataSpectrumBeforePhotonicSubstraction = (AliCFDataGrid *) ((AliCFDataGrid *)GetSpectrum(GetContainer(kDataContainer),fStepData))->Clone();
385 dataSpectrumBeforePhotonicSubstraction->SetName("dataSpectrumBeforePhotonicSubstraction");
388 // Background Estimate
389 AliCFContainer *photonicContainer = GetContainer(kPhotonicBackground);
390 if(!photonicContainer){
391 AliError("Photonic background container not found");
395 Int_t stepbackground = 0;
396 AliCFDataGrid *photonicGrid = new AliCFDataGrid("ContaminationGrid","ContaminationGrid",*photonicContainer,stepbackground);
399 spectrumPhotonicSubtracted->Add(photonicGrid,-1.0);
402 TH1D *photonicsubtractedspectrum = (TH1D *) spectrumPhotonicSubtracted->Project(0);
403 AliHFEtools::NormaliseBinWidth(photonicsubtractedspectrum);
404 TH1D *newrawspectrum = (TH1D *) dataSpectrumBeforePhotonicSubstraction->Project(0);
405 AliHFEtools::NormaliseBinWidth(newrawspectrum);
406 fQA->AddResultAt(photonicsubtractedspectrum,AliHFEInclusiveSpectrumQA::kAfterSPB);
407 fQA->AddResultAt(newrawspectrum,AliHFEInclusiveSpectrumQA::kBeforeSPB);
408 fQA->DrawSubtractPhotonicBackground();
410 return spectrumPhotonicSubtracted;
414 //____________________________________________________________
415 AliCFDataGrid *AliHFEInclusiveSpectrum::CorrectParametrizedEfficiency(AliCFDataGrid* const bgsubpectrum){
418 // Apply TPC pid efficiency correction from parametrisation
421 // Data in the right format
422 AliCFDataGrid *dataGrid = 0x0;
424 dataGrid = bgsubpectrum;
428 AliCFContainer *dataContainer = GetContainer(kDataContainer);
430 AliError("Data Container not available");
433 dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData);
435 AliCFDataGrid *result = (AliCFDataGrid *) dataGrid->Clone();
436 result->SetName("ParametrizedEfficiencyBefore");
437 THnSparse *h = result->GetGrid();
438 Int_t nbdimensions = h->GetNdimensions();
439 //printf("CorrectParametrizedEfficiency::We have dimensions %d\n",nbdimensions);
440 AliCFContainer *dataContainer = GetContainer(kDataContainer);
442 AliError("Data Container not available");
445 AliCFContainer *dataContainerbis = (AliCFContainer *) dataContainer->Clone();
446 dataContainerbis->Add(dataContainerbis,-1.0);
449 Int_t* coord = new Int_t[nbdimensions];
450 memset(coord, 0, sizeof(Int_t) * nbdimensions);
451 Double_t* points = new Double_t[nbdimensions];
453 ULong64_t nEntries = h->GetNbins();
454 for (ULong64_t i = 0; i < nEntries; ++i) {
456 Double_t value = h->GetBinContent(i, coord);
457 //Double_t valuecontainer = dataContainerbis->GetBinContent(coord,fStepData);
458 //printf("Value %f, and valuecontainer %f\n",value,valuecontainer);
460 // Get the bin co-ordinates given an coord
461 for (Int_t j = 0; j < nbdimensions; ++j)
462 points[j] = h->GetAxis(j)->GetBinCenter(coord[j]);
464 if (!fEfficiencyFunction->IsInside(points))
466 TF1::RejectPoint(kFALSE);
468 // Evaulate function at points
469 Double_t valueEfficiency = fEfficiencyFunction->EvalPar(points, NULL);
470 //printf("Value efficiency is %f\n",valueEfficiency);
472 if(valueEfficiency > 0.0) {
473 h->SetBinContent(coord,value/valueEfficiency);
474 dataContainerbis->SetBinContent(coord,fStepData,value/valueEfficiency);
475 Double_t error = h->GetBinError(i);
476 h->SetBinError(coord,error/valueEfficiency);
477 dataContainerbis->SetBinError(coord,fStepData,error/valueEfficiency);
485 AliCFDataGrid *resultt = new AliCFDataGrid("spectrumEfficiencyParametrized", "Data Grid for spectrum after Efficiency parametrized", *dataContainerbis,fStepData);
488 TH1D *afterE = (TH1D *) resultt->Project(0);
489 AliHFEtools::NormaliseBinWidth(afterE);
490 TH1D *beforeE = (TH1D *) dataGrid->Project(0);
491 AliHFEtools::NormaliseBinWidth(beforeE);
492 fQA->AddResultAt(afterE,AliHFEInclusiveSpectrumQA::kAfterPE);
493 fQA->AddResultAt(beforeE,AliHFEInclusiveSpectrumQA::kBeforePE);
494 fQA->AddResultAt(fEfficiencyFunction,AliHFEInclusiveSpectrumQA::kPEfficiency);
495 fQA->DrawCorrectWithEfficiency(AliHFEInclusiveSpectrumQA::kParametrized);
497 if(fNbDimensions>=2) {
498 fQA->AddResultAt((TObject *) resultt,AliHFEInclusiveSpectrumQA::kAfterPEND);
499 fQA->AddResultAt((TObject *) dataGrid,AliHFEInclusiveSpectrumQA::kBeforePEND);
500 fQA->AddResultAt((TObject *) fEfficiencyFunction,AliHFEInclusiveSpectrumQA::kPEfficiencyND);
501 fQA->DrawCorrectWithEfficiencyND(AliHFEInclusiveSpectrumQA::kParametrized);
507 //____________________________________________________________
508 AliCFDataGrid *AliHFEInclusiveSpectrum::CorrectV0Efficiency(AliCFDataGrid* const bgsubpectrum){
511 // Apply TPC pid efficiency correction from V0
514 AliCFContainer *v0Container = GetContainer(kDataContainerV0);
516 AliError("V0 Container not available");
521 AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*v0Container);
522 efficiencyD->CalculateEfficiency(fStepAfterCutsV0,fStepBeforeCutsV0);
524 // Data in the right format
525 AliCFDataGrid *dataGrid = 0x0;
527 dataGrid = bgsubpectrum;
530 AliCFContainer *dataContainer = GetContainer(kDataContainer);
532 AliError("Data Container not available");
535 dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData);
539 AliCFDataGrid *result = (AliCFDataGrid *) dataGrid->Clone();
540 result->ApplyEffCorrection(*efficiencyD);
543 TH1D *afterE = (TH1D *) result->Project(0);
544 AliHFEtools::NormaliseBinWidth(afterE);
545 TH1D *beforeE = (TH1D *) dataGrid->Project(0);
546 AliHFEtools::NormaliseBinWidth(beforeE);
547 TH1D* efficiencyDproj = (TH1D *) efficiencyD->Project(0);
548 fQA->AddResultAt(afterE,AliHFEInclusiveSpectrumQA::kAfterV0);
549 fQA->AddResultAt(beforeE,AliHFEInclusiveSpectrumQA::kBeforeV0);
550 fQA->AddResultAt(efficiencyDproj,AliHFEInclusiveSpectrumQA::kV0Efficiency);
551 fQA->DrawCorrectWithEfficiency(AliHFEInclusiveSpectrumQA::kV0);
553 if(fNbDimensions>=2) {
554 fQA->AddResultAt((TObject *) result,AliHFEInclusiveSpectrumQA::kAfterV0ND);
555 fQA->AddResultAt((TObject *) dataGrid,AliHFEInclusiveSpectrumQA::kBeforeV0ND);
556 fQA->AddResultAt((TObject *) efficiencyD,AliHFEInclusiveSpectrumQA::kV0EfficiencyND);
557 fQA->DrawCorrectWithEfficiencyND(AliHFEInclusiveSpectrumQA::kV0);
564 //____________________________________________________________
565 THnSparse *AliHFEInclusiveSpectrum::Unfold(AliCFDataGrid* const bgsubpectrum){
568 // Return the unfolded spectrum
571 AliCFContainer *mcContainer = GetContainer(kMCContainerMC);
573 AliError("MC Container not available");
578 AliError("No Correlation map available");
583 AliCFDataGrid *dataGrid = 0x0;
585 dataGrid = bgsubpectrum;
589 AliCFContainer *dataContainer = GetContainer(kDataContainer);
591 AliError("Data Container not available");
595 dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData);
599 AliCFDataGrid* guessedGrid = new AliCFDataGrid("guessed","",*mcContainer, fStepGuessedUnfolding);
600 THnSparse* guessedTHnSparse = ((AliCFGridSparse*)guessedGrid->GetData())->GetGrid();
603 AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*mcContainer);
604 efficiencyD->CalculateEfficiency(fStepMC,fStepTrue);
608 AliCFUnfolding unfolding("unfolding","",fNbDimensions,fCorrelation,efficiencyD->GetGrid(),dataGrid->GetGrid(),guessedTHnSparse,1.e-06,0,fNumberOfIterations);
609 if(fSetSmoothing) unfolding.UseSmoothing();
613 THnSparse* result = unfolding.GetUnfolded();
614 THnSparse* residual = unfolding.GetEstMeasured();
617 TH1D *residualh = (TH1D *) residual->Projection(0);
618 TH1D *beforeE = (TH1D *) dataGrid->Project(0);
619 TH1D* efficiencyDproj = (TH1D *) efficiencyD->Project(0);
620 TH1D *afterE = (TH1D *) result->Projection(0);
621 AliHFEtools::NormaliseBinWidth(residualh);
622 AliHFEtools::NormaliseBinWidth(beforeE);
623 AliHFEtools::NormaliseBinWidth(afterE);
624 fQA->AddResultAt(residualh,AliHFEInclusiveSpectrumQA::kResidualU);
625 fQA->AddResultAt(afterE,AliHFEInclusiveSpectrumQA::kAfterU);
626 fQA->AddResultAt(beforeE,AliHFEInclusiveSpectrumQA::kBeforeU);
627 fQA->AddResultAt(efficiencyDproj,AliHFEInclusiveSpectrumQA::kUEfficiency);
628 fQA->DrawUnfolding();
630 return (THnSparse *) result->Clone();
633 //____________________________________________________________
634 AliCFDataGrid *AliHFEInclusiveSpectrum::CorrectForEfficiency(AliCFDataGrid* const bgsubpectrum){
637 // Apply unfolding and efficiency correction together to bgsubspectrum
640 AliCFContainer *mcContainer = GetContainer(kMCContainerESD);
642 AliError("MC Container not available");
647 AliCFEffGrid* efficiencyD = new AliCFEffGrid("efficiency","",*mcContainer);
648 efficiencyD->CalculateEfficiency(fStepMC,fStepTrue);
650 // Data in the right format
651 AliCFDataGrid *dataGrid = 0x0;
653 dataGrid = bgsubpectrum;
657 AliCFContainer *dataContainer = GetContainer(kDataContainer);
659 AliError("Data Container not available");
663 dataGrid = new AliCFDataGrid("dataGrid","dataGrid",*dataContainer, fStepData);
667 AliCFDataGrid *result = (AliCFDataGrid *) dataGrid->Clone();
668 result->ApplyEffCorrection(*efficiencyD);
671 TH1D *afterE = (TH1D *) result->Project(0);
672 AliHFEtools::NormaliseBinWidth(afterE);
673 TH1D *beforeE = (TH1D *) dataGrid->Project(0);
674 AliHFEtools::NormaliseBinWidth(beforeE);
675 TH1D* efficiencyDproj = (TH1D *) efficiencyD->Project(0);
676 fQA->AddResultAt(afterE,AliHFEInclusiveSpectrumQA::kAfterMCE);
677 fQA->AddResultAt(beforeE,AliHFEInclusiveSpectrumQA::kBeforeMCE);
678 fQA->AddResultAt(efficiencyDproj,AliHFEInclusiveSpectrumQA::kMCEfficiency);
679 fQA->DrawCorrectWithEfficiency(AliHFEInclusiveSpectrumQA::kMC);
681 if(fNbDimensions>=2) {
682 fQA->AddResultAt((TObject *) result,AliHFEInclusiveSpectrumQA::kAfterMCEND);
683 fQA->AddResultAt((TObject *) dataGrid,AliHFEInclusiveSpectrumQA::kBeforeMCEND);
684 fQA->AddResultAt((TObject *) efficiencyD,AliHFEInclusiveSpectrumQA::kMCEfficiencyND);
685 fQA->DrawCorrectWithEfficiencyND(AliHFEInclusiveSpectrumQA::kMC);
691 //____________________________________________________________
692 void AliHFEInclusiveSpectrum::WriteResults(const char *filename)
698 AliCFContainer *dataContainer = GetContainer(kDataContainer);
699 AliCFContainer *mcContainer = GetContainer(kMCContainerMC);
700 TObject *unfolded = 0x0;
701 TObject *correctedspectrum = 0x0;
703 unfolded = fQA->GetResult(AliHFEInclusiveSpectrumQA::kFinalResultUnfolded);
704 correctedspectrum = fQA->GetResult(AliHFEInclusiveSpectrumQA::kFinalResultDirectEfficiency);
707 TFile *file = TFile::Open(filename,"recreate");
708 if(dataContainer) dataContainer->Write("data");
709 if(mcContainer) mcContainer->Write("mcefficiency");
710 if(fCorrelation) fCorrelation->Write("correlationmatrix");
711 if(unfolded) unfolded->Write("unfoldedspectrum");
712 if(correctedspectrum) correctedspectrum->Write("correctedspectrum");
713 if(fQA) fQA->Write("QAResults");