fix for bug #70582 (change from L. Molnar)
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnCutBetheBloch.cxx
CommitLineData
5eb970a4 1//
2// Class AliRsnCutBetheBloch
3//
4// General implementation of a single cut strategy, which can be:
5// - a value contained in a given interval [--> IsBetween() ]
6// - a value equal to a given reference [--> MatchesValue()]
7//
8// In all cases, the reference value(s) is (are) given as data members
9// and each kind of cut requires a given value type (Int, UInt, Double),
10// but the cut check procedure is then automatized and chosen thanks to
11// an enumeration of the implemented cut types.
12// At the end, the user (or any other point which uses this object) has
13// to use the method IsSelected() to check if this cut has been passed.
14//
15// authors: Martin Vala (martin.vala@cern.ch)
16// Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
17//
18#include "TMath.h"
19
a59975b3 20#include "AliExternalTrackParam.h"
21
5eb970a4 22#include "AliRsnDaughter.h"
23#include "AliRsnCutBetheBloch.h"
24
25ClassImp(AliRsnCutBetheBloch)
26
27//_________________________________________________________________________________________________
28AliRsnCutBetheBloch::AliRsnCutBetheBloch() :
4fbb2459 29 AliRsnCut(),
30 fCorrect(kTRUE),
31 fMIP(50.0),
32 fType(AliPID::kUnknown)
5eb970a4 33{
34//
35// Default constructor.
36//
37
38 fConst[0] = fConst[1] = fConst[2] = fConst[3] = fConst[4] = 0.0;
39}
40
41//_________________________________________________________________________________________________
42AliRsnCutBetheBloch::AliRsnCutBetheBloch
5eb970a4 43(const char *name, Double_t fractionRange, AliPID::EParticleType type, Double_t mip, Bool_t correct) :
4fbb2459 44 AliRsnCut(name, 0.0, fractionRange),
45 fCorrect(correct),
46 fMIP(mip),
47 fType(type)
5eb970a4 48{
49//
50// Main constructor.
51// the cut range is the relative fraction of the value:
52// BB*(1-fraction) < TPC < BB*(1+fraction)
53// which means:
54// -fraction < (TPC - BB)/BB < fraction
55//
56
5eb970a4 57 fConst[0] = fConst[1] = fConst[2] = fConst[3] = fConst[4] = 0.0;
58}
59
60//_____________________________________________________________________________
a59975b3 61Double_t AliRsnCutBetheBloch::BetheBloch(AliRsnDaughter * const trackRef)
5eb970a4 62{
63//
64// Computes the theoretical dE/dx according to
65// a given mass hypothesis, from which betaGamma is computed
66//
67// This is the empirical ALEPH parameterization of the Bethe-Bloch formula.
68// It is normalized to 1 at the minimum.
69//
70// The default values for the kp* parameters are for ALICE TPC.
71// The value is computed in MIP units, multiplied by 50 to have it in energy.
72//
73
cff22d54 74 AliPID pid;
75 Double_t mass = pid.ParticleMass(fType);
76
867852c5 77 // get the track momentum at the inner wall of TPC: if absent cut is not passed
a59975b3 78 AliExternalTrackParam track(*trackRef->GetRefESD()->GetInnerParam());
79
80 Double_t betaGamma = track.P() / mass;
5eb970a4 81 Double_t beta = betaGamma / TMath::Sqrt(1.0 + betaGamma * betaGamma);
82 Double_t aa = TMath::Power(beta, fConst[3]);
83 Double_t bb = TMath::Power(1.0 / betaGamma, fConst[4]);
84
85 bb = TMath::Log(fConst[2] + bb);
86
87 Double_t out = (fConst[1] - aa - bb) * fConst[0] / aa;
88
4fbb2459 89 if (fCorrect) {
5eb970a4 90 Double_t kMeanCorr = 0.1;
91 Double_t meanCorr = (1 + (out - 1) * kMeanCorr);
92 out *= meanCorr;
93 }
94
95 return out * fMIP;
96}
97
6f4a992c 98//_____________________________________________________________________________
99Double_t AliRsnCutBetheBloch::RelDiff(AliRsnDaughter *track)
100{
101//
102// Relative difference between BB value and TPC signal
103//
104
105 if (!track->GetRefESD()) return -99999.9;
106
107 // compute Bethe-Bloch with the given mass hypothesis
108 Double_t bb = BetheBloch(track);
109 return TMath::Abs((track->GetRefESD()->GetTPCsignal() - bb) / bb);
110}
111
5eb970a4 112//_________________________________________________________________________________________________
113Bool_t AliRsnCutBetheBloch::IsSelected(ETarget tgt, AliRsnDaughter *track)
114{
115//
116// Cut checker.
117//
118
119 // coherence check
4fbb2459 120 if (tgt != AliRsnCut::kParticle) {
1d8bfdee 121 AliError(Form("[%s] Wrong target. Skipping cut", GetName()));
5eb970a4 122 return kTRUE;
123 }
124
cff22d54 125 // if the required PID of the track is not the same as the
126 // reference of the cut, the cut is automatically skipped
127 if (track->RequiredPID() != fType) return kTRUE;
128
5eb970a4 129 // retrieve the TPC signal
4fbb2459 130 AliESDtrack *esd = track->GetRefESD();
5eb970a4 131 if (!esd) {
132 AliError("ESD information unavailable");
133 return kTRUE;
134 }
ceaa78d3 135 if (!track->GetRefESD()->GetInnerParam()) {
136 AliDebug(AliLog::kDebug+2, "Rejecting a track with no info at the TPC inner wall");
137 return kFALSE;
138 }
5eb970a4 139
5eb970a4 140 // the cut range is the relative fraction of the value:
141 // BB*(1-fraction) < TPC < BB*(1+fraction)
142 // which means:
143 // -fraction < (TPC - BB)/BB < fraction
144 // so we must compute the cut value accordingly
6f4a992c 145 fCutValueD = RelDiff(track);
5eb970a4 146
147 // then, this cut is checked inside the range
148 return OkRange();
149}
150
151//_________________________________________________________________________________________________
e79f56bd 152Bool_t AliRsnCutBetheBloch::IsSelected(ETarget /*tgt*/, AliRsnPairParticle* /*pair*/)
5eb970a4 153{
154//
155// Cut checker
156//
157
158 AliWarning("Cannot apply this cut to pairs");
159 return kTRUE;
160}
161
162//_________________________________________________________________________________________________
e79f56bd 163Bool_t AliRsnCutBetheBloch::IsSelected(ETarget /*tgt*/, AliRsnEvent* /*event*/)
5eb970a4 164{
165//
166// Cut checker
167//
168
169 AliWarning("Cannot apply this cut to events");
170 return kTRUE;
171}
172
173//_________________________________________________________________________________________________
e79f56bd 174Bool_t AliRsnCutBetheBloch::IsSelected(ETarget /*tgt*/, AliRsnEvent* /*ev1*/, AliRsnEvent* /*ev2*/)
5eb970a4 175{
176//
177// Cut checker
178//
179
180 AliWarning("Cannot apply this cut to event mixing");
181 return kTRUE;
182}