PWGHF/hfe/AliHFEbayesPIDqa.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 // QA class for Bayes PID
  17 // Plot likelihoods vs various observables
  18 // combined likelihood as well as indiv. detector response
  19 // also plot mass spectrum using TOF for prior determination
  20 //
  21 // Author:
  22 //   Yvonne Pachmayer <Y.Pachmayer@gsi.de>
  23 //
  24 #include <TBrowser.h>
  25 #include <TClass.h>
  26 #include <TH2.h>
  27 #include <THnSparse.h>
  28 #include <TString.h>
  29
  30 #include "AliAODTrack.h"
  31 #include "AliESDtrack.h"
  32 #include "AliLog.h"
  33 #include "AliPID.h"
  34 #include "AliPIDResponse.h"
  35
  36 #include "AliHFEcollection.h"
  37 #include "AliHFEpidBase.h"
  38 #include "AliHFEpidQAmanager.h"
  39 #include "AliHFEpidTPC.h"
  40 #include "AliHFEpidTOF.h"
  41 #include "AliHFEpidTRD.h"
  42 #include "AliHFEpidBayes.h"
  43 #include "AliHFEtools.h"
  44 #include "AliHFEbayesPIDqa.h"
  45 #include "AliHFEdetPIDqa.h"
  46
  47 ClassImp(AliHFEbayesPIDqa)
  48
  49 //_________________________________________________________
  50 AliHFEbayesPIDqa::AliHFEbayesPIDqa():
  51     AliHFEdetPIDqa()
  52   , fHistos(NULL)
  53 {
  54   //
  55   // Dummy constructor
  56   //
  57 }
  58
  59 //_________________________________________________________
  60 AliHFEbayesPIDqa::AliHFEbayesPIDqa(const char* name):
  61     AliHFEdetPIDqa(name, "QA for Bayes")
  62   , fHistos(NULL)
  63 {
  64   //
  65   // Default constructor
  66   //
  67 }
  68
  69 //____________________________________________________________
  70 AliHFEbayesPIDqa::AliHFEbayesPIDqa(const AliHFEbayesPIDqa &o):
  71     AliHFEdetPIDqa(o),
  72     fHistos(NULL)
  73 {
  74   //
  75   // Copy constructor
  76   //
  77 }
  78
  79 //____________________________________________________________
  80 AliHFEbayesPIDqa &AliHFEbayesPIDqa::operator=(const AliHFEbayesPIDqa &o){
  81   //
  82   // Make assignment
  83   //
  84   if(this != &o){
  85     AliHFEdetPIDqa::operator=(o);
  86     fHistos = o.fHistos;
  87   }
  88   return *this;
  89 }
  90
  91
  92 //_________________________________________________________
  93 AliHFEbayesPIDqa::~AliHFEbayesPIDqa(){
  94   //
  95   // Destructor
  96   //
  97   if(fHistos) delete fHistos;
  98 }
  99 //_________________________________________________________
 100 void AliHFEbayesPIDqa::Copy(TObject &o) const {
 101   //
 102   // Make copy
 103   //
 104   AliHFEbayesPIDqa &target = dynamic_cast<AliHFEbayesPIDqa &>(o);
 105   if(target.fHistos){
 106     delete target.fHistos;
 107     target.fHistos = NULL;
 108   }
 109   if(fHistos) target.fHistos = new AliHFEcollection(*fHistos);
 110
 111 }
 112
 113 //_________________________________________________________
 114 Long64_t AliHFEbayesPIDqa::Merge(TCollection *coll){
 115   //
 116   // Merge with other objects
 117   //
 118   if(!coll) return 0;
 119   if(coll->IsEmpty()) return 1;
 120
 121   TIter it(coll);
 122   AliHFEbayesPIDqa *refQA = NULL;
 123   TObject *o = NULL;
 124   Long64_t count = 0;
 125   TList listHistos;
 126   while((o = it())){
 127     refQA = dynamic_cast<AliHFEbayesPIDqa *>(o);
 128     if(!refQA) continue;
 129
 130     listHistos.Add(refQA->fHistos);
 131     count++;
 132   }
 133   fHistos->Merge(&listHistos);
 134   return count + 1;
 135 }
 136 //_________________________________________________________
 137 void AliHFEbayesPIDqa::Initialize(){
 138   //
 139   // Define Histograms
 140   //
 141
 142   fHistos = new AliHFEcollection("BAYESqahistos", "Collection of Bayes QA histograms");
 143
 144   CreateDetectorSignalHistograms();
 145
 146 //  fBAYESpid = new AliHFEpidBayes("QAbayesPID");
 147
 148 }
 149
 150
 151 void AliHFEbayesPIDqa::CreateDetectorSignalHistograms(){
 152   //
 153   // Create Histogram for Probability Studies
 154     //
 155
 156     Int_t kPbins = 1000;  //fQAmanager->HasHighResolutionHistos() ? 1000 : 100;
 157     Int_t kSigmaBins = 300; //fQAmanager->HasHighResolutionHistos() ? 600 : 150;
 158     Int_t trdLikelihoodBins = 100; // fQAmanager->HasHighResolutionHistos() ? 200 : 100;
 159     const Int_t kPIDbins = AliPID::kSPECIES + 1;
 160     const Int_t kProbbins = 100;
 161     const Int_t kSteps = 2;
 162     const Int_t kCentralityBins = 11;
 163     const Double_t kMinPID = -1;
 164     const Double_t kMinP = 0.;
 165     const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
 166     const Double_t kMaxP = 10.;
 167
 168     Int_t nBinsContr[8] =  {kPIDbins, kPbins, kProbbins, kSigmaBins, kSigmaBins, trdLikelihoodBins, kSteps, kCentralityBins};
 169     Double_t minContr[8] = {kMinPID, kMinP, 0,  -10., -7., 0., 0., 0.};
 170     Double_t maxContr[8] = {kMaxPID, kMaxP, 1,    5.,   8., 1., 2., 11.};
 171     fHistos->CreateTHnSparse("control", "Control; species; p [GeV/c]; Comb prob; TPC sigma; TOF sigma; TRD electron Likelihood; selection step; centrality", 8, nBinsContr, minContr,maxContr);
 172
 173
 174
 175     Int_t nBinsTOFmass[5] = {kPIDbins, kPbins, 150, kSteps, kCentralityBins};
 176     Double_t minTOFmass[5] = {kMinPID, kMinP, 0., 0., 0.};
 177     Double_t maxTOFmass[5] = {kMaxPID, kMaxP, 1.5, 2., 11.};
 178     fHistos->CreateTHnSparse("tofmass", "TOF mass; species; p [GeV/c]; TOF mass; selection step", 5, nBinsTOFmass, minTOFmass, maxTOFmass);
 179
 180
 181 }
 182 //_________________________________________________________
 183 void AliHFEbayesPIDqa::ProcessTrack(const AliHFEpidObject *track, AliHFEdetPIDqa::EStep_t step){
 184   //
 185   // Fill TPC histograms
 186   //
 187   AliDebug(1, Form("QA started for BAYES PID for step %d", (Int_t)step));
 188   AliHFEpidObject::AnalysisType_t anatype = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
 189   Float_t centrality = track->GetCentrality();
 190   Int_t species = track->GetAbInitioPID();
 191   if(species >= AliPID::kSPECIES) species = -1;
 192
 193   Double_t probComb[AliPID::kSPECIES]={0.};
 194   AliHFEpidBayes *bayespid = dynamic_cast<AliHFEpidBayes *>(fQAmanager->GetDetectorPID(AliHFEpid::kBAYESpid));
 195   const AliPIDResponse *pidResponseBayes = bayespid ? bayespid->GetPIDResponse() : NULL;
 196   if(!pidResponseBayes){
 197     AliError("No PID Response available");
 198     return;
 199   }
 200   bayespid->CalcCombProb(track,pidResponseBayes, probComb);
 201
 202   Double_t contentSignal[8];
 203   contentSignal[0] = species;
 204   contentSignal[1] = track->GetRecTrack()->P();
 205   contentSignal[2] = probComb[AliPID::kElectron];
 206 //  contentSignal[2] = contentSignal[2] = pidResponse ? pidResponse->NumberOfSigmasTPC(track->GetRecTrack(), AliPID::kElectron) : 0.;
 207
 208   AliHFEpidTOF *tofpid = dynamic_cast<AliHFEpidTOF *>(fQAmanager->GetDetectorPID(AliHFEpid::kTOFpid));
 209   const AliPIDResponse *pidResponseTOF = tofpid ? tofpid->GetPIDResponse() : NULL;
 210   if(!pidResponseTOF){
 211     AliError("No PID response available");
 212     return;
 213   }
 214   contentSignal[4] = pidResponseTOF->NumberOfSigmasTOF(track->GetRecTrack(), AliPID::kElectron);
 215
 216   AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fQAmanager->GetDetectorPID(AliHFEpid::kTPCpid));
 217   const AliPIDResponse *pidResponse = tpcpid ? tpcpid->GetPIDResponse() : NULL;
 218   if(!pidResponse){
 219     AliError("No PID response available");
 220     return;
 221   }
 222   contentSignal[3] = pidResponse->NumberOfSigmasTPC(track->GetRecTrack(), AliPID::kElectron);
 223
 224   AliHFEpidTRD *trdpid = dynamic_cast<AliHFEpidTRD *>(fQAmanager->GetDetectorPID(AliHFEpid::kTRDpid));
 225   contentSignal[5] = trdpid ? trdpid->GetElectronLikelihood(static_cast<const AliVTrack*>(track->GetRecTrack()), anatype) : -10;
 226
 227   contentSignal[6] = step;
 228   contentSignal[7] = centrality;
 229   fHistos->Fill("control", contentSignal);
 230
 231   Double_t contentFill[5];
 232   contentFill[0]=contentSignal[0];
 233   contentFill[1]=contentSignal[1];
 234   contentFill[3]=contentSignal[6];
 235   contentFill[4]=contentSignal[7];
 236
 237   Double_t masscalcfromtof=CalcTOFMass(track);
 238   contentFill[2]=masscalcfromtof;
 239   fHistos->Fill("tofmass", contentFill);
 240 }
 241
 242 Double_t AliHFEbayesPIDqa::CalcTOFMass(const AliHFEpidObject *track){
 243     //
 244     // Calc TOF mass
 245     //
 246
 247    Double_t mass=-99; //GeV
 248    Double_t length=((AliESDtrack*)track->GetRecTrack())->GetIntegratedLength();
 249    if (length<=0) return mass;
 250    Double_t tofTime=((AliESDtrack*)track->GetRecTrack())->GetTOFsignal();//in ps
 251    Double_t tof= tofTime*1E-3; // ns, average T0 fill subtracted, no info from T0detector
 252    if (tof<=0) return mass;
 253    Double_t c=TMath::C()*1.E-9;// m/ns
 254    length =length*0.01; // in meters
 255    tof=tof*c;
 256    Double_t fact= (tof/length)*(tof/length) -1.;
 257    Double_t mom=track->GetRecTrack()->P();
 258    if(mom==0) return mass;
 259
 260    if(fact<=0) {
 261        mass = -mom*TMath::Sqrt(-fact);
 262    }else{
 263        mass = mom*TMath::Sqrt(fact);
 264    }
 265
 266    return mass;
 267 }
 268