laura is hosted by Hepforge, IPPP Durham
Laura++  v3r3
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauLASSRes.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 "LauLASSRes.hh"
19 #include "LauResonanceInfo.hh"
20 
22 
23 
24 LauLASSRes::LauLASSRes(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) :
25  LauAbsResonance(resInfo, resPairAmpInt, daughters),
26  q0_(0.0),
27  mDaugSum_(0.0),
28  mDaugSumSq_(0.0),
29  mDaugDiff_(0.0),
30  mDaugDiffSq_(0.0),
31  resMassSq_(0.0),
32  r_(0),
33  a_(0),
34  cutOff_(0.0)
35 {
36  // Default values for LASS parameters
37  cutOff_ = 1.8;
38  const Double_t rVal = 3.32;
39  const Double_t aVal = 2.07;
40 
41  const TString& parNameBase = this->getSanitisedName();
42 
43  TString rName(parNameBase);
44  rName += "_r";
45  r_ = resInfo->getExtraParameter( rName );
46  if ( r_ == 0 ) {
47  r_ = new LauParameter( rName, rVal, 0.0, 10.0, kTRUE );
48  r_->secondStage(kTRUE);
49  resInfo->addExtraParameter( r_ );
50  }
51 
52  TString aName(parNameBase);
53  aName += "_a";
54  a_ = resInfo->getExtraParameter( aName );
55  if ( a_ == 0 ) {
56  a_ = new LauParameter( aName, aVal, 0.0, 10.0, kTRUE );
57  a_->secondStage(kTRUE);
58  resInfo->addExtraParameter( a_ );
59  }
60 }
61 
63 {
64 }
65 
67 {
68  // Create the mass sums and differences
69  Double_t massDaug1 = this->getMassDaug1();
70  Double_t massDaug2 = this->getMassDaug2();
71 
72  mDaugSum_ = massDaug1 + massDaug2;
74 
75  mDaugDiff_ = massDaug1 - massDaug2;
77 
78  Int_t resSpin = this->getSpin();
79  if (resSpin != 0) {
80  std::cerr << "WARNING in LauLASSRes::amplitude : Resonance spin is " << resSpin << "." << std::endl;
81  std::cerr << " : LASS amplitude is only for scalers, resetting spin to 0." << std::endl;
82  this->changeResonance( -1.0, -1.0, 0 );
83  }
84 
85  this->calcQ0();
86 }
87 
89 {
90  // Decay momentum of either daughter in the resonance rest frame
91  // when resonance mass = rest-mass value, m_0 (PDG value)
92 
93  resMass_ = this->getMass();
95 
96  q0_ = TMath::Sqrt((resMassSq_ - mDaugSumSq_)*(resMassSq_ - mDaugDiffSq_))/(2.0*resMass_);
97 }
98 
99 LauComplex LauLASSRes::resAmp(Double_t mass, Double_t /*spinTerm*/)
100 {
101  LauComplex resAmplitude(0.0, 0.0);
102  LauComplex bkgAmplitude(0.0, 0.0);
103  LauComplex totAmplitude(0.0, 0.0);
104 
105  if (mass < 1e-10) {
106  std::cerr << "WARNING in LauLASSRes::amplitude : mass < 1e-10." << std::endl;
107  return LauComplex(0.0, 0.0);
108  }
109 
110  //---------------------------
111  // First do the resonant part
112  //---------------------------
113 
114  // Calculate the width of the resonance (as a function of mass)
115  // q is the momentum of either daughter in the resonance rest-frame
116  const Double_t q = this->getQ();
117  const Double_t resMass = this->getMass();
118  const Double_t resWidth = this->getWidth();
119 
120  // If the mass is floating and their value have changed
121  // we need to recalculate everything that assumes this value
122  if ( (!this->fixMass()) && resMass != resMass_ ) {
123  this->calcQ0();
124  }
125 
126  Double_t qRatio = q/q0_;
127  Double_t totWidth = resWidth*qRatio*(resMass/mass);
128 
129  Double_t massSq = mass*mass;
130  Double_t massSqTerm = resMassSq_ - massSq;
131 
132  // Compute the complex amplitude
133  resAmplitude = LauComplex(massSqTerm, resMass*totWidth);
134 
135  // Scale by the denominator factor
136  resAmplitude.rescale((resMassSq_*resWidth/q0_)/(massSqTerm*massSqTerm + resMassSq_*totWidth*totWidth));
137 
138  // Calculate the phase shift term
139  const Double_t rVal = this->getEffectiveRange();
140  const Double_t aVal = this->getScatteringLength();
141  const Double_t tandeltaB = (2.0*aVal*q)/(2.0 + aVal*rVal*q*q);
142  const Double_t tanSq = tandeltaB*tandeltaB;
143  const Double_t cos2PhaseShift = (1.0 - tanSq) / (1.0 + tanSq);
144  const Double_t sin2PhaseShift = 2.0*tandeltaB / (1.0 + tanSq);
145  LauComplex phaseShift(cos2PhaseShift, sin2PhaseShift);
146 
147  // Multiply by the phase shift
148  resAmplitude = resAmplitude * phaseShift;
149 
150 
151  //--------------------------------
152  // Now do the effective range part
153  //--------------------------------
154 
155  // Form the real and imaginary parts
156  const Double_t qcotdeltaB = 1.0/aVal + (rVal*q*q)/2.0;
157 
158  // Compute the complex amplitude
159  bkgAmplitude = LauComplex(qcotdeltaB, q);
160 
161  // Scale by the numerator and denominator factors
162  bkgAmplitude.rescale(mass/(qcotdeltaB*qcotdeltaB + q*q));
163 
164 
165  //------------------
166  // Add them together
167  //------------------
168 
169  if (mass > cutOff_) {
170  totAmplitude = resAmplitude;
171  } else {
172  totAmplitude = bkgAmplitude + resAmplitude;
173  }
174 
175  return totAmplitude;
176 }
177 
178 const std::vector<LauParameter*>& LauLASSRes::getFloatingParameters()
179 {
180  this->clearFloatingParameters();
181 
182  if ( ! this->fixMass() ) {
183  this->addFloatingParameter( this->getMassPar() );
184  }
185 
186  if ( ! this->fixWidth() ) {
187  this->addFloatingParameter( this->getWidthPar() );
188  }
189 
190  if ( ! this->fixEffectiveRange() ) {
191  this->addFloatingParameter( r_ );
192  }
193 
194  if ( ! this->fixScatteringLength() ) {
195  this->addFloatingParameter( a_ );
196  }
197 
198  return this->getParameters();
199 }
200 
201 void LauLASSRes::setResonanceParameter(const TString& name, const Double_t value)
202 {
203  // Set various parameters for the LASS lineshape dynamics
204  if (name == "a") {
205  this->setScatteringLength(value);
206  std::cout << "INFO in LauLASSRes::setResonanceParameter : Setting LASS Scattering Length = " << this->getScatteringLength() << std::endl;
207  } else if (name == "r") {
208  this->setEffectiveRange(value);
209  std::cout << "INFO in LauLASSRes::setResonanceParameter : Setting LASS Effective Range = " << this->getEffectiveRange() << std::endl;
210  } else {
211  std::cerr << "WARNING in LauLASSRes::setResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
212  }
213 }
214 
216 {
217  if (name == "a") {
218  if ( a_->fixed() ) {
219  a_->fixed( kFALSE );
220  this->addFloatingParameter( a_ );
221  } else {
222  std::cerr << "WARNING in LauLASSRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
223  }
224  } else if (name == "r") {
225  if ( r_->fixed() ) {
226  r_->fixed( kFALSE );
227  this->addFloatingParameter( r_ );
228  } else {
229  std::cerr << "WARNING in LauLASSRes::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
230  }
231  } else {
232  std::cerr << "WARNING in LauLASSRes::fixResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
233  }
234 }
235 
237 {
238  if (name == "a") {
239  return a_;
240  } else if (name == "r") {
241  return r_;
242  } else {
243  std::cerr << "WARNING in LauLASSRes::getResonanceParameter: Parameter name not reconised." << std::endl;
244  return 0;
245  }
246 }
247 
248 void LauLASSRes::setEffectiveRange(const Double_t r)
249 {
250  r_->value( r );
251  r_->genValue( r );
252  r_->initValue( r );
253 }
254 
255 void LauLASSRes::setScatteringLength(const Double_t a)
256 {
257  a_->value( a );
258  a_->genValue( a );
259  a_->initValue( a );
260 }
261 
Double_t getQ() const
Get the current value of the daughter momentum in the resonance rest frame.
LauParameter * getMassPar()
Get the mass parameter of the resonance.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
Class for defining the LASS resonance model.
Definition: LauLASSRes.hh:31
virtual void floatResonanceParameter(const TString &name)
Allow the various parameters to float in the fit.
Definition: LauLASSRes.cc:215
Bool_t fixWidth() const
Get the status of resonance width (fixed or released)
Double_t getMass() const
Get the mass of the resonance.
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)
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 getScatteringLength() const
Get the scattering length range parameter.
Definition: LauLASSRes.hh:115
Double_t mDaugDiffSq_
Square of mDaugDiff.
Definition: LauLASSRes.hh:155
Bool_t fixScatteringLength() const
See if the scattering length parameter is fixed or floating.
Definition: LauLASSRes.hh:127
Double_t getMassDaug1() const
Get the mass of daughter 1.
Bool_t fixMass() const
Get the status of resonance mass (fixed or released)
Double_t getMassDaug2() const
Get the mass of daughter 2.
Double_t q0_
Decay momentum of either daughter in the resonance rest frame.
Definition: LauLASSRes.hh:147
Double_t resMassSq_
Square of the resonance mass.
Definition: LauLASSRes.hh:159
virtual void initialise()
Initialise the model.
Definition: LauLASSRes.cc:66
void addFloatingParameter(LauParameter *param)
Add parameter to the list of floating parameters.
LauParameter * a_
LASS scattering length parameter.
Definition: LauLASSRes.hh:164
Bool_t fixEffectiveRange() const
See if the effective range parameter is fixed or floating.
Definition: LauLASSRes.hh:121
Double_t mDaugDiff_
Difference between the daughter masses.
Definition: LauLASSRes.hh:153
Double_t getEffectiveRange() const
Get the effective range parameter.
Definition: LauLASSRes.hh:109
std::vector< LauParameter * > & getParameters()
Access the list of floating parameters.
LauParameter * getWidthPar()
Get the width parameter of the resonance.
virtual LauComplex resAmp(Double_t mass, Double_t spinTerm)
Complex resonant amplitude.
Definition: LauLASSRes.cc:99
Class for defining the fit parameter objects.
Definition: LauParameter.hh:35
virtual void setResonanceParameter(const TString &name, const Double_t value)
Set value of a resonance parameter.
Definition: LauLASSRes.cc:201
File containing declaration of LauLASSRes class.
Double_t mDaugSum_
Sum of the daughter masses.
Definition: LauLASSRes.hh:149
Double_t getWidth() const
Get the width of the resonance.
void setScatteringLength(const Double_t a)
Set the scattering length parameter value.
Definition: LauLASSRes.cc:255
LauParameter * r_
LASS effective range parameter.
Definition: LauLASSRes.hh:162
Abstract class for defining type for resonance amplitude models (Breit-Wigner, Flatte etc...
Double_t cutOff_
LASS cut off.
Definition: LauLASSRes.hh:167
void rescale(Double_t scaleVal)
Scale this by a factor.
Definition: LauComplex.hh:285
Double_t resMass_
The resonance mass.
Definition: LauLASSRes.hh:157
Double_t initValue() const
The initial value of the parameter.
File containing LauConstants namespace.
Class for defining a complex number.
Definition: LauComplex.hh:47
void setEffectiveRange(const Double_t r)
Set the effective range parameter value.
Definition: LauLASSRes.cc:248
Double_t value() const
The value of the parameter.
Double_t mDaugSumSq_
Square of mDaugSum.
Definition: LauLASSRes.hh:151
Int_t getSpin() const
Get the spin of the resonance.
virtual const std::vector< LauParameter * > & getFloatingParameters()
Retrieve the resonance parameters, e.g. so that they can be loaded into a fit.
Definition: LauLASSRes.cc:178
Double_t genValue() const
The value generated for the parameter.
void calcQ0()
Utility function to calculate the q0 value.
Definition: LauLASSRes.cc:88
virtual ~LauLASSRes()
Destructor.
Definition: LauLASSRes.cc:62
virtual LauParameter * getResonanceParameter(const TString &name)
Access the given resonance parameter.
Definition: LauLASSRes.cc:236
void clearFloatingParameters()
Clear list of floating parameters.