PWGHF/hfe/AliHFEitsPIDqa.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 // Class AliHFEitsPIDqa
  17 // Monitoring ITS PID in the HFE PID montioring framework. The following
  18 // quantities are monitored:
  19 //   ITS dE/dx (Number of sigmas)
  20 // (Always as function of momentum, particle species and centrality
  21 // before and after cut)
  22 // More information about the PID monitoring framework can be found in
  23 // AliHFEpidQAmanager.cxx and AliHFEdetPIDqa.cxx
  24 //
  25 // Author:
  26 //    Yvonne Pachmayer <pachmay@physi.uni-heidelberg.de>
  27 //
  28 #include <TBrowser.h>
  29 #include <TClass.h>
  30 #include <TH2.h>
  31 #include <THnSparse.h>
  32 #include <TString.h>
  33
  34 #include "AliAODTrack.h"
  35 #include "AliESDtrack.h"
  36 #include "AliLog.h"
  37 #include "AliPID.h"
  38 #include "AliPIDResponse.h"
  39
  40 #include "AliHFEcollection.h"
  41 #include "AliHFEpidBase.h"
  42 #include "AliHFEpidQAmanager.h"
  43 #include "AliHFEpidITS.h"
  44 #include "AliHFEtools.h"
  45 #include "AliHFEitsPIDqa.h"
  46
  47 ClassImp(AliHFEitsPIDqa)
  48
  49 //_________________________________________________________
  50 AliHFEitsPIDqa::AliHFEitsPIDqa():
  51     AliHFEdetPIDqa()
  52   , fHistos(NULL)
  53   , fBrowseCentrality(-1)
  54 {
  55   //
  56   // Dummy constructor
  57   //
  58 }
  59
  60 //_________________________________________________________
  61 AliHFEitsPIDqa::AliHFEitsPIDqa(const char* name):
  62     AliHFEdetPIDqa(name, "QA for ITS")
  63   , fHistos(NULL)
  64   , fBrowseCentrality(-1)
  65 {
  66   //
  67   // Default constructor
  68   //
  69 }
  70
  71 //_________________________________________________________
  72 AliHFEitsPIDqa::AliHFEitsPIDqa(const AliHFEitsPIDqa &o):
  73     AliHFEdetPIDqa(o)
  74   , fHistos(NULL)
  75   , fBrowseCentrality(o.fBrowseCentrality)
  76 {
  77   //
  78   // Copy constructor
  79   //
  80   o.Copy(*this);
  81 }
  82
  83 //_________________________________________________________
  84 AliHFEitsPIDqa &AliHFEitsPIDqa::operator=(const AliHFEitsPIDqa &o){
  85   //
  86   // Do assignment
  87   //
  88   AliHFEdetPIDqa::operator=(o);
  89   if(&o != this) o.Copy(*this);
  90   return *this;
  91 }
  92
  93 //_________________________________________________________
  94 AliHFEitsPIDqa::~AliHFEitsPIDqa(){
  95   //
  96   // Destructor
  97   //
  98   if(fHistos) delete fHistos;
  99 }
 100
 101 //_________________________________________________________
 102 void AliHFEitsPIDqa::Copy(TObject &o) const {
 103   //
 104   // Make copy
 105   //
 106   AliHFEitsPIDqa &target = dynamic_cast<AliHFEitsPIDqa &>(o);
 107   if(target.fHistos){
 108     delete target.fHistos;
 109     target.fHistos = NULL;
 110   }
 111   if(fHistos) target.fHistos = new AliHFEcollection(*fHistos);
 112   target.fBrowseCentrality = fBrowseCentrality;
 113 }
 114
 115 //_________________________________________________________
 116 Long64_t AliHFEitsPIDqa::Merge(TCollection *coll){
 117   //
 118   // Merge with other objects
 119   //
 120   if(!coll) return 0;
 121   if(coll->IsEmpty()) return 1;
 122
 123   TIter it(coll);
 124   AliHFEitsPIDqa *refQA = NULL;
 125   TObject *o = NULL;
 126   Long64_t count = 0;
 127   TList listHistos;
 128   while((o = it())){
 129     refQA = dynamic_cast<AliHFEitsPIDqa *>(o);
 130     if(!refQA) continue;
 131
 132     listHistos.Add(refQA->fHistos);
 133     count++;
 134   }
 135   fHistos->Merge(&listHistos);
 136   return count + 1;
 137 }
 138
 139 //_________________________________________________________
 140 void AliHFEitsPIDqa::Browse(TBrowser *b){
 141   //
 142   // Browse the PID QA
 143   //
 144   if(b){
 145     if(fHistos){
 146       b->Add(fHistos, fHistos->GetName());
 147
 148       // Make Projections of the nsigma Spectra and add them to a new Folder
 149       TString specnames[AliPID::kSPECIES+1] = {"All", "Electrons", "Muon", "Pions", "Kaon", "Protons"};
 150       Int_t specind[AliPID::kSPECIES+1] = {-1, AliPID::kElectron, AliPID::kMuon, AliPID::kPion, AliPID::kKaon, AliPID::kProton};
 151       TList *listNsigma = new TList;
 152       listNsigma->SetOwner();
 153
 154       TH2 *hptr = NULL;
 155       for(Int_t ispec = 0; ispec < AliPID::kSPECIES+1; ispec++){
 156         for(Int_t istep = 0; istep < 2; istep++){
 157           hptr = MakeSpectrumNSigma(static_cast<AliHFEdetPIDqa::EStep_t>(istep), specind[ispec], fBrowseCentrality);
 158           hptr->SetName(Form("hITSnsigma%s%s%s", specnames[ispec].Data(), istep == 0 ? "Before" : "After", fBrowseCentrality == -1 ? "MinBias" : Form("Cent%d", fBrowseCentrality)));
 159           listNsigma->Add(hptr);
 160         }
 161       }
 162
 163       b->Add(listNsigma, "Projections NSigma");
 164     }
 165   }
 166 }
 167
 168 //_________________________________________________________
 169 void AliHFEitsPIDqa::Initialize(){
 170   //
 171   // Define Histograms
 172   //
 173
 174   fHistos = new AliHFEcollection("itsqahistos", "Collection of ITS QA histograms");
 175
 176   // Make common binning
 177   const Int_t kNdim = 5;
 178   const Int_t kPIDbins = AliPID::kSPECIES + 1;
 179   const Int_t kSteps = 2;
 180   const Int_t kCentralityBins = 11;
 181   const Double_t kMinPID = -1;
 182   const Double_t kMinP = 0.;
 183   const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
 184   const Double_t kMaxP = 2.;
 185
 186   // Quantities where one can switch between low and high resolution
 187   Int_t kPbins = fQAmanager->HasHighResolutionHistos() ? 100 : 10;
 188   Int_t kSigmaBins = fQAmanager->HasHighResolutionHistos() ? 1400 : 240;
 189
 190   // 1st histogram: ITS sigmas: (species, p, nsigma, step, centrality)
 191   Int_t nBinsSigma[kNdim] = {kPIDbins, kPbins, kSigmaBins, kSteps, kCentralityBins};
 192   Double_t minSigma[kNdim] = {kMinPID, kMinP, -12., 0., 0.};
 193   Double_t maxSigma[kNdim] = {kMaxPID, kMaxP, 12., 2., 11.};
 194   fHistos->CreateTHnSparse("itsnSigma", "ITS signal; species; p [GeV/c]; ITS signal [a.u.]; Selection Step; Centrality; Eta", kNdim, nBinsSigma, minSigma, maxSigma);
 195
 196   // General ITS QA
 197 }
 198
 199 //_________________________________________________________
 200 void AliHFEitsPIDqa::ProcessTrack(const AliHFEpidObject *track, AliHFEdetPIDqa::EStep_t step){
 201   //
 202   // Fill ITS histograms
 203   //
 204   AliDebug(1, Form("QA started for ITS PID for step %d", (Int_t)step));
 205   //AliHFEpidObject::AnalysisType_t anatype = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
 206   Float_t centrality = track->GetCentrality();
 207   Int_t species = track->GetAbInitioPID();
 208   if(species >= AliPID::kSPECIES) species = -1;
 209
 210   AliHFEpidITS *itspid = dynamic_cast<AliHFEpidITS *>(fQAmanager->GetDetectorPID(AliHFEpid::kITSpid));
 211   Double_t contentSignal[5];
 212   contentSignal[0] = species;
 213   contentSignal[1] = track->GetRecTrack()->P();
 214   contentSignal[2] = itspid ? itspid->GetITSNsigmaCorrected(track->GetRecTrack()) : -100.;
 215   contentSignal[3] = step;
 216   contentSignal[4] = centrality;
 217
 218   fHistos->Fill("itsnSigma", contentSignal);
 219 }
 220
 221 //_________________________________________________________
 222 TH2 *AliHFEitsPIDqa::MakeSpectrumNSigma(AliHFEdetPIDqa::EStep_t istep, Int_t species, Int_t centralityClass){
 223   //
 224   // Plot the Spectrum
 225   //
 226   THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("itsnSigma"));
 227   if(!hSignal) return NULL;
 228   hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
 229   if(species >= 0 && species < AliPID::kSPECIES)
 230     hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
 231   TString hname = Form("hITSsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
 232           htitle = Form("ITS dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
 233   if(species > -1){
 234     hname += AliPID::ParticleName(species);
 235     htitle += Form(" for %ss", AliPID::ParticleName(species));
 236   }
 237   TString centname, centtitle;
 238   Bool_t hasCentralityInfo = kTRUE;
 239   if(centralityClass > -1){
 240     if(hSignal->GetNdimensions() < 5){
 241       AliError("Centrality Information not available");
 242       centname = centtitle = "MinBias";
 243       hasCentralityInfo= kFALSE;
 244     } else {
 245       // Project centrality classes
 246       // -1 is Min. Bias
 247       hSignal->GetAxis(4)->SetRange(centralityClass+1, centralityClass+1);
 248       centname = Form("Cent%d", centralityClass);
 249       centtitle = Form(" Centrality %d", centralityClass);
 250     }
 251   } else {
 252     centname = centtitle = "MinBias";
 253     hasCentralityInfo= kFALSE;
 254   }
 255   hname += centtitle;
 256   htitle += centtitle;
 257
 258   TH2 *hTmp = hSignal->Projection(2,1);
 259   hTmp->SetName(hname.Data());
 260   hTmp->SetTitle(htitle.Data());
 261   hTmp->SetStats(kFALSE);
 262   hTmp->GetXaxis()->SetTitle("p [GeV/c]");
 263   hTmp->GetYaxis()->SetTitle("ITS dE/dx - <dE/dx>|_{el} [#sigma]");
 264   hSignal->GetAxis(3)->SetRange(0, hSignal->GetAxis(3)->GetNbins());
 265   hSignal->GetAxis(0)->SetRange(0, hSignal->GetAxis(0)->GetNbins());
 266   if(hasCentralityInfo) hSignal->GetAxis(4)->SetRange(0, hSignal->GetAxis(4)->GetNbins());
 267   return hTmp;
 268 }
 269