ComPWA
Common Partial-Wave-Analysis Framework
Flatte.hpp
Go to the documentation of this file.
1 // Copyright (c) 2013, 2017 The ComPWA Team.
2 // This file is part of the ComPWA framework, check
3 // https://github.com/ComPWA/ComPWA/license.txt for details.
4 
5 #ifndef COMPWA_PHYSICS_DYNAMICS_FLATTE_HPP_
6 #define COMPWA_PHYSICS_DYNAMICS_FLATTE_HPP_
7 
8 #include "Coupling.hpp"
9 #include "FormFactor.hpp"
11 
12 namespace ComPWA {
13 namespace Physics {
14 namespace Dynamics {
15 
16 namespace Flatte {
19  std::shared_ptr<ComPWA::FunctionTree::FitParameter> G;
21  std::vector<Coupling> HiddenCouplings;
22 };
23 
25 inline std::complex<double>
26 flatteCouplingTerm(double sqrtS, double mR, double coupling, double massA,
27  double massB, unsigned int J, double mesonRadius,
28  std::shared_ptr<FormFactor> FormFactorFunctor) {
29  auto qR = std::abs(qValueAC(mR, massA, massB));
30  auto phspR = phspFactorAC(sqrtS, massA, massB);
31  auto ffR = FormFactorFunctor->operator()(qR *qR, J, mesonRadius);
32  auto qS = std::abs(qValueAC(sqrtS, massA, massB));
33  auto barrierA = FormFactorFunctor->operator()(qS *qS, J, mesonRadius) / ffR;
34 
35  // Calculate normalized vertex functions vtxA(s_R)
36  std::complex<double> vtxA(1, 0); // spin==0
37  if (J > 0 || FormFactorFunctor->getName() == "CrystalBarrel") {
38  vtxA = ffR * std::pow(qR, J);
39  }
40  auto width = couplingToWidth(mR, coupling, vtxA, phspR);
41  // Including the factor qTermA, as suggested by PDG 2014, Chapter 47.2,
42  // leads to an amplitude that doesn't converge.
43  // qTermA = qValue(sqrtS,massA1,massA2) / qValue(mR,massA1,massA2);
44  // termA = gammaA * barrierA * barrierA * std::pow(qTermA, (double)2 * J +
45  // 1);
46 
47  return (width * barrierA * barrierA);
48 }
49 
60 inline std::complex<double>
61 dynamicalFunction(double mSq, double mR, double gA, std::complex<double> termA,
62  std::complex<double> termB,
63  std::complex<double> termC = std::complex<double>(0, 0)) {
64  std::complex<double> i(0, 1);
65  double sqrtS = sqrt(mSq);
66 
67  std::complex<double> denom = std::complex<double>(mR * mR - mSq, 0) +
68  (-1.0) * i * sqrtS * (termA + termB + termC);
69 
70  std::complex<double> result = std::complex<double>(gA, 0) / denom;
71 
72 #ifndef NDEBUG
73  if (std::isnan(result.real()) || std::isnan(result.imag())) {
74  std::cout << "AmpFlatteRes::dynamicalFunction() | " << mR << " " << mSq
75  << " " << termA << " " << termB << " " << termC << std::endl;
76  return 0;
77  }
78 #endif
79 
80  return result;
81 }
82 
101 inline std::complex<double>
102 dynamicalFunction(double mSq, double mR, double massA1, double massA2,
103  double gA, double massB1, double massB2, double couplingB,
104  double massC1, double massC2, double couplingC,
105  unsigned int L, double mesonRadius,
106  std::shared_ptr<FormFactor> FormFactorFunctor) {
107  double sqrtS = sqrt(mSq);
108 
109  // channel A - signal channel
110  auto termA = flatteCouplingTerm(sqrtS, mR, gA, massA1, massA2, L, mesonRadius,
111  FormFactorFunctor);
112  // channel B - hidden channel
113  auto termB = flatteCouplingTerm(sqrtS, mR, couplingB, massB1, massB2, L,
114  mesonRadius, FormFactorFunctor);
115 
116  // channel C - hidden channel
117  std::complex<double> termC;
118  if (couplingC != 0.0) {
119  termC = flatteCouplingTerm(sqrtS, mR, couplingC, massC1, massC2, L,
120  mesonRadius, FormFactorFunctor);
121  }
122  return dynamicalFunction(mSq, mR, gA, termA, termB, termC);
123 }
124 
125 std::shared_ptr<ComPWA::FunctionTree::TreeNode> createFunctionTree(
126  InputInfo Params,
127  std::shared_ptr<ComPWA::FunctionTree::Value<std::vector<double>>>
128  InvMassSquared);
129 
130 } // namespace Flatte
131 
133 public:
134  FlatteStrategy(std::shared_ptr<FormFactor> FF)
135  : Strategy(ComPWA::FunctionTree::ParType::MCOMPLEX, "Flatte"),
136  FormFactorFunctor(FF) {}
137 
138  virtual void execute(ComPWA::FunctionTree::ParameterList &paras,
139  std::shared_ptr<ComPWA::FunctionTree::Parameter> &out);
140 
141 private:
142  std::shared_ptr<FormFactor> FormFactorFunctor;
143 };
144 
145 } // namespace Dynamics
146 } // namespace Physics
147 } // namespace ComPWA
148 
149 #endif
std::complex< double > phspFactorAC(double sqrtS, double ma, double mb)
Two body phsp factor.
Definition: FormFactor.hpp:42
std::shared_ptr< FormFactor > FormFactorFunctor
Form factor function object.
std::shared_ptr< FormFactor > FormFactorFunctor
Definition: Flatte.hpp:142
std::complex< double > dynamicalFunction(double mSq, double mR, double gA, std::complex< double > termA, std::complex< double > termB, std::complex< double > termC=std::complex< double >(0, 0))
Dynamical function for two coupled channel approach.
Definition: Flatte.hpp:61
FlatteStrategy(std::shared_ptr< FormFactor > FF)
Definition: Flatte.hpp:134
std::shared_ptr< ComPWA::FunctionTree::TreeNode > createFunctionTree(InputInfo Params, std::shared_ptr< ComPWA::FunctionTree::Value< std::vector< double >>> InvMassSquared)
Definition: Flatte.cpp:17
std::complex< double > couplingToWidth(double mR, double g, std::complex< double > gamma, std::complex< double > phspFactor)
Convert width to complex coupling.
Definition: Coupling.hpp:18
std::vector< Coupling > HiddenCouplings
Coupling parameters and final state masses for multiple hidden channels.
Definition: Flatte.hpp:21
unsigned int L
Orbital Angular Momentum between two daughters in Resonance decay.
std::complex< double > qValueAC(double sqrtS, double ma, double mb)
Calculate Break-up momentum.
Definition: FormFactor.hpp:95
std::complex< double > flatteCouplingTerm(double sqrtS, double mR, double coupling, double massA, double massB, unsigned int J, double mesonRadius, std::shared_ptr< FormFactor > FormFactorFunctor)
Helper function to calculate the coupling terms for the Flatte formular.
Definition: Flatte.hpp:26
ParType
Enums for the type of the parameter, should be extended if an new parameter type is added...
Definition: Parameter.hpp:34
Virtual base class for operations of FunctionTree nodes.
Definition: Functions.hpp:30
std::shared_ptr< ComPWA::FunctionTree::FitParameter > G
Coupling to signal channel.
Definition: Flatte.hpp:19
This class provides a list of parameters and values of different types.