New class for dE/dx analysis (comparison with Bethe-Bloch, check of response function...
authorprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 10 Jun 2009 15:58:53 +0000 (15:58 +0000)
committerprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 10 Jun 2009 15:58:53 +0000 (15:58 +0000)
ITS/AliITSdEdxAnalyzer.cxx [new file with mode: 0644]
ITS/AliITSdEdxAnalyzer.h [new file with mode: 0644]
ITS/CMake_libITSrec.txt
ITS/ITSrecLinkDef.h
ITS/MakeITSdEdxAnalysis.C [new file with mode: 0644]
ITS/libITSrec.pkg

diff --git a/ITS/AliITSdEdxAnalyzer.cxx b/ITS/AliITSdEdxAnalyzer.cxx
new file mode 100644 (file)
index 0000000..c5c981f
--- /dev/null
@@ -0,0 +1,250 @@
+/**************************************************************************
+ * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id$ */
+
+///////////////////////////////////////////////////////////////////
+//                                                               //
+// Implementation of the base class for dEdx analysis            //
+// Origin: F.Prino, Torino, prino@to.infn.it                     //
+//                                                               //
+///////////////////////////////////////////////////////////////////
+
+#include <TString.h>
+#include <TFile.h>
+#include <TParticlePDG.h>
+#include <TDatabasePDG.h>
+#include "AliITSdEdxAnalyzer.h"
+#include "AliExternalTrackParam.h"
+#include "AliITSpidESD.h"
+#include "AliPID.h"
+
+ClassImp(AliITSdEdxAnalyzer)
+
+const Int_t AliITSdEdxAnalyzer::fgkPdgCode[kNParticles]={211,321,2212};
+const Int_t AliITSdEdxAnalyzer::fgkLayerCode[kNValuesdEdx]={3,4,5,6,0};
+
+//______________________________________________________________________
+AliITSdEdxAnalyzer::AliITSdEdxAnalyzer() :  
+  TObject(),
+  fNPBins(10),
+  fPMin(0.1),
+  fPMax(1.1),
+  fHistodEdx(0),
+  fHistoDeltadEdx(0),
+  fHistodEdxVsP(0),
+  fThickness(0.03),
+  fDensity(2.33),
+  fBBmodel(0),
+  fMIP(82.),
+  fTPCpidCut(0.5)
+{
+  // default constructor
+  BookHistos();
+}
+//______________________________________________________________________
+AliITSdEdxAnalyzer::AliITSdEdxAnalyzer(const Int_t npbins, const Float_t pmin, const Float_t pmax) :  
+  TObject(),
+  fNPBins(npbins),
+  fPMin(pmin),
+  fPMax(pmax),
+  fHistodEdx(0),
+  fHistoDeltadEdx(0),
+  fHistodEdxVsP(0),
+  fThickness(0.03),
+  fDensity(2.33),
+  fBBmodel(0),
+  fMIP(82.),
+  fTPCpidCut(0.5)
+{
+  // standard constructor
+  BookHistos();
+}
+//______________________________________________________________________
+AliITSdEdxAnalyzer::~AliITSdEdxAnalyzer(){
+  // destructor
+  DeleteHistos();
+}
+//______________________________________________________________________
+void AliITSdEdxAnalyzer::SetMomentumBins(const Int_t npbins, const Float_t pmin, const Float_t pmax){
+  // Kill exisiting histos, set new binning, book new histos
+  DeleteHistos();
+  fNPBins=npbins;
+  fPMin=pmin;
+  fPMax=pmax;
+  BookHistos();
+}
+//______________________________________________________________________
+void AliITSdEdxAnalyzer::DeleteHistos(){
+  //
+  if(fHistodEdx){
+    for(Int_t i=0; i<GetNHistos();i++) delete fHistodEdx[i];
+    delete [] fHistodEdx;
+    fHistodEdx=0;
+  }
+  if(fHistoDeltadEdx){
+    for(Int_t i=0; i<GetNHistos();i++) delete fHistoDeltadEdx[i];
+    delete [] fHistoDeltadEdx;
+    fHistoDeltadEdx=0;
+  }
+  if(fHistodEdxVsP){
+    for(Int_t i=0; i<GetNHistos2();i++) delete fHistodEdxVsP[i];
+    delete [] fHistodEdxVsP;
+    fHistodEdxVsP=0;
+  }
+}
+//______________________________________________________________________
+void AliITSdEdxAnalyzer::BookHistos(){
+  //
+  fHistodEdx=new TH1F*[GetNHistos()];
+  fHistoDeltadEdx=new TH1F*[GetNHistos()];
+  for(Int_t iSpecie=0; iSpecie<kNParticles; iSpecie++){
+    for(Int_t iPBin=0; iPBin<fNPBins; iPBin++){
+      for(Int_t iVal=0; iVal<kNValuesdEdx; iVal++){
+       TString hisnam1=Form("hdEdx%dpbin%dpart%d",fgkLayerCode[iVal],iPBin,fgkPdgCode[iSpecie]);
+       TString histit1=Form("dEdx layer %d (keV) pbin%d part%d",fgkLayerCode[iVal],iPBin,fgkPdgCode[iSpecie]);
+       if(iVal==kNValuesdEdx-1) histit1=Form("dEdx trunc. mean (keV) pbin%d part%d",iPBin,fgkPdgCode[iSpecie]);
+       fHistodEdx[GetIndex(iSpecie,iPBin,iVal)]=new TH1F(hisnam1.Data(),histit1.Data(),100,0.,500.);
+       fHistodEdx[GetIndex(iSpecie,iPBin,iVal)]->GetXaxis()->SetTitle(histit1.Data());
+       TString hisnam2=Form("hdeltadEdx%dpbin%dpart%d",fgkLayerCode[iVal],iPBin,fgkPdgCode[iSpecie]);
+       TString histit2=Form("(dEdx_l%d-BetheBloch)/BetheBloch pbin%d part%d",fgkLayerCode[iVal],iPBin,fgkPdgCode[iSpecie]);
+       if(iVal==kNValuesdEdx-1) histit1=Form("(dEdx_truncmean-BetheBloch)/BetheBloch pbin%d part%d",iPBin,fgkPdgCode[iSpecie]);
+       fHistoDeltadEdx[GetIndex(iSpecie,iPBin,iVal)]=new TH1F(hisnam2.Data(),histit2.Data(),100,-0.5,0.5);
+       fHistoDeltadEdx[GetIndex(iSpecie,iPBin,iVal)]->GetXaxis()->SetTitle(histit2.Data());
+      }
+    }
+  }
+
+  fHistodEdxVsP=new TH2F*[GetNHistos2()];
+  for(Int_t iSpecie=0; iSpecie<kNParticles; iSpecie++){
+    for(Int_t iVal=0; iVal<kNValuesdEdx; iVal++){
+      TString hisnam=Form("hdEdx%dVsPpart%d",fgkLayerCode[iVal],fgkPdgCode[iSpecie]);
+      TString histit=Form("dEdx layer %d (keV) vs P part%d",fgkLayerCode[iVal],fgkPdgCode[iSpecie]);
+      if(iVal==kNValuesdEdx-1) histit=Form("dEdx trunc. mean (keV) vs P part%d",fgkPdgCode[iSpecie]);      
+      fHistodEdxVsP[GetIndex2(iSpecie,iVal)]=new TH2F(hisnam.Data(),histit.Data(),50,fPMin,fPMax,50,0.,500.);
+      histit.ReplaceAll(" vs P "," ");
+      fHistodEdxVsP[GetIndex2(iSpecie,iVal)]->GetYaxis()->SetTitle(histit.Data());
+      fHistodEdxVsP[GetIndex2(iSpecie,iVal)]->GetXaxis()->SetTitle("Momentum (GeV/c)");
+    }
+  }
+}
+//______________________________________________________________________
+void AliITSdEdxAnalyzer:: WriteHistos(TString filename) const {
+  TFile *outfile=new TFile(filename.Data(),"recreate");
+  outfile->cd();
+  for(Int_t i=0; i<GetNHistos();i++){ 
+    fHistodEdx[i]->Write();
+    fHistoDeltadEdx[i]->Write();
+  }
+  for(Int_t i=0; i<GetNHistos2();i++) fHistodEdxVsP[i]->Write();
+  outfile->Close();
+  delete outfile;
+}
+//______________________________________________________________________
+void AliITSdEdxAnalyzer::ReadEvent(AliESDEvent* esd, AliStack* stack){
+  // Fill histos 
+  // if stack is !=0 MC truth is used to define particle specie
+  // else TPC pid is used to define particle specie
+
+  if(!esd) AliFatal("Bad ESD event");
+  for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
+    AliESDtrack* track = esd->GetTrack(iTrack);
+    Double_t trmean=track->GetITSsignal();
+    Double_t sig[4];
+    track->GetITSdEdxSamples(sig);
+    Double_t preco=track->P();
+    Int_t iPBin=GetMomentumBin(preco);
+    if(iPBin==-1) continue;
+    Int_t iSpecie=-1;
+    Double_t dedxbb=0;
+    if(stack){
+      Int_t lab=track->GetLabel();
+      if(lab<0) continue;
+      TParticle* part=(TParticle*)stack->Particle(lab);
+      Int_t absPdgCode=TMath::Abs(part->GetPdgCode());
+      iSpecie=GetSpecieBin(absPdgCode);
+      if(iSpecie==-1) continue;
+      dedxbb=BetheBloch(part);
+    }else{   
+      iSpecie=GetPaticleIdFromTPC(track);
+      if(iSpecie==-1) continue;
+      TParticlePDG* p=TDatabasePDG::Instance()->GetParticle(fgkPdgCode[iSpecie]);    
+      dedxbb=BetheBloch(preco,p->Mass());
+    }
+    for(Int_t ilay=0;ilay<4;ilay++){    
+      if(sig[ilay]>0.){
+       fHistodEdx[GetIndex(iSpecie,iPBin,ilay)]->Fill(sig[ilay]);
+       fHistoDeltadEdx[GetIndex(iSpecie,iPBin,ilay)]->Fill((sig[ilay]-dedxbb)/dedxbb);
+       fHistodEdxVsP[GetIndex2(iSpecie,ilay)]->Fill(preco,sig[ilay]);
+      }
+    }
+    if(trmean>0.){
+      fHistodEdx[GetIndex(iSpecie,iPBin,4)]->Fill(trmean);
+      fHistoDeltadEdx[GetIndex(iSpecie,iPBin,4)]->Fill((trmean-dedxbb)/dedxbb);
+      fHistodEdxVsP[GetIndex2(iSpecie,4)]->Fill(preco,trmean);
+    }
+  }
+}
+//______________________________________________________________________
+Int_t AliITSdEdxAnalyzer::GetPaticleIdFromTPC(const AliESDtrack* track) const {
+  //
+  Double_t tpcpid[AliPID::kSPECIES];
+  track->GetTPCpid(tpcpid);
+  Int_t maxPos=-1;
+  Float_t maxProb=0.;
+  for(Int_t is=0; is<AliPID::kSPECIES; is++){
+    if(tpcpid[is]>maxProb){
+      maxProb=tpcpid[is];
+      maxPos=is;
+    }
+  }
+  Int_t iSpecie=-1;
+  if(maxProb>fTPCpidCut){
+    if(maxPos==AliPID::kPion) iSpecie=0;
+    else if(maxPos==AliPID::kKaon) iSpecie=1;
+    else if(maxPos==AliPID::kProton) iSpecie=2;    
+  }
+  return iSpecie;
+}
+//______________________________________________________________________
+Double_t AliITSdEdxAnalyzer::BetheBloch(const Float_t p, const Float_t m) const {
+  //
+  Double_t dedxbb=0.;
+  switch(fBBmodel){
+  case 0:
+    Double_t betagamma=p/m;
+    Double_t conv=fDensity*1E6*fThickness/116.31*fMIP;
+    dedxbb=conv*AliExternalTrackParam::BetheBlochSolid(betagamma);
+    break;
+  case 1:
+    dedxbb=fMIP*AliITSpidESD::Bethe(p,m);
+    break;
+  }
+  return dedxbb;
+}
+//______________________________________________________________________
+TGraph* AliITSdEdxAnalyzer::GetBetheBlochGraph(const Int_t pdgcode) const {
+  // Fills a TGraph with Bethe-Bloch expe
+  TGraph* g=new TGraph(0);
+  TParticlePDG* p=TDatabasePDG::Instance()->GetParticle(pdgcode);
+  Float_t mass=p->Mass();
+  for(Int_t ip=0; ip<100;ip++){
+    Float_t p=fPMin+(ip+0.5)*(fPMax-fPMin)/100.;
+    g->SetPoint(ip,p,BetheBloch(p,mass));
+  }
+  g->GetXaxis()->SetTitle("Momentum (GeV/c)");
+  g->GetYaxis()->SetTitle("dE/dx (keV/300 #mum)");
+  return g;
+}
diff --git a/ITS/AliITSdEdxAnalyzer.h b/ITS/AliITSdEdxAnalyzer.h
new file mode 100644 (file)
index 0000000..64a7199
--- /dev/null
@@ -0,0 +1,132 @@
+#ifndef ALIITSDEDXANALYZER_H
+#define ALIITSDEDXANALYZER_H
+/* Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+///////////////////////////////////////////////////////////////////
+//                                                               //
+// Class for:                                                    //
+// - building histograms of dE/dx in bins of pt/layer/particle   //
+// - comparing dEdx signal with Bethe-Bloch expected value       //
+// Origin: F.Prino, Torino, prino@to.infn.it                     //
+//                                                               //
+///////////////////////////////////////////////////////////////////
+
+#include <TObject.h>
+#include <TParticle.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TGraph.h>
+#include "AliStack.h"
+#include "AliLog.h"
+#include "AliESDEvent.h"
+
+class AliITSdEdxAnalyzer : public TObject {
+
+ public:
+  AliITSdEdxAnalyzer();
+  AliITSdEdxAnalyzer(const Int_t npbins, const Float_t pmin, const Float_t pmax);
+  virtual ~AliITSdEdxAnalyzer();
+
+  void SetMaterial(Float_t dens, Float_t thick){
+    fDensity=dens;
+    fThickness=thick;
+  }
+  void SetMomentumBins(const Int_t npbins, const Float_t pmin, const Float_t pmax);
+  void SetUseBBFromAliExternalTrackParam() {fBBmodel=0;}
+  void SetUseBBFromAliITSpidESD() {fBBmodel=1;}
+  void SetMIPdEdx(Float_t val){fMIP=val;}
+  void SetTPCMinimumPIDProb(Float_t min){fTPCpidCut=min;}
+
+  void ReadEvent(AliESDEvent* ev, AliStack* stack=0);
+
+
+  Double_t BetheBloch(TParticle* part) const{
+    return BetheBloch(part->P(),part->GetMass());
+  }
+  Double_t BetheBloch(const Float_t p, const Float_t m) const;
+
+
+  TH1F* GetSingleLayerdEdxHisto(const Int_t lay, const Int_t pdgcode, const Int_t pbin) const {
+    if(lay<3 || lay>6) AliFatal("Wrong LayerNumber");
+    return fHistodEdx[GetIndex(GetSpecieBin(pdgcode),pbin,lay-3)];
+  }
+  TH1F* GetTruncatedMeandEdxHisto(const Int_t pdgcode, const Int_t pbin) const {
+    return fHistodEdx[GetIndex(GetSpecieBin(pdgcode),pbin,4)];
+  }
+
+  TH1F* GetSingleLayerDeltadEdxHisto(const Int_t lay, const Int_t pdgcode, const Int_t pbin) const {
+    if(lay<3 || lay>6) AliFatal("Wrong LayerNumber");
+    return fHistoDeltadEdx[GetIndex(GetSpecieBin(pdgcode),pbin,lay-3)];
+  }
+  TH1F* GetTruncatedMeanDeltadEdxHisto(const Int_t pdgcode, const Int_t pbin) const {
+    return fHistoDeltadEdx[GetIndex(GetSpecieBin(pdgcode),pbin,4)];
+  }
+
+  TH2F* GetSingleLayerdEdxVsPHisto(const Int_t lay, const Int_t pdgcode) const {
+    if(lay<3 || lay>6) AliFatal("Wrong LayerNumber");
+    return fHistodEdxVsP[GetIndex2(GetSpecieBin(pdgcode),lay-3)];
+  }
+  TH2F* GetTruncatedMeandEdxVsPHisto(const Int_t pdgcode) const {
+    return fHistodEdxVsP[GetIndex2(GetSpecieBin(pdgcode),4)];
+  }
+
+  TGraph* GetBetheBlochGraph(const Int_t pdgcode) const;
+
+  void WriteHistos(TString filename="ITS.dEdx.root") const;
+
+ protected:
+
+  Int_t GetNHistos() const {return kNParticles*kNValuesdEdx*fNPBins;}
+  Int_t GetNHistos2() const {return kNParticles*kNValuesdEdx;}
+  Int_t GetIndex(const Int_t iSp, const Int_t iP, const Int_t iVal) const {
+    return iVal+kNValuesdEdx*iP+fNPBins*kNValuesdEdx*iSp;
+  }
+  Int_t GetIndex2(const Int_t iSp, const Int_t iVal) const {
+    return iVal+kNValuesdEdx*iSp;
+  }
+  Int_t GetMomentumBin(Float_t p) const{
+    Int_t iBin=(Int_t)((p-fPMin)/(fPMax-fPMin)*fNPBins);
+    if(iBin>=0 && iBin<fNPBins)return iBin; 
+    return -1;
+  }
+  Int_t GetSpecieBin(Int_t absPdgCode) const {
+    for(Int_t iS=0; iS<kNParticles; iS++){
+      if(absPdgCode==fgkPdgCode[iS]) return iS;
+    }
+    return -1;
+  }
+
+  Int_t GetPaticleIdFromTPC(const AliESDtrack* track) const;
+  void BookHistos();
+  void DeleteHistos();
+  
+ private:
+
+  AliITSdEdxAnalyzer(const AliITSdEdxAnalyzer& dum);
+  AliITSdEdxAnalyzer& operator=(const AliITSdEdxAnalyzer& dum);
+
+  enum {kNParticles=3};  
+  enum {kNValuesdEdx=5};
+
+  static const Int_t fgkPdgCode[kNParticles]; // initialized in the cxx
+  static const Int_t fgkLayerCode[kNValuesdEdx]; // initialized in the cxx
+  
+
+  Int_t   fNPBins;         // Number of Momentum bins
+  Float_t fPMin;           // Minimum P
+  Float_t fPMax;           // Maximum P
+  TH1F**  fHistodEdx;      // Array of histograms of dEdx_meas in bins of pt
+  TH1F**  fHistoDeltadEdx; // Array of histograms of dEdx_meas-dEdx_expected in bins of pt
+  TH2F**  fHistodEdxVsP;   // Array of histograms of dEdx_meas vs. pt
+  Float_t fThickness;      // detector thickness (cm)
+  Float_t fDensity;        // detector density (g/cm3)
+  Int_t   fBBmodel;        // is there MC truth info?
+  Float_t fMIP;            // normalization for MIP
+  Float_t fTPCpidCut;      // minimum probability for PID in TPC   
+
+  ClassDef(AliITSdEdxAnalyzer,0);
+};
+#endif
index a2c1e74..c38d460 100644 (file)
@@ -38,6 +38,7 @@ set(SRCS
                AliITSpidESD.cxx 
                AliITSpidESD1.cxx 
                AliITSpidESD2.cxx 
+               AliITSdEdxAnalyzer.cxx
                 AliITSPident.cxx 
                 AliITSSteerPid.cxx 
                 AliITSPidParItem.cxx 
index f8cf030..56cea4c 100644 (file)
@@ -68,6 +68,7 @@
 #pragma link C++ class AliITSpidESD+;
 #pragma link C++ class AliITSpidESD1+;
 #pragma link C++ class AliITSpidESD2+;
+#pragma link C++ class AliITSdEdxAnalyzer+;
 //beam test classes
 #pragma link C++ class AliITSBeamTestDig+;
 #pragma link C++ class AliITSBeamTestDigSPD+;
diff --git a/ITS/MakeITSdEdxAnalysis.C b/ITS/MakeITSdEdxAnalysis.C
new file mode 100644 (file)
index 0000000..f3719bb
--- /dev/null
@@ -0,0 +1,119 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include <TCanvas.h>
+#include <TClassTable.h>
+#include <TGraph.h>
+#include <TTree.h>
+#include <TH1.h>
+#include <TFile.h>
+#include <TH2.h>
+#include <TInterpreter.h>
+#include <TStyle.h>
+#include "AliHeader.h"
+#include "AliITSdEdxAnalyzer.h"
+#include "AliRunLoader.h"
+#include "AliESDEvent.h"
+#endif
+
+Bool_t MakeITSdEdxAnalysis(TString path=".",Int_t nmombins=16){
+
+  if (gClassTable->GetID("AliRun") < 0) {
+    gInterpreter->ExecuteMacro("loadlibs.C");
+  }
+  TString galicefile = path+"/galice.root";
+  TString esdfile = path+"/AliESDs.root";
+
+  AliRunLoader* runLoader = AliRunLoader::Open(galicefile.Data());
+  if (!runLoader) {
+    printf("Error in getting run loader");
+    return kFALSE;
+  }
+  runLoader->LoadHeader();
+  runLoader->LoadKinematics();
+
+  TFile* esdFile = TFile::Open(esdfile.Data());
+  if (!esdFile || !esdFile->IsOpen()) {
+    printf("Error in opening ESD file");
+    return kFALSE;
+  }
+
+  AliESDEvent * esd = new AliESDEvent;
+  TTree* tree = (TTree*) esdFile->Get("esdTree");
+  if (!tree) {
+    printf("Error: no ESD tree found");
+    return kFALSE;
+  }
+  esd->ReadFromTree(tree);
+  AliITSdEdxAnalyzer* enan=new AliITSdEdxAnalyzer(nmombins,0.1,3.1);
+  enan->SetMIPdEdx(79.);
+
+  for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
+    runLoader->GetEvent(iEvent);
+    tree->GetEvent(iEvent);
+    if (!esd) {
+      printf("Error: no ESD object found for event %d", iEvent);
+      return kFALSE;
+    }
+    AliStack* stack = runLoader->Stack();
+    enan->ReadEvent(esd,stack);
+  }
+  enan->WriteHistos();
+
+
+
+  TCanvas* c1=new TCanvas("c1","dEdx - pions - layer 4");
+  c1->Divide(4,4);
+  for(Int_t i=0; i<nmombins;i++){
+    TH1F* h=enan->GetSingleLayerdEdxHisto(4,211,i);
+    c1->cd(i+1);
+    h->Draw();
+  }
+
+
+  TCanvas* c1b=new TCanvas("c1b","dEdx - pions - Truncated mean");
+  c1b->Divide(4,4);
+  for(Int_t i=0; i<nmombins;i++){
+    TH1F* h=enan->GetTruncatedMeandEdxHisto(211,i);
+    c1b->cd(i+1);
+    h->Draw();
+  }
+
+  TCanvas* c2=new TCanvas("c2","Delta dEdx - pions - layer 4");
+  c2->Divide(4,4);
+  for(Int_t i=0; i<nmombins;i++){
+    TH1F* h=enan->GetSingleLayerDeltadEdxHisto(4,211,i);
+    c2->cd(i+1);
+    h->Draw();
+  }
+
+  TCanvas* c2b=new TCanvas("c2b","Delta dEdx - pions - Truncated Mean");
+  c2b->Divide(4,4);
+  for(Int_t i=0; i<nmombins;i++){
+    TH1F* h=enan->GetTruncatedMeanDeltadEdxHisto(211,i);
+    c2b->cd(i+1);
+    h->Draw();
+  }
+
+  gStyle->SetPalette(1);
+
+  TCanvas* c3=0;
+  c3=new TCanvas("c3","dEdx vs. p - pions - layer 4");
+  TH2F* h2=enan->GetSingleLayerdEdxVsPHisto(4,211);
+  h2->Draw("colz");
+  enan->SetUseBBFromAliExternalTrackParam();
+  TGraph* gpion1=enan->GetBetheBlochGraph(211);
+  gpion1->Draw("LSAME");
+  enan->SetUseBBFromAliITSpidESD();
+  TGraph* gpion2=enan->GetBetheBlochGraph(211);
+  gpion2->SetLineColor(2);
+  gpion2->Draw("LSAME");
+
+  TCanvas* c3b=0;
+  c3b=new TCanvas("c3b","dEdx vs. p - pions - Truncated Mean");
+  TH2F* h2b=enan->GetTruncatedMeandEdxVsPHisto(211);
+  h2b->Draw("colz");
+  gpion1->Draw("LSAME");
+  gpion2->Draw("LSAME");
+
+  return kTRUE;
+}
index abc6f09..fd513ec 100644 (file)
@@ -37,6 +37,7 @@ SRCS =        AliITSDetTypeRec.cxx \
                AliITSpidESD.cxx \
                AliITSpidESD1.cxx \
                AliITSpidESD2.cxx \
+               AliITSdEdxAnalyzer.cxx \
                 AliITSPident.cxx \
                 AliITSSteerPid.cxx \
                 AliITSPidParItem.cxx \