laura is hosted by Hepforge, IPPP Durham
Laura++  v3r5
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauPolarFormFactorSymNR.cc
Go to the documentation of this file.
1 
2 /*
3 Copyright 2018 University of Warwick
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 */
17 
18 /*
19 Laura++ package authors:
20 John Back
21 Paul Harrison
22 Thomas Latham
23 */
24 
29 #include <iostream>
30 
31 #include "TMath.h"
32 
33 #include "LauConstants.hh"
35 #include "LauDaughters.hh"
36 #include "LauParameter.hh"
37 #include "LauResonanceInfo.hh"
38 
40 
41 
42 LauPolarFormFactorSymNR::LauPolarFormFactorSymNR(LauResonanceInfo* resInfo, const LauAbsResonance::LauResonanceModel resType, const Int_t resPairAmpInt, const LauDaughters* daughters) :
43  LauAbsResonance(resInfo, resPairAmpInt, daughters),
44  lambda_(0),
45  model_(resType)
46 
47 {
48  TString parName = this->getSanitisedName();
49  parName += "_lambda";
50  lambda_ = resInfo->getExtraParameter( parName );
51  if ( lambda_ == 0 ) {
52  lambda_ = new LauParameter( parName, 1.0, 0.0, 10.0, kTRUE );
53  lambda_->secondStage(kTRUE);
54  resInfo->addExtraParameter( lambda_ );
55  }
56 }
57 
59 {
60 }
61 
63 {
64  const LauDaughters* daughters = this->getDaughters();
65  if ( ! daughters->gotSymmetricalDP() ) {
66  std::cerr << "WARNING in LauPolarFormFactorSymNR::initialise : Dalitz plot is symmetric - this lineshape is not appropriate." << std::endl;
67  }
68 
69  Int_t resPairAmpInt = this->getPairInt();
70  if ( resPairAmpInt == 3 ) {
71  std::cerr << "WARNING in LauPolarFormFactorSymNR::initialise : This lineshape is intended to be on the symmetrised axes of the DP." << std::endl;
72  }
73 
75  std::cerr << "WARNING in LauPolarFormFactorSymNR::initialise : Unknown model requested, defaulting to Polar Form Factor." << std::endl;
77  }
78 
79  if ( (model_ != LauAbsResonance::PolarFFSymNR) && (this->getSpin() != 0) ) {
80  std::cerr << "WARNING in LauPolarFormFactorSymNR::initialise : Non-zero spin will be ignored for this model - perhaps you should use LauAbsResonance::PolarFFSymNRNoInter instead" << std::endl;
81  }
82 
83  // NB we do not need to call setSpinType(LauAbsResonance::Legendre) here (as is done in LauPolarFormFactorNR) since override the amplitude method and explicitly use calcLegendrePoly
84 }
85 
86 
87 LauComplex LauPolarFormFactorSymNR::resAmp(Double_t mass, Double_t spinTerm)
88 {
89  std::cerr << "ERROR in LauPolarFormFactorSymNR : This method should never be called." << std::endl;
90  std::cerr << " : Returning zero amplitude for mass = " << mass << " and spinTerm = " << spinTerm << "." << std::endl;
91  return LauComplex(0.0, 0.0);
92 }
93 
95 {
96  // This function returns the complex dynamical amplitude for a Polar Form Factor Non-Resonant distribution
97 
98  // Calculate for symmetric DPs, e.g. 3pi or 3K, by using shapeNo = 1 or 2
99  // Have s<->t symmetry already done in Dynamics flip function.
100  // For Kpipi or similar plots, one can use the separate terms
101  // and consider them as two separate components with their own mag and phase.
102  // For this shapeNo = 3 and shapeNo = 4 need to be used to create the two
103  // individual amplitudes (with the same value of lambda).
104 
105  // Calculate Mandelstam variables.
106  // s = m_13^2, t = m_23^2
107  const Double_t s = kinematics->getm13Sq();
108  const Double_t t = kinematics->getm23Sq();
109 
110  Double_t magnitude(1.0);
111 
112  const Double_t lambda = this->getLambda();
113 
115 
116  magnitude = 1.0/(1.0 + s /(lambda*lambda)) + 1.0/(1.0 + t /(lambda*lambda));
117 
119 
120  magnitude = (s <= t) ? 1.0/(1.0 + s /(lambda*lambda)) : 1.0/(1.0 + t /(lambda*lambda));
121  }
122 
123  LauComplex resAmplitude(magnitude, 0.0);
124 
125  return resAmplitude;
126 
127 }
128 
129 const std::vector<LauParameter*>& LauPolarFormFactorSymNR::getFloatingParameters()
130 {
131  this->clearFloatingParameters();
132 
133  if ( ! this->fixLambda() ) {
134  this->addFloatingParameter( lambda_ );
135  }
136 
137  return this->getParameters();
138 }
139 
140 void LauPolarFormFactorSymNR::setResonanceParameter(const TString& name, const Double_t value)
141 {
142  // Set various parameters for the lineshape
143  if (name == "lambda") {
144  this->setLambda(value);
145  std::cout << "INFO in LauPolarFormFactorSymNR::setResonanceParameter : Setting parameter lambda = " << this->getLambda() << std::endl;
146  }
147  else {
148  std::cerr << "WARNING in LauPolarFormFactorSymNR::setResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
149  }
150 }
151 
153 {
154  if (name == "lambda") {
155  if ( lambda_->fixed() ) {
156  lambda_->fixed( kFALSE );
157  this->addFloatingParameter( lambda_ );
158  } else {
159  std::cerr << "WARNING in LauPolarFormFactorSymNR::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl;
160  }
161  }
162  else {
163  std::cerr << "WARNING in LauPolarFormFactorSymNR::fixResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl;
164  }
165 }
166 
168 {
169  if (name == "lambda") {
170  return lambda_;
171  }
172  else {
173  std::cerr << "WARNING in LauPolarFormFactorSymNR::getResonanceParameter: Parameter name not reconised." << std::endl;
174  return 0;
175  }
176 }
177 
178 void LauPolarFormFactorSymNR::setLambda(const Double_t lambda)
179 {
180  lambda_->value( lambda );
181  lambda_->genValue( lambda );
182  lambda_->initValue( lambda );
183 }
184 
185 
Bool_t fixed() const
Check whether the parameter is fixed or floated.
Bool_t fixLambda() const
See if the lambda parameter is fixed or floating.
virtual ~LauPolarFormFactorSymNR()
Destructor.
File containing declaration of LauResonanceInfo class.
Double_t getLambda() const
Get the parameter lambda, the NR shape parameter.
ClassImp(LauAbsCoeffSet)
LauParameter()
Default constructor.
Definition: LauParameter.cc:44
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:47
virtual LauParameter * getResonanceParameter(const TString &name)
Access the given resonance parameter.
File containing declaration of LauDaughters class.
void setLambda(const Double_t lambda)
Set the parameter lambda, the NR shape parameter.
Double_t getm23Sq() const
Get the m23 invariant mass square.
const LauDaughters * getDaughters() const
Access the daughters object.
Int_t getPairInt() const
Get the integer to identify which DP axis the resonance belongs to.
virtual const std::vector< LauParameter * > & getFloatingParameters()
Retrieve the resonance parameters, e.g. so that they can be loaded into a fit.
void addFloatingParameter(LauParameter *param)
Add parameter to the list of floating parameters.
LauAbsResonance::LauResonanceModel model_
The model to use.
Bool_t gotSymmetricalDP() const
Is Dalitz plot symmetric, i.e. 2 identical particles.
Definition: LauDaughters.hh:80
File containing declaration of LauPolarFormFactorSymNR class.
File containing declaration of LauParameter class.
std::vector< LauParameter * > & getParameters()
Access the list of floating parameters.
virtual LauComplex amplitude(const LauKinematics *kinematics)
Get the complex dynamical amplitude.
Class for defining the fit parameter objects.
Definition: LauParameter.hh:49
LauResonanceModel
Define the allowed resonance types.
Abstract class for defining type for resonance amplitude models (Breit-Wigner, Flatte etc...
Double_t getm13Sq() const
Get the m13 invariant mass square.
virtual LauComplex resAmp(Double_t mass, Double_t spinTerm)
Complex resonant amplitude.
Double_t initValue() const
The initial value of the parameter.
File containing LauConstants namespace.
Class for defining a complex number.
Definition: LauComplex.hh:61
Class for calculating 3-body kinematic quantities.
Double_t value() const
The value of the parameter.
virtual void setResonanceParameter(const TString &name, const Double_t value)
Set value of the various parameters.
virtual void initialise()
Initialise the model.
Int_t getSpin() const
Get the spin of the resonance.
Double_t genValue() const
The value generated for the parameter.
Class for defining a nonresonant form factor model.
virtual void floatResonanceParameter(const TString &name)
Allow the various parameters to float in the fit.
void clearFloatingParameters()
Clear list of floating parameters.