1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 ////////////////////////////////////////////////////////////////////////////////
18 // This class contains all code which is used to compute any of the values
19 // which can be of interest within a resonance analysis. Besides the obvious
20 // invariant mass, it allows to compute other utility values on all possible
21 // targets, in order to allow a wide spectrum of binning and checks.
22 // When needed, this object can also define a binning in the variable which
23 // it is required to compute, which is used for initializing axes of output
24 // histograms (see AliRsnFunction).
25 // The value computation requires this object to be passed the object whose
26 // informations will be used. This object can be of any allowed input type
27 // (track, Daughter, event), then this class must inherit from AliRsnTarget.
28 // Then, when value computation is attempted, a check on target type is done
29 // and computation is successful only if expected target matches that of the
31 // In some cases, the value computation can require a support external object,
32 // which must then be passed to this class. It can be of any type inheriting
35 // authors: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
36 // M. Vala (martin.vala@cern.ch)
38 ////////////////////////////////////////////////////////////////////////////////
40 #include <Riostream.h>
41 #include "AliVVertex.h"
42 #include "AliMultiplicity.h"
43 #include "AliESDtrackCuts.h"
44 #include "AliESDpid.h"
45 #include "AliAODPid.h"
46 #include "AliCentrality.h"
47 #include "AliESDUtils.h"
48 #include "AliPIDResponse.h"
50 #include "AliRsnEvent.h"
51 #include "AliRsnDaughter.h"
52 #include "AliRsnMother.h"
53 #include "AliRsnDaughterDef.h"
54 #include "AliRsnDaughterDef.h"
56 #include "AliRsnValueDaughter.h"
58 ClassImp(AliRsnValueDaughter)
60 //_____________________________________________________________________________
61 AliRsnValueDaughter::AliRsnValueDaughter(const char *name, EType type) :
62 AliRsnValue(name, AliRsnTarget::kDaughter),
70 //_____________________________________________________________________________
71 AliRsnValueDaughter::AliRsnValueDaughter(const AliRsnValueDaughter ©) :
80 //_____________________________________________________________________________
81 AliRsnValueDaughter &AliRsnValueDaughter::operator=(const AliRsnValueDaughter ©)
84 // Assignment operator.
85 // Works like copy constructor.
88 AliRsnValue::operator=(copy);
96 //_____________________________________________________________________________
97 const char *AliRsnValueDaughter::GetTypeName() const
100 // This method returns a string to give a name to each possible
101 // computation value.
105 case kP: return "SingleTrackPtot";
106 case kPt: return "SingleTrackPt";
107 case kPtpc: return "SingleTrackPtpc";
108 case kEta: return "SingleTrackEta";
109 case kITSsignal: return "SingleTrackITSsignal";
110 case kTPCsignal: return "SingleTrackTPCsignal";
111 case kTOFsignal: return "SingleTrackTOFsignal";
112 case kTPCnsigmaPi: return "SingleTrackTPCnsigmaPion";
113 case kTPCnsigmaK: return "SingleTrackTPCnsigmaKaon";
114 case kTPCnsigmaP: return "SingleTrackTPCnsigmaProton";
115 case kTOFnsigmaPi: return "SingleTrackTOFnsigmaPion";
116 case kTOFnsigmaK: return "SingleTrackTOFnsigmaKaon";
117 case kTOFnsigmaP: return "SingleTrackTOFnsigmaProton";
118 default: return "Undefined";
122 //_____________________________________________________________________________
123 Bool_t AliRsnValueDaughter::Eval(TObject *object)
126 // Evaluation of the required value.
127 // Checks that the passed object is of the right type
128 // and if this check is successful, computes the required value.
129 // The output of the function tells if computing was successful,
130 // and the values must be taken with GetValue().
133 // coherence check, which also casts object
134 // to AliRsnTarget data members and returns kFALSE
135 // in case the object is NULL
136 if (!TargetOK(object)) return kFALSE;
138 // set a reference to the mother momentum
139 AliVParticle *ref = fDaughter->GetRef();
140 AliVParticle *refMC = fDaughter->GetRefMC();
141 AliVTrack *track = fDaughter->Ref2Vtrack();
142 if (fUseMCInfo && !refMC) {
146 if (!fUseMCInfo && !ref) {
151 // compute value depending on types in the enumeration
152 // if the type does not match any available choice, or if
153 // the computation is not doable due to any problem
154 // (not initialized support object, wrong values, risk of floating point errors)
155 // the method returng kFALSE and sets the computed value to a meaningless number
158 fComputedValue = (fUseMCInfo ? refMC->P() : ref->P());
161 fComputedValue = (fUseMCInfo ? refMC->Pt() : ref->Pt());
164 fComputedValue = (fUseMCInfo ? refMC->Eta() : ref->Eta());
168 fComputedValue = track->GetTPCmomentum();
171 AliWarning("Cannot get TPC momentum for non-track object");
172 fComputedValue = 0.0;
177 fComputedValue = track->GetITSsignal();
180 AliWarning("Cannot get ITS signal for non-track object");
181 fComputedValue = 0.0;
186 fComputedValue = track->GetTPCsignal();
189 AliWarning("Cannot get TPC signal for non-track object");
190 fComputedValue = 0.0;
195 fComputedValue = track->GetTOFsignal();
198 AliWarning("Cannot get TOF signal for non-track object");
199 fComputedValue = 0.0;
204 AliPIDResponse *pid = fEvent->GetPIDResponse();
205 fComputedValue = pid->NumberOfSigmasTPC(track, AliPID::kPion);
208 AliWarning("Cannot get TOF signal for non-track object");
209 fComputedValue = 0.0;
214 AliPIDResponse *pid = fEvent->GetPIDResponse();
215 fComputedValue = pid->NumberOfSigmasTPC(track, AliPID::kKaon);
218 AliWarning("Cannot get TOF signal for non-track object");
219 fComputedValue = 0.0;
224 AliPIDResponse *pid = fEvent->GetPIDResponse();
225 fComputedValue = pid->NumberOfSigmasTPC(track, AliPID::kProton);
228 AliWarning("Cannot get TOF signal for non-track object");
229 fComputedValue = 0.0;
234 AliPIDResponse *pid = fEvent->GetPIDResponse();
235 fComputedValue = pid->NumberOfSigmasTOF(track, AliPID::kPion);
238 AliWarning("Cannot get TOF signal for non-track object");
239 fComputedValue = 0.0;
244 AliPIDResponse *pid = fEvent->GetPIDResponse();
245 fComputedValue = pid->NumberOfSigmasTOF(track, AliPID::kKaon);
248 AliWarning("Cannot get TOF signal for non-track object");
249 fComputedValue = 0.0;
254 AliPIDResponse *pid = fEvent->GetPIDResponse();
255 fComputedValue = pid->NumberOfSigmasTOF(track, AliPID::kProton);
258 AliWarning("Cannot get TOF signal for non-track object");
259 fComputedValue = 0.0;
263 AliError(Form("[%s] Invalid value type for this computation", GetName()));