PWGLF/SPECTRA/Nuclei/B2/AliLnSpectra.cxx

   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  **************************************************************************/
  15
  16 // (invariant) differential yield
  17 // author: Eulogio Serradilla <eulogio.serradilla@cern.ch>
  18
  19 #include <Riostream.h>
  20 #include <TFile.h>
  21 #include <TH1D.h>
  22 #include <TString.h>
  23 #include <TMath.h>
  24 #include <TGraphErrors.h>
  25 #include <TF1.h>
  26
  27 #include "AliLnSpectra.h"
  28 #include "B2.h"
  29
  30 ClassImp(AliLnSpectra)
  31
  32 AliLnSpectra::AliLnSpectra(const TString& particle, const TString& ptFilename, const TString& tag, const TString& outputFilename, const TString& otag)
  33 : TObject()
  34 , fParticle(particle)
  35 , fPtFilename(ptFilename)
  36 , fTag(tag)
  37 , fOutputFilename(outputFilename)
  38 , fOutputTag(otag)
  39 , fYMin(-0.5)
  40 , fYMax(0.5)
  41 , fINEL(1)
  42 {
  43 //
  44 // constructor
  45 //
  46 }
  47
  48 AliLnSpectra::~AliLnSpectra()
  49 {
  50 //
  51 // destructor
  52 //
  53 }
  54
  55 Int_t AliLnSpectra::Exec()
  56 {
  57 //
  58 // (invariant) differential yield projection in pt for |y| < ymax-ymin
  59 //
  60         using namespace std;
  61
  62         TFile* finput = new TFile(fPtFilename.Data(), "read");
  63         if(finput->IsZombie()) exit(1);
  64
  65         // pt and number of events
  66
  67         TH1D* hPt = FindObj<TH1D>(finput, fParticle + "_Pt");
  68         TH1D* hStats = FindObj<TH1D>(finput, fParticle + "_Stats");
  69
  70         Double_t nEvent = hStats->Integral(3,3);
  71         if(fINEL) nEvent = hStats->Integral(4,4);
  72
  73         // ouputfile
  74
  75         TFile* foutput = new TFile(fOutputFilename.Data(),"recreate");
  76
  77         foutput->mkdir(fOutputTag.Data());
  78         foutput->cd(fOutputTag.Data());
  79
  80         // differential yield
  81
  82         TGraphErrors* grDYieldPt = this->GetDiffYieldPt(hPt, nEvent, fParticle + "_DiffYield_Pt");
  83         grDYieldPt->Write();
  84
  85         // invariant differential yield
  86
  87         TGraphErrors* grInvDYieldPt = this->GetInvDiffYieldPt(grDYieldPt, fParticle + "_InvDiffYield_Pt");
  88         grInvDYieldPt->Write();
  89
  90         // clean
  91
  92         delete grDYieldPt;
  93         delete grInvDYieldPt;
  94
  95         delete foutput;
  96         delete finput;
  97
  98         return 0;
  99 }
 100
 101 TGraphErrors* AliLnSpectra::GetDiffYieldPt(const TH1D* hPt, Int_t nEvent, const TString& name) const
 102 {
 103 //
 104 // projection of the differential yield in pt
 105 // only statistical error of N is taken into account
 106 //
 107         TGraphErrors* grDYieldPt = new TGraphErrors();
 108         grDYieldPt->SetName(name.Data());
 109
 110         Double_t dy = fYMax - fYMin;
 111
 112         TAxis* xaxis = hPt->GetXaxis();
 113
 114         for(Int_t i=1, j=0; i < xaxis->GetNbins()+1; ++i)
 115         {
 116                 Double_t pt  = xaxis->GetBinCenter(i);
 117                 Double_t dpt = xaxis->GetBinWidth(i);
 118
 119                 Double_t n = hPt->GetBinContent(i);
 120                 Double_t eN = hPt->GetBinError(i);
 121
 122                 if(n==0) continue;
 123
 124                 Double_t yield = n/(nEvent*dpt*dy);
 125                 Double_t yError = yield*eN/n;
 126
 127                 grDYieldPt->SetPoint(j, pt, yield);
 128                 grDYieldPt->SetPointError(j++, dpt/2., yError);
 129         }
 130
 131         return grDYieldPt;
 132 }
 133
 134 TGraphErrors* AliLnSpectra::GetInvDiffYieldPt(const TGraphErrors* grDYieldPt, const TString& name) const
 135 {
 136 //
 137 // projection of the invariant differential yield in pt
 138 // only statistical error of N is taken into account
 139 //
 140         TGraphErrors* grInvDYieldPt = new TGraphErrors();
 141         grInvDYieldPt->SetName(name.Data());
 142
 143         for(Int_t i=0; i < grDYieldPt->GetN(); ++i)
 144         {
 145                 Double_t pt, y;
 146                 grDYieldPt->GetPoint(i,pt,y);
 147
 148                 Double_t errpt = grDYieldPt->GetErrorX(i);
 149                 Double_t erry = grDYieldPt->GetErrorY(i);
 150
 151                 y = y/(2.*TMath::Pi()*pt);
 152                 erry = erry / (2.*TMath::Pi()*pt);
 153
 154                 grInvDYieldPt->SetPoint(i, pt, y);
 155                 grInvDYieldPt->SetPointError(i, errpt, erry);
 156         }
 157
 158         return grInvDYieldPt;
 159 }
 160
 161 TGraphErrors* AliLnSpectra::GetInvDiffYieldPt(const TH1D* hPt, Int_t nEvent, const TString& name) const
 162 {
 163 //
 164 // projection of the invariant differential yield in pt
 165 // only statistical error of N is taken into account
 166 //
 167         TGraphErrors* grDYieldPt = this->GetDiffYieldPt(hPt, nEvent, "_tmp_diff_yield_");
 168         TGraphErrors* grInvDYieldPt = this->GetInvDiffYieldPt(grDYieldPt, name);
 169
 170         delete grDYieldPt;
 171         return grInvDYieldPt;
 172 }