Ad-hoc implementation of cuts for pions/kaons/protons for 2010 analysis
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnValuePair.cxx
CommitLineData
2895972e 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////////////////////////////////////////////////////////////////////////////////
17//
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, pair, 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
30// passed object.
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
33// from TObject.
34//
35// authors: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
36// M. Vala (martin.vala@cern.ch)
37//
38////////////////////////////////////////////////////////////////////////////////
39
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
49#include "AliRsnEvent.h"
50#include "AliRsnDaughter.h"
51#include "AliRsnMother.h"
52#include "AliRsnPairDef.h"
53#include "AliRsnDaughterDef.h"
54
55#include "AliRsnValuePair.h"
56
57ClassImp(AliRsnValuePair)
58
59//_____________________________________________________________________________
60AliRsnValuePair::AliRsnValuePair(const char *name, EType type) :
61 AliRsnValue(name, AliRsnTarget::kMother),
62 fType(type)
63{
64//
65// Constructor
66//
67}
68
69//_____________________________________________________________________________
70AliRsnValuePair::AliRsnValuePair(const AliRsnValuePair& copy) :
71 AliRsnValue(copy),
72 fType(copy.fType)
73{
74//
75// Copy constructor
76//
77}
78
79//_____________________________________________________________________________
80AliRsnValuePair& AliRsnValuePair::operator=(const AliRsnValuePair& copy)
81{
82//
83// Assignment operator.
84// Works like copy constructor.
85//
86
87 AliRsnValue::operator=(copy);
88 fType = copy.fType;
89
90 return (*this);
91}
92
93//_____________________________________________________________________________
94const char* AliRsnValuePair::GetTypeName() const
95{
96//
97// This method returns a string to give a name to each possible
98// computation value.
99//
100
101 switch (fType) {
102 case kPt: return "PairPt";
103 case kPz: return "PairPz";
104 case kInvMass: return "PairInvMass";
105 case kInvMassRes: return "PairInvMassResolution";
106 case kEta: return "PairEta";
107 case kMt: return "PairMt";
108 case kY: return "PairY";
109 case kPtRatio: return "PairPtRatio";
110 case kDipAngle: return "PairDipAngle";
111 case kCosThetaStar: return "PairCosThetaStar";
112 default: return "Undefined";
113 }
114}
115
116//_____________________________________________________________________________
117Bool_t AliRsnValuePair::Eval(TObject *object)
118{
119//
120// Evaluation of the required value.
121// In this implementation, fills the member 4-vectors with data
122// coming from the object passed as argument, and then returns the value
123//
124
125 // coherence check, which also casts object
126 // to AliRsnTarget data members and returns kFALSE
127 // in case the object is NULL
128 if (!TargetOK(object)) return kFALSE;
129
130 // set a reference to the mother momentum
131 TLorentzVector &sum = fMother->Sum(fUseMCInfo);
132 TLorentzVector &ref = fMother->Ref(fUseMCInfo);
133 TLorentzVector &p1 = fMother->GetDaughter(0)->P(fUseMCInfo);
134 TLorentzVector &p2 = fMother->GetDaughter(1)->P(fUseMCInfo);
135
136 // compute value depending on types in the enumeration
137 // if the type does not match any available choice, or if
138 // the computation is not doable due to any problem
139 // (not initialized support object, wrong values, risk of floating point errors)
140 // the method returng kFALSE and sets the computed value to a meaningless number
141 switch (fType) {
142 case kPt:
143 fComputedValue = sum.Perp();
144 return kTRUE;
145 case kInvMass:
146 fComputedValue = sum.M();
147 return kTRUE;
148 case kEta:
149 fComputedValue = sum.Eta();
150 return kTRUE;
151 case kInvMassRes:
152 fComputedValue = fMother->Sum(kFALSE).M() - fMother->Sum(kTRUE).M();
153 fComputedValue /= fMother->Sum(kTRUE).M();
154 return kTRUE;
155 case kMt:
156 fComputedValue = ref.Mt();
157 return kTRUE;
158 case kY:
159 fComputedValue = ref.Rapidity();
160 return kTRUE;
161 case kPtRatio:
162 fComputedValue = TMath::Abs(p1.Perp() - p2.Perp());
163 fComputedValue /= TMath::Abs(p1.Perp() + p2.Perp());
164 return kTRUE;
165 case kDipAngle:
166 fComputedValue = p1.Perp() * p2.Perp() + p1.Z() * p2.Z();
167 fComputedValue /= p1.Mag() * p2.Mag();
168 return kTRUE;
169 case kCosThetaStar:
170 fComputedValue = fMother->CosThetaStar();
171 return kTRUE;
172 default:
173 AliError(Form("[%s] Invalid value type for this computation", GetName()));
174 return kFALSE;
175 }
176}