Page MenuHomeHEPForge
Differential D41 Diff 206 src/LauRelBreitWignerRes.cc

Changeset View

Standalone View

View Options

src/LauRelBreitWignerRes.cc

Show All 31 Lines
#include "LauRelBreitWignerRes.hh" #include "LauRelBreitWignerRes.hh"
ClassImp(LauRelBreitWignerRes) ClassImp(LauRelBreitWignerRes)
LauRelBreitWignerRes::LauRelBreitWignerRes(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) : LauRelBreitWignerRes::LauRelBreitWignerRes(LauResonanceInfo* resInfo, const Int_t resPairAmpInt, const LauDaughters* daughters) :
LauAbsResonance(resInfo, resPairAmpInt, daughters), LauAbsResonance(resInfo, resPairAmpInt, daughters),
q0_(0.0), q0_(0.0),
p0_(0.0),
pstar0_(0.0),
erm0_(0.0),
resMass_(0.0), resMass_(0.0),
resMassSq_(0.0), resMassSq_(0.0),
resWidth_(0.0), resWidth_(0.0),
resRadius_(0.0), resRadius_(0.0),
parRadius_(0.0), parRadius_(0.0),
mDaugSum_(0.0), mDaugSum_(0.0),
mDaugSumSq_(0.0), mDaugSumSq_(0.0),
mDaugDiff_(0.0), mDaugDiff_(0.0),
mDaugDiffSq_(0.0), mDaugDiffSq_(0.0),
mParentSq_(0.0), mParentSq_(0.0),
mBachSq_(0.0), mBachSq_(0.0),
FR0_(1.0), FR0_(1.0)
FP0_(1.0)
{ {
} }
LauRelBreitWignerRes::~LauRelBreitWignerRes() LauRelBreitWignerRes::~LauRelBreitWignerRes()
{ {
} }
void LauRelBreitWignerRes::initialise() void LauRelBreitWignerRes::initialise()
Show All 19 Linesvoid LauRelBreitWignerRes::initialise()
mDaugDiffSq_ = mDaugDiff_*mDaugDiff_; mDaugDiffSq_ = mDaugDiff_*mDaugDiff_;
mParentSq_ = massParent*massParent; mParentSq_ = massParent*massParent;
mBachSq_ = massBachelor*massBachelor; mBachSq_ = massBachelor*massBachelor;
// Create an effective resonance pole mass to protect against resonances // Create an effective resonance pole mass to protect against resonances
// that are below threshold // that are below threshold
Double_t effResMass = resMass_; Double_t effResMass = resMass_;
Double_t effResMassSq = resMassSq_; Double_t effResMassSq = resMassSq_;
if (resMassSq_ - mDaugSumSq_ < 0.0 || resMass_ > massParent - massBachelor){ if ( resMassSq_ - mDaugSumSq_ < 0.0 ) {
Double_t minMass = mDaugSum_; Double_t minMass = mDaugSum_;
Double_t maxMass = massParent - massBachelor; Double_t maxMass = massParent - massBachelor;
Double_t tanhTerm = std::tanh( (resMass_ - ((minMass + maxMass)/2))/(maxMass-minMass)); Double_t tanhTerm = std::tanh( (resMass_ - ((minMass + maxMass)/2))/(maxMass-minMass));
effResMass = minMass + (maxMass-minMass)*(1+tanhTerm)/2; effResMass = minMass + (maxMass-minMass)*(1+tanhTerm)/2;
effResMassSq = effResMass*effResMass; effResMassSq = effResMass*effResMass;
} }
// Decay momentum of either daughter in the resonance rest frame // Decay momentum of either daughter in the resonance rest frame
// when resonance mass = rest-mass value, m_0 (PDG value) // when resonance mass = rest-mass value, m_0 (PDG value)
Double_t term1 = effResMassSq - mDaugSumSq_; Double_t term1 = effResMassSq - mDaugSumSq_;
Double_t term2 = effResMassSq - mDaugDiffSq_; Double_t term2 = effResMassSq - mDaugDiffSq_;
Double_t term12 = term1*term2; Double_t term12 = term1*term2;
if (term12 > 0.0) { if (term12 > 0.0) {
q0_ = TMath::Sqrt(term12)/(2.0*effResMass); q0_ = TMath::Sqrt(term12)/(2.0*effResMass);
} else { } else {
q0_ = 0.0; q0_ = 0.0;
} }
// Momentum of the bachelor particle in the resonance rest frame // Calculate the resonance Blatt-Weisskopf form factor for the case when m = m_0
// when resonance mass = rest-mass value, m_0 (PDG value)
Double_t eBach = (mParentSq_ - effResMassSq - mBachSq_)/(2.0*effResMass);
Double_t termBach = eBach*eBach - mBachSq_;
if ( eBach<0.0 || termBach<0.0 ) {
p0_ = 0.0;
} else {
p0_ = TMath::Sqrt( termBach );
}
// Momentum of the bachelor particle in the parent rest frame
// when resonance mass = rest-mass value, m_0 (PDG value)
Double_t eStarBach = (mParentSq_ + mBachSq_ - effResMassSq)/(2.0*massParent);
Double_t termStarBach = eStarBach*eStarBach - mBachSq_;
if ( eStarBach<0.0 || termStarBach<0.0 ) {
pstar0_ = 0.0;
} else {
pstar0_ = TMath::Sqrt( termStarBach );
}
// Covariant factor when resonance mass = rest-mass value, m_0 (PDF value)
erm0_ = (mParentSq_ + effResMassSq - mBachSq_)/(2.0*massParent*effResMass);
this->calcCovFactor( erm0_ );
// Calculate the Blatt-Weisskopf form factor for the case when m = m_0
FR0_ = 1.0; FR0_ = 1.0;
FP0_ = 1.0; if ( this->getSpin() > 0 ) {
const Int_t resSpin = this->getSpin();
if ( resSpin > 0 ) {
const LauBlattWeisskopfFactor* resBWFactor = this->getResBWFactor(); const LauBlattWeisskopfFactor* resBWFactor = this->getResBWFactor();
const LauBlattWeisskopfFactor* parBWFactor = this->getParBWFactor(); if ( resBWFactor != nullptr ) {
FR0_ = (resBWFactor!=0) ? resBWFactor->calcFormFactor(q0_) : 1.0; FR0_ = resBWFactor->calcFormFactor(q0_);
switch ( parBWFactor->getRestFrame() ) {
case LauBlattWeisskopfFactor::ResonanceFrame:
FP0_ = (parBWFactor!=0) ? parBWFactor->calcFormFactor(p0_) : 1.0;
break;
case LauBlattWeisskopfFactor::ParentFrame:
FP0_ = (parBWFactor!=0) ? parBWFactor->calcFormFactor(pstar0_) : 1.0;
break;
case LauBlattWeisskopfFactor::Covariant:
{
Double_t covFactor = this->getCovFactor();
if ( resSpin > 2 ) {
covFactor = TMath::Power( covFactor, 1.0/resSpin );
} else if ( resSpin == 2 ) {
covFactor = TMath::Sqrt( covFactor );
}
FP0_ = (parBWFactor!=0) ? parBWFactor->calcFormFactor(pstar0_*covFactor) : 1.0;
break;
}
} }
} }
} }
LauComplex LauRelBreitWignerRes::resAmp(Double_t mass, Double_t spinTerm) LauComplex LauRelBreitWignerRes::resAmp(Double_t mass, Double_t spinTerm)
{ {
// This function returns the complex dynamical amplitude for a Breit-Wigner resonance, // This function returns the complex dynamical amplitude for a Breit-Wigner resonance,
// given the invariant mass and cos(helicity) values. // given the invariant mass and cos(helicity) values.
Show All 35 LinesLauComplex LauRelBreitWignerRes::resAmp(Double_t mass, Double_t spinTerm)
const Double_t p = this->getP(); const Double_t p = this->getP();
const Double_t pstar = this->getPstar(); const Double_t pstar = this->getPstar();
// Get barrier scaling factors // Get barrier scaling factors
Double_t fFactorR(1.0); Double_t fFactorR(1.0);
Double_t fFactorB(1.0); Double_t fFactorB(1.0);
if ( resSpin > 0 ) { if ( resSpin > 0 ) {
const LauBlattWeisskopfFactor* resBWFactor = this->getResBWFactor(); const LauBlattWeisskopfFactor* resBWFactor = this->getResBWFactor();
if ( resBWFactor != nullptr ) {
fFactorR = resBWFactor->calcFormFactor(q);
}
const LauBlattWeisskopfFactor* parBWFactor = this->getParBWFactor(); const LauBlattWeisskopfFactor* parBWFactor = this->getParBWFactor();
fFactorR = (resBWFactor!=0) ? resBWFactor->calcFormFactor(q) : 1.0; if ( parBWFactor != nullptr ) {
switch ( parBWFactor->getRestFrame() ) { switch ( parBWFactor->getRestFrame() ) {
case LauBlattWeisskopfFactor::ResonanceFrame: case LauBlattWeisskopfFactor::ResonanceFrame:
fFactorB = (parBWFactor!=0) ? parBWFactor->calcFormFactor(p) : 1.0; fFactorB = parBWFactor->calcFormFactor(p);
break; break;
case LauBlattWeisskopfFactor::ParentFrame: case LauBlattWeisskopfFactor::ParentFrame:
fFactorB = (parBWFactor!=0) ? parBWFactor->calcFormFactor(pstar) : 1.0; fFactorB = parBWFactor->calcFormFactor(pstar);
break; break;
case LauBlattWeisskopfFactor::Covariant: case LauBlattWeisskopfFactor::Covariant:
{ {
Double_t covFactor = this->getCovFactor(); Double_t covFactor = this->getCovFactor();
if ( resSpin > 2 ) { if ( resSpin > 2 ) {
covFactor = TMath::Power( covFactor, 1.0/resSpin ); covFactor = TMath::Power( covFactor, 1.0/resSpin );
} else if ( resSpin == 2 ) { } else if ( resSpin == 2 ) {
covFactor = TMath::Sqrt( covFactor ); covFactor = TMath::Sqrt( covFactor );
} }
fFactorB = (parBWFactor!=0) ? parBWFactor->calcFormFactor(pstar*covFactor) : 1.0; fFactorB = parBWFactor->calcFormFactor(pstar*covFactor);
break; break;
} }
} }
} }
const Double_t fFactorRRatio = fFactorR/FR0_; }
const Double_t fFactorBRatio = fFactorB/FP0_;
// If ignoreMomenta is set, set the total width simply as the pole width, and do // If ignoreMomenta is set, set the total width simply as the pole width, and do
// not include any momentum-dependent barrier factors (set them to unity) // not include any momentum-dependent barrier factors (set them to unity)
Double_t totWidth(resWidth); Double_t totWidth(resWidth);
if (!this->ignoreMomenta()) { if (!this->ignoreMomenta()) {
const Double_t qRatio = q/q0_; const Double_t qRatio = q/q0_;
Double_t qTerm(0.0); Double_t qTerm(0.0);
if (resSpin == 0) { if (resSpin == 0) {
qTerm = qRatio; qTerm = qRatio;
} else if (resSpin == 1) { } else if (resSpin == 1) {
qTerm = qRatio*qRatio*qRatio; qTerm = qRatio*qRatio*qRatio;
} else if (resSpin == 2) { } else if (resSpin == 2) {
qTerm = qRatio*qRatio*qRatio*qRatio*qRatio; qTerm = qRatio*qRatio*qRatio*qRatio*qRatio;
} else { } else {
qTerm = TMath::Power(qRatio, 2.0*resSpin + 1.0); qTerm = TMath::Power(qRatio, 2.0*resSpin + 1.0);
} }
const Double_t fFactorRRatio = fFactorR/FR0_;
totWidth = resWidth*qTerm*(resMass/mass)*fFactorRRatio*fFactorRRatio; totWidth = resWidth*qTerm*(resMass/mass)*fFactorRRatio*fFactorRRatio;
} }
const Double_t massSq = mass*mass; const Double_t massSq = mass*mass;
const Double_t massSqTerm = resMassSq_ - massSq; const Double_t massSqTerm = resMassSq_ - massSq;
// Compute the complex amplitude // Compute the complex amplitude
resAmplitude = LauComplex(massSqTerm, resMass*totWidth); resAmplitude = LauComplex(massSqTerm, resMass*totWidth);
// Scale by the denominator factor, as well as the spin term // Scale by the denominator factor, as well as the spin term
Double_t scale = spinTerm/(massSqTerm*massSqTerm + resMassSq_*totWidth*totWidth); Double_t scale = spinTerm/(massSqTerm*massSqTerm + resMassSq_*totWidth*totWidth);
// Include Blatt-Weisskopf barrier factors // Include Blatt-Weisskopf barrier factors
if (!this->ignoreBarrierScaling()) { if (!this->ignoreBarrierScaling()) {
scale *= fFactorRRatio*fFactorBRatio; scale *= fFactorR * fFactorB;
} }
resAmplitude.rescale(scale); resAmplitude.rescale(scale);
return resAmplitude; return resAmplitude;
} }
const std::vector<LauParameter*>& LauRelBreitWignerRes::getFloatingParameters() const std::vector<LauParameter*>& LauRelBreitWignerRes::getFloatingParameters()
Show All 19 Lines