diff --git a/src/LauAbsResonance.cc b/src/LauAbsResonance.cc index 4b654e7..37c9be7 100644 --- a/src/LauAbsResonance.cc +++ b/src/LauAbsResonance.cc @@ -1,719 +1,720 @@ /* Copyright 2004 University of Warwick Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0 Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License. */ /* Laura++ package authors: John Back Paul Harrison Thomas Latham */ /*! \file LauAbsResonance.cc \brief File containing implementation of LauAbsResonance class. */ #include #include "TSystem.h" #include "LauAbsResonance.hh" #include "LauConstants.hh" #include "LauDaughters.hh" #include "LauKinematics.hh" #include "LauParameter.hh" #include "LauResonanceInfo.hh" ClassImp(LauAbsResonance) bool LauAbsResonance::isIncoherentModel(LauResonanceModel model) { switch(model) { case BW: case RelBW: case GS: case Flatte: case Sigma: case Kappa: case Dabba: case LASS: case LASS_BW: case LASS_NR: case EFKLLM: case KMatrix: case FlatNR: case NRModel: case BelleNR: case PowerLawNR: case BelleSymNR: case BelleSymNRNoInter: case TaylorNR: case PolNR: case POLE: case PolarFFNR: case PolarFFSymNR: case PolarFFSymNRNoInter: case Rescattering: + case Rescattering2: case RescatteringNoInter: case MIPW_MagPhase: case MIPW_RealImag: case RhoOmegaMix_GS: case RhoOmegaMix_RBW: case RhoOmegaMix_GS_1: case RhoOmegaMix_RBW_1: break; case GaussIncoh: return true; } return false; } // Constructor 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) { if ( resInfo == 0 ) { std::cerr << "ERROR in LauAbsResonance constructor : null LauResonanceInfo object provided" << std::endl; gSystem->Exit(EXIT_FAILURE); } if ( daughters_ == 0 ) { std::cerr << "ERROR in LauAbsResonance constructor : null LauDaughters object provided" << std::endl; gSystem->Exit(EXIT_FAILURE); } nameParent_ = this->getNameParent(); nameDaug1_ = this->getNameDaug1(); nameDaug2_ = this->getNameDaug2(); nameBachelor_ = this->getNameBachelor(); massParent_ = this->getMassParent(); massDaug1_ = this->getMassDaug1(); massDaug2_ = this->getMassDaug2(); massBachelor_ = this->getMassBachelor(); chargeParent_ = this->getChargeParent(); chargeDaug1_ = this->getChargeDaug1(); chargeDaug2_ = this->getChargeDaug2(); chargeBachelor_ = this->getChargeBachelor(); // check that the total charge adds up to that of the resonance Int_t totalCharge = chargeDaug1_ + chargeDaug2_; if ( (totalCharge != resCharge_) && (resPairAmpInt_ != 0) ) { std::cerr << "ERROR in LauAbsResonance : Total charge of daughters = " << totalCharge << ". Resonance charge = " << resCharge_ << "." << std::endl; gSystem->Exit(EXIT_FAILURE); } } // 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) { if ( daughters_ == 0 ) { std::cerr << "ERROR in LauAbsResonance constructor : null LauDaughters object provided" << std::endl; gSystem->Exit(EXIT_FAILURE); } nameParent_ = this->getNameParent(); nameDaug1_ = this->getNameDaug1(); nameDaug2_ = this->getNameDaug2(); nameBachelor_ = this->getNameBachelor(); massParent_ = this->getMassParent(); massDaug1_ = this->getMassDaug1(); massDaug2_ = this->getMassDaug2(); massBachelor_ = this->getMassBachelor(); chargeParent_ = this->getChargeParent(); chargeDaug1_ = this->getChargeDaug1(); chargeDaug2_ = this->getChargeDaug2(); chargeBachelor_ = this->getChargeBachelor(); // Since we haven't been provided with a LauResonanceInfo object we can just // set the change of the resonance to be the sum of the daughter charges resCharge_ = chargeDaug1_ + chargeDaug2_; } // Destructor LauAbsResonance::~LauAbsResonance() { } LauComplex LauAbsResonance::amplitude(const LauKinematics* kinematics) { // Use LauKinematics interface for amplitude // For resonance made from tracks i, j, we need the momenta // of tracks i and k in the i-j rest frame for spin helicity calculations // in the Zemach tensor formalism. // Also need the momentum of track k in the parent rest-frame for // calculation of the Blatt-Weisskopf factors. mass_ = 0.0; cosHel_ = 0.0; q_ = 0.0; p_ = 0.0; pstar_ = 0.0; erm_ = 1.0; covFactor_ = 1.0; if (resPairAmpInt_ == 1) { mass_ = kinematics->getm23(); cosHel_ = kinematics->getc23(); q_ = kinematics->getp2_23(); p_ = kinematics->getp1_23(); pstar_ = kinematics->getp1_Parent(); erm_ = kinematics->getcov23(); } else if (resPairAmpInt_ == 2) { mass_ = kinematics->getm13(); cosHel_ = kinematics->getc13(); q_ = kinematics->getp1_13(); p_ = kinematics->getp2_13(); pstar_ = kinematics->getp2_Parent(); erm_ = kinematics->getcov13(); } else if (resPairAmpInt_ == 3) { mass_ = kinematics->getm12(); cosHel_ = kinematics->getc12(); q_ = kinematics->getp1_12(); p_ = kinematics->getp3_12(); pstar_ = kinematics->getp3_Parent(); erm_ = kinematics->getcov12(); } else { std::cerr << "ERROR in LauAbsResonance::amplitude : Nonsense setup of resPairAmp array." << std::endl; gSystem->Exit(EXIT_FAILURE); } if (this->flipHelicity()) { cosHel_ *= -1.0; } if (this->ignoreMomenta()) { q_ = 1.0; p_ = 1.0; pstar_ = 1.0; erm_ = 1.0; } // Calculate the spin factors Double_t spinTerm(1.0); Double_t pProd(1.0); if (!this->ignoreSpin()) { switch ( this->getSpinType() ) { case Zemach_P: pProd = q_*p_; spinTerm = this->calcZemachSpinFactor( pProd ); break; case Zemach_Pstar: pProd = q_*pstar_; spinTerm = this->calcZemachSpinFactor( pProd ); break; case Covariant: pProd = q_*pstar_; spinTerm = this->calcCovSpinFactor( pProd ); break; case Legendre: spinTerm = this->calcLegendrePoly(); break; } } // Calculate the full amplitude LauComplex resAmplitude = this->resAmp(mass_, spinTerm); return resAmplitude; } void LauAbsResonance::calcCovFactor( const Double_t erm ) { if (resSpin_ == 0) { covFactor_ = 1.0; } else if (resSpin_ == 1) { covFactor_ = erm; } else if (resSpin_ == 2) { covFactor_ = erm*erm + 0.5; } else if (resSpin_ == 3) { covFactor_ = erm*(erm*erm + 1.5); } else if (resSpin_ == 4) { covFactor_ = (8.*erm*erm*erm*erm + 24.*erm*erm + 3.)/35.; } else if (resSpin_ > 4) { std::cerr << "WARNING in LauAbsResonance::calcCovFactor : covariant spin factor cannot (yet) be fully calculated for spin >= 5" << std::endl; std::cerr << " : the function of sqrt(1 + (p/mParent)^2) part will be missing" << std::endl; covFactor_ = 1.0; } } Double_t LauAbsResonance::calcCovSpinFactor( const Double_t pProd ) { if (resSpin_ == 0) { covFactor_ = 1.0; return 1.0; } // Covariant spin factor is (p* q)^L * f_L(erm) * P_L(cosHel) Double_t spinFactor(pProd); for ( Int_t i(1); i < resSpin_; ++i ) { spinFactor *= pProd; } this->calcCovFactor( erm_ ); spinFactor *= covFactor_; spinFactor *= this->calcLegendrePoly(); return spinFactor; } Double_t LauAbsResonance::calcZemachSpinFactor( const Double_t pProd ) const { // Calculate the spin factors // // These are calculated as follows // // -2^j * (q*p)^j * cj * Pj(cosHel) // // where Pj(coshHel) is the jth order Legendre polynomial and // // cj = j! / (2j-1)!! if (resSpin_ == 0) { return 1.0; } Double_t spinFactor(pProd); for ( Int_t i(1); i < resSpin_; ++i ) { spinFactor *= pProd; } spinFactor *= this->calcLegendrePoly(); return spinFactor; } Double_t LauAbsResonance::calcLegendrePoly( const Double_t cosHel ) { cosHel_ = cosHel; return this->calcLegendrePoly(); } Double_t LauAbsResonance::calcLegendrePoly() const { Double_t legPol = 1.0; if (resSpin_ == 1) { // Calculate vector resonance Legendre polynomial legPol = -2.0*cosHel_; } else if (resSpin_ == 2) { // Calculate tensor resonance Legendre polynomial legPol = 4.0*(3.0*cosHel_*cosHel_ - 1.0)/3.0; } else if (resSpin_ == 3) { // Calculate spin 3 resonance Legendre polynomial legPol = -8.0*(5.0*cosHel_*cosHel_*cosHel_ - 3.0*cosHel_)/5.0; } else if (resSpin_ == 4) { // Calculate spin 4 resonance Legendre polynomial legPol = 16.0*(35.0*cosHel_*cosHel_*cosHel_*cosHel_ - 30.0*cosHel_*cosHel_ + 3.0)/35.0; } else if (resSpin_ == 5) { // Calculate spin 5 resonance Legendre polynomial legPol = -32.0*(63.0*cosHel_*cosHel_*cosHel_*cosHel_*cosHel_ - 70.0*cosHel_*cosHel_*cosHel_ + 15.0*cosHel_)/63.0; } else if (resSpin_ > 5) { std::cerr << "WARNING in LauAbsResonance::calcLegendrePoly : Legendre polynomials not (yet) implemented for spin > 5" << std::endl; } return legPol; } void LauAbsResonance::changeResonance(const Double_t newMass, const Double_t newWidth, const Int_t newSpin) { if (newMass > 0.0) { resMass_->valueAndRange(newMass,0.0,3.0*newMass); resMass_->initValue(newMass); resMass_->genValue(newMass); std::cout << "INFO in LauAbsResonance::changeResonance : Setting mass to " << resMass_->value() << std::endl; } if (newWidth > 0.0) { resWidth_->valueAndRange(newWidth,0.0,3.0*newWidth); resWidth_->initValue(newWidth); resWidth_->genValue(newWidth); std::cout << "INFO in LauAbsResonance::changeResonance : Setting width to " << resWidth_->value() << std::endl; } if (newSpin > -1) { resSpin_ = newSpin; std::cout << "INFO in LauAbsResonance::changeResonance : Setting spin to " << resSpin_ << std::endl; } } void LauAbsResonance::changeBWBarrierRadii(const Double_t resRadius, const Double_t parRadius) { if ( resRadius >= 0.0 && resBWFactor_ != 0 ) { LauParameter* resBWRadius = resBWFactor_->getRadiusParameter(); resBWRadius->value(resRadius); resBWRadius->initValue(resRadius); resBWRadius->genValue(resRadius); std::cout << "INFO in LauAbsResonance::changeBWBarrierRadii : Setting resonance factor radius to " << resBWRadius->value() << std::endl; } if ( parRadius >= 0.0 && parBWFactor_ != 0 ) { LauParameter* parBWRadius = parBWFactor_->getRadiusParameter(); parBWRadius->value(parRadius); parBWRadius->initValue(parRadius); parBWRadius->genValue(parRadius); std::cout << "INFO in LauAbsResonance::changeBWBarrierRadii : Setting parent factor radius to " << parBWRadius->value() << std::endl; } } void LauAbsResonance::setResonanceParameter(const TString& name, const Double_t value) { //This function should always be overwritten if needed in classes inheriting from LauAbsResonance. std::cerr << "WARNING in LauAbsResonance::setResonanceParameter : Unable to set parameter \"" << name << "\" to value: " << value << "." << std::endl; } void LauAbsResonance::floatResonanceParameter(const TString& name) { //This function should always be overwritten if needed in classes inheriting from LauAbsResonance. std::cerr << "WARNING in LauAbsResonance::floatResonanceParameter : Unable to release parameter \"" << name << "\"." << std::endl; } LauParameter* LauAbsResonance::getResonanceParameter(const TString& name) { //This function should always be overwritten if needed in classes inheriting from LauAbsResonance. std::cerr << "WARNING in LauAbsResonance::getResonanceParameter : Unable to get parameter \"" << name << "\"." << std::endl; return 0; } void LauAbsResonance::addFloatingParameter( LauParameter* param ) { if ( param == 0 ) { return; } if ( param->clone() ) { resParameters_.push_back( param->parent() ); } else { resParameters_.push_back( param ); } } void LauAbsResonance::fixBarrierRadii(const Bool_t fixResRad, const Bool_t fixParRad) { if ( resBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::fixBarrierRadii : resonance barrier factor not present, cannot fix/float it" << std::endl; return; } if ( parBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::fixBarrierRadii : parent barrier factor not present, cannot fix/float it" << std::endl; return; } LauParameter* resBWRadius = resBWFactor_->getRadiusParameter(); resBWRadius->fixed(fixResRad); LauParameter* parBWRadius = parBWFactor_->getRadiusParameter(); parBWRadius->fixed(fixParRad); } Bool_t LauAbsResonance::fixResRadius() const { if ( resBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::fixResRadius : resonance barrier factor not present" << std::endl; return kTRUE; } LauParameter* bwRadius = resBWFactor_->getRadiusParameter(); return bwRadius->fixed(); } Bool_t LauAbsResonance::fixParRadius() const { if ( parBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::fixParRadius : parent barrier factor not present" << std::endl; return kTRUE; } LauParameter* bwRadius = parBWFactor_->getRadiusParameter(); return bwRadius->fixed(); } Double_t LauAbsResonance::getResRadius() const { if ( resBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::getResRadius : resonance barrier factor not present" << std::endl; return -1.0; } LauParameter* bwRadius = resBWFactor_->getRadiusParameter(); return bwRadius->unblindValue(); } Double_t LauAbsResonance::getParRadius() const { if ( parBWFactor_ == 0 ) { std::cerr << "WARNING in LauAbsResonance::getParRadius : parent barrier factor not present" << std::endl; return -1.0; } LauParameter* bwRadius = parBWFactor_->getRadiusParameter(); return bwRadius->unblindValue(); } Double_t LauAbsResonance::getMassParent() const { // Get the parent mass Double_t mass(LauConstants::mB); if (daughters_) { mass = daughters_->getMassParent(); } return mass; } Double_t LauAbsResonance::getMassDaug1() const { // Get the daughter mass Double_t mass(LauConstants::mPi); if (daughters_) { if (resPairAmpInt_ == 1) { mass = daughters_->getMassDaug2(); } else if (resPairAmpInt_ == 2) { mass = daughters_->getMassDaug1(); } else if (resPairAmpInt_ == 3) { mass = daughters_->getMassDaug1(); } } return mass; } Double_t LauAbsResonance::getMassDaug2() const { // Get the daughter mass Double_t mass(LauConstants::mPi); if (daughters_) { if (resPairAmpInt_ == 1) { mass = daughters_->getMassDaug3(); } else if (resPairAmpInt_ == 2) { mass = daughters_->getMassDaug3(); } else if (resPairAmpInt_ == 3) { mass = daughters_->getMassDaug2(); } } return mass; } Double_t LauAbsResonance::getMassBachelor() const { // Get the bachelor mass Double_t mass(LauConstants::mPi); if (daughters_) { if (resPairAmpInt_ == 1) { mass = daughters_->getMassDaug1(); } else if (resPairAmpInt_ == 2) { mass = daughters_->getMassDaug2(); } else if (resPairAmpInt_ == 3) { mass = daughters_->getMassDaug3(); } } return mass; } Int_t LauAbsResonance::getChargeParent() const { // Get the parent charge Int_t charge(0); if (daughters_) { charge = daughters_->getChargeParent(); } return charge; } Int_t LauAbsResonance::getChargeDaug1() const { // Get the daughter charge Int_t charge(0); if (daughters_) { if (resPairAmpInt_ == 1) { charge = daughters_->getChargeDaug2(); } else if (resPairAmpInt_ == 2) { charge = daughters_->getChargeDaug1(); } else if (resPairAmpInt_ == 3) { charge = daughters_->getChargeDaug1(); } } return charge; } Int_t LauAbsResonance::getChargeDaug2() const { // Get the daughter charge Int_t charge(0); if (daughters_) { if (resPairAmpInt_ == 1) { charge = daughters_->getChargeDaug3(); } else if (resPairAmpInt_ == 2) { charge = daughters_->getChargeDaug3(); } else if (resPairAmpInt_ == 3) { charge = daughters_->getChargeDaug2(); } } return charge; } Int_t LauAbsResonance::getChargeBachelor() const { // Get the bachelor charge Int_t charge(0); if (daughters_) { if (resPairAmpInt_ == 1) { charge = daughters_->getChargeDaug1(); } else if (resPairAmpInt_ == 2) { charge = daughters_->getChargeDaug2(); } else if (resPairAmpInt_ == 3) { charge = daughters_->getChargeDaug3(); } } return charge; } TString LauAbsResonance::getNameParent() const { // Get the parent name TString name(""); if (daughters_) { name = daughters_->getNameParent(); } return name; } TString LauAbsResonance::getNameDaug1() const { // Get the daughter name TString name(""); if (daughters_) { if (resPairAmpInt_ == 1) { name = daughters_->getNameDaug2(); } else if (resPairAmpInt_ == 2) { name = daughters_->getNameDaug1(); } else if (resPairAmpInt_ == 3) { name = daughters_->getNameDaug1(); } } return name; } TString LauAbsResonance::getNameDaug2() const { // Get the daughter name TString name(""); if (daughters_) { if (resPairAmpInt_ == 1) { name = daughters_->getNameDaug3(); } else if (resPairAmpInt_ == 2) { name = daughters_->getNameDaug3(); } else if (resPairAmpInt_ == 3) { name = daughters_->getNameDaug2(); } } return name; } TString LauAbsResonance::getNameBachelor() const { // Get the bachelor name TString name(""); if (daughters_) { if (resPairAmpInt_ == 1) { name = daughters_->getNameDaug1(); } else if (resPairAmpInt_ == 2) { name = daughters_->getNameDaug2(); } else if (resPairAmpInt_ == 3) { name = daughters_->getNameDaug3(); } } return name; } diff --git a/src/LauPolarFormFactorNR.cc b/src/LauPolarFormFactorNR.cc index 94dc5f3..b7567b0 100644 --- a/src/LauPolarFormFactorNR.cc +++ b/src/LauPolarFormFactorNR.cc @@ -1,148 +1,148 @@ /* Copyright 2004 University of Warwick Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0 Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License. */ /* Laura++ package authors: John Back Paul Harrison Thomas Latham */ /*! \file LauBelleNR.cc \brief File containing implementation of LauBelleNR class. */ #include #include "TMath.h" #include "LauConstants.hh" #include "LauPolarFormFactorNR.hh" #include "LauDaughters.hh" #include "LauParameter.hh" #include "LauResonanceInfo.hh" ClassImp(LauPolarFormFactorNR) LauPolarFormFactorNR::LauPolarFormFactorNR(LauResonanceInfo* resInfo, const LauAbsResonance::LauResonanceModel resType, const Int_t resPairAmpInt, const LauDaughters* daughters) : LauAbsResonance(resInfo, resPairAmpInt, daughters), lambda_(0), model_(resType) { TString parName = this->getSanitisedName(); parName += "_lambda"; lambda_ = resInfo->getExtraParameter( parName ); if ( lambda_ == 0 ) { lambda_ = new LauParameter( parName, 1.0, 0.0, 10.0, kTRUE ); lambda_->secondStage(kTRUE); resInfo->addExtraParameter( lambda_ ); } } LauPolarFormFactorNR::~LauPolarFormFactorNR() { } void LauPolarFormFactorNR::initialise() { const LauDaughters* daughters = this->getDaughters(); Int_t resPairAmpInt = this->getPairInt(); if ( daughters->gotSymmetricalDP() && resPairAmpInt != 3 ) { std::cerr << "WARNING in LauPolarFormFactorNR::initialise : Dalitz plot is symmetric - this lineshape is not appropriate." << std::endl; } if ( model_ != LauAbsResonance::PolarFFNR) { std::cerr << "WARNING in LauPolarFormFactorNR::initialise : Unknown model requested, defaulting to Polar Form Factor." << std::endl; model_ = LauAbsResonance::PolarFFNR; } } -LauComplex LauPolarFormFactorNR::resAmp(Double_t mass, Double_t spinTerm) +LauComplex LauPolarFormFactorNR::resAmp(Double_t mass, Double_t) { Double_t magnitude(1.0); Double_t lambda = this->getLambda(); magnitude = 1.0/(1.0 + mass*mass /(lambda*lambda)); LauComplex resAmplitude(magnitude, 0.0); return resAmplitude; } const std::vector& LauPolarFormFactorNR::getFloatingParameters() { this->clearFloatingParameters(); if ( ! this->fixLambda() ) { this->addFloatingParameter( lambda_ ); } return this->getParameters(); } void LauPolarFormFactorNR::setResonanceParameter(const TString& name, const Double_t value) { // Set various parameters for the lineshape if (name == "lambda") { this->setLambda(value); std::cout << "INFO in LauPolarFormFactorNR::setResonanceParameter : Setting parameter lambda = " << this->getLambda() << std::endl; } else { std::cerr << "WARNING in LauPolarFormFactorNR::setResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl; } } void LauPolarFormFactorNR::floatResonanceParameter(const TString& name) { if (name == "lambda") { if ( lambda_->fixed() ) { lambda_->fixed( kFALSE ); this->addFloatingParameter( lambda_ ); } else { std::cerr << "WARNING in LauPolarFormFactorNR::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else { std::cerr << "WARNING in LauPolarFormFactorNR::fixResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl; } } LauParameter* LauPolarFormFactorNR::getResonanceParameter(const TString& name) { if (name == "lambda") { return lambda_; } else { std::cerr << "WARNING in LauPolarFormFactorNR::getResonanceParameter: Parameter name not reconised." << std::endl; return 0; } } void LauPolarFormFactorNR::setLambda(const Double_t lambda) { lambda_->value( lambda ); lambda_->genValue( lambda ); lambda_->initValue( lambda ); } diff --git a/src/LauRescattering2Res.cc b/src/LauRescattering2Res.cc index 0fb5cc8..da5f30c 100644 --- a/src/LauRescattering2Res.cc +++ b/src/LauRescattering2Res.cc @@ -1,658 +1,658 @@ // Copyright University of Warwick 2004 - 2014. // Distributed under the Boost Software License, Version 1.0. // (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) // Authors: // Thomas Latham // John Back // Paul Harrison /*! \file LauRescattering2Res.cc \brief File containing implementation of LauRescattering2Res class. */ #include #include "LauConstants.hh" #include "LauRescattering2Res.hh" #include "LauResonanceInfo.hh" ClassImp(LauRescattering2Res) LauRescattering2Res::LauRescattering2Res(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) : LauAbsResonance(resInfo, resPairAmpInt, daughters), B1_(0),B2_(0),B3_(0), C1_(0),C2_(0),C3_(0),C4_(0),C5_(0), D0_(0),D1_(0),D2_(0),D3_(0), F1_(0),F2_(0),F3_(0),F4_(0) { // Default values for parameters, taken from const Double_t B1Val(23.6); const Double_t B2Val(29.4); const Double_t B3Val(0.6); const Double_t C1Val(34.39); const Double_t C2Val(4.4); const Double_t C3Val(-32.9); const Double_t C4Val(-16.); const Double_t C5Val(7.4); const Double_t D0Val(0.59); const Double_t D1Val(-0.38); const Double_t D2Val(0.12); const Double_t D3Val(-0.09); const Double_t F1Val(-0.043); const Double_t F2Val(-0.008); const Double_t F3Val(-0.28); const Double_t F4Val(0.026); const TString& parNameBase = this->getSanitisedName(); TString B1Name(parNameBase); B1Name += "_B1"; B1_ = resInfo->getExtraParameter( B1Name ); if ( B1_ == 0 ) { B1_ = new LauParameter( B1Name, B1Val, -100.,100., kTRUE ); B1_->secondStage(kTRUE); resInfo->addExtraParameter( B1_ ); } TString B2Name(parNameBase); B2Name += "_B2"; B2_ = resInfo->getExtraParameter( B2Name ); if ( B2_ == 0 ) { B2_ = new LauParameter( B2Name, B2Val, -100.,100., kTRUE ); B2_->secondStage(kTRUE); resInfo->addExtraParameter( B2_ ); } TString B3Name(parNameBase); B3Name += "_B3"; B3_ = resInfo->getExtraParameter( B3Name ); if ( B3_ == 0 ) { B3_ = new LauParameter( B3Name, B3Val, -100.,100., kTRUE ); B3_->secondStage(kTRUE); resInfo->addExtraParameter( B3_ ); } TString C1Name(parNameBase); C1Name += "_C1"; C1_ = resInfo->getExtraParameter( C1Name ); if ( C1_ == 0 ) { C1_ = new LauParameter( C1Name, C1Val, -100.,100., kTRUE ); C1_->secondStage(kTRUE); resInfo->addExtraParameter( C1_ ); } TString C2Name(parNameBase); C2Name += "_C2"; C2_ = resInfo->getExtraParameter( C2Name ); if ( C2_ == 0 ) { C2_ = new LauParameter( C2Name, C2Val, -100.,100., kTRUE ); C2_->secondStage(kTRUE); resInfo->addExtraParameter( C2_ ); } TString C3Name(parNameBase); C3Name += "_C3"; C3_ = resInfo->getExtraParameter( C3Name ); if ( C3_ == 0 ) { C3_ = new LauParameter( C3Name, C3Val, -100.,100., kTRUE ); C3_->secondStage(kTRUE); resInfo->addExtraParameter( C3_ ); } TString C4Name(parNameBase); C4Name += "_C4"; C4_ = resInfo->getExtraParameter( C4Name ); if ( C4_ == 0 ) { C4_ = new LauParameter( C4Name, C4Val, -100.,100., kTRUE ); C4_->secondStage(kTRUE); resInfo->addExtraParameter( C4_ ); } TString C5Name(parNameBase); C5Name += "_C5"; C5_ = resInfo->getExtraParameter( C5Name ); if ( C5_ == 0 ) { C5_ = new LauParameter( C5Name, C5Val, -100.,100., kTRUE ); C5_->secondStage(kTRUE); resInfo->addExtraParameter( C5_ ); } TString D0Name(parNameBase); D0Name += "_D0"; D0_ = resInfo->getExtraParameter( D0Name ); if ( D0_ == 0 ) { D0_ = new LauParameter( D0Name, D0Val, -100.,100., kTRUE ); D0_->secondStage(kTRUE); resInfo->addExtraParameter( D0_ ); } TString D1Name(parNameBase); D1Name += "_D1"; D1_ = resInfo->getExtraParameter( D1Name ); if ( D1_ == 0 ) { D1_ = new LauParameter( D1Name, D1Val, -100.,100., kTRUE ); D1_->secondStage(kTRUE); resInfo->addExtraParameter( D1_ ); } TString D2Name(parNameBase); D2Name += "_D2"; D2_ = resInfo->getExtraParameter( D2Name ); if ( D2_ == 0 ) { D2_ = new LauParameter( D2Name, D2Val,-100.,100., kTRUE ); D2_->secondStage(kTRUE); resInfo->addExtraParameter( D2_ ); } TString D3Name(parNameBase); D3Name += "_D3"; D3_ = resInfo->getExtraParameter( D3Name ); if ( D3_ == 0 ) { D3_ = new LauParameter( D3Name, D3Val, -100.,100., kTRUE ); D3_->secondStage(kTRUE); resInfo->addExtraParameter( D3_ ); } TString F1Name(parNameBase); F1Name += "_F1"; F1_ = resInfo->getExtraParameter( F1Name ); if ( F1_ == 0 ) { F1_ = new LauParameter( F1Name, F1Val, -100.,100., kTRUE ); F1_->secondStage(kTRUE); resInfo->addExtraParameter( F1_ ); } TString F2Name(parNameBase); F2Name += "_F2"; F2_ = resInfo->getExtraParameter( F2Name ); if ( F2_ == 0 ) { F2_ = new LauParameter( F2Name, F2Val, -100.,100., kTRUE ); F2_->secondStage(kTRUE); resInfo->addExtraParameter( F2_ ); } TString F3Name(parNameBase); F3Name += "_F3"; F3_ = resInfo->getExtraParameter( F3Name ); if ( F3_ == 0 ) { F3_ = new LauParameter( F3Name, F3Val, -100.,100., kTRUE ); F3_->secondStage(kTRUE); resInfo->addExtraParameter( F3_ ); } TString F4Name(parNameBase); F4Name += "_F4"; F4_ = resInfo->getExtraParameter( F4Name ); if ( F4_ == 0 ) { F4_ = new LauParameter( F4Name, F4Val, -100.,100., kTRUE ); F4_->secondStage(kTRUE); resInfo->addExtraParameter( F4_ ); } } LauRescattering2Res::~LauRescattering2Res() { delete B1_;delete B2_;delete B3_; delete C1_;delete C2_;delete C3_;delete C4_;delete C5_; delete D0_;delete D1_;delete D2_;delete D3_; delete F1_;delete F2_;delete F3_;delete F4_; } void LauRescattering2Res::initialise() { sqr_tmin[1] = 2.00*LauConstants::mK; sqr_tmin[2] = 1.47; sqr_tmax[1] = 1.47; sqr_tmax[2] = 2.00; B0_ = 226.5*TMath::Pi()/180.0 + B1_->value() - B2_->value() + B3_->value(); C0_ = phi00(sqr_tmax[1]*sqr_tmax[1], 1) + C1_->value() - C2_->value() + C3_->value() - C4_->value() + C5_->value(); F0_ = g00(sqr_tmax[1]*sqr_tmax[1], 1) + F1_->value() - F2_->value() + F3_->value() - F4_->value(); std::cout << "##### B = " << this->B1_->value()<<", "<< this->B2_->value()<<", "<< this->B3_->value()<C2_->value()<<", "<< this->C3_->value()<<", "<< this->C4_->value()<<", "<< this->C5_->value()<D1_->value()<<", "<< this->D2_->value()<<", "<< this->D3_->value()<< std::endl; std::cout << "##### F = " << this->F1_->value()<<", "<< this->F2_->value()<<", "<< this->F3_->value()<<", "<< this->F4_->value()< sqr_tmax[1] && mass < sqr_tmax[2]) i = 2; if (i == 0) { std::cout << " ERROR MASS = " << mass <<" out of limmits mI, mII" << std::endl; return LauComplex(0,0); } Double_t sqr_t = mass; Double_t mag = this->g00(sqr_t, i); Double_t phase = this->phi00(sqr_t, i); return LauComplex(mag*TMath::Cos(phase), mag*TMath::Sin(phase)); } -Double_t LauRescattering2Res::pn(const Double_t x, const Double_t n) const +Double_t LauRescattering2Res::pn(const Double_t x_, const Double_t n) const { if (n==0) return 1.0; - if (n==1) return x; - return 2*x*pn(x,n-1) - pn(x,n-2); + if (n==1) return x_; + return 2*x_*pn(x_,n-1) - pn(x_,n-2); } Double_t LauRescattering2Res::x(const Double_t sqr_t, const Int_t i) const { return 2.0*(sqr_t-sqr_tmin[i])/(sqr_tmax[i]-sqr_tmin[i]) -1.0; } Double_t LauRescattering2Res::phi00(const Double_t sqr_t, const Int_t i) const { Double_t x_t = x(sqr_t, i); if (i==1) return B0_*this->pn(x_t,0)+ B1_->value()*this->pn(x_t,1)+ B2_->value()*this->pn(x_t,2)+ B3_->value()*this->pn(x_t,3); if (i==2) return C0_*this->pn(x_t,0)+ C1_->value()*this->pn(x_t,1)+ C2_->value()*this->pn(x_t,2)+ C3_->value()*this->pn(x_t,3)+ C4_->value()*this->pn(x_t,4)+ C5_->value()*this->pn(x_t,5); return 0; } Double_t LauRescattering2Res::g00(const Double_t sqr_t, const Int_t i) const { Double_t x_t = x(sqr_t, i); if (i==1) return D0_->value()*this->pn(x_t,0)+ D1_->value()*this->pn(x_t,1)+ D2_->value()*this->pn(x_t,2)+ D3_->value()*this->pn(x_t,3); if (i==2) return F0_*this->pn(x_t,0)+ F1_->value()*this->pn(x_t,1)+ F2_->value()*this->pn(x_t,2)+ F3_->value()*this->pn(x_t,3)+ F4_->value()*this->pn(x_t,4); return 0; } // TODO up to here! const std::vector& LauRescattering2Res::getFloatingParameters() { this->clearFloatingParameters(); if ( ! this->fixB1Parameter() ) { this->addFloatingParameter( B1_ ); } if ( ! this->fixB2Parameter() ) { this->addFloatingParameter( B2_ ); } if ( ! this->fixB3Parameter() ) { this->addFloatingParameter( B3_ ); } if ( ! this->fixC1Parameter() ) { this->addFloatingParameter( C1_ ); } if ( ! this->fixC2Parameter() ) { this->addFloatingParameter( C2_ ); } if ( ! this->fixC3Parameter() ) { this->addFloatingParameter( C3_ ); } if ( ! this->fixC4Parameter() ) { this->addFloatingParameter( C4_ ); } if ( ! this->fixC5Parameter() ) { this->addFloatingParameter( C5_ ); } if ( ! this->fixD0Parameter() ) { this->addFloatingParameter( D0_ ); } if ( ! this->fixD1Parameter() ) { this->addFloatingParameter( D1_ ); } if ( ! this->fixD2Parameter() ) { this->addFloatingParameter( D2_ ); } if ( ! this->fixD3Parameter() ) { this->addFloatingParameter( D3_ ); } if ( ! this->fixF1Parameter() ) { this->addFloatingParameter( F1_ ); } if ( ! this->fixF2Parameter() ) { this->addFloatingParameter( F2_ ); } if ( ! this->fixF3Parameter() ) { this->addFloatingParameter( F3_ ); } if ( ! this->fixF4Parameter() ) { this->addFloatingParameter( F4_ ); } return this->getParameters(); } void LauRescattering2Res::setResonanceParameter(const TString& name, const Double_t value) { // Set various parameters for the NRAmplitude lineshape dynamics if (name == "B1") { this->setB1Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude B1 = " << this->getB1Parameter() << std::endl; } else if (name == "B2") { this->setB2Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude B2 = " << this->getB2Parameter() << std::endl; } else if (name == "B3") { this->setB3Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude B3 = " << this->getB3Parameter() << std::endl; } else if (name == "C1") { this->setC1Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude C1 = " << this->getC1Parameter() << std::endl; } else if (name == "C2") { this->setC2Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude C2 = " << this->getC2Parameter() << std::endl; } else if (name == "C3") { this->setC3Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude C3 = " << this->getC3Parameter() << std::endl; } else if (name == "C4") { this->setC4Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude C4 = " << this->getC4Parameter() << std::endl; } else if (name == "C5") { this->setC5Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude C5 = " << this->getC5Parameter() << std::endl; } else if (name == "D0") { this->setD0Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude D0 = " << this->getD0Parameter() << std::endl; } else if (name == "D1") { this->setD1Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude D1 = " << this->getD1Parameter() << std::endl; } else if (name == "D2") { this->setD2Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude D2 = " << this->getD2Parameter() << std::endl; } else if (name == "D3") { this->setD3Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude D3 = " << this->getD3Parameter() << std::endl; } else if (name == "F1") { this->setF1Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude F1 = " << this->getF1Parameter() << std::endl; } else if (name == "F2") { this->setF2Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude F2 = " << this->getF2Parameter() << std::endl; } else if (name == "F3") { this->setF3Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude F3 = " << this->getF3Parameter() << std::endl; } else if (name == "F4") { this->setF4Parameter(value); std::cout << "INFO in LauRescattering2Res::setResonanceParameter : Setting NRAmplitude F4 = " << this->getF4Parameter() << std::endl; } else { std::cerr << "WARNING in LauRescattering2Res::setResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl; } } void LauRescattering2Res::floatResonanceParameter(const TString& name) { if (name == "B1") { if ( B1_->fixed() ) { B1_->fixed( kFALSE ); this->addFloatingParameter( B1_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "B2") { if ( B2_->fixed() ) { B2_->fixed( kFALSE ); this->addFloatingParameter( B2_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "B3") { if ( B3_->fixed() ) { B3_->fixed( kFALSE ); this->addFloatingParameter( B3_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "C1") { if ( C1_->fixed() ) { C1_->fixed( kFALSE ); this->addFloatingParameter( C1_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "C2") { if ( C2_->fixed() ) { C2_->fixed( kFALSE ); this->addFloatingParameter( C2_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "C3") { if ( C3_->fixed() ) { C3_->fixed( kFALSE ); this->addFloatingParameter( C3_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "C4") { if ( C4_->fixed() ) { C4_->fixed( kFALSE ); this->addFloatingParameter( C4_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "C5") { if ( C5_->fixed() ) { C5_->fixed( kFALSE ); this->addFloatingParameter( C5_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "D0") { if ( D0_->fixed() ) { D0_->fixed( kFALSE ); this->addFloatingParameter( D0_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "D1") { if ( D1_->fixed() ) { D1_->fixed( kFALSE ); this->addFloatingParameter( D1_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "D2") { if ( D2_->fixed() ) { D2_->fixed( kFALSE ); this->addFloatingParameter( D2_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "D3") { if ( D3_->fixed() ) { D3_->fixed( kFALSE ); this->addFloatingParameter( D3_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "F1") { if ( F1_->fixed() ) { F1_->fixed( kFALSE ); this->addFloatingParameter( F1_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "F2") { if ( F2_->fixed() ) { F2_->fixed( kFALSE ); this->addFloatingParameter( F2_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "F3") { if ( F3_->fixed() ) { F3_->fixed( kFALSE ); this->addFloatingParameter( F3_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else if (name == "4F") { if ( F4_->fixed() ) { F4_->fixed( kFALSE ); this->addFloatingParameter( F4_ ); } else { std::cerr << "WARNING in LauRescattering2Res::floatResonanceParameter: Parameter already floating. No parameter changes made." << std::endl; } } else { std::cerr << "WARNING in LauRescattering2Res::fixResonanceParameter: Parameter name not reconised. No parameter changes made." << std::endl; } } LauParameter* LauRescattering2Res::getResonanceParameter(const TString& name) { if (name == "B1") { return B1_; } else if (name == "B2") { return B2_; } else if (name == "B3") { return B3_; } else if (name == "C1") { return C1_; } else if (name == "C2") { return C2_; } else if (name == "C3") { return C3_; } else if (name == "C4") { return C4_; } else if (name == "C5") { return C5_; } else if (name == "D0") { return D0_; } else if (name == "D1") { return D1_; } else if (name == "D2") { return D2_; } else if (name == "D3") { return D3_; } else if (name == "F1") { return F1_; } else if (name == "F2") { return F2_; } else if (name == "F3") { return F3_; } else if (name == "F4") { return F4_; } else { std::cerr << "WARNING in LauRescattering2Res::getResonanceParameter: Parameter name not reconised." << std::endl; return 0; } } void LauRescattering2Res::setB1Parameter(const Double_t B1) { B1_->value( B1 ); B1_->genValue( B1 ); B1_->initValue( B1 ); } void LauRescattering2Res::setB2Parameter(const Double_t B2) { B2_->value( B2 ); B2_->genValue( B2 ); B2_->initValue( B2 ); } void LauRescattering2Res::setB3Parameter(const Double_t B3) { B3_->value( B3 ); B3_->genValue( B3 ); B3_->initValue( B3 ); } void LauRescattering2Res::setC1Parameter(const Double_t C1) { C1_->value( C1 ); C1_->genValue( C1 ); C1_->initValue( C1 ); } void LauRescattering2Res::setC2Parameter(const Double_t C2) { C2_->value( C2 ); C2_->genValue( C2 ); C2_->initValue( C2 ); } void LauRescattering2Res::setC3Parameter(const Double_t C3) { C3_->value( C3 ); C3_->genValue( C3 ); C3_->initValue( C3 ); } void LauRescattering2Res::setC4Parameter(const Double_t C4) { C4_->value( C4 ); C4_->genValue( C4 ); C4_->initValue( C4 ); } void LauRescattering2Res::setC5Parameter(const Double_t C5) { C5_->value( C5 ); C5_->genValue( C5 ); C5_->initValue( C5 ); } void LauRescattering2Res::setD0Parameter(const Double_t D0) { D0_->value( D0 ); D0_->genValue( D0 ); D0_->initValue( D0 ); } void LauRescattering2Res::setD1Parameter(const Double_t D1) { D1_->value( D1 ); D1_->genValue( D1 ); D1_->initValue( D1 ); } void LauRescattering2Res::setD2Parameter(const Double_t D2) { D2_->value( D2 ); D2_->genValue( D2 ); D2_->initValue( D2 ); } void LauRescattering2Res::setD3Parameter(const Double_t D3) { D3_->value( D3 ); D3_->genValue( D3 ); D3_->initValue( D3 ); } void LauRescattering2Res::setF1Parameter(const Double_t F1) { F1_->value( F1 ); F1_->genValue( F1 ); F1_->initValue( F1 ); } void LauRescattering2Res::setF2Parameter(const Double_t F2) { F2_->value( F2 ); F2_->genValue( F2 ); F2_->initValue( F2 ); } void LauRescattering2Res::setF3Parameter(const Double_t F3) { F3_->value( F3 ); F3_->genValue( F3 ); F3_->initValue( F3 ); } void LauRescattering2Res::setF4Parameter(const Double_t F4) { F4_->value( F4 ); F4_->genValue( F4 ); F4_->initValue( F4 ); } diff --git a/src/LauResonanceMaker.cc b/src/LauResonanceMaker.cc index 03caa26..d2e49b7 100644 --- a/src/LauResonanceMaker.cc +++ b/src/LauResonanceMaker.cc @@ -1,983 +1,990 @@ /* Copyright 2004 University of Warwick Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0 Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License. */ /* Laura++ package authors: John Back Paul Harrison Thomas Latham */ /*! \file LauResonanceMaker.cc \brief File containing implementation of LauResonanceMaker class. */ #include #include "LauAbsResonance.hh" #include "LauBelleNR.hh" #include "LauBelleSymNR.hh" #include "LauBreitWignerRes.hh" #include "LauDabbaRes.hh" #include "LauDaughters.hh" #include "LauEFKLLMRes.hh" #include "LauFlatteRes.hh" #include "LauFlatNR.hh" #include "LauGaussIncohRes.hh" #include "LauGounarisSakuraiRes.hh" #include "LauKappaRes.hh" #include "LauLASSRes.hh" #include "LauLASSBWRes.hh" #include "LauLASSNRRes.hh" #include "LauModIndPartWaveMagPhase.hh" #include "LauModIndPartWaveRealImag.hh" #include "LauNRAmplitude.hh" #include "LauRescatteringRes.hh" +#include "LauRescattering2Res.hh" #include "LauPolNR.hh" #include "LauPoleRes.hh" #include "LauPolarFormFactorNR.hh" #include "LauPolarFormFactorSymNR.hh" #include "LauRelBreitWignerRes.hh" #include "LauResonanceInfo.hh" #include "LauResonanceMaker.hh" #include "LauRhoOmegaMix.hh" #include "LauSigmaRes.hh" ClassImp(LauResonanceMaker); LauResonanceMaker* LauResonanceMaker::resonanceMaker_ = 0; LauResonanceMaker::LauResonanceMaker() : nResDefMax_(0), bwBarrierType_(LauBlattWeisskopfFactor::BWPrimeBarrier), bwRestFrame_(LauBlattWeisskopfFactor::ResonanceFrame), spinFormalism_(LauAbsResonance::Zemach_P), summaryPrinted_(kFALSE) { this->createResonanceVector(); this->setDefaultBWRadius( LauBlattWeisskopfFactor::Parent, 4.0 ); } LauResonanceMaker::~LauResonanceMaker() { for ( std::vector::iterator iter = bwIndepFactors_.begin(); iter != bwIndepFactors_.end(); ++iter ) { delete *iter; } bwIndepFactors_.clear(); for ( BWFactorCategoryMap::iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) { delete iter->second.bwFactor_; } bwFactors_.clear(); } LauResonanceMaker& LauResonanceMaker::get() { if ( resonanceMaker_ == 0 ) { resonanceMaker_ = new LauResonanceMaker(); } return *resonanceMaker_; } void LauResonanceMaker::createResonanceVector() { // Function to create all possible resonances that this class supports. // Also add in the sigma and kappa - but a special paramterisation is used // instead of the PDG "pole mass and width" values. std::cout << "INFO in LauResonanceMaker::createResonanceVector : Setting up possible resonance states..." << std::endl; LauResonanceInfo* neutral(0); LauResonanceInfo* positve(0); LauResonanceInfo* negatve(0); // Define the resonance names and store them in the array resInfo_.clear(); resInfo_.reserve(100); // rho resonances name, mass, width, spin, charge, default BW category, BW radius parameter (defaults to 4.0) // rho(770) neutral = new LauResonanceInfo("rho0(770)", 0.77526, 0.1478, 1, 0, LauBlattWeisskopfFactor::Light, 5.3); positve = new LauResonanceInfo("rho+(770)", 0.77511, 0.1491, 1, 1, LauBlattWeisskopfFactor::Light, 5.3); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // The following two lines of code are placed here in order to allow the following, rather niche, scenario: // The LauRhoOmegaMix code permits (through the use of the optional independentPar argument of LauResonanceInfo::addExtraParameter) the magnitude and phase of the rho/omega mixing to potentially differ between the decay of the parent particle to rho0 X and the parent antiparticle to rho0 Xbar. // This can be acheived by using the rho0(770) record in one case and the rho0(770)_COPY record in the other. neutral = neutral->createSharedParameterRecord("rho0(770)_COPY"); resInfo_.push_back( neutral ); // rho(1450) neutral = new LauResonanceInfo("rho0(1450)", 1.465, 0.400, 1, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("rho+(1450)", 1.465, 0.400, 1, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // rho_3(1690) neutral = new LauResonanceInfo("rho0_3(1690)", 1.686, 0.186, 3, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("rho+_3(1690)", 1.686, 0.186, 3, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // rho(1700) neutral = new LauResonanceInfo("rho0(1700)", 1.720, 0.250, 1, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("rho+(1700)", 1.720, 0.250, 1, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // rho(1900) neutral = new LauResonanceInfo("rho0(1900)", 1.909, 0.130, 1, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("rho+(1900)", 1.909, 0.130, 1, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // rho_3(1990) neutral = new LauResonanceInfo("rho0_3(1990)", 1.982, 0.188, 3, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("rho+_3(1990)", 1.982, 0.188, 3, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K* resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // K*(892) neutral = new LauResonanceInfo("K*0(892)", 0.89581, 0.0474, 1, 0, LauBlattWeisskopfFactor::Kstar, 3.0); positve = new LauResonanceInfo("K*+(892)", 0.89166, 0.0508, 1, 1, LauBlattWeisskopfFactor::Kstar, 3.0); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K*(1410) neutral = new LauResonanceInfo("K*0(1410)", 1.414, 0.232, 1, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("K*+(1410)", 1.414, 0.232, 1, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K*_0(1430) neutral = new LauResonanceInfo("K*0_0(1430)", 1.425, 0.270, 0, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("K*+_0(1430)", 1.425, 0.270, 0, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // LASS nonresonant model neutral = neutral->createSharedParameterRecord("LASSNR0"); positve = positve->createSharedParameterRecord("LASSNR+"); negatve = negatve->createSharedParameterRecord("LASSNR-"); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K*_2(1430) neutral = new LauResonanceInfo("K*0_2(1430)", 1.4324, 0.109, 2, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("K*+_2(1430)", 1.4256, 0.0985, 2, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K*(1680) neutral = new LauResonanceInfo("K*0(1680)", 1.717, 0.322, 1, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("K*+(1680)", 1.717, 0.322, 1, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // K*(1950) neutral = new LauResonanceInfo("K*0_0(1950)", 1.945, 0.201, 0, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("K*+_0(1950)", 1.945, 0.201, 0, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // phi resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // phi(1020) neutral = new LauResonanceInfo("phi(1020)", 1.019461, 0.004266, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // phi(1680) neutral = new LauResonanceInfo("phi(1680)", 1.680, 0.150, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // f_0(980) neutral = new LauResonanceInfo("f_0(980)", 0.990, 0.070, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1270) neutral = new LauResonanceInfo("f_2(1270)", 1.2751, 0.1851, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_0(1370) neutral = new LauResonanceInfo("f_0(1370)", 1.370, 0.350, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f'_0(1300), from Belle's Kspipi paper neutral = new LauResonanceInfo("f'_0(1300)", 1.449, 0.126, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1430) neutral = new LauResonanceInfo("f_2(1430)", 1.430, 0.150, 2, 0, LauBlattWeisskopfFactor::Light ); // PDG width in the range 13 - 150 resInfo_.push_back( neutral ); // f_0(1500) neutral = new LauResonanceInfo("f_0(1500)", 1.505, 0.109, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f'_2(1525) neutral = new LauResonanceInfo("f'_2(1525)", 1.525, 0.073, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1565) neutral = new LauResonanceInfo("f_2(1565)", 1.562, 0.134, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1640) neutral = new LauResonanceInfo("f_2(1640)", 1.639, 0.099, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_0(1710) neutral = new LauResonanceInfo("f_0(1710)", 1.722, 0.135, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1810) neutral = new LauResonanceInfo("f_2(1810)", 1.816, 0.197, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1910) neutral = new LauResonanceInfo("f_2(1910)", 1.903, 0.196, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(1950) neutral = new LauResonanceInfo("f_2(1950)", 1.944, 0.472, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_2(2010) neutral = new LauResonanceInfo("f_2(2010)", 2.011, 0.202, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_0(2020) neutral = new LauResonanceInfo("f_0(2020)", 1.992, 0.442, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_4(2050) neutral = new LauResonanceInfo("f_4(2050)", 2.018, 0.237, 4, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // f_0(2100) neutral = new LauResonanceInfo("f_0(2100)", 2.101, 0.224, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // omega resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // omega(782) neutral = new LauResonanceInfo("omega(782)", 0.78265, 0.00849, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // a resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // a_0(980) neutral = new LauResonanceInfo("a0_0(980)", 0.980, 0.092, 0, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("a+_0(980)", 0.980, 0.092, 0, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // a_0(1450) neutral = new LauResonanceInfo("a0_0(1450)", 1.474, 0.265, 0, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("a+_0(1450)", 1.474, 0.265, 0, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // a_2(1320) neutral = new LauResonanceInfo("a0_2(1320)", 1.3190, 0.1050, 2, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("a+_2(1320)", 1.3190, 0.1050, 2, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // charmonium resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // chi_c0 neutral = new LauResonanceInfo("chi_c0", 3.41475, 0.0105, 0, 0, LauBlattWeisskopfFactor::Charmonium ); resInfo_.push_back( neutral ); // chi_c1 neutral = new LauResonanceInfo("chi_c1", 3.51066, 0.00084, 0, 0, LauBlattWeisskopfFactor::Charmonium ); resInfo_.push_back( neutral ); // chi_c2 neutral = new LauResonanceInfo("chi_c2", 3.55620, 0.00193, 2, 0, LauBlattWeisskopfFactor::Charmonium ); resInfo_.push_back( neutral ); // X(3872) neutral = new LauResonanceInfo("X(3872)", 3.87169, 0.0012, 1, 0, LauBlattWeisskopfFactor::Charmonium ); resInfo_.push_back( neutral ); // unknown scalars name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // sigma neutral = new LauResonanceInfo("sigma0", 0.475, 0.550, 0, 0, LauBlattWeisskopfFactor::Light ); positve = new LauResonanceInfo("sigma+", 0.475, 0.550, 0, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // kappa neutral = new LauResonanceInfo("kappa0", 0.682, 0.547, 0, 0, LauBlattWeisskopfFactor::Kstar ); positve = new LauResonanceInfo("kappa+", 0.682, 0.547, 0, 1, LauBlattWeisskopfFactor::Kstar ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // dabba neutral = new LauResonanceInfo("dabba0", 2.098, 0.520, 0, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("dabba+", 2.098, 0.520, 0, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // excited charm states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // D* neutral = new LauResonanceInfo("D*0", 2.00696, 0.0021, 1, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D*+", 2.01026, 83.4e-6, 1, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D*_0 neutral = new LauResonanceInfo("D*0_0", 2.318, 0.267, 0, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D*+_0", 2.403, 0.283, 0, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D*_2 //AVERAGE--neutral = new LauResonanceInfo("D*0_2", 2.4618, 0.049, 2, 0 ); neutral = new LauResonanceInfo("D*0_2", 2.4626, 0.049, 2, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D*+_2", 2.4643, 0.037, 2, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D1(2420) neutral = new LauResonanceInfo("D0_1(2420)", 2.4214, 0.0274, 1, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D+_1(2420)", 2.4232, 0.025, 1, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D(2600) //OLD--neutral = new LauResonanceInfo("D0(2600)", 2.6087, 0.093, 0, 0 ); //OLD--positve = new LauResonanceInfo("D+(2600)", 2.6213, 0.093, 0, 1 ); neutral = new LauResonanceInfo("D0(2600)", 2.612, 0.093, 0, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D+(2600)", 2.612, 0.093, 0, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D(2760) //OLD-- neutral = new LauResonanceInfo("D0(2760)", 2.7633, 0.061, 1, 0 ); //OLD-- positve = new LauResonanceInfo("D+(2760)", 2.7697, 0.061, 1, 1 ); neutral = new LauResonanceInfo("D0(2760)", 2.761, 0.063, 1, 0, LauBlattWeisskopfFactor::Charm ); positve = new LauResonanceInfo("D+(2760)", 2.761, 0.063, 1, 1, LauBlattWeisskopfFactor::Charm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // D(2900) neutral = new LauResonanceInfo("D0(3000)", 3.0, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm ); resInfo_.push_back( neutral ); // D(3400) neutral = new LauResonanceInfo("D0(3400)", 3.4, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm ); resInfo_.push_back( neutral ); // excited strange charm name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // Ds* positve = new LauResonanceInfo("Ds*+", 2.1121, 0.0019, 1, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Ds0*(2317) positve = new LauResonanceInfo("Ds*+_0(2317)", 2.3177, 0.0038, 0, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Ds2*(2573) positve = new LauResonanceInfo("Ds*+_2(2573)", 2.5719, 0.017, 2, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Ds1*(2700) positve = new LauResonanceInfo("Ds*+_1(2700)", 2.709, 0.117, 1, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Ds1*(2860) positve = new LauResonanceInfo("Ds*+_1(2860)", 2.862, 0.180, 1, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Ds3*(2860) positve = new LauResonanceInfo("Ds*+_3(2860)", 2.862, 0.058, 3, 1, LauBlattWeisskopfFactor::StrangeCharm ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // excited bottom states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // B* neutral = new LauResonanceInfo("B*0", 5.3252, 0.00, 1, 0, LauBlattWeisskopfFactor::Beauty, 6.0); positve = new LauResonanceInfo("B*+", 5.3252, 0.00, 1, 1, LauBlattWeisskopfFactor::Beauty, 6.0); negatve = positve->createChargeConjugate(); resInfo_.push_back( neutral ); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // excited strange bottom name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // Bs* neutral = new LauResonanceInfo("Bs*0", 5.4154, 0.00, 1, 0, LauBlattWeisskopfFactor::StrangeBeauty, 6.0); resInfo_.push_back( neutral ); // nonresonant models name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0) // Phase-space nonresonant model neutral = new LauResonanceInfo("NonReson", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // Theory-based nonresonant model neutral = new LauResonanceInfo("NRModel", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // Belle nonresonant models neutral = new LauResonanceInfo("BelleSymNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("BelleNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); positve = new LauResonanceInfo("BelleNR+", 0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); neutral = new LauResonanceInfo("BelleNR_Swave", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); positve = new LauResonanceInfo("BelleNR_Swave+",0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); neutral = new LauResonanceInfo("BelleNR_Pwave", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); positve = new LauResonanceInfo("BelleNR_Pwave+",0.0, 0.0, 1, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); neutral = new LauResonanceInfo("BelleNR_Dwave", 0.0, 0.0, 2, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); positve = new LauResonanceInfo("BelleNR_Dwave+",0.0, 0.0, 2, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); neutral = new LauResonanceInfo("BelleNR_Fwave", 0.0, 0.0, 3, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); positve = new LauResonanceInfo("BelleNR_Fwave+",0.0, 0.0, 3, 1, LauBlattWeisskopfFactor::Light ); negatve = positve->createChargeConjugate(); resInfo_.push_back( positve ); resInfo_.push_back( negatve ); // Taylor expansion nonresonant model neutral = new LauResonanceInfo("NRTaylor", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // Polynomial nonresonant models neutral = new LauResonanceInfo("PolNR_S0", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolNR_S1", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolNR_S2", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolNR_P0", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolNR_P1", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolNR_P2", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // Fake resonances for S-Wave splines neutral = new LauResonanceInfo("Spline_S0", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("Spline_S0_Bar", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // Polar Form Factor nonresonant model neutral = new LauResonanceInfo("PolarFFSymNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); neutral = new LauResonanceInfo("PolarFFNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); // PiPi-KK Inelastic Scattering neutral = new LauResonanceInfo("Rescattering", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light ); resInfo_.push_back( neutral ); nResDefMax_ = resInfo_.size(); } void LauResonanceMaker::setBWType(const LauBlattWeisskopfFactor::BarrierType bwType) { // Check whether any BW factors have been created and bail out if so if ( ! bwIndepFactors_.empty() ) { std::cerr << "ERROR in LauResonanceMaker::setBWType : some barrier factors have already been created - cannot change the barrier type now!" << std::endl; return; } for ( BWFactorCategoryMap::const_iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) { if ( iter->second.bwFactor_ != 0 ) { std::cerr << "ERROR in LauResonanceMaker::setBWType : some barrier factors have already been created - cannot change the barrier type now!" << std::endl; return; } } bwBarrierType_ = bwType; } void LauResonanceMaker::setBWBachelorRestFrame(const LauBlattWeisskopfFactor::RestFrame restFrame) { // Check whether any BW factors have been created and bail out if so if ( ! bwIndepFactors_.empty() ) { std::cerr << "ERROR in LauResonanceMaker::setBWBachelorRestFrame : some barrier factors have already been created - cannot change the rest frame now!" << std::endl; return; } for ( BWFactorCategoryMap::const_iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) { if ( iter->second.bwFactor_ != 0 ) { std::cerr << "ERROR in LauResonanceMaker::setBWBachelorRestFrame : some barrier factors have already been created - cannot change the rest frame now!" << std::endl; return; } } bwRestFrame_ = restFrame; } void LauResonanceMaker::setSpinFormalism(const LauAbsResonance::LauSpinType spinType) { if ( summaryPrinted_ ) { std::cerr << "ERROR in LauResonanceMaker::setSpinFormalism : cannot redefine the spin formalism after creating one or more resonances" << std::endl; return; } spinFormalism_ = spinType; } void LauResonanceMaker::setDefaultBWRadius(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const Double_t bwRadius) { if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) { std::cerr << "WARNING in LauResonanceMaker::setDefaultBWRadius : cannot set radius values for Default or Indep categories" << std::endl; return; } // Check if we have an information object for this category BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory ); if ( factor_iter != bwFactors_.end() ) { // If so, we can set the value in the information object BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second; categoryInfo.defaultRadius_ = bwRadius; // Then we can check if a LauBlattWeisskopfFactor object has been created for this category LauBlattWeisskopfFactor* bwFactor = categoryInfo.bwFactor_; if ( bwFactor != 0 ) { // If it has then we can also set its radius value directly LauParameter* radius = bwFactor->getRadiusParameter(); radius->value(bwRadius); radius->initValue(bwRadius); radius->genValue(bwRadius); } } else { // If not then we just store the value to be used later BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory]; categoryInfo.bwFactor_ = 0; categoryInfo.defaultRadius_ = bwRadius; categoryInfo.radiusFixed_ = kTRUE; } } void LauResonanceMaker::fixBWRadius(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const Bool_t fixRadius) { if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) { std::cerr << "WARNING in LauResonanceMaker::fixBWRadius : cannot fix/float radius values for Default or Indep categories" << std::endl; return; } // Check if we have an information object for this category BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory ); if ( factor_iter != bwFactors_.end() ) { // If so, we can set the value in the information object BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second; categoryInfo.radiusFixed_ = fixRadius; // Then we can check if a LauBlattWeisskopfFactor object has been created for this category LauBlattWeisskopfFactor* bwFactor = categoryInfo.bwFactor_; if ( bwFactor != 0 ) { // If it has then we can also fix/float its radius value directly LauParameter* radius = bwFactor->getRadiusParameter(); radius->fixed(fixRadius); } } else { // If not then we just store the value to be used later BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory]; categoryInfo.bwFactor_ = 0; categoryInfo.defaultRadius_ = -1.0; categoryInfo.radiusFixed_ = fixRadius; } } LauBlattWeisskopfFactor* LauResonanceMaker::getBWFactor( const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const LauResonanceInfo* resInfo ) { LauBlattWeisskopfFactor* bwFactor(0); // If this is an independent factor, create it and add it to the list of independent factors, then return it if ( bwCategory == LauBlattWeisskopfFactor::Indep ) { bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory ); bwIndepFactors_.push_back(bwFactor); return bwFactor; } // Otherwise, look up the category in the category information map BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory ); if ( factor_iter == bwFactors_.end() ) { // If the category is currently undefined we need to create it bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory ); BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory]; categoryInfo.bwFactor_ = bwFactor; categoryInfo.defaultRadius_ = bwFactor->getRadiusParameter()->value(); categoryInfo.radiusFixed_ = kTRUE; } else { // If it exists, we can check if the factor object has been created BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second; if ( categoryInfo.bwFactor_ != 0 ) { // If so, simply clone it const UInt_t resSpin = resInfo->getSpin(); bwFactor = categoryInfo.bwFactor_->createClone( resSpin ); } else { // Otherwise we need to create it, using the default value if it has been set if ( categoryInfo.defaultRadius_ >= 0.0 ) { bwFactor = new LauBlattWeisskopfFactor( *resInfo, categoryInfo.defaultRadius_, bwBarrierType_, bwRestFrame_, bwCategory ); } else { bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory ); } categoryInfo.bwFactor_ = bwFactor; // Set whether the radius should be fixed/floated LauParameter* radius = bwFactor->getRadiusParameter(); radius->fixed( categoryInfo.radiusFixed_ ); } } return bwFactor; } LauAbsResonance* LauResonanceMaker::getResonance(const LauDaughters* daughters, const TString& resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory) { // Routine to return the appropriate LauAbsResonance object given the resonance // name (resName), which daughter is the bachelor track (resPairAmpInt = 1,2 or 3), // and the resonance type ("BW" = Breit-Wigner, "Flatte" = Flatte distribution). // If this is the first resonance we are making, first print a summary of the formalism if ( ! summaryPrinted_ ) { std::cout << "INFO in LauResonanceMaker::getResonance : Freezing amplitude formalism:" << std::endl; switch ( spinFormalism_ ) { case LauAbsResonance::Zemach_P : std::cout << " : Spin factors use Zemach spin tensors, with bachelor momentum in resonance rest frame" << std::endl; break; case LauAbsResonance::Zemach_Pstar : std::cout << " : Spin factors use Zemach spin tensors, with bachelor momentum in parent rest frame" << std::endl; break; case LauAbsResonance::Covariant : std::cout << " : Spin factors use Covariant spin tensors" << std::endl; break; case LauAbsResonance::Legendre : std::cout << " : Spin factors are just Legendre polynomials" << std::endl; break; } switch ( bwBarrierType_ ) { case LauBlattWeisskopfFactor::BWBarrier : std::cout << " : Blatt-Weisskopf barrier factors are the 'non-primed' form" << std::endl; break; case LauBlattWeisskopfFactor::BWPrimeBarrier : std::cout << " : Blatt-Weisskopf barrier factors are the 'primed' form" << std::endl; break; case LauBlattWeisskopfFactor::ExpBarrier : std::cout << " : Blatt-Weisskopf barrier factors are the exponential form" << std::endl; break; } switch ( bwRestFrame_ ) { case LauBlattWeisskopfFactor::ParentFrame : std::cout << " : Blatt-Weisskopf barrier factors use bachelor momentum in parent rest frame" << std::endl; break; case LauBlattWeisskopfFactor::ResonanceFrame : std::cout << " : Blatt-Weisskopf barrier factors use bachelor momentum in resonance rest frame" << std::endl; break; case LauBlattWeisskopfFactor::Covariant : std::cout << " : Blatt-Weisskopf barrier factors use covariant expression" << std::endl; break; } summaryPrinted_ = kTRUE; } // Loop over all possible resonance states we have defined in // createResonanceVector() until we get a match with the name of the resonance LauResonanceInfo* resInfo(0); for (std::vector::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) { if (resName == (*iter)->getName()) { // We have recognised the resonance name. std::cout<<"INFO in LauResonanceMaker::getResonance : Creating resonance: "<getBWCategory(); } LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo ); LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo ); theResonance->setBarrierRadii( resBWFactor, parBWFactor ); break; } case LauAbsResonance::GS : { // Gounaris-Sakurai function to try and model the rho(770) better std::cout<<" : Using Gounaris-Sakurai lineshape. "<getBWCategory(); } LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo ); LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo ); theResonance->setBarrierRadii( resBWFactor, parBWFactor ); break; } case LauAbsResonance::Flatte : // Flatte distribution - coupled channel Breit-Wigner std::cout<<" : Using Flatte lineshape. "<getBWCategory(); } LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo ); LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo ); theResonance->setBarrierRadii( resBWFactor, parBWFactor ); break; } // Set the spin formalism choice theResonance->setSpinType( spinFormalism_ ); return theResonance; } Int_t LauResonanceMaker::resTypeInt(const TString& name) const { // Internal function that returns the resonance integer, specified by the // order of the resonance vector defined in createResonanceVector(), // for a given resonance name. Int_t resTypInt(-99); Int_t i(0); for (std::vector::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) { if (name.BeginsWith((*iter)->getName(), TString::kExact) == kTRUE) { // We have recognised the resonance from those that are available resTypInt = i; return resTypInt; } ++i; } return resTypInt; } void LauResonanceMaker::printAll( std::ostream& stream ) const { for ( std::vector::const_iterator iter = resInfo_.begin(); iter != resInfo_.end(); ++iter ) { stream << (**iter) << std::endl; } } LauResonanceInfo* LauResonanceMaker::getResInfo(const TString& resName) const { LauResonanceInfo* resInfo(0); for (std::vector::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) { if (resName == (*iter)->getName()) { // We have recognised the resonance name. resInfo = (*iter); // stop looping break; } } return resInfo; }