STEER/AliAODpidUtil.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 /* $Id: AliAODpidUtil.cxx 38329 2010-01-17 19:17:24Z hristov $ */
  17
  18 //-----------------------------------------------------------------
  19 //           Implementation of the combined PID class
  20 //           For the AOD Class
  21 //           containing information on the particle identification
  22 //      Origin: Rosa Romita, GSI, r.romita@gsi.de
  23 //-----------------------------------------------------------------
  24
  25 #include "AliLog.h"
  26 #include "AliPID.h"
  27 #include "AliAODpidUtil.h"
  28 #include "AliAODEvent.h"
  29 #include "AliAODTrack.h"
  30 #include "AliAODPid.h"
  31 #include "AliTRDPIDResponse.h"
  32
  33 ClassImp(AliAODpidUtil)
  34
  35   Int_t AliAODpidUtil::MakePID(AliAODTrack *track,Float_t TimeZeroTOF,Double_t *p) const {
  36   //
  37   //  Calculate probabilities for all detectors, except if TPConly==kTRUE
  38   //  and combine PID
  39   //
  40   //   Option TPConly==kTRUE is used during reconstruction,
  41   //  because ITS tracking uses TPC pid
  42   //  HMPID and TRD pid are done in detector reconstructors
  43   //
  44
  45   /*
  46     Float_t TimeZeroTOF = 0;
  47     if (subtractT0)
  48     TimeZeroTOF = event->GetT0();
  49   */
  50   Int_t ns=AliPID::kSPECIES;
  51   Double_t tpcPid[AliPID::kSPECIES];
  52   MakeTPCPID(track,tpcPid);
  53   Double_t itsPid[AliPID::kSPECIES];
  54   Double_t tofPid[AliPID::kSPECIES];
  55   Double_t trdPid[AliPID::kSPECIES];
  56   MakeITSPID(track,itsPid);
  57   MakeTOFPID(track, TimeZeroTOF,tofPid);
  58   //MakeHMPIDPID(track);
  59   MakeTRDPID(track,trdPid);
  60   for (Int_t j=0; j<ns; j++) {
  61     p[j]=tpcPid[j]*itsPid[j]*tofPid[j]*trdPid[j];
  62   }
  63
  64   return 0;
  65 }
  66 //_________________________________________________________________________
  67 void AliAODpidUtil::MakeTPCPID(AliAODTrack *track,Double_t *p) const
  68 {
  69   //
  70   //  TPC pid using bethe-bloch and gaussian response
  71   //
  72
  73   Double_t mom = track->P();
  74   AliAODPid *pidObj = track->GetDetPid();
  75   if (pidObj) mom = pidObj->GetTPCmomentum();
  76
  77   Double_t dedx=pidObj->GetTPCsignal();
  78   Bool_t mismatch=kTRUE;
  79
  80   for (Int_t j=0; j<AliPID::kSPECIES; j++) {
  81     AliPID::EParticleType type=AliPID::EParticleType(j);
  82     Double_t bethe=fTPCResponse.GetExpectedSignal(mom,type);
  83     Double_t sigma=fTPCResponse.GetExpectedSigma(mom,50,type);
  84     if (TMath::Abs(dedx-bethe) > fRange*sigma) {
  85       p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
  86     } else {
  87       p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
  88       mismatch=kFALSE;
  89     }
  90
  91   }
  92
  93   if (mismatch)
  94     for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
  95
  96
  97   return;
  98 }
  99 //_________________________________________________________________________
 100 void AliAODpidUtil::MakeITSPID(AliAODTrack *track,Double_t *p) const
 101 {
 102   //
 103   // ITS PID
 104   //  1) Truncated mean method
 105   //
 106
 107
 108   Double_t mom=track->P();
 109   AliAODPid *pidObj = track->GetDetPid();
 110
 111   Double_t dedx=pidObj->GetITSsignal();
 112   Bool_t mismatch=kTRUE;
 113   for (Int_t j=0; j<AliPID::kSPECIES; j++) {
 114     Double_t mass=AliPID::ParticleMass(j);//GeV/c^2
 115     Double_t bethe=fITSResponse.Bethe(mom,mass);
 116     Double_t sigma=fITSResponse.GetResolution(bethe);
 117     if (TMath::Abs(dedx-bethe) > fRange*sigma) {
 118       p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
 119     } else {
 120       p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
 121       mismatch=kFALSE;
 122     }
 123
 124     // Check for particles heavier than (AliPID::kSPECIES - 1)
 125
 126   }
 127
 128   if (mismatch)
 129     for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1./AliPID::kSPECIES;
 130
 131   return;
 132
 133 }
 134 //_________________________________________________________________________
 135 void AliAODpidUtil::MakeTOFPID(AliAODTrack *track, Float_t TimeZeroTOF,Double_t *p) const
 136 {
 137   //
 138   //   TOF PID using gaussian response
 139   //
 140
 141   Double_t time[AliPID::kSPECIESN];
 142   Double_t sigma[AliPID::kSPECIESN];
 143   AliAODPid *pidObj = track->GetDetPid();
 144   pidObj->GetIntegratedTimes(time);
 145
 146   for (Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++) {
 147     sigma[iPart] = fTOFResponse.GetExpectedSigma(track->P(),time[iPart],AliPID::ParticleMass(iPart));
 148   }
 149
 150   AliDebugGeneral("AliESDpid::MakeTOFPID",2,
 151                   Form("Expected TOF signals [ps]: %f %f %f %f %f",
 152                        time[AliPID::kElectron],
 153                        time[AliPID::kMuon],
 154                        time[AliPID::kPion],
 155                        time[AliPID::kKaon],
 156                        time[AliPID::kProton]));
 157
 158   AliDebugGeneral("AliESDpid::MakeTOFPID",2,
 159                   Form("Expected TOF std deviations [ps]: %f %f %f %f %f",
 160                        sigma[AliPID::kElectron],
 161                        sigma[AliPID::kMuon],
 162                        sigma[AliPID::kPion],
 163                        sigma[AliPID::kKaon],
 164                        sigma[AliPID::kProton]
 165                        ));
 166
 167   Double_t tof = pidObj->GetTOFsignal() - TimeZeroTOF;
 168
 169   Bool_t mismatch = kTRUE;
 170   for (Int_t j=0; j<AliPID::kSPECIES; j++) {
 171     Double_t sig = sigma[j];
 172     if (TMath::Abs(tof-time[j]) > fRange*sig) {
 173       p[j] = TMath::Exp(-0.5*fRange*fRange)/sig;
 174     } else
 175       p[j] = TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sig*sig))/sig;
 176
 177     // Check the mismatching
 178     Double_t mass = AliPID::ParticleMass(j);
 179     Double_t pm = fTOFResponse.GetMismatchProbability(track->P(),mass);
 180     if (p[j]>pm) mismatch = kFALSE;
 181
 182     // Check for particles heavier than (AliPID::kSPECIES - 1)
 183
 184   }
 185
 186   if (mismatch)
 187     for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
 188
 189   return;
 190 }
 191 //_________________________________________________________________________
 192 void AliAODpidUtil::MakeTRDPID(AliAODTrack *track,Double_t *p) const
 193 {
 194
 195   // Method to recalculate the TRD PID probabilities
 196   AliAODPid *pidObj = track->GetDetPid();
 197   Float_t *mom=pidObj->GetTRDmomentum();
 198   Double_t *dedx=pidObj->GetTRDsignal();
 199   Bool_t norm=kTRUE;
 200   fTRDResponse.GetResponse(pidObj->GetTRDnSlices(),dedx,mom,p,norm);
 201   return;
 202 }