laura is hosted by Hepforge, IPPP Durham
Laura++  v3r0p1
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauSigmaRes.cc
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2004 - 2014.
3 // Distributed under the Boost Software License, Version 1.0.
4 // (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5 
6 // Authors:
7 // Thomas Latham
8 // John Back
9 // Paul Harrison
10 
15 #include <iostream>
16 
17 #include "LauConstants.hh"
18 #include "LauSigmaRes.hh"
19 #include "LauResonanceInfo.hh"
20 
22 
23 
24 LauSigmaRes::LauSigmaRes(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) :
25  LauAbsResonance(resInfo, resPairAmpInt, daughters),
26  mPiSq4_(4.0*LauConstants::mPiSq),
27  sAdler_(LauConstants::mPiSq*0.5),
28  b1_(0),
29  b2_(0),
30  a_(0),
31  m0_(0)
32 {
33  // Initialise various constants
34  mPiSq4_ = 4.0*LauConstants::mPiSq;
35  sAdler_ = LauConstants::mPiSq*0.5; // Adler zero at 0.5*(mpi)^2
36 
37  // constant factors from BES data
38  const Double_t b1Val = 0.5843;
39  const Double_t b2Val = 1.6663;
40  const Double_t aVal = 1.082;
41  const Double_t m0Val = 0.9264;
42 
43  const TString& parNameBase = this->getSanitisedName();
44 
45  TString b1Name(parNameBase);
46  b1Name += "_b1";
47  b1_ = resInfo->getExtraParameter( b1Name );
48  if ( b1_ == 0 ) {
49  b1_ = new LauParameter( b1Name, b1Val, 0.0, 100.0, kTRUE );
50  b1_->secondStage(kTRUE);
51  resInfo->addExtraParameter( b1_ );
52  }
53 
54  TString b2Name(parNameBase);
55  b2Name += "_b2";
56  b2_ = resInfo->getExtraParameter( b2Name );
57  if ( b2_ == 0 ) {
58  b2_ = new LauParameter( b2Name, b2Val, 0.0, 100.0, kTRUE );
59  b2_->secondStage(kTRUE);
60  resInfo->addExtraParameter( b2_ );
61  }
62 
63  TString aName(parNameBase);
64  aName += "_A";
65  a_ = resInfo->getExtraParameter( aName );
66  if ( a_ == 0 ) {
67  a_ = new LauParameter( aName, aVal, 0.0, 10.0, kTRUE );
68  a_->secondStage(kTRUE);
69  resInfo->addExtraParameter( a_ );
70  }
71 
72  TString m0Name(parNameBase);
73  m0Name += "_m0";
74  m0_ = resInfo->getExtraParameter( m0Name );
75  if ( m0_ == 0 ) {
76  m0_ = new LauParameter( m0Name, m0Val, 0.0, 10.0, kTRUE );
77  m0_->secondStage(kTRUE);
78  resInfo->addExtraParameter( m0_ );
79  }
80 }
81 
83 {
84 }
85 
87 {
88  this->checkDaughterTypes();
89 
90  Double_t resSpin = this->getSpin();
91  if (resSpin != 0) {
92  std::cerr << "WARNING in LauSigmaRes::initialise : Resonance spin is " << resSpin << "." << std::endl;
93  std::cerr << " : Sigma amplitude is only for scalers, resetting spin to 0." << std::endl;
94  this->changeResonance( -1.0, -1.0, 0 );
95  }
96 }
97 
99 {
100  // Check that the daughter tracks are the same type. Otherwise issue a warning
101  // and set the type to be pion for the Sigma distribution.
102  Int_t resPairAmpInt = this->getPairInt();
103  if (resPairAmpInt < 1 || resPairAmpInt > 3) {
104  std::cerr << "WARNING in LauSigmaRes::checkDaughterTypes : resPairAmpInt = " << resPairAmpInt << " is out of the range [1,2,3]." << std::endl;
105  return;
106  }
107 
108  const TString& nameDaug1 = this->getNameDaug1();
109  const TString& nameDaug2 = this->getNameDaug2();
110  if (!nameDaug1.CompareTo(nameDaug2, TString::kExact)) {
111  // Daughter types agree. Find out if we have pion or kaon system
112  if (!nameDaug1.Contains("pi")) {
113  std::cerr << "ERROR in LauSigmaRes::checkDaughterTypes : Sigma model is using daughters \""<<nameDaug1<<"\" and \""<<nameDaug2<<"\", which are not pions." << std::endl;
114  }
115  }
116 }
117 
118 LauComplex LauSigmaRes::resAmp(Double_t mass, Double_t spinTerm)
119 {
120  // This function returns the complex dynamical amplitude for a Sigma distribution
121  // given the invariant mass and cos(helicity) values.
122 
123  // First check that the appropriate daughters are either pi+pi- or K+K-
124  // Check that the daughter tracks are the same type. Otherwise issue a warning
125  // and set the type to be pion for the Sigma distribution. Returns the
126  // integer defined by the enum LauSigmaRes::SigmaPartType.
127 
128  Double_t s = mass*mass; // Invariant mass squared combination for the system
129  Double_t rho(0.0); // Phase-space factor
130  if (s > mPiSq4_) {rho = TMath::Sqrt(1.0 - mPiSq4_/s);}
131 
132  const Double_t m0Val = this->getM0Value();
133  const Double_t m0Sq = m0Val * m0Val;
134 
135  const Double_t aVal = this->getAValue();
136  const Double_t b1Val = this->getB1Value();
137  const Double_t b2Val = this->getB2Value();
138 
139  Double_t f = b2Val*s + b1Val; // f(s) function
140  Double_t numerator = s - sAdler_;
141  Double_t denom = m0Sq - sAdler_;
142  Double_t gamma(0.0);
143  if (TMath::Abs(denom) > 1e-10 && TMath::Abs(aVal) > 1e-10) {
144  // Decay width of the system
145  gamma = rho*(numerator/denom)*f*TMath::Exp(-(s - m0Sq)/aVal);
146  }
147 
148  // Now form the complex amplitude - use relativistic BW form (without barrier factors)
149  // Note that the M factor in the denominator is not the "pole" at ~500 MeV, but is
150  // m0_ = 0.9264, the mass when the phase shift goes through 90 degrees.
151 
152  Double_t dMSq = m0Sq - s;
153  Double_t widthTerm = gamma*m0Val;
154  LauComplex resAmplitude(dMSq, widthTerm);
155 
156  Double_t denomFactor = dMSq*dMSq + widthTerm*widthTerm;
157 
158  Double_t invDenomFactor = 0.0;
159  if (denomFactor > 1e-10) {invDenomFactor = 1.0/denomFactor;}
160 
161  resAmplitude.rescale(spinTerm*invDenomFactor);
162 
163  return resAmplitude;
164 
165 }
166 
167 const std::vector<LauParameter*>& LauSigmaRes::getFloatingParameters()
168 {
169  this->clearFloatingParameters();
170 
171  if ( ! this->fixB1Value() ) {
172  this->addFloatingParameter( b1_ );
173  }
174 
175  if ( ! this->fixB2Value() ) {
176  this->addFloatingParameter( b2_ );
177  }
178 
179  if ( ! this->fixAValue() ) {
180  this->addFloatingParameter( a_ );
181  }
182 
183  if ( ! this->fixM0Value() ) {
184  this->addFloatingParameter( m0_ );
185  }
186 
187  return this->getParameters();
188 }
189 
190 void LauSigmaRes::setResonanceParameter(const TString& name, const Double_t value)
191 {
192  // Set various parameters for the lineshape
193  if (name == "b1") {
194  this->setB1Value(value);
195  std::cout << "INFO in LauSigmaRes::setResonanceParameter : Setting parameter b1 = " << this->getB1Value() << std::endl;
196  }
197  else if (name == "b2") {
198  this->setB2Value(value);
199  std::cout << "INFO in LauSigmaRes::setResonanceParameter : Setting parameter b2 = " << this->getB2Value() << std::endl;
200  }
201  else if (name == "A") {
202  this->setAValue(value);
203  std::cout << "INFO in LauSigmaRes::setResonanceParameter : Setting parameter A = " << this->getAValue() << std::endl;
204  }
205  else if (name == "m0") {
206  this->setM0Value(value);
207  std::cout << "INFO in LauSigmaRes::setResonanceParameter : Setting parameter m0 = " << this->getM0Value() << std::endl;
208  }
209  else {
210  std::cerr << "WARNING in LauSigmaRes::setResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
211  }
212 }
213 
215 {
216  if (name == "b1") {
217  if ( b1_->fixed() ) {
218  b1_->fixed( kFALSE );
219  this->addFloatingParameter( b1_ );
220  } else {
221  std::cerr << "WARNING in LauSigmaRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
222  }
223  } else if (name == "b2") {
224  if ( b2_->fixed() ) {
225  b2_->fixed( kFALSE );
226  this->addFloatingParameter( b2_ );
227  } else {
228  std::cerr << "WARNING in LauSigmaRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
229  }
230  } else if (name == "A") {
231  if ( a_->fixed() ) {
232  a_->fixed( kFALSE );
233  this->addFloatingParameter( a_ );
234  } else {
235  std::cerr << "WARNING in LauSigmaRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
236  }
237  } else if (name == "m0") {
238  if ( m0_->fixed() ) {
239  m0_->fixed( kFALSE );
240  this->addFloatingParameter( m0_ );
241  } else {
242  std::cerr << "WARNING in LauSigmaRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
243  }
244  } else {
245  std::cerr << "WARNING in LauSigmaRes::fixResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
246  }
247 }
248 
250 {
251  if (name == "b1") {
252  return b1_;
253  } else if (name == "b2") {
254  return b2_;
255  } else if (name == "A") {
256  return a_;
257  } else if (name == "m0") {
258  return m0_;
259  } else {
260  std::cerr << "WARNING in LauSigmaRes::getResonanceParameter: Parameter name not reconised." << std::endl;
261  return 0;
262  }
263 }
264 
265 void LauSigmaRes::setB1Value(const Double_t b1)
266 {
267  b1_->value( b1 );
268  b1_->genValue( b1 );
269  b1_->initValue( b1 );
270 }
271 
272 void LauSigmaRes::setB2Value(const Double_t b2)
273 {
274  b2_->value( b2 );
275  b2_->genValue( b2 );
276  b2_->initValue( b2 );
277 }
278 
279 void LauSigmaRes::setAValue(const Double_t A)
280 {
281  a_->value( A );
282  a_->genValue( A );
283  a_->initValue( A );
284 }
285 
286 void LauSigmaRes::setM0Value(const Double_t m0)
287 {
288  m0_->value( m0 );
289  m0_->genValue( m0 );
290  m0_->initValue( m0 );
291 }
292 
virtual ~LauSigmaRes()
Destructor.
Definition: LauSigmaRes.cc:82
Bool_t fixed() const
Check whether the parameter is fixed or floated.
void setAValue(const Double_t A)
Set the A parameter.
Definition: LauSigmaRes.cc:279
Bool_t fixB2Value() const
Fix the b2 parameter value.
Definition: LauSigmaRes.hh:139
void changeResonance(const Double_t newMass, const Double_t newWidth, const Int_t newSpin)
Allow the mass, width and spin of the resonance to be changed.
File containing declaration of LauResonanceInfo class.
ClassImp(LauAbsCoeffSet)
Bool_t fixB1Value() const
Fix the b1 parameter value.
Definition: LauSigmaRes.hh:133
LauParameter()
Default constructor.
Definition: LauParameter.cc:30
Class for defining the properties of a resonant particle.
const TString & name() const
The parameter name.
Class that defines the particular 3-body decay under study.
Definition: LauDaughters.hh:33
Double_t getAValue() const
Get the A parameter value.
Definition: LauSigmaRes.hh:121
Double_t getB2Value() const
Get the b2 parameter value.
Definition: LauSigmaRes.hh:115
TString getNameDaug2() const
Get the name of the second daughter of the resonance.
Bool_t fixM0Value() const
Fix the m0 parameter value.
Definition: LauSigmaRes.hh:151
virtual void initialise()
Initialise the model.
Definition: LauSigmaRes.cc:86
virtual void floatResonanceParameter(const TString &name)
Allow the various parameters to float in the fit.
Definition: LauSigmaRes.cc:214
Int_t getPairInt() const
Get the integer to identify which DP axis the resonance belongs to.
Double_t sAdler_
Defined as 0.5*mPi*mPi.
Definition: LauSigmaRes.hh:173
Class for defining the Sigma resonance model.
Definition: LauSigmaRes.hh:31
void addFloatingParameter(LauParameter *param)
Add parameter to the list of floating parameters.
Double_t getB1Value() const
Get the b1 parameter value.
Definition: LauSigmaRes.hh:109
File containing declaration of LauSigmaRes class.
virtual void setResonanceParameter(const TString &name, const Double_t value)
Set value of the various parameters.
Definition: LauSigmaRes.cc:190
virtual LauParameter * getResonanceParameter(const TString &name)
Access the given resonance parameter.
Definition: LauSigmaRes.cc:249
LauParameter * b2_
Factor from BES data.
Definition: LauSigmaRes.hh:178
const Double_t mPiSq
Square of pi+- mass.
Definition: LauConstants.hh:65
std::vector< LauParameter * > & getParameters()
Access the list of floating parameters.
LauParameter * a_
Factor from BES data.
Definition: LauSigmaRes.hh:180
Class for defining the fit parameter objects.
Definition: LauParameter.hh:34
void checkDaughterTypes() const
Check that both daughters are the same type of particle.
Definition: LauSigmaRes.cc:98
void setB2Value(const Double_t b2)
Set the b2 parameter.
Definition: LauSigmaRes.cc:272
void setB1Value(const Double_t b1)
Set the b1 parameter.
Definition: LauSigmaRes.cc:265
virtual LauComplex resAmp(Double_t mass, Double_t spinTerm)
Complex resonant ampltiude.
Definition: LauSigmaRes.cc:118
Abstract class for defining type for resonance amplitude models (Breit-Wigner, Flatte etc...
void rescale(Double_t scaleVal)
Scale this by a factor.
Definition: LauComplex.hh:285
LauParameter * m0_
Factor from BES data.
Definition: LauSigmaRes.hh:182
Double_t initValue() const
The initial value of the parameter.
File containing LauConstants namespace.
Class for defining a complex number.
Definition: LauComplex.hh:47
TString getNameDaug1() const
Get the name of the first daughter of the resonance.
Double_t value() const
The value of the parameter.
Double_t getM0Value() const
Get the m0 parameter value.
Definition: LauSigmaRes.hh:127
virtual const std::vector< LauParameter * > & getFloatingParameters()
Retrieve the resonance parameters, e.g. so that they can be loaded into a fit.
Definition: LauSigmaRes.cc:167
void setM0Value(const Double_t m0)
Set the m0 parameter.
Definition: LauSigmaRes.cc:286
Bool_t fixAValue() const
Fix the A parameter value.
Definition: LauSigmaRes.hh:145
Double_t mPiSq4_
Defined as 4*mPi*mPi.
Definition: LauSigmaRes.hh:171
Int_t getSpin() const
Get the spin of the resonance.
Double_t genValue() const
The value generated for the parameter.
LauParameter * b1_
Factor from BES data.
Definition: LauSigmaRes.hh:176
void clearFloatingParameters()
Clear list of floating parameters.