diff --git a/doc/ReleaseNotes.md b/doc/ReleaseNotes.md --- a/doc/ReleaseNotes.md +++ b/doc/ReleaseNotes.md @@ -1,5 +1,9 @@ # Laura++ release notes +1st February 2021 Dan Johnson +* Allow floating of parameters in the K-matrix + - see https://phab.hepforge.org/T59 + 2nd December 2020 Thomas Latham * Fix LauFormulaPar to follow change in behaviour of TFormula - see https://phab.hepforge.org/T129 diff --git a/inc/LauAbsResonance.hh b/inc/LauAbsResonance.hh --- a/inc/LauAbsResonance.hh +++ b/inc/LauAbsResonance.hh @@ -473,37 +473,37 @@ LauAbsResonance& operator=(const LauAbsResonance& rhs); //! Information on the resonance - LauResonanceInfo* resInfo_; + LauResonanceInfo* resInfo_{0}; //! Information on the particles const LauDaughters* daughters_; //! Parent name - TString nameParent_; + TString nameParent_{""}; //! Daughter 1 name - TString nameDaug1_; + TString nameDaug1_{""}; //! Daughter 2 name - TString nameDaug2_; + TString nameDaug2_{""}; //! Bachelor name - TString nameBachelor_; + TString nameBachelor_{""}; //! Parent charge - Int_t chargeParent_; + Int_t chargeParent_{0}; //! Daughter 1 charge - Int_t chargeDaug1_; + Int_t chargeDaug1_{0}; //! Daughter 2 charge - Int_t chargeDaug2_; + Int_t chargeDaug2_{0}; //! Bachelor charge - Int_t chargeBachelor_; + Int_t chargeBachelor_{0}; //! Parent mass - Double_t massParent_; + Double_t massParent_{0.0}; //! Daughter 1 mass - Double_t massDaug1_; + Double_t massDaug1_{0.0}; //! Daughter 2 mass - Double_t massDaug2_; + Double_t massDaug2_{0.0}; // Bachelor mass - Double_t massBachelor_; + Double_t massBachelor_{0.0}; //! Resonance name TString resName_; @@ -512,49 +512,49 @@ TString sanitisedName_; //! Resonance mass - LauParameter* resMass_; + LauParameter* resMass_{0}; //! Resonance width - LauParameter* resWidth_; + LauParameter* resWidth_{0}; //! All parameters of the resonance std::vector resParameters_; //! Resonance spin - Int_t resSpin_; + Int_t resSpin_{0}; //! Resonance charge - Int_t resCharge_; + Int_t resCharge_{0}; //! DP axis identifier Int_t resPairAmpInt_; //! Blatt Weisskopf barrier for parent decay - LauBlattWeisskopfFactor* parBWFactor_; + LauBlattWeisskopfFactor* parBWFactor_{0}; //! Blatt Weisskopf barrier for resonance decay - LauBlattWeisskopfFactor* resBWFactor_; + LauBlattWeisskopfFactor* resBWFactor_{0}; //! Spin formalism - LauSpinType spinType_; + LauSpinType spinType_{Zemach_P}; //! Boolean to flip helicity - Bool_t flipHelicity_; + Bool_t flipHelicity_{kFALSE}; //! Boolean to ignore the momentum factors in both the spin factor and the mass-dependent width - Bool_t ignoreMomenta_; + Bool_t ignoreMomenta_{kFALSE}; //! Boolean to set the spinTerm to unity always - Bool_t ignoreSpin_; + Bool_t ignoreSpin_{kFALSE}; //! Boolean to ignore barrier factor scaling in the amplitude numerator, they are still used for the mass-dependent width - Bool_t ignoreBarrierScaling_; + Bool_t ignoreBarrierScaling_{kFALSE}; // Event kinematics information //! Invariant mass - Double_t mass_; + Double_t mass_{0.0}; //! Helicity angle cosine - Double_t cosHel_; + Double_t cosHel_{0.0}; //! Daughter momentum in resonance rest frame - Double_t q_; + Double_t q_{0.0}; //! Bachelor momentum in resonance rest frame - Double_t p_; + Double_t p_{0.0}; //! Bachelor momentum in parent rest frame - Double_t pstar_; + Double_t pstar_{0.0}; //! Covariant factor /*! @@ -566,10 +566,10 @@ \see LauKinematics::getcov13 \see LauKinematics::getcov23 */ - Double_t erm_; + Double_t erm_{1.0}; //! Covariant factor (full spin-dependent expression) - Double_t covFactor_; + Double_t covFactor_{1.0}; ClassDef(LauAbsResonance,0) // Abstract resonance class diff --git a/inc/LauKMatrixProdPole.hh b/inc/LauKMatrixProdPole.hh --- a/inc/LauKMatrixProdPole.hh +++ b/inc/LauKMatrixProdPole.hh @@ -44,7 +44,7 @@ class LauKMatrixProdPole : public LauAbsResonance { - public: + public: //! Constructor /*! \param [in] poleName name of the pole @@ -54,34 +54,40 @@ \param [in] daughters the daughter particles \param [in] useProdAdler boolean to turn on/off the production Adler zero factor */ - LauKMatrixProdPole(const TString& poleName, Int_t poleIndex, Int_t resPairAmpInt, - LauKMatrixPropagator* propagator, const LauDaughters* daughters, - Bool_t useProdAdler = kFALSE); + LauKMatrixProdPole( const TString& poleName, Int_t poleIndex, Int_t resPairAmpInt, + LauKMatrixPropagator* propagator, const LauDaughters* daughters, + Bool_t useProdAdler = kFALSE); //! Destructor - virtual ~LauKMatrixProdPole(); + virtual ~LauKMatrixProdPole(); // Initialise the model - virtual void initialise() {return;} + virtual void initialise() {return;} //! The amplitude calculation /*! \param [in] kinematics the kinematic variables of the current event \return the complex amplitude */ - virtual LauComplex amplitude(const LauKinematics* kinematics); + virtual LauComplex amplitude(const LauKinematics* kinematics); //! Get the resonance model type - /*! - \return the resonance model type - */ + /*! + \return the resonance model type + */ virtual LauAbsResonance::LauResonanceModel getResonanceModel() const {return LauAbsResonance::KMatrix;} + //! Retrieve the resonance parameters, e.g. so that they can be loaded into a fit + /*! + \return floating parameters of the resonance + */ + virtual const std::vector& getFloatingParameters(); + protected: //! Function not meant to be called, amplitude is called directly in this case - virtual LauComplex resAmp(Double_t mass, Double_t spinTerm); + virtual LauComplex resAmp(Double_t mass, Double_t spinTerm); - private: + private: //! Copy constructor (not implemented) LauKMatrixProdPole(const LauKMatrixProdPole& rhs); @@ -89,14 +95,14 @@ LauKMatrixProdPole& operator=(const LauKMatrixProdPole& rhs); //! The K-matrix propagator - LauKMatrixPropagator* thePropagator_; + LauKMatrixPropagator* thePropagator_; //! The number of the pole Int_t poleIndex_; - //! Boolean to turn on/off the production Adler zero factor - Bool_t useProdAdler_; + //! Boolean to turn on/off the production Adler zero factor + Bool_t useProdAdler_; - ClassDef(LauKMatrixProdPole, 0) // K-matrix production pole + ClassDef(LauKMatrixProdPole, 0) // K-matrix production pole }; diff --git a/inc/LauKMatrixProdSVP.hh b/inc/LauKMatrixProdSVP.hh --- a/inc/LauKMatrixProdSVP.hh +++ b/inc/LauKMatrixProdSVP.hh @@ -76,6 +76,12 @@ \return the resonance model type */ virtual LauAbsResonance::LauResonanceModel getResonanceModel() const {return LauAbsResonance::KMatrix;} + + //! Retrieve the resonance parameters, e.g. so that they can be loaded into a fit + /*! + \return floating parameters of the resonance + */ + const std::vector& getFloatingParameters(); protected: //! Function not meant to be called diff --git a/inc/LauKMatrixPropagator.hh b/inc/LauKMatrixPropagator.hh --- a/inc/LauKMatrixPropagator.hh +++ b/inc/LauKMatrixPropagator.hh @@ -49,7 +49,7 @@ class LauKMatrixPropagator { - public: + public: //! Constructor /*! \param [in] name name of the propagator @@ -59,12 +59,12 @@ \param [in] nPoles the number of poles \param [in] rowIndex this specifies which row of the propagator should be used when summing over the amplitude channels */ - LauKMatrixPropagator(const TString& name, const TString& paramFileName, - Int_t resPairAmpInt, Int_t nChannels, Int_t nPoles, - Int_t rowIndex = 1); + LauKMatrixPropagator( const TString& name, const TString& paramFileName, + Int_t resPairAmpInt, Int_t nChannels, Int_t nPoles, + Int_t rowIndex = 1 ); //! Destructor - virtual ~LauKMatrixPropagator(); + virtual ~LauKMatrixPropagator(); //! Calculate the invariant mass squared s /*! @@ -82,25 +82,25 @@ /*! \param [in] inputFile name of the input file */ - void setParameters(const TString& inputFile); + void setParameters(const TString& inputFile); - //! Get the scattering K matrix - /*! - \return the real, symmetric scattering K matrix + //! Get the scattering K matrix + /*! + \return the real, symmetric scattering K matrix */ - TMatrixD getKMatrix() const {return ScattKMatrix_;} + TMatrixD getKMatrix() const {return ScattKMatrix_;} - //! Get the real part of the propagator full matrix - /*! - \return the real part of the propagator full matrix - */ - TMatrixD getRealPropMatrix() const {return realProp_;} + //! Get the real part of the propagator full matrix + /*! + \return the real part of the propagator full matrix + */ + TMatrixD getRealPropMatrix() const {return realProp_;} - //! Get the negative imaginary part of the full propagator matrix - /*! - \return the negative imaginary part of the full propagator matrix - */ - TMatrixD getNegImagPropMatrix() const {return negImagProp_;} + //! Get the negative imaginary part of the full propagator matrix + /*! + \return the negative imaginary part of the full propagator matrix + */ + TMatrixD getNegImagPropMatrix() const {return negImagProp_;} //! Get the real part of the term of the propagator /*! @@ -123,6 +123,14 @@ */ Double_t getPoleDenomTerm(Int_t poleIndex) const; + //! Get coupling parameters, set according to the input file + /*! + \param [in] poleIndex number of the required pole + \param [in] channelIndex number of the required channel + \return the parameter of the coupling constant + */ + LauParameter& getPoleMassSqParameter(Int_t poleIndex); + //! Get coupling constants that were loaded from the input file /*! \param [in] poleIndex number of the required pole @@ -131,6 +139,14 @@ */ Double_t getCouplingConstant(Int_t poleIndex, Int_t channelIndex) const; + //! Get coupling parameters, set according to the input file + /*! + \param [in] poleIndex number of the required pole + \param [in] channelIndex number of the required channel + \return the parameter of the coupling constant + */ + LauParameter& getCouplingParameter(Int_t poleIndex, Int_t channelIndex); + //! Get scattering constants that were loaded from the input file /*! \param [in] channel1Index number of the first channel index @@ -139,6 +155,44 @@ */ Double_t getScatteringConstant(Int_t channel1Index, Int_t channel2Index) const; + //! Get scattering parameters, set according to the input file + /*! + \param [in] channel1Index number of the first channel index + \param [in] channel2Index number of the second channel index + \return the parameter of the scattering constant + */ + LauParameter& getScatteringParameter(Int_t channel1Index, Int_t channel2Index); + + //! Get mSq0 production parameter + /*! + \return the mSq0 parameter + */ + LauParameter& getmSq0() {return mSq0_;} + + //! Get s0Scatt production parameter + /*! + \return the s0Scatt parameter + */ + LauParameter& gets0Scatt() {return s0Scatt_;} + + //! Get s0 production parameter + /*! + \return the s0Prod parameter + */ + LauParameter& gets0Prod() {return s0Prod_;} + + //! Get sA production parameter + /*! + \return the sA parameter + */ + LauParameter& getsA() {return sA_;} + + //! Get sA0 production parameter + /*! + \return the sA0 parameter + */ + LauParameter& getsA0() {return sA0_;} + //! Get the "slowly-varying part" term of the amplitude /*! \return the svp term @@ -154,62 +208,60 @@ //! Get the DP axis identifier /*! - /return the value to identify the DP axis in question + \return the value to identify the DP axis in question */ Int_t getResPairAmpInt() const {return resPairAmpInt_;} //! Get the number of channels /*! - /return the number of channels + \return the number of channels */ Int_t getNChannels() const {return nChannels_;} //! Get the number of poles /*! - /return the number of poles + \return the number of poles */ Int_t getNPoles() const {return nPoles_;} //! Get the propagator name /*! - /return the name of the propagator + \return the name of the propagator */ TString getName() const {return name_;} - //! Get the unitary transition amplitude for the given channel - /*! - \param [in] s The invariant mass squared - \param [in] channel The index number of the channel process - \return the complex amplitude T + //! Get the unitary transition amplitude for the given channel + /*! + \param [in] s The invariant mass squared + \param [in] channel The index number of the channel process + \return the complex amplitude T */ - LauComplex getTransitionAmp(Double_t s, Int_t channel); + LauComplex getTransitionAmp(Double_t s, Int_t channel); - - //! Get the complex phase space term for the given channel and invariant mass squared - /*! - \param [in] s The invariant mass squared - \param [in] channel The index number of the channel process - \return the complex phase space term rho(channel, channel) + //! Get the complex phase space term for the given channel and invariant mass squared + /*! + \param [in] s The invariant mass squared + \param [in] channel The index number of the channel process + \return the complex phase space term rho(channel, channel) */ - LauComplex getPhaseSpaceTerm(Double_t s, Int_t channel); + LauComplex getPhaseSpaceTerm(Double_t s, Int_t channel); - //! Get the Adler zero factor, which is set when updatePropagator is called - /*! - \return the Adler zero factor + //! Get the Adler zero factor, which is set when updatePropagator is called + /*! + \return the Adler zero factor */ - Double_t getAdlerZero() const {return adlerZeroFactor_;} + Double_t getAdlerZero() const {return adlerZeroFactor_;} - //! Get the THat amplitude for the given s and channel number - /*! - \param [in] s The invariant mass squared + //! Get the THat amplitude for the given s and channel number + /*! + \param [in] s The invariant mass squared \param [in] channel The index number of the channel process \return the complex THat amplitude */ - LauComplex getTHat(Double_t s, Int_t channel); - + LauComplex getTHat(Double_t s, Int_t channel); - protected: + protected: //! Calculate the scattering K-matrix for the given value of s /*! \param [in] s the invariant mass squared @@ -318,22 +370,22 @@ Bool_t checkPhaseSpaceType(Int_t phaseSpaceInt) const; - //! Get the unitary transition amplitude matrix for the given kinematics - /*! - \param [in] kinematics The pointer to the constant kinematics + //! Get the unitary transition amplitude matrix for the given kinematics + /*! + \param [in] kinematics The pointer to the constant kinematics */ - void getTMatrix(const LauKinematics* kinematics); + void getTMatrix(const LauKinematics* kinematics); - //! Get the unitary transition amplitude matrix for the given kinematics - /*! - \param [in] s The invariant mass squared of the system + //! Get the unitary transition amplitude matrix for the given kinematics + /*! + \param [in] s The invariant mass squared of the system */ - void getTMatrix(Double_t s); + void getTMatrix(Double_t s); - //! Get the square root of the phase space matrix - void getSqrtRhoMatrix(); - - private: + //! Get the square root of the phase space matrix + void getSqrtRhoMatrix(); + + private: //! Copy constructor (not implemented) LauKMatrixPropagator(const LauKMatrixPropagator& rhs); @@ -341,39 +393,39 @@ LauKMatrixPropagator& operator=(const LauKMatrixPropagator& rhs); //! Create a map for the K-matrix parameters - typedef std::map > KMatrixParamMap; + typedef std::map > KMatrixParamMap; - //! Initialise and set the dimensions for the internal matrices and parameter arrays - void initialiseMatrices(); + //! Initialise and set the dimensions for the internal matrices and parameter arrays + void initialiseMatrices(); - //! Store the (phase space) channel indices from a line in the parameter file - /*! - \param [in] theLine Vector of strings corresponding to the line from the parameter file + //! Store the (phase space) channel indices from a line in the parameter file + /*! + \param [in] theLine Vector of strings corresponding to the line from the parameter file */ - void storeChannels(const std::vector& theLine); + void storeChannels(const std::vector& theLine); - //! Store the pole mass and couplings from a line in the parameter file - /*! - \param [in] theLine Vector of strings corresponding to the line from the parameter file + //! Store the pole mass and couplings from a line in the parameter file + /*! + \param [in] theLine Vector of strings corresponding to the line from the parameter file */ - void storePole(const std::vector& theLine); + void storePole(const std::vector& theLine); - //! Store the scattering coefficients from a line in the parameter file - /*! - \param [in] theLine Vector of strings corresponding to the line from the parameter file + //! Store the scattering coefficients from a line in the parameter file + /*! + \param [in] theLine Vector of strings corresponding to the line from the parameter file */ - void storeScattering(const std::vector& theLine); + void storeScattering(const std::vector& theLine); - //! Store miscelleanous parameters from a line in the parameter file - /*! - \param [in] keyword the name of the parameter to be set - \param [in] parString the string containing the value of the parameter + //! Store miscelleanous parameters from a line in the parameter file + /*! + \param [in] keyword the name of the parameter to be set + \param [in] parString the string containing the value of the parameter */ - void storeParameter(const TString& keyword, const TString& parString); + void storeParameter(const TString& keyword, const TString& parString); - //! String to store the propagator name + //! String to store the propagator name TString name_; //! Name of the input parameter file TString paramFileName_; @@ -415,10 +467,10 @@ TMatrixD ReSqrtRhoMatrix_; //! Imaginary part of the square root of the phase space density diagonal matrix TMatrixD ImSqrtRhoMatrix_; - //! Real part of the unitary T matrix - TMatrixD ReTMatrix_; - //! Imaginary part of the unitary T matrix - TMatrixD ImTMatrix_; + //! Real part of the unitary T matrix + TMatrixD ReTMatrix_; + //! Imaginary part of the unitary T matrix + TMatrixD ImTMatrix_; //! Number of channels Int_t nChannels_; @@ -426,12 +478,12 @@ Int_t nPoles_; //! Vector of squared pole masses - std::vector mSqPoles_; + std::vector mSqPoles_; //! Array of coupling constants - LauParArray gCouplings_; + LauParArray gCouplings_; //! Array of scattering SVP values - LauParArray fScattering_; + LauParArray fScattering_; //! Vector of phase space types std::vector phaseSpaceTypes_; @@ -490,8 +542,8 @@ //! Control the output of the functions Bool_t verbose_; - //! Control if scattering constants are channel symmetric: f_ji = f_ij - Bool_t scattSymmetry_; + //! Control if scattering constants are channel symmetric: f_ji = f_ij + Bool_t scattSymmetry_; ClassDef(LauKMatrixPropagator,0) // K-matrix amplitude model }; diff --git a/src/LauAbsResonance.cc b/src/LauAbsResonance.cc --- a/src/LauAbsResonance.cc +++ b/src/LauAbsResonance.cc @@ -85,30 +85,13 @@ LauAbsResonance::LauAbsResonance(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) : resInfo_(resInfo), daughters_(daughters), - nameParent_(""), nameDaug1_(""), nameDaug2_(""), nameBachelor_(""), - chargeParent_(0), chargeDaug1_(0), chargeDaug2_(0), chargeBachelor_(0), - massParent_(0.0), massDaug1_(0.0), massDaug2_(0.0), massBachelor_(0.0), resName_( (resInfo!=0) ? resInfo->getName() : "" ), sanitisedName_( (resInfo!=0) ? resInfo->getSanitisedName() : "" ), resMass_( (resInfo!=0) ? resInfo->getMass() : 0 ), resWidth_( (resInfo!=0) ? resInfo->getWidth() : 0 ), resSpin_( (resInfo!=0) ? resInfo->getSpin() : 0 ), resCharge_( (resInfo!=0) ? resInfo->getCharge() : 0 ), - resPairAmpInt_(resPairAmpInt), - parBWFactor_(0), - resBWFactor_(0), - spinType_(Zemach_P), - flipHelicity_(kFALSE), - ignoreMomenta_(kFALSE), - ignoreSpin_(kFALSE), - ignoreBarrierScaling_(kFALSE), - mass_(0.0), - cosHel_(0.0), - q_(0.0), - p_(0.0), - pstar_(0.0), - erm_(1.0), - covFactor_(1.0) + resPairAmpInt_(resPairAmpInt) { if ( resInfo == 0 ) { std::cerr << "ERROR in LauAbsResonance constructor : null LauResonanceInfo object provided" << std::endl; @@ -143,32 +126,10 @@ // Constructor LauAbsResonance::LauAbsResonance(const TString& resName, const Int_t resPairAmpInt, const LauDaughters* daughters) : - resInfo_(0), daughters_(daughters), - nameParent_(""), nameDaug1_(""), nameDaug2_(""), nameBachelor_(""), - chargeParent_(0), chargeDaug1_(0), chargeDaug2_(0), chargeBachelor_(0), - massParent_(0.0), massDaug1_(0.0), massDaug2_(0.0), massBachelor_(0.0), resName_(resName), sanitisedName_(resName), - resMass_(0), - resWidth_(0), - resSpin_(0), - resCharge_(0), - resPairAmpInt_(resPairAmpInt), - parBWFactor_(0), - resBWFactor_(0), - spinType_(Zemach_P), - flipHelicity_(kFALSE), - ignoreMomenta_(kFALSE), - ignoreSpin_(kFALSE), - ignoreBarrierScaling_(kFALSE), - mass_(0.0), - cosHel_(0.0), - q_(0.0), - p_(0.0), - pstar_(0.0), - erm_(1.0), - covFactor_(1.0) + resPairAmpInt_(resPairAmpInt) { if ( daughters_ == 0 ) { std::cerr << "ERROR in LauAbsResonance constructor : null LauDaughters object provided" << std::endl; diff --git a/src/LauKMatrixProdPole.cc b/src/LauKMatrixProdPole.cc --- a/src/LauKMatrixProdPole.cc +++ b/src/LauKMatrixProdPole.cc @@ -67,8 +67,6 @@ { // Calculate the amplitude for the K-matrix production pole. - // First, make sure the K-matrix propagator is up-to-date for - // the given centre-of-mass squared value ("s") from the kinematics. LauComplex amp(0.0, 0.0); if (thePropagator_ == 0) { @@ -106,3 +104,29 @@ return amp; } + +const std::vector& LauKMatrixProdPole::getFloatingParameters() +{ + + this->clearFloatingParameters(); + + Int_t nChannels = thePropagator_->getNChannels(); + + for (int jChannel = 0 ; jChannel < nChannels ; jChannel++) + { + LauParameter& par_gj_ = thePropagator_->getCouplingParameter(poleIndex_, jChannel); + if ( !par_gj_.fixed() ) + { + this->addFloatingParameter( &par_gj_ ); + } + } + + LauParameter& par_polemasssq_ = thePropagator_->getPoleMassSqParameter(poleIndex_); + if ( !par_polemasssq_.fixed() ) + { + this->addFloatingParameter( &par_polemasssq_ ); + } + + return this->getParameters(); + +} \ No newline at end of file diff --git a/src/LauKMatrixProdSVP.cc b/src/LauKMatrixProdSVP.cc --- a/src/LauKMatrixProdSVP.cc +++ b/src/LauKMatrixProdSVP.cc @@ -94,3 +94,53 @@ } +const std::vector& LauKMatrixProdSVP::getFloatingParameters() +{ + + this->clearFloatingParameters(); + + Int_t nChannels = thePropagator_->getNChannels(); + + for (int jChannel = 0 ; jChannel < nChannels ; jChannel++) + { + LauParameter& par_f_ = thePropagator_->getScatteringParameter(channelIndex_, jChannel); + if ( !par_f_.fixed() ) + { + this->addFloatingParameter( &par_f_ ); + } + + } + + LauParameter& par_mSq0_ = thePropagator_->getmSq0(); + if ( !par_mSq0_.fixed() ) + { + this->addFloatingParameter( &par_mSq0_ ); + } + + LauParameter& par_s0Scatt_ = thePropagator_->gets0Scatt(); + if ( !par_s0Scatt_.fixed() ) + { + this->addFloatingParameter( &par_s0Scatt_ ); + } + + LauParameter& par_s0Prod_ = thePropagator_->gets0Prod(); + if ( !par_s0Prod_.fixed() ) + { + this->addFloatingParameter( &par_s0Prod_ ); + } + + LauParameter& par_sA_ = thePropagator_->getsA(); + if ( !par_sA_.fixed() ) + { + this->addFloatingParameter( &par_sA_ ); + } + + LauParameter& par_sA0_ = thePropagator_->getsA0(); + if ( !par_sA0_.fixed() ) + { + this->addFloatingParameter( &par_sA0_ ); + } + + return this->getParameters(); + +} \ No newline at end of file diff --git a/src/LauKMatrixPropagator.cc b/src/LauKMatrixPropagator.cc --- a/src/LauKMatrixPropagator.cc +++ b/src/LauKMatrixPropagator.cc @@ -45,7 +45,7 @@ ClassImp(LauKMatrixPropagator) LauKMatrixPropagator::LauKMatrixPropagator(const TString& name, const TString& paramFile, - Int_t resPairAmpInt, Int_t nChannels, + Int_t resPairAmpInt, Int_t nChannels, Int_t nPoles, Int_t rowIndex) : name_(name), paramFileName_(paramFile), @@ -452,7 +452,7 @@ Double_t poleMass = std::atof(theLine[2].c_str()); Double_t poleMassSq = poleMass*poleMass; - LauParameter mPoleParam(poleMassSq); + LauParameter mPoleParam(Form("KM_%s_poleMassSq_%i",name_.Data(),poleIndex),poleMassSq); mSqPoles_[poleIndex] = mPoleParam; cout<<"Added bare pole mass "<= nPoles_) ) { + std::cerr << "ERROR from LauKMatrixPropagator::getPoleMassSqParameter(). Invalid pole." << std::endl; + gSystem->Exit(EXIT_FAILURE); + } + + return mSqPoles_[poleIndex]; +} + Double_t LauKMatrixPropagator::getCouplingConstant(Int_t poleIndex, Int_t channelIndex) const { @@ -663,6 +673,18 @@ } +LauParameter& LauKMatrixPropagator::getCouplingParameter(Int_t poleIndex, Int_t channelIndex) +{ + + if ( (parametersSet_ == kFALSE) || (poleIndex < 0 || poleIndex >= nPoles_) || (channelIndex < 0 || channelIndex >= nChannels_) ) { + std::cerr << "ERROR from LauKMatrixPropagator::getCouplingParameter(). Invalid coupling." << std::endl; + gSystem->Exit(EXIT_FAILURE); + } + + //std::cout << "Minvalue + range for " << poleIndex << ", " << channelIndex << ": " << gCouplings_[poleIndex][channelIndex].minValue() << " => + " << gCouplings_[poleIndex][channelIndex].range() << + // " and init value: " << gCouplings_[poleIndex][channelIndex].initValue() << std::endl; + return gCouplings_[poleIndex][channelIndex]; +} Double_t LauKMatrixPropagator::getScatteringConstant(Int_t channel1Index, Int_t channel2Index) const { @@ -676,6 +698,17 @@ } +LauParameter& LauKMatrixPropagator::getScatteringParameter(Int_t channel1Index, Int_t channel2Index) +{ + + if ( (parametersSet_ == kFALSE) || (channel1Index < 0 || channel1Index >= nChannels_) || (channel2Index < 0 || channel2Index >= nChannels_) ) { + std::cerr << "ERROR from LauKMatrixPropagator::getScatteringParameter(). Invalid chanel index." << std::endl; + gSystem->Exit(EXIT_FAILURE); + } + + return fScattering_[channel1Index][channel2Index]; +} + Double_t LauKMatrixPropagator::calcSVPTerm(Double_t s, Double_t s0) const {