Page MenuHomeHEPForge
Files F11222447

No OneTemporary
Actions

View Options

diff --git a/src/LauAbsCoeffSet.cc b/src/LauAbsCoeffSet.cc
index 093f552..da9acbe 100644
--- a/src/LauAbsCoeffSet.cc
+++ b/src/LauAbsCoeffSet.cc
@@ -1,366 +1,368 @@
/*
Copyright 2006 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 LauAbsCoeffSet.cc
\brief File containing implementation of LauAbsCoeffSet class.
*/
#include <fstream>
#include <iostream>
#include "TString.h"
#include "LauAbsCoeffSet.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauRandom.hh"
ClassImp(LauAbsCoeffSet);
TRandom* LauAbsCoeffSet::randomiser_ = nullptr;
Double_t LauAbsCoeffSet::minMagnitude_ = -10.0;
Double_t LauAbsCoeffSet::maxMagnitude_ = 10.0;
Double_t LauAbsCoeffSet::minPhase_ = -LauConstants::threePi;
Double_t LauAbsCoeffSet::maxPhase_ = LauConstants::threePi;
Double_t LauAbsCoeffSet::minRealImagPart_ = -10.0;
Double_t LauAbsCoeffSet::maxRealImagPart_ = 10.0;
Double_t LauAbsCoeffSet::minDelta_ = -2.0;
Double_t LauAbsCoeffSet::maxDelta_ = 2.0;
LauAbsCoeffSet::LauAbsCoeffSet(const TString& theName, const TString& theBaseName) :
name_{theName},
basename_{theBaseName}
{
}
TRandom* LauAbsCoeffSet::getRandomiser()
{
if ( randomiser_ == nullptr ) {
randomiser_ = LauRandom::zeroSeedRandom();
}
return randomiser_;
}
void LauAbsCoeffSet::index(const UInt_t newIndex)
{
index_ = newIndex;
const TString oldBaseName{ this->baseName() };
TString basename{ oldBaseName };
basename += newIndex;
basename += "_";
this->baseName(basename);
std::vector<LauParameter*> pars { this->getParameters() };
for ( LauParameter* par : pars ) {
this->adjustName( *par, oldBaseName );
}
}
void LauAbsCoeffSet::adjustName(LauParameter& par, const TString& oldBaseName)
{
TString theName{ par.name() };
if ( theName.BeginsWith( oldBaseName ) && theName != oldBaseName ) {
theName.Remove(0,oldBaseName.Length());
}
theName.Prepend(this->baseName());
par.name(theName);
}
void LauAbsCoeffSet::setParameterValue(const TString& parName, const Double_t value, const Bool_t init)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::setParameterValue : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->value( value );
if ( init ) {
par->genValue( value );
par->initValue( value );
}
}
void LauAbsCoeffSet::setParameterError(const TString& parName, const Double_t error)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::setParameterError : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->error( error );
}
void LauAbsCoeffSet::fixParameter(const TString& parName)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::fixParameter : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->fixed( kTRUE );
}
void LauAbsCoeffSet::floatParameter(const TString& parName)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::floatParameter : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->fixed( kFALSE );
}
void LauAbsCoeffSet::blindParameter(const TString& parName, const TString& blindingString, const Double_t width)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::blindParameter : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->blindParameter( blindingString, width );
}
void LauAbsCoeffSet::addGaussianConstraint(const TString& parName, const Double_t mean, const Double_t width)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::addGaussianConstraint : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
par->addGaussianConstraint( mean, width );
}
void LauAbsCoeffSet::addSuffixToParameterName(const TString& parName, const TString& suffix)
{
LauParameter* par { this->findParameter( parName ) };
if ( par == nullptr ) {
std::cerr << "ERROR in LauAbsCoeffSet::addSuffixToParameterName : Unable to find parameter \"" << parName << "\"" << std::endl;
return;
}
TString newName{ par->name() };
if ( ! suffix.BeginsWith('_') ) {
newName += "_";
}
newName += suffix;
par->name( newName );
}
LauParameter* LauAbsCoeffSet::findParameter(const TString& parName)
{
std::vector<LauParameter*> pars { this->getParameters() };
for ( LauParameter* par : pars ) {
const TString& iName { par->name() };
if ( iName.EndsWith( parName ) ) {
return par;
}
}
return nullptr;
}
#include "LauBelleCPCoeffSet.hh"
#include "LauCartesianCPCoeffSet.hh"
#include "LauCartesianGammaCPCoeffSet.hh"
#include "LauCleoCPCoeffSet.hh"
#include "LauMagPhaseCoeffSet.hh"
#include "LauMagPhaseCPCoeffSet.hh"
#include "LauNSCCartesianCPCoeffSet.hh"
#include "LauPolarGammaCPCoeffSet.hh"
#include "LauRealImagCoeffSet.hh"
#include "LauRealImagCPCoeffSet.hh"
#include "LauRealImagGammaCPCoeffSet.hh"
+//! \cond DOXYGEN_IGNORE
void nlohmann::adl_serializer<LauAbsCoeffSet>::to_json( json& j, const LauAbsCoeffSet& t )
{
switch ( t.type() ) {
case LauCoeffType::MagPhase :
nlohmann::adl_serializer<LauMagPhaseCoeffSet>::to_json( j, static_cast<const LauMagPhaseCoeffSet&>(t));
break;
case LauCoeffType::RealImag :
nlohmann::adl_serializer<LauRealImagCoeffSet>::to_json( j, static_cast<const LauRealImagCoeffSet&>(t));
break;
case LauCoeffType::BelleCP :
nlohmann::adl_serializer<LauBelleCPCoeffSet>::to_json( j, static_cast<const LauBelleCPCoeffSet&>(t));
break;
case LauCoeffType::CartesianCP :
nlohmann::adl_serializer<LauCartesianCPCoeffSet>::to_json( j, static_cast<const LauCartesianCPCoeffSet&>(t));
break;
case LauCoeffType::CartesianGammaCP :
nlohmann::adl_serializer<LauCartesianGammaCPCoeffSet>::to_json( j, static_cast<const LauCartesianGammaCPCoeffSet&>(t));
break;
case LauCoeffType::CleoCP :
nlohmann::adl_serializer<LauCleoCPCoeffSet>::to_json( j, static_cast<const LauCleoCPCoeffSet&>(t));
break;
case LauCoeffType::MagPhaseCP :
nlohmann::adl_serializer<LauMagPhaseCPCoeffSet>::to_json( j, static_cast<const LauMagPhaseCPCoeffSet&>(t));
break;
case LauCoeffType::NSCCartesianCP :
nlohmann::adl_serializer<LauNSCCartesianCPCoeffSet>::to_json( j, static_cast<const LauNSCCartesianCPCoeffSet&>(t));
break;
case LauCoeffType::PolarGammaCP :
nlohmann::adl_serializer<LauPolarGammaCPCoeffSet>::to_json( j, static_cast<const LauPolarGammaCPCoeffSet&>(t));
break;
case LauCoeffType::RealImagCP :
nlohmann::adl_serializer<LauRealImagCPCoeffSet>::to_json( j, static_cast<const LauRealImagCPCoeffSet&>(t));
break;
case LauCoeffType::RealImagGammaCP :
nlohmann::adl_serializer<LauRealImagGammaCPCoeffSet>::to_json( j, static_cast<const LauRealImagGammaCPCoeffSet&>(t));
break;
}
}
+//! \endcond DOXYGEN_IGNORE
std::vector<std::unique_ptr<LauAbsCoeffSet>> LauAbsCoeffSet::readFromJson( const TString& fileName )
{
using nlohmann::json;
std::ifstream in(fileName, std::ios_base::in);
if ( ! in ) {
std::cerr << "ERROR in LauAbsCoeffSet::readFromJson : couldn't open file \"" << fileName << "\"" << std::endl;
return {};
}
json j;
in >> j;
const auto nCoeffs { j.at("nCoeffs").get<std::size_t>() };
std::vector<std::unique_ptr<LauAbsCoeffSet>> coeffs;
coeffs.reserve( nCoeffs );
for ( auto& coeff : j.at("coeffs") ) {
const LauCoeffType type { coeff.at("type").get<LauCoeffType>() };
switch ( type ) {
case LauCoeffType::MagPhase :
coeffs.emplace_back( std::make_unique<LauMagPhaseCoeffSet>( coeff.get<LauMagPhaseCoeffSet>() ) );
break;
case LauCoeffType::RealImag :
coeffs.emplace_back( std::make_unique<LauRealImagCoeffSet>( coeff.get<LauRealImagCoeffSet>() ) );
break;
case LauCoeffType::BelleCP :
coeffs.emplace_back( std::make_unique<LauBelleCPCoeffSet>( coeff.get<LauBelleCPCoeffSet>() ) );
break;
case LauCoeffType::CartesianCP :
coeffs.emplace_back( std::make_unique<LauCartesianCPCoeffSet>( coeff.get<LauCartesianCPCoeffSet>() ) );
break;
case LauCoeffType::CartesianGammaCP :
coeffs.emplace_back( std::make_unique<LauCartesianGammaCPCoeffSet>( coeff.get<LauCartesianGammaCPCoeffSet>() ) );
break;
case LauCoeffType::CleoCP :
coeffs.emplace_back( std::make_unique<LauCleoCPCoeffSet>( coeff.get<LauCleoCPCoeffSet>() ) );
break;
case LauCoeffType::MagPhaseCP :
coeffs.emplace_back( std::make_unique<LauMagPhaseCPCoeffSet>( coeff.get<LauMagPhaseCPCoeffSet>() ) );
break;
case LauCoeffType::NSCCartesianCP :
coeffs.emplace_back( std::make_unique<LauNSCCartesianCPCoeffSet>( coeff.get<LauNSCCartesianCPCoeffSet>() ) );
break;
case LauCoeffType::PolarGammaCP :
coeffs.emplace_back( std::make_unique<LauPolarGammaCPCoeffSet>( coeff.get<LauPolarGammaCPCoeffSet>() ) );
break;
case LauCoeffType::RealImagCP :
coeffs.emplace_back( std::make_unique<LauRealImagCPCoeffSet>( coeff.get<LauRealImagCPCoeffSet>() ) );
break;
case LauCoeffType::RealImagGammaCP :
coeffs.emplace_back( std::make_unique<LauRealImagGammaCPCoeffSet>( coeff.get<LauRealImagGammaCPCoeffSet>() ) );
break;
}
}
return coeffs;
}
void LauAbsCoeffSet::writeToJson( const TString& fileName, const std::vector<std::unique_ptr<LauAbsCoeffSet>>& coeffs )
{
using nlohmann::json;
json j;
j["nCoeffs"] = coeffs.size();
j["coeffs"] = json::array();
for ( auto& coeffset : coeffs ) {
j["coeffs"].push_back( *coeffset );
}
std::ofstream out{fileName, std::ios_base::out};
if ( ! out ) {
std::cerr << "ERROR in LauAbsCoeffSet::writeToJson : couldn't open file \"" << fileName << "\" for writing. No file will be written!" << std::endl;
return;
}
out << j.dump(4);
out << std::endl;
}
std::ostream& operator<<( std::ostream& os, const LauCoeffType type )
{
switch ( type ) {
case LauCoeffType::MagPhase :
os << "MagPhase";
break;
case LauCoeffType::RealImag :
os << "RealImag";
break;
case LauCoeffType::BelleCP :
os << "BelleCP";
break;
case LauCoeffType::CartesianCP :
os << "CartesianCP";
break;
case LauCoeffType::CartesianGammaCP :
os << "CartesianGammaCP";
break;
case LauCoeffType::CleoCP :
os << "CleoCP";
break;
case LauCoeffType::MagPhaseCP :
os << "MagPhaseCP";
break;
case LauCoeffType::NSCCartesianCP :
os << "NSCCartesianCP";
break;
case LauCoeffType::PolarGammaCP :
os << "PolarGammaCP";
break;
case LauCoeffType::RealImagCP :
os << "RealImagCP";
break;
case LauCoeffType::RealImagGammaCP :
os << "RealImagGammaCP";
break;
}
return os;
}
diff --git a/src/LauBelleCPCoeffSet.cc b/src/LauBelleCPCoeffSet.cc
index 8b0740f..45ffc23 100644
--- a/src/LauBelleCPCoeffSet.cc
+++ b/src/LauBelleCPCoeffSet.cc
@@ -1,386 +1,388 @@
/*
Copyright 2006 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 LauBelleCPCoeffSet.cc
\brief File containing implementation of LauBelleCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauBelleCPCoeffSet.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauBelleCPCoeffSet)
LauBelleCPCoeffSet::LauBelleCPCoeffSet(const TString& compName, const Double_t a, const Double_t delta, const Double_t b, const Double_t phi,
const Bool_t aFixed, const Bool_t deltaFixed, const Bool_t bFixed, const Bool_t phiFixed,
const Bool_t bSecondStage, const Bool_t phiSecondStage) :
LauAbsCoeffSet{ compName },
a_{ std::make_unique<LauParameter>("A", a, minMagnitude_, maxMagnitude_, aFixed) },
b_{ std::make_unique<LauParameter>("B", b, minMagnitude_, maxMagnitude_, bFixed) },
delta_{ std::make_unique<LauParameter>("Delta", delta, minPhase_, maxPhase_, deltaFixed) },
phi_{ std::make_unique<LauParameter>("Phi", phi, minPhase_, maxPhase_, phiFixed) },
acp_{ "ACP", (-2.0*b*TMath::Cos(phi))/(1.0+b*b), -1.0, 1.0, bFixed&&phiFixed }
{
if (bSecondStage && !bFixed) {
b_->secondStage(kTRUE);
b_->initValue(0.0);
}
if (phiSecondStage && !phiFixed) {
phi_->secondStage(kTRUE);
phi_->initValue(0.0);
}
}
LauBelleCPCoeffSet::LauBelleCPCoeffSet(const LauBelleCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) {
a_.reset( rhs.a_->createClone(constFactor) );
} else {
a_ = std::make_unique<LauParameter>("A", rhs.a_->value(), minMagnitude_, maxMagnitude_, rhs.a_->fixed());
if ( rhs.a_->blind() ) {
const LauBlind* blinder { rhs.a_->blinder() };
a_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
b_.reset( rhs.b_->createClone(constFactor) );
} else {
b_ = std::make_unique<LauParameter>("B", rhs.b_->value(), minMagnitude_, maxMagnitude_, rhs.b_->fixed());
if ( rhs.b_->blind() ) {
const LauBlind* blinder { rhs.b_->blinder() };
b_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase ) {
delta_.reset( rhs.delta_->createClone(constFactor) );
} else {
delta_ = std::make_unique<LauParameter>("Delta", rhs.delta_->value(), minPhase_, maxPhase_, rhs.delta_->fixed());
if ( rhs.delta_->blind() ) {
const LauBlind* blinder { rhs.delta_->blinder() };
delta_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
phi_.reset( rhs.phi_->createClone(constFactor) );
} else {
phi_ = std::make_unique<LauParameter>("Phi", rhs.phi_->value(), minPhase_, maxPhase_, rhs.phi_->fixed());
if ( rhs.phi_->blind() ) {
const LauBlind* blinder { rhs.phi_->blinder() };
phi_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
void LauBelleCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauBelleCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"a-magnitude = "<<a_->value()<<",\t";
std::cout<<"delta = "<<delta_->value()<<",\t";
std::cout<<"b-magnitude = "<<b_->value()<<",\t";
std::cout<<"phi = "<<phi_->value()<<"."<<std::endl;
}
void LauBelleCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & a-Magnitude & delta & b-Magnitude & phi \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauBelleCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, a_->value());
stream<<" \\pm ";
print.printFormat(stream, a_->error());
stream<<"$ & $";
print.printFormat(stream, delta_->value());
stream<<" \\pm ";
print.printFormat(stream, delta_->error());
stream<<"$ & $";
print.printFormat(stream, b_->value());
stream<<" \\pm ";
print.printFormat(stream, b_->error());
stream<<"$ & $";
print.printFormat(stream, phi_->value());
stream<<" \\pm ";
print.printFormat(stream, phi_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauBelleCPCoeffSet::randomiseInitValues()
{
if (a_->fixed() == kFALSE) {
// Choose an a-magnitude between 0.0 and 2.0
const Double_t mag { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
a_->initValue(mag); a_->value(mag);
}
if (b_->fixed() == kFALSE && b_->secondStage() == kFALSE) {
// Choose a b-magnitude between 0.0 and 0.1
const Double_t mag { LauAbsCoeffSet::getRandomiser()->Rndm()*0.1 };
b_->initValue(mag); b_->value(mag);
}
if (delta_->fixed() == kFALSE) {
// Choose a phase between +- pi
const Double_t phase { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
delta_->initValue(phase); delta_->value(phase);
}
if (phi_->fixed() == kFALSE && phi_->secondStage() == kFALSE) {
// Choose a phase between +- pi
const Double_t phase { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
phi_->initValue(phase); phi_->value(phase);
}
}
void LauBelleCPCoeffSet::finaliseValues()
{
// retrieve the current values from the parameters
Double_t aVal { a_->value() };
Double_t bVal { b_->value() };
Double_t deltaVal { delta_->value() };
Double_t phiVal { phi_->value() };
// Check whether we have a negative "a" magnitude.
// If so make it positive and add pi to the "delta" phase.
if (aVal < 0.0) {
aVal *= -1.0;
deltaVal += LauConstants::pi;
}
// Check whether we have a negative "b" magnitude.
// If so make it positive and add pi to the "phi" phase.
if (bVal < 0.0) {
bVal *= -1.0;
phiVal += LauConstants::pi;
}
// Check now whether the phases lies in the right range (-pi to pi).
Bool_t deltaWithinRange{kFALSE};
Bool_t phiWithinRange{kFALSE};
while (deltaWithinRange == kFALSE && phiWithinRange == kFALSE) {
if (deltaVal > -LauConstants::pi && deltaVal < LauConstants::pi) {
deltaWithinRange = kTRUE;
} else {
// Not within the specified range
if (deltaVal > LauConstants::pi) {
deltaVal -= LauConstants::twoPi;
} else if (deltaVal < -LauConstants::pi) {
deltaVal += LauConstants::twoPi;
}
}
if (phiVal > -LauConstants::pi && phiVal < LauConstants::pi) {
phiWithinRange = kTRUE;
} else {
// Not within the specified range
if (phiVal > LauConstants::pi) {
phiVal -= LauConstants::twoPi;
} else if (phiVal < -LauConstants::pi) {
phiVal += LauConstants::twoPi;
}
}
}
// A further problem can occur when the generated phase is close to -pi or pi.
// The phase can wrap over to the other end of the scale -
// this leads to artificially large pulls so we wrap it back.
const Double_t genDelta { delta_->genValue() };
const Double_t genPhi { phi_->genValue() };
Double_t diff { deltaVal - genDelta };
if (diff > LauConstants::pi) {
deltaVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
deltaVal += LauConstants::twoPi;
}
diff = phiVal - genPhi;
if (diff > LauConstants::pi) {
phiVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
phiVal += LauConstants::twoPi;
}
// finally store the new values in the parameters
// and update the pulls
a_->value(aVal); a_->updatePull();
b_->value(bVal); b_->updatePull();
delta_->value(deltaVal); delta_->updatePull();
phi_->value(phiVal); phi_->updatePull();
}
const LauComplex& LauBelleCPCoeffSet::particleCoeff()
{
const LauComplex aTerm{a_->unblindValue()*TMath::Cos(delta_->unblindValue()), a_->unblindValue()*TMath::Sin(delta_->unblindValue())};
const LauComplex bTerm{b_->unblindValue()*TMath::Cos(phi_->unblindValue()), b_->unblindValue()*TMath::Sin(phi_->unblindValue())};
particleCoeff_.setRealImagPart(1.0,0.0);
particleCoeff_ += bTerm;
particleCoeff_ *= aTerm;
return particleCoeff_;
}
const LauComplex& LauBelleCPCoeffSet::antiparticleCoeff()
{
const LauComplex aTerm{a_->unblindValue()*TMath::Cos(delta_->unblindValue()), a_->unblindValue()*TMath::Sin(delta_->unblindValue())};
const LauComplex bTerm{b_->unblindValue()*TMath::Cos(phi_->unblindValue()), b_->unblindValue()*TMath::Sin(phi_->unblindValue())};
antiparticleCoeff_.setRealImagPart(1.0,0.0);
antiparticleCoeff_ -= bTerm;
antiparticleCoeff_ *= aTerm;
return antiparticleCoeff_;
}
void LauBelleCPCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, const Bool_t init )
{
const LauComplex sum { coeff + coeffBar };
const LauComplex diff { coeff - coeffBar };
const LauComplex ratio { diff / sum };
const Double_t aVal{ 0.5 * sum.abs() };
const Double_t deltaVal{ sum.arg() };
const Double_t bVal{ ratio.abs() };
const Double_t phiVal{ ratio.arg() };
a_->value( aVal );
delta_->value( deltaVal );
b_->value( bVal );
phi_->value( phiVal );
if ( init ) {
a_->genValue( aVal );
delta_->genValue( deltaVal );
b_->genValue( bVal );
phi_->genValue( phiVal );
a_->initValue( aVal );
delta_->initValue( deltaVal );
b_->initValue( bVal );
phi_->initValue( phiVal );
}
}
LauParameter LauBelleCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const Double_t value { (-2.0*b_->value()*TMath::Cos(phi_->value()))/(1.0+b_->value()*b_->value()) };
// is it fixed?
const Bool_t fixed { b_->fixed() && phi_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauBelleCPCoeffSet* LauBelleCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauBelleCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauBelleCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauBelleCPCoeffSet nlohmann::adl_serializer<LauBelleCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::BelleCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t a { j.at("a").get<Double_t>() };
const Double_t delta { j.at("delta").get<Double_t>() };
const Double_t b { j.at("b").get<Double_t>() };
const Double_t phi { j.at("phi").get<Double_t>() };
const Bool_t aFixed { j.at("aFixed").get<Bool_t>() };
const Bool_t deltaFixed { j.at("deltaFixed").get<Bool_t>() };
const Bool_t bFixed { j.at("bFixed").get<Bool_t>() };
const Bool_t phiFixed { j.at("phiFixed").get<Bool_t>() };
// TODO - these should be optional?
const Bool_t bSecondStage { j.at("bSecondStage").get<Bool_t>() };
const Bool_t phiSecondStage { j.at("phiSecondStage").get<Bool_t>() };
return { name, a, delta, b, phi, aFixed, deltaFixed, bFixed, phiFixed, bSecondStage, phiSecondStage };
}
void nlohmann::adl_serializer<LauBelleCPCoeffSet>::to_json(json& j, const LauBelleCPCoeffSet& t)
{
j["type"] = LauCoeffType::BelleCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["a"] = pars[0]->value();
j["delta"] = pars[1]->value();
j["b"] = pars[2]->value();
j["phi"] = pars[3]->value();
j["aFixed"] = pars[0]->fixed();
j["deltaFixed"] = pars[1]->fixed();
j["bFixed"] = pars[2]->fixed();
j["phiFixed"] = pars[3]->fixed();
j["bSecondStage"] = pars[2]->secondStage();
j["phiSecondStage"] = pars[3]->secondStage();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauCartesianCPCoeffSet.cc b/src/LauCartesianCPCoeffSet.cc
index f4acdb4..0c30dec 100644
--- a/src/LauCartesianCPCoeffSet.cc
+++ b/src/LauCartesianCPCoeffSet.cc
@@ -1,325 +1,327 @@
/*
Copyright 2006 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 LauCartesianCPCoeffSet.cc
\brief File containing implementation of LauCartesianCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauCartesianCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauCartesianCPCoeffSet)
LauCartesianCPCoeffSet::LauCartesianCPCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Double_t deltaX, const Double_t deltaY,
const Bool_t xFixed, const Bool_t yFixed, const Bool_t deltaXFixed, const Bool_t deltaYFixed,
const Bool_t deltaXSecondStage, const Bool_t deltaYSecondStage) :
LauAbsCoeffSet{ compName },
x_{ std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed) },
y_{ std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) },
deltaX_{ std::make_unique<LauParameter>("DeltaX", deltaX, minDelta_, maxDelta_, deltaXFixed) },
deltaY_{ std::make_unique<LauParameter>("DeltaY", deltaY, minDelta_, maxDelta_, deltaYFixed) },
acp_{ "ACP", -2.0*(x*deltaX + y*deltaY)/(x*x + deltaX*deltaX + y*y + deltaY*deltaY), -1.0, 1.0, deltaXFixed&&deltaYFixed }
{
if (deltaXSecondStage && !deltaXFixed) {
deltaX_->secondStage(kTRUE);
deltaX_->initValue(0.0);
}
if (deltaYSecondStage && !deltaYFixed) {
deltaY_->secondStage(kTRUE);
deltaY_->initValue(0.0);
}
}
LauCartesianCPCoeffSet::LauCartesianCPCoeffSet(const LauCartesianCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
deltaX_.reset( rhs.deltaX_->createClone(constFactor) );
deltaY_.reset( rhs.deltaY_->createClone(constFactor) );
} else {
deltaX_ = std::make_unique<LauParameter>("DeltaX", rhs.deltaX_->value(), minDelta_, maxDelta_, rhs.deltaX_->fixed());
deltaY_ = std::make_unique<LauParameter>("DeltaY", rhs.deltaY_->value(), minDelta_, maxDelta_, rhs.deltaY_->fixed());
if ( rhs.deltaX_->secondStage() && !rhs.deltaX_->fixed() ) {
deltaX_->secondStage(kTRUE);
deltaX_->initValue(0.0);
}
if ( rhs.deltaY_->secondStage() && !rhs.deltaY_->fixed() ) {
deltaY_->secondStage(kTRUE);
deltaY_->initValue(0.0);
}
if ( rhs.deltaX_->blind() ) {
const LauBlind* blinder { rhs.deltaX_->blinder() };
deltaX_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
if ( rhs.deltaY_->blind() ) {
const LauBlind* blinder { rhs.deltaY_->blinder() };
deltaY_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauCartesianCPCoeffSet::getParameters()
{
return { x_.get(), y_.get(), deltaX_.get(), deltaY_.get() };
}
std::vector<const LauParameter*> LauCartesianCPCoeffSet::getParameters() const
{
return { x_.get(), y_.get(), deltaX_.get(), deltaY_.get() };
}
void LauCartesianCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauCartesianCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"x = "<<x_->value()<<",\t";
std::cout<<"y = "<<y_->value()<<",\t";
std::cout<<"Delta x = "<<deltaX_->value()<<",\t";
std::cout<<"Delta y = "<<deltaY_->value()<<"."<<std::endl;
}
void LauCartesianCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $\\Delta$ Real Part & $\\Delta$ Imaginary Part \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauCartesianCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
print.printFormat(stream, deltaX_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaX_->error());
stream<<"$ & $";
print.printFormat(stream, deltaY_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaY_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauCartesianCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
y_->initValue(value); y_->value(value);
}
if (deltaX_->fixed() == kFALSE && deltaX_->secondStage() == kFALSE) {
// Choose a value for "Delta X" between -0.5 and 0.5
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*1.0 - 0.5 };
deltaX_->initValue(value); deltaX_->value(value);
}
if (deltaY_->fixed() == kFALSE && deltaY_->secondStage() == kFALSE) {
// Choose a value for "Delta Y" between -0.5 and 0.5
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*1.0 - 0.5 };
deltaY_->initValue(value); deltaY_->value(value);
}
}
void LauCartesianCPCoeffSet::finaliseValues()
{
// update the pulls
x_->updatePull();
y_->updatePull();
deltaX_->updatePull();
deltaY_->updatePull();
}
const LauComplex& LauCartesianCPCoeffSet::particleCoeff()
{
particleCoeff_.setRealImagPart( x_->unblindValue() + deltaX_->unblindValue(), y_->unblindValue() + deltaY_->unblindValue() );
return particleCoeff_;
}
const LauComplex& LauCartesianCPCoeffSet::antiparticleCoeff()
{
antiparticleCoeff_.setRealImagPart( x_->unblindValue() - deltaX_->unblindValue(), y_->unblindValue() - deltaY_->unblindValue() );
return antiparticleCoeff_;
}
void LauCartesianCPCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, const Bool_t init )
{
LauComplex average{ coeff };
average += coeffBar;
average.rescale( 0.5 );
const Double_t xVal{ average.re() };
const Double_t yVal{ average.im() };
const Double_t deltaXVal{ coeff.re() - average.re() };
const Double_t deltaYVal{ coeff.im() - average.im() };
x_->value( xVal );
y_->value( yVal );
deltaX_->value( deltaXVal );
deltaY_->value( deltaYVal );
if ( init ) {
x_->genValue( xVal );
y_->genValue( yVal );
deltaX_->genValue( deltaXVal );
deltaY_->genValue( deltaYVal );
x_->initValue( xVal );
y_->initValue( yVal );
deltaX_->initValue( deltaXVal );
deltaY_->initValue( deltaYVal );
}
}
LauParameter LauCartesianCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const Double_t numer { x_->value()*deltaX_->value() + y_->value()*deltaY_->value() };
const Double_t denom { x_->value()*x_->value() + deltaX_->value()*deltaX_->value() + y_->value()*y_->value() + deltaY_->value()*deltaY_->value() };
const Double_t value { -2.0*numer/denom };
// is it fixed?
const Bool_t fixed { deltaX_->fixed() && deltaY_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauCartesianCPCoeffSet* LauCartesianCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauCartesianCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauCartesianCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauCartesianCPCoeffSet nlohmann::adl_serializer<LauCartesianCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::CartesianCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauCartesianCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t deltaX { j.at("deltaX").get<Double_t>() };
const Double_t deltaY { j.at("deltaY").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t deltaXFixed { j.at("deltaXFixed").get<Bool_t>() };
const Bool_t deltaYFixed { j.at("deltaYFixed").get<Bool_t>() };
const Bool_t deltaXSecondStage { j.at("deltaXSecondStage").get<Bool_t>() };
const Bool_t deltaYSecondStage { j.at("deltaYSecondStage").get<Bool_t>() };
return { name, x, y, deltaX, deltaY, xFixed, yFixed, deltaXFixed, deltaYFixed, deltaXSecondStage, deltaYSecondStage };
}
void nlohmann::adl_serializer<LauCartesianCPCoeffSet>::to_json(json& j, const LauCartesianCPCoeffSet& t)
{
j["type"] = LauCoeffType::CartesianCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["deltaX"] = pars[2]->value();
j["deltaY"] = pars[3]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["deltaXFixed"] = pars[2]->fixed();
j["deltaYFixed"] = pars[3]->fixed();
j["deltaXSecondStage"] = pars[2]->secondStage();
j["deltaYSecondStage"] = pars[3]->secondStage();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauCartesianGammaCPCoeffSet.cc b/src/LauCartesianGammaCPCoeffSet.cc
index 086e8e0..8e5d182 100644
--- a/src/LauCartesianGammaCPCoeffSet.cc
+++ b/src/LauCartesianGammaCPCoeffSet.cc
@@ -1,377 +1,379 @@
/*
Copyright 2014 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 LauCartesianGammaCPCoeffSet.cc
\brief File containing implementation of LauCartesianGammaCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauCartesianGammaCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauCartesianGammaCPCoeffSet)
LauCartesianGammaCPCoeffSet::LauCartesianGammaCPCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Double_t xCP, const Double_t yCP, const Double_t deltaXCP, const Double_t deltaYCP,
const Bool_t xFixed, const Bool_t yFixed, const Bool_t xCPFixed, const Bool_t yCPFixed, const Bool_t deltaXCPFixed, const Bool_t deltaYCPFixed,
const Bool_t deltaXCPSecondStage, const Bool_t deltaYCPSecondStage) :
LauAbsCoeffSet{ compName },
x_{ std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed) },
y_{ std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) },
xCP_{ std::make_unique<LauParameter>("XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed) },
yCP_{ std::make_unique<LauParameter>("YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed) },
deltaXCP_{ std::make_unique<LauParameter>("DeltaXCP", deltaXCP, minDelta_, maxDelta_, deltaXCPFixed) },
deltaYCP_{ std::make_unique<LauParameter>("DeltaYCP", deltaYCP, minDelta_, maxDelta_, deltaYCPFixed) },
nonCPPart_{ x, y },
cpPart_{ 1+xCP+deltaXCP, yCP+deltaYCP },
cpAntiPart_{ 1+xCP-deltaXCP, yCP-deltaYCP },
particleCoeff_{ nonCPPart_ * cpPart_ },
antiparticleCoeff_{ nonCPPart_ * cpAntiPart_ },
acp_{ "ACP", (antiparticleCoeff_.abs2()-particleCoeff_.abs2())/(antiparticleCoeff_.abs2()+particleCoeff_.abs2()), -1.0, 1.0, deltaXCPFixed&&deltaYCPFixed }
{
if (deltaXCPSecondStage && !deltaXCPFixed) {
deltaXCP_->secondStage(kTRUE);
deltaXCP_->initValue(0.0);
}
if (deltaYCPSecondStage && !deltaYCPFixed) {
deltaYCP_->secondStage(kTRUE);
deltaYCP_->initValue(0.0);
}
}
LauCartesianGammaCPCoeffSet::LauCartesianGammaCPCoeffSet(const LauCartesianGammaCPCoeffSet& rhs, CloneOption cloneOption, Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
nonCPPart_{ rhs.nonCPPart_ },
cpPart_{ rhs.cpPart_ },
cpAntiPart_{ rhs.cpAntiPart_ },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
xCP_.reset( rhs.xCP_->createClone(constFactor) );
yCP_.reset( rhs.yCP_->createClone(constFactor) );
deltaXCP_.reset( rhs.deltaXCP_->createClone(constFactor) );
deltaYCP_.reset( rhs.deltaYCP_->createClone(constFactor) );
} else {
xCP_ = std::make_unique<LauParameter>("XCP", rhs.xCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.xCP_->fixed());
if ( rhs.xCP_->blind() ) {
const LauBlind* blinder { rhs.xCP_->blinder() };
xCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
yCP_ = std::make_unique<LauParameter>("YCP", rhs.yCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.yCP_->fixed());
if ( rhs.yCP_->blind() ) {
const LauBlind* blinder { rhs.yCP_->blinder() };
yCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
deltaXCP_ = std::make_unique<LauParameter>("DeltaXCP", rhs.deltaXCP_->value(), minDelta_, maxDelta_, rhs.deltaXCP_->fixed());
deltaYCP_ = std::make_unique<LauParameter>("DeltaYCP", rhs.deltaYCP_->value(), minDelta_, maxDelta_, rhs.deltaYCP_->fixed());
if ( rhs.deltaXCP_->secondStage() && !rhs.deltaXCP_->fixed() ) {
deltaXCP_->secondStage(kTRUE);
deltaXCP_->initValue(0.0);
}
if ( rhs.deltaYCP_->secondStage() && !rhs.deltaYCP_->fixed() ) {
deltaYCP_->secondStage(kTRUE);
deltaYCP_->initValue(0.0);
}
if ( rhs.deltaXCP_->blind() ) {
const LauBlind* blinder { rhs.deltaXCP_->blinder() };
deltaXCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
if ( rhs.deltaYCP_->blind() ) {
const LauBlind* blinder { rhs.deltaYCP_->blinder() };
deltaYCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauCartesianGammaCPCoeffSet::getParameters()
{
return {
x_.get(),
y_.get(),
xCP_.get(),
yCP_.get(),
deltaXCP_.get(),
deltaYCP_.get()
};
}
std::vector<const LauParameter*> LauCartesianGammaCPCoeffSet::getParameters() const
{
return {
x_.get(),
y_.get(),
xCP_.get(),
yCP_.get(),
deltaXCP_.get(),
deltaYCP_.get()
};
}
void LauCartesianGammaCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauCartesianGammaCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"x = "<<x_->value()<<",\t";
std::cout<<"y = "<<y_->value()<<",\t";
std::cout<<"xCP = "<<xCP_->value()<<",\t";
std::cout<<"yCP = "<<yCP_->value()<<",\t";
std::cout<<"Delta xCP = "<<deltaXCP_->value()<<",\t";
std::cout<<"Delta yCP = "<<deltaYCP_->value()<<"."<<std::endl;
}
void LauCartesianGammaCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & Real CP Part & Imaginary CP Part & $\\Delta$ Real CP Part & $\\Delta$ Imaginary CP Part \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauCartesianGammaCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
print.printFormat(stream, xCP_->value());
stream<<" \\pm ";
print.printFormat(stream, xCP_->error());
stream<<"$ & $";
print.printFormat(stream, yCP_->value());
stream<<" \\pm ";
print.printFormat(stream, yCP_->error());
stream<<"$ & $";
print.printFormat(stream, deltaXCP_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaXCP_->error());
stream<<"$ & $";
print.printFormat(stream, deltaYCP_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaYCP_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauCartesianGammaCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
y_->initValue(value); y_->value(value);
}
if (xCP_->fixed() == kFALSE) {
// Choose a value for "XCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
xCP_->initValue(value); xCP_->value(value);
}
if (yCP_->fixed() == kFALSE) {
// Choose a value for "YCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
yCP_->initValue(value); yCP_->value(value);
}
if (deltaXCP_->fixed() == kFALSE && deltaXCP_->secondStage() == kFALSE) {
// Choose a value for "Delta XCP" between -0.5 and 0.5
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*1.0 - 0.5 };
deltaXCP_->initValue(value); deltaXCP_->value(value);
}
if (deltaYCP_->fixed() == kFALSE && deltaYCP_->secondStage() == kFALSE) {
// Choose a value for "Delta YCP" between -0.5 and 0.5
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*1.0 - 0.5 };
deltaYCP_->initValue(value); deltaYCP_->value(value);
}
}
void LauCartesianGammaCPCoeffSet::finaliseValues()
{
// update the pulls
x_->updatePull();
y_->updatePull();
xCP_->updatePull();
yCP_->updatePull();
deltaXCP_->updatePull();
deltaYCP_->updatePull();
}
const LauComplex& LauCartesianGammaCPCoeffSet::particleCoeff()
{
nonCPPart_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
cpPart_.setRealImagPart( 1.0 + xCP_->unblindValue() + deltaXCP_->unblindValue(), yCP_->unblindValue() + deltaYCP_->unblindValue() );
particleCoeff_ = nonCPPart_ * cpPart_;
return particleCoeff_;
}
const LauComplex& LauCartesianGammaCPCoeffSet::antiparticleCoeff()
{
nonCPPart_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
cpAntiPart_.setRealImagPart( 1.0 + xCP_->unblindValue() - deltaXCP_->unblindValue(), yCP_->unblindValue() - deltaYCP_->unblindValue() );
antiparticleCoeff_ = nonCPPart_ * cpAntiPart_;
return antiparticleCoeff_;
}
void LauCartesianGammaCPCoeffSet::setCoeffValues( const LauComplex&, const LauComplex&, const Bool_t )
{
std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters" << std::endl;
}
LauParameter LauCartesianGammaCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const LauComplex nonCPPart{ x_->value(), y_->value() };
const LauComplex cpPart{ 1.0 + xCP_->value() + deltaXCP_->value(), yCP_->value() + deltaYCP_->value() };
const LauComplex cpAntiPart{ 1.0 + xCP_->value() - deltaXCP_->value(), yCP_->value() - deltaYCP_->value() };
const LauComplex partCoeff { nonCPPart * cpPart };
const LauComplex antiCoeff { nonCPPart * cpAntiPart };
const Double_t numer { antiCoeff.abs2() - partCoeff.abs2() };
const Double_t denom { antiCoeff.abs2() + partCoeff.abs2() };
const Double_t value { numer/denom };
// is it fixed?
const Bool_t fixed { deltaXCP_->fixed() && deltaYCP_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauCartesianGammaCPCoeffSet* LauCartesianGammaCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauCartesianGammaCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauCartesianGammaCPCoeffSet nlohmann::adl_serializer<LauCartesianGammaCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::CartesianGammaCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauCartesianGammaCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t xCP { j.at("xCP").get<Double_t>() };
const Double_t yCP { j.at("yCP").get<Double_t>() };
const Double_t deltaXCP { j.at("deltaXCP").get<Double_t>() };
const Double_t deltaYCP { j.at("deltaYCP").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t xCPFixed { j.at("xCPFixed").get<Bool_t>() };
const Bool_t yCPFixed { j.at("yCPFixed").get<Bool_t>() };
const Bool_t deltaXCPFixed { j.at("deltaXCPFixed").get<Bool_t>() };
const Bool_t deltaYCPFixed { j.at("deltaYCPFixed").get<Bool_t>() };
const Bool_t deltaXCPSecondStage { j.at("deltaXCPSecondStage").get<Bool_t>() };
const Bool_t deltaYCPSecondStage { j.at("deltaYCPSecondStage").get<Bool_t>() };
return { name, x, y, xCP, yCP, deltaXCP, deltaYCP, xFixed, yFixed, xCPFixed, yCPFixed, deltaXCPFixed, deltaYCPFixed, deltaXCPSecondStage, deltaYCPSecondStage };
}
void nlohmann::adl_serializer<LauCartesianGammaCPCoeffSet>::to_json(json& j, const LauCartesianGammaCPCoeffSet& t)
{
j["type"] = LauCoeffType::CartesianGammaCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["xCP"] = pars[2]->value();
j["yCP"] = pars[3]->value();
j["deltaXCP"] = pars[4]->value();
j["deltaYCP"] = pars[5]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["xCPFixed"] = pars[2]->fixed();
j["yCPFixed"] = pars[3]->fixed();
j["deltaXCPFixed"] = pars[4]->fixed();
j["deltaYCPFixed"] = pars[5]->fixed();
j["deltaXCPSecondStage"] = pars[4]->secondStage();
j["deltaYCPSecondStage"] = pars[5]->secondStage();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauCleoCPCoeffSet.cc b/src/LauCleoCPCoeffSet.cc
index 9699313..5436911 100644
--- a/src/LauCleoCPCoeffSet.cc
+++ b/src/LauCleoCPCoeffSet.cc
@@ -1,390 +1,391 @@
/*
Copyright 2006 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 LauCleoCPCoeffSet.cc
\brief File containing implementation of LauCleoCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauCleoCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauCleoCPCoeffSet)
LauCleoCPCoeffSet::LauCleoCPCoeffSet(const TString& compName, const Double_t a, const Double_t delta, const Double_t b, const Double_t phi,
const Bool_t aFixed, const Bool_t deltaFixed, const Bool_t bFixed, const Bool_t phiFixed,
const Bool_t bSecondStage, const Bool_t phiSecondStage) :
LauAbsCoeffSet{ compName },
a_{ std::make_unique<LauParameter>("A", a, minMagnitude_, maxMagnitude_, aFixed) },
b_{ std::make_unique<LauParameter>("B", b, minMagnitude_, maxMagnitude_, bFixed) },
delta_{ std::make_unique<LauParameter>("Delta", delta, minPhase_, maxPhase_, deltaFixed) },
phi_{ std::make_unique<LauParameter>("Phi", phi, minPhase_, maxPhase_, phiFixed) },
acp_{ "ACP", (-2.0*a*b)/(a*a + b*b), -1.0, 1.0, bFixed&&phiFixed }
{
if (bSecondStage && !bFixed) {
b_->secondStage(kTRUE);
b_->initValue(0.0);
}
if (phiSecondStage && !phiFixed) {
phi_->secondStage(kTRUE);
phi_->initValue(0.0);
}
}
LauCleoCPCoeffSet::LauCleoCPCoeffSet(const LauCleoCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) {
a_.reset( rhs.a_->createClone(constFactor) );
} else {
a_ = std::make_unique<LauParameter>("A", rhs.a_->value(), minMagnitude_, maxMagnitude_, rhs.a_->fixed());
if ( rhs.a_->blind() ) {
const LauBlind* blinder { rhs.a_->blinder() };
a_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
b_.reset( rhs.b_->createClone(constFactor) );
} else {
b_ = std::make_unique<LauParameter>("B", rhs.b_->value(), minMagnitude_, maxMagnitude_, rhs.b_->fixed());
if ( rhs.b_->blind() ) {
const LauBlind* blinder { rhs.b_->blinder() };
b_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase ) {
delta_.reset( rhs.delta_->createClone(constFactor) );
} else {
delta_ = std::make_unique<LauParameter>("Delta", rhs.delta_->value(), minPhase_, maxPhase_, rhs.delta_->fixed());
if ( rhs.delta_->blind() ) {
const LauBlind* blinder { rhs.delta_->blinder() };
delta_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
phi_.reset( rhs.phi_->createClone(constFactor) );
} else {
phi_ = std::make_unique<LauParameter>("Phi", rhs.phi_->value(), minPhase_, maxPhase_, rhs.phi_->fixed());
if ( rhs.phi_->blind() ) {
const LauBlind* blinder { rhs.phi_->blinder() };
phi_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauCleoCPCoeffSet::getParameters()
{
return { a_.get(), delta_.get(), b_.get(), phi_.get() };
}
std::vector<const LauParameter*> LauCleoCPCoeffSet::getParameters() const
{
return { a_.get(), delta_.get(), b_.get(), phi_.get() };
}
void LauCleoCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauCleoCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"a-magnitude = "<<a_->value()<<",\t";
std::cout<<"delta = "<<delta_->value()<<",\t";
std::cout<<"b-magnitude = "<<b_->value()<<",\t";
std::cout<<"phi = "<<phi_->value()<<"."<<std::endl;
}
void LauCleoCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & a-Magnitude & delta & b-Magnitude & phi \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauCleoCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, a_->value());
stream<<" \\pm ";
print.printFormat(stream, a_->error());
stream<<"$ & $";
print.printFormat(stream, delta_->value());
stream<<" \\pm ";
print.printFormat(stream, delta_->error());
stream<<"$ & $";
print.printFormat(stream, b_->value());
stream<<" \\pm ";
print.printFormat(stream, b_->error());
stream<<"$ & $";
print.printFormat(stream, phi_->value());
stream<<" \\pm ";
print.printFormat(stream, phi_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauCleoCPCoeffSet::randomiseInitValues()
{
if (a_->fixed() == kFALSE) {
// Choose an a-magnitude between 0.0 and 2.0
const Double_t mag { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
a_->initValue(mag); a_->value(mag);
}
if (b_->fixed() == kFALSE && b_->secondStage() == kFALSE) {
// Choose a b-magnitude between 0.0 and 0.1
const Double_t mag { LauAbsCoeffSet::getRandomiser()->Rndm()*0.1 };
b_->initValue(mag); b_->value(mag);
}
if (delta_->fixed() == kFALSE) {
// Choose a phase between +- pi
const Double_t phase { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
delta_->initValue(phase); delta_->value(phase);
}
if (phi_->fixed() == kFALSE && phi_->secondStage() == kFALSE) {
// Choose a phase between +- pi
const Double_t phase { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
phi_->initValue(phase); phi_->value(phase);
}
}
void LauCleoCPCoeffSet::finaliseValues()
{
// retrieve the current values from the parameters
Double_t aVal { a_->value() };
Double_t bVal { b_->value() };
Double_t deltaVal { delta_->value() };
Double_t phiVal { phi_->value() };
// Check whether we have a negative "a" magnitude.
// If so make it positive and add pi to the "delta" phase.
if (aVal < 0.0) {
aVal *= -1.0;
bVal *= -1.0;
deltaVal += LauConstants::pi;
}
// Check now whether the phases lies in the right range (-pi to pi).
Bool_t deltaWithinRange{kFALSE};
Bool_t phiWithinRange{kFALSE};
while (deltaWithinRange == kFALSE && phiWithinRange == kFALSE) {
if (deltaVal > -LauConstants::pi && deltaVal < LauConstants::pi) {
deltaWithinRange = kTRUE;
} else {
// Not within the specified range
if (deltaVal > LauConstants::pi) {
deltaVal -= LauConstants::twoPi;
} else if (deltaVal < -LauConstants::pi) {
deltaVal += LauConstants::twoPi;
}
}
if (phiVal > -LauConstants::pi && phiVal < LauConstants::pi) {
phiWithinRange = kTRUE;
} else {
// Not within the specified range
if (phiVal > LauConstants::pi) {
phiVal -= LauConstants::twoPi;
} else if (phiVal < -LauConstants::pi) {
phiVal += LauConstants::twoPi;
}
}
}
// A further problem can occur when the generated phase is close to -pi or pi.
// The phase can wrap over to the other end of the scale -
// this leads to artificially large pulls so we wrap it back.
const Double_t genDelta { delta_->genValue() };
const Double_t genPhi { phi_->genValue() };
Double_t diff { deltaVal - genDelta };
if (diff > LauConstants::pi) {
deltaVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
deltaVal += LauConstants::twoPi;
}
diff = phiVal - genPhi;
if (diff > LauConstants::pi) {
phiVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
phiVal += LauConstants::twoPi;
}
// finally store the new values in the parameters
// and update the pulls
a_->value(aVal); a_->updatePull();
b_->value(bVal); b_->updatePull();
delta_->value(deltaVal); delta_->updatePull();
phi_->value(phiVal); phi_->updatePull();
}
const LauComplex& LauCleoCPCoeffSet::particleCoeff()
{
const Double_t magnitude { a_->unblindValue() + b_->unblindValue() };
const Double_t phase { delta_->unblindValue() + phi_->unblindValue() };
particleCoeff_.setRealImagPart(magnitude*TMath::Cos(phase), magnitude*TMath::Sin(phase));
return particleCoeff_;
}
const LauComplex& LauCleoCPCoeffSet::antiparticleCoeff()
{
const Double_t magnitude { a_->unblindValue() - b_->unblindValue() };
const Double_t phase { delta_->unblindValue() - phi_->unblindValue() };
antiparticleCoeff_.setRealImagPart(magnitude*TMath::Cos(phase), magnitude*TMath::Sin(phase));
return antiparticleCoeff_;
}
void LauCleoCPCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, const Bool_t init )
{
const Double_t mag { coeff.abs() };
const Double_t magBar { coeffBar.abs() };
const Double_t phase { coeff.arg() };
const Double_t phaseBar { coeffBar.arg() };
const Double_t aVal{ 0.5 * ( mag + magBar ) };
const Double_t deltaVal{ 0.5 * ( phase + phaseBar ) };
const Double_t bVal{ 0.5 * ( mag - magBar ) };
const Double_t phiVal{ 0.5 * ( phase - phaseBar ) };
a_->value( aVal );
delta_->value( deltaVal );
b_->value( bVal );
phi_->value( phiVal );
if ( init ) {
a_->genValue( aVal );
delta_->genValue( deltaVal );
b_->genValue( bVal );
phi_->genValue( phiVal );
a_->initValue( aVal );
delta_->initValue( deltaVal );
b_->initValue( bVal );
phi_->initValue( phiVal );
}
}
LauParameter LauCleoCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const Double_t numer { -2.0*a_->value()*b_->value() };
const Double_t denom { a_->value()*a_->value() + b_->value()*b_->value() };
const Double_t value { numer/denom };
// is it fixed?
const Bool_t fixed { a_->fixed() && b_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauCleoCPCoeffSet* LauCleoCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauCleoCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauCleoCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
-
+//! \cond DOXYGEN_IGNORE
LauCleoCPCoeffSet nlohmann::adl_serializer<LauCleoCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::CleoCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauCleoCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t a { j.at("a").get<Double_t>() };
const Double_t delta { j.at("delta").get<Double_t>() };
const Double_t b { j.at("b").get<Double_t>() };
const Double_t phi { j.at("phi").get<Double_t>() };
const Bool_t aFixed { j.at("aFixed").get<Bool_t>() };
const Bool_t deltaFixed { j.at("deltaFixed").get<Bool_t>() };
const Bool_t bFixed { j.at("bFixed").get<Bool_t>() };
const Bool_t phiFixed { j.at("phiFixed").get<Bool_t>() };
// TODO - these should be optional?
const Bool_t bSecondStage { j.at("bSecondStage").get<Bool_t>() };
const Bool_t phiSecondStage { j.at("phiSecondStage").get<Bool_t>() };
return { name, a, delta, b, phi, aFixed, deltaFixed, bFixed, phiFixed, bSecondStage, phiSecondStage };
}
void nlohmann::adl_serializer<LauCleoCPCoeffSet>::to_json(json& j, const LauCleoCPCoeffSet& t)
{
j["type"] = LauCoeffType::CleoCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["a"] = pars[0]->value();
j["delta"] = pars[1]->value();
j["b"] = pars[2]->value();
j["phi"] = pars[3]->value();
j["aFixed"] = pars[0]->fixed();
j["deltaFixed"] = pars[1]->fixed();
j["bFixed"] = pars[2]->fixed();
j["phiFixed"] = pars[3]->fixed();
j["bSecondStage"] = pars[2]->secondStage();
j["phiSecondStage"] = pars[3]->secondStage();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauMagPhaseCPCoeffSet.cc b/src/LauMagPhaseCPCoeffSet.cc
index 61f568f..9ad75c0 100644
--- a/src/LauMagPhaseCPCoeffSet.cc
+++ b/src/LauMagPhaseCPCoeffSet.cc
@@ -1,358 +1,360 @@
/*
Copyright 2011 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 LauMagPhaseCPCoeffSet.cc
\brief File containing implementation of LauMagPhaseCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauMagPhaseCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauMagPhaseCPCoeffSet)
LauMagPhaseCPCoeffSet::LauMagPhaseCPCoeffSet(const TString& compName, const Double_t mag, const Double_t phase, const Double_t magBar, const Double_t phaseBar,
const Bool_t magFixed, const Bool_t phaseFixed, const Bool_t magBarFixed, const Bool_t phaseBarFixed) :
LauAbsCoeffSet{ compName },
mag_{ std::make_unique<LauParameter>("Mag", mag, minMagnitude_, maxMagnitude_, magFixed) },
phase_{ std::make_unique<LauParameter>("Phase", phase, minPhase_, maxPhase_, phaseFixed) },
magBar_{ std::make_unique<LauParameter>("MagBar", magBar, minMagnitude_, maxMagnitude_, magBarFixed) },
phaseBar_{ std::make_unique<LauParameter>("PhaseBar", phaseBar, minPhase_, maxPhase_, phaseBarFixed) },
acp_{ "ACP", (magBar*magBar - mag*mag)/(magBar*magBar + mag*mag), -1.0, 1.0 }
{
}
LauMagPhaseCPCoeffSet::LauMagPhaseCPCoeffSet(const LauMagPhaseCPCoeffSet& rhs, CloneOption cloneOption, Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) {
mag_.reset( rhs.mag_->createClone(constFactor) );
magBar_.reset( rhs.magBar_->createClone(constFactor) );
} else {
mag_ = std::make_unique<LauParameter>("Mag", rhs.mag_->value(), minMagnitude_, maxMagnitude_, rhs.mag_->fixed());
if ( rhs.mag_->blind() ) {
const LauBlind* blinder { rhs.mag_->blinder() };
mag_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
magBar_ = std::make_unique<LauParameter>("MagBar", rhs.magBar_->value(), minMagnitude_, maxMagnitude_, rhs.magBar_->fixed());
if ( rhs.magBar_->blind() ) {
const LauBlind* blinder { rhs.magBar_->blinder() };
magBar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase ) {
phase_.reset( rhs.phase_->createClone(constFactor) );
phaseBar_.reset( rhs.phaseBar_->createClone(constFactor) );
} else {
phase_ = std::make_unique<LauParameter>("Phase", rhs.phase_->value(), minPhase_, maxPhase_, rhs.phase_->fixed());
if ( rhs.phase_->blind() ) {
const LauBlind* blinder { rhs.phase_->blinder() };
phase_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
phaseBar_ = std::make_unique<LauParameter>("PhaseBar", rhs.phaseBar_->value(), minPhase_, maxPhase_, rhs.phaseBar_->fixed());
if ( rhs.phaseBar_->blind() ) {
const LauBlind* blinder { rhs.phaseBar_->blinder() };
phaseBar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauMagPhaseCPCoeffSet::getParameters()
{
return { mag_.get(), phase_.get(), magBar_.get(), phaseBar_.get() };
}
std::vector<const LauParameter*> LauMagPhaseCPCoeffSet::getParameters() const
{
return { mag_.get(), phase_.get(), magBar_.get(), phaseBar_.get() };
}
void LauMagPhaseCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauMagPhaseCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"mag = "<<mag_->value()<<",\t";
std::cout<<"phase = "<<phase_->value()<<",\t";
std::cout<<"magBar = "<<magBar_->value()<<",\t";
std::cout<<"phaseBar = "<<phaseBar_->value()<<"."<<std::endl;
}
void LauMagPhaseCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Magnitude & Phase & Magnitude_bar & Phase_bar \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauMagPhaseCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, mag_->value());
stream<<" \\pm ";
print.printFormat(stream, mag_->error());
stream<<"$ & $";
print.printFormat(stream, phase_->value());
stream<<" \\pm ";
print.printFormat(stream, phase_->error());
stream<<"$ & $";
print.printFormat(stream, magBar_->value());
stream<<" \\pm ";
print.printFormat(stream, magBar_->error());
stream<<"$ & $";
print.printFormat(stream, phaseBar_->value());
stream<<" \\pm ";
print.printFormat(stream, phaseBar_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauMagPhaseCPCoeffSet::randomiseInitValues()
{
if (mag_->fixed() == kFALSE) {
// Choose a value for "magnitude" between 0.0 and 2.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
mag_->initValue(value); mag_->value(value);
}
if (phase_->fixed() == kFALSE) {
// Choose a phase between +- pi
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
phase_->initValue(value); phase_->value(value);
}
if (magBar_->fixed() == kFALSE) {
// Choose a value for "magnitude" between 0.0 and 2.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
magBar_->initValue(value); magBar_->value(value);
}
if (phaseBar_->fixed() == kFALSE) {
// Choose a phase between +- pi
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
phaseBar_->initValue(value); phaseBar_->value(value);
}
}
void LauMagPhaseCPCoeffSet::finaliseValues()
{
// retrieve the current values from the parameters
Double_t mVal{ mag_->value() };
Double_t pVal{ phase_->value() };
Double_t mBarVal{ magBar_->value() };
Double_t pBarVal{ phaseBar_->value() };
// Check whether we have a negative magnitude.
// If so make it positive and add pi to the phase.
if (mVal < 0.0) {
mVal *= -1.0;
pVal += LauConstants::pi;
}
if (mBarVal < 0.0) {
mBarVal *= -1.0;
pBarVal += LauConstants::pi;
}
// Check now whether the phases lies in the right range (-pi to pi).
Bool_t pWithinRange{kFALSE};
Bool_t pBarWithinRange{kFALSE};
while (pWithinRange == kFALSE && pBarWithinRange == kFALSE) {
if (pVal > -LauConstants::pi && pVal < LauConstants::pi) {
pWithinRange = kTRUE;
} else {
// Not within the specified range
if (pVal > LauConstants::pi) {
pVal -= LauConstants::twoPi;
} else if (pVal < -LauConstants::pi) {
pVal += LauConstants::twoPi;
}
}
if (pBarVal > -LauConstants::pi && pBarVal < LauConstants::pi) {
pBarWithinRange = kTRUE;
} else {
// Not within the specified range
if (pBarVal > LauConstants::pi) {
pBarVal -= LauConstants::twoPi;
} else if (pBarVal < -LauConstants::pi) {
pBarVal += LauConstants::twoPi;
}
}
}
// A further problem can occur when the generated phase is close to -pi or pi.
// The phase can wrap over to the other end of the scale -
// this leads to artificially large pulls so we wrap it back.
const Double_t genPhase { phase_->genValue() };
const Double_t genPhaseBar { phaseBar_->genValue() };
Double_t diff { pVal - genPhase };
if (diff > LauConstants::pi) {
pVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
pVal += LauConstants::twoPi;
}
diff = pBarVal - genPhaseBar;
if (diff > LauConstants::pi) {
pBarVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
pBarVal += LauConstants::twoPi;
}
// finally store the new values in the parameters
// and update the pulls
mag_->value(mVal); mag_->updatePull();
phase_->value(pVal); phase_->updatePull();
magBar_->value(mBarVal); magBar_->updatePull();
phaseBar_->value(pBarVal); phaseBar_->updatePull();
}
const LauComplex& LauMagPhaseCPCoeffSet::particleCoeff()
{
particleCoeff_.setRealImagPart( mag_->unblindValue()*TMath::Cos(phase_->unblindValue()), mag_->unblindValue()*TMath::Sin(phase_->unblindValue()) );
return particleCoeff_;
}
const LauComplex& LauMagPhaseCPCoeffSet::antiparticleCoeff()
{
antiparticleCoeff_.setRealImagPart( magBar_->unblindValue()*TMath::Cos(phaseBar_->unblindValue()), magBar_->unblindValue()*TMath::Sin(phaseBar_->unblindValue()) );
return antiparticleCoeff_;
}
void LauMagPhaseCPCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, Bool_t init )
{
const Double_t magVal{ coeff.abs() };
const Double_t phaseVal{ coeff.arg() };
const Double_t magBarVal{ coeffBar.abs() };
const Double_t phaseBarVal{ coeffBar.arg() };
mag_->value( magVal );
phase_->value( phaseVal );
magBar_->value( magBarVal );
phaseBar_->value( phaseBarVal );
if ( init ) {
mag_->genValue( magVal );
phase_->genValue( phaseVal );
magBar_->genValue( magBarVal );
phaseBar_->genValue( phaseBarVal );
mag_->initValue( magVal );
phase_->initValue( phaseVal );
magBar_->initValue( magBarVal );
phaseBar_->initValue( phaseBarVal );
}
}
LauParameter LauMagPhaseCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const Double_t value { (magBar_->value()*magBar_->value() - mag_->value()*mag_->value())/(magBar_->value()*magBar_->value() + mag_->value()*mag_->value()) };
// is it fixed?
const Bool_t fixed { magBar_->fixed() && mag_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauMagPhaseCPCoeffSet* LauMagPhaseCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude ) ) {
std::cerr << "ERROR in LauMagPhaseCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauMagPhaseCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauMagPhaseCPCoeffSet nlohmann::adl_serializer<LauMagPhaseCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::MagPhaseCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauMagPhaseCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t mag { j.at("mag").get<Double_t>() };
const Double_t phase { j.at("phase").get<Double_t>() };
const Double_t magBar { j.at("magBar").get<Double_t>() };
const Double_t phaseBar { j.at("phaseBar").get<Double_t>() };
const Bool_t magFixed { j.at("magFixed").get<Bool_t>() };
const Bool_t phaseFixed { j.at("phaseFixed").get<Bool_t>() };
const Bool_t magBarFixed { j.at("magBarFixed").get<Bool_t>() };
const Bool_t phaseBarFixed { j.at("phaseBarFixed").get<Bool_t>() };
return { name, mag, phase, magBar, phaseBar, magFixed, phaseFixed, magBarFixed, phaseBarFixed };
}
void nlohmann::adl_serializer<LauMagPhaseCPCoeffSet>::to_json(json& j, const LauMagPhaseCPCoeffSet& t)
{
j["type"] = LauCoeffType::MagPhaseCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["mag"] = pars[0]->value();
j["phase"] = pars[1]->value();
j["magBar"] = pars[2]->value();
j["phaseBar"] = pars[3]->value();
j["magFixed"] = pars[0]->fixed();
j["phaseFixed"] = pars[1]->fixed();
j["magBarFixed"] = pars[2]->fixed();
j["phaseBarFixed"] = pars[3]->fixed();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauMagPhaseCoeffSet.cc b/src/LauMagPhaseCoeffSet.cc
index 443a833..1d1e4ee 100644
--- a/src/LauMagPhaseCoeffSet.cc
+++ b/src/LauMagPhaseCoeffSet.cc
@@ -1,257 +1,259 @@
/*
Copyright 2006 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 LauMagPhaseCoeffSet.cc
\brief File containing implementation of LauMagPhaseCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauMagPhaseCoeffSet.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauMagPhaseCoeffSet)
LauMagPhaseCoeffSet::LauMagPhaseCoeffSet(const TString& compName, const Double_t magnitude, const Double_t phase, const Bool_t magFixed, const Bool_t phaseFixed) :
LauAbsCoeffSet{ compName },
magnitude_{ std::make_unique<LauParameter>("A",magnitude,minMagnitude_,maxMagnitude_,magFixed) },
phase_{ std::make_unique<LauParameter>("Delta",phase,minPhase_,maxPhase_,phaseFixed) }
{
}
LauMagPhaseCoeffSet::LauMagPhaseCoeffSet(const LauMagPhaseCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
coeff_{ rhs.coeff_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) {
magnitude_.reset( rhs.magnitude_->createClone(constFactor) );
} else {
magnitude_ = std::make_unique<LauParameter>("A", rhs.magnitude_->value(), minMagnitude_, maxMagnitude_, rhs.magnitude_->fixed());
if ( rhs.magnitude_->blind() ) {
const LauBlind* blinder { rhs.magnitude_->blinder() };
magnitude_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase ) {
phase_.reset( rhs.phase_->createClone(constFactor) );
} else {
phase_ = std::make_unique<LauParameter>("Delta", rhs.phase_->value(), minPhase_, maxPhase_, rhs.phase_->fixed());
if ( rhs.phase_->blind() ) {
const LauBlind* blinder { rhs.phase_->blinder() };
phase_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauMagPhaseCoeffSet::getParameters()
{
return { magnitude_.get(), phase_.get() };
}
std::vector<const LauParameter*> LauMagPhaseCoeffSet::getParameters() const
{
return { magnitude_.get(), phase_.get() };
}
void LauMagPhaseCoeffSet::printParValues() const
{
std::cout<<"INFO in LauMagPhaseCoeffSet::printParValues : Component \""<<this->name()<<"\" has magnitude = "<<magnitude_->value()<<" and phase = "<<phase_->value()<<"."<<std::endl;
}
void LauMagPhaseCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Magnitude & Phase \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauMagPhaseCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, magnitude_->value());
stream<<" \\pm ";
print.printFormat(stream, magnitude_->error());
stream<<"$ & $";
print.printFormat(stream, phase_->value());
stream<<" \\pm ";
print.printFormat(stream, phase_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauMagPhaseCoeffSet::randomiseInitValues()
{
if (magnitude_->fixed() == kFALSE) {
// Choose a magnitude between 0.0 and 2.0
const Double_t mag { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
magnitude_->initValue(mag); magnitude_->value(mag);
}
if (phase_->fixed() == kFALSE) {
// Choose a phase between +- pi
const Double_t phase { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
phase_->initValue(phase); phase_->value(phase);
}
}
void LauMagPhaseCoeffSet::finaliseValues()
{
// retrieve the current values from the parameters
Double_t mag { magnitude_->value() };
Double_t phase { phase_->value() };
// Check whether we have a negative magnitude.
// If so make it positive and add pi to the phase.
if (mag < 0.0) {
mag *= -1.0;
phase += LauConstants::pi;
}
// Check now whether the phase lies in the right range (-pi to pi).
Bool_t withinRange{kFALSE};
while (withinRange == kFALSE) {
if (phase > -LauConstants::pi && phase <= LauConstants::pi) {
withinRange = kTRUE;
} else {
// Not within the specified range
if (phase > LauConstants::pi) {
phase -= LauConstants::twoPi;
} else if (phase <= -LauConstants::pi) {
phase += LauConstants::twoPi;
}
}
}
// A further problem can occur when the generated phase is close to -pi or pi.
// The phase can wrap over to the other end of the scale -
// this leads to artificially large pulls so we wrap it back.
const Double_t genPhase { phase_->genValue() };
const Double_t diff { phase - genPhase };
if (diff > LauConstants::pi) {
phase -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
phase += LauConstants::twoPi;
}
// finally store the new values in the parameters
// and update the pulls
magnitude_->value(mag); magnitude_->updatePull();
phase_->value(phase); phase_->updatePull();
}
const LauComplex& LauMagPhaseCoeffSet::particleCoeff()
{
coeff_.setRealImagPart(magnitude_->unblindValue()*TMath::Cos(phase_->unblindValue()), magnitude_->unblindValue()*TMath::Sin(phase_->unblindValue()));
return coeff_;
}
const LauComplex& LauMagPhaseCoeffSet::antiparticleCoeff()
{
return this->particleCoeff();
}
void LauMagPhaseCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, const Bool_t init )
{
LauComplex average{ coeff };
average += coeffBar;
average.rescale( 0.5 );
const Double_t magVal{ average.abs() };
const Double_t phaseVal{ average.arg() };
magnitude_->value( magVal );
phase_->value( phaseVal );
if ( init ) {
magnitude_->genValue( magVal );
phase_->genValue( phaseVal );
magnitude_->initValue( magVal );
phase_->initValue( phaseVal );
}
}
LauParameter LauMagPhaseCoeffSet::acp()
{
const TString parName{ this->baseName() + "_ACP" };
return LauParameter{ parName, 0.0 };
}
LauMagPhaseCoeffSet* LauMagPhaseCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude ) ) {
std::cerr << "ERROR in LauMagPhaseCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauMagPhaseCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauMagPhaseCoeffSet nlohmann::adl_serializer<LauMagPhaseCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::MagPhase ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauMagPhaseCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t mag { j.at("mag").get<Double_t>() };
const Double_t phase { j.at("phase").get<Double_t>() };
const Bool_t magFixed { j.at("magFixed").get<Bool_t>() };
const Bool_t phaseFixed { j.at("phaseFixed").get<Bool_t>() };
return { name, mag, phase, magFixed, phaseFixed };
}
void nlohmann::adl_serializer<LauMagPhaseCoeffSet>::to_json(json& j, const LauMagPhaseCoeffSet& t)
{
j["type"] = LauCoeffType::MagPhase;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["mag"] = pars[0]->value();
j["phase"] = pars[1]->value();
j["magFixed"] = pars[0]->fixed();
j["phaseFixed"] = pars[1]->fixed();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauNSCCartesianCPCoeffSet.cc b/src/LauNSCCartesianCPCoeffSet.cc
index 1fcf04d..0338f9a 100644
--- a/src/LauNSCCartesianCPCoeffSet.cc
+++ b/src/LauNSCCartesianCPCoeffSet.cc
@@ -1,406 +1,407 @@
/*
Copyright 2015 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 LauNSCCartesianCPCoeffSet.cc
\brief File containing implementation of LauNSCCartesianCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauNSCCartesianCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
#include "LauRandom.hh"
ClassImp(LauNSCCartesianCPCoeffSet)
LauNSCCartesianCPCoeffSet::LauNSCCartesianCPCoeffSet(const TString& compName,
const Bool_t finalStateIsF,
const Double_t x, const Double_t y, const Double_t deltaX, const Double_t deltaY,
const Bool_t xFixed, const Bool_t yFixed, const Bool_t deltaXFixed, const Bool_t deltaYFixed,
const Bool_t deltaXSecondStage, const Bool_t deltaYSecondStage,
const Double_t xPrime, const Double_t yPrime, const Double_t deltaXPrime, const Double_t deltaYPrime,
const Bool_t xPrimeFixed, const Bool_t yPrimeFixed, const Bool_t deltaXPrimeFixed, const Bool_t deltaYPrimeFixed,
const Bool_t deltaXPrimeSecondStage, const Bool_t deltaYPrimeSecondStage ) :
LauAbsCoeffSet{ compName },
finalStateIsF_{ finalStateIsF },
x_{ std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed) },
y_{ std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) },
deltaX_{ std::make_unique<LauParameter>("DeltaX", deltaX, minDelta_, maxDelta_, deltaXFixed) },
deltaY_{ std::make_unique<LauParameter>("DeltaY", deltaY, minDelta_, maxDelta_, deltaYFixed) },
xPrime_{ std::make_unique<LauParameter>("XPrime", xPrime, minRealImagPart_, maxRealImagPart_, xPrimeFixed) },
yPrime_{ std::make_unique<LauParameter>("YPrime", yPrime, minRealImagPart_, maxRealImagPart_, yPrimeFixed) },
deltaXPrime_{ std::make_unique<LauParameter>("DeltaXPrime", deltaXPrime, minDelta_, maxDelta_, deltaXPrimeFixed) },
deltaYPrime_{ std::make_unique<LauParameter>("DeltaYPrime", deltaYPrime, minDelta_, maxDelta_, deltaYPrimeFixed) },
acp_{ "ACP", 0.0, -1.0, 1.0, kTRUE }
{
if (deltaXSecondStage && !deltaXFixed) {
deltaX_->secondStage(kTRUE);
deltaX_->initValue(0.0);
}
if (deltaYSecondStage && !deltaYFixed) {
deltaY_->secondStage(kTRUE);
deltaY_->initValue(0.0);
}
if (deltaXPrimeSecondStage && !deltaXPrimeFixed) {
deltaXPrime_->secondStage(kTRUE);
deltaXPrime_->initValue(0.0);
}
if (deltaYPrimeSecondStage && !deltaYPrimeFixed) {
deltaYPrime_->secondStage(kTRUE);
deltaYPrime_->initValue(0.0);
}
}
LauNSCCartesianCPCoeffSet::LauNSCCartesianCPCoeffSet(const LauNSCCartesianCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
finalStateIsF_{ rhs.finalStateIsF_ },
coeffAf_{ rhs.coeffAf_ },
coeffAfbar_{ rhs.coeffAfbar_ },
coeffAbarf_{ rhs.coeffAbarf_ },
coeffAbarfbar_{ rhs.coeffAbarfbar_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
xPrime_.reset( rhs.xPrime_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
xPrime_ = std::make_unique<LauParameter>("XPrime", rhs.xPrime_->value(), minRealImagPart_, maxRealImagPart_, rhs.xPrime_->fixed());
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
yPrime_.reset( rhs.yPrime_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
yPrime_ = std::make_unique<LauParameter>("YPrime", rhs.yPrime_->value(), minRealImagPart_, maxRealImagPart_, rhs.yPrime_->fixed());
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
deltaX_.reset( rhs.deltaX_->createClone(constFactor) );
deltaY_.reset( rhs.deltaY_->createClone(constFactor) );
deltaXPrime_.reset( rhs.deltaXPrime_->createClone(constFactor) );
deltaYPrime_.reset( rhs.deltaYPrime_->createClone(constFactor) );
} else {
deltaX_ = std::make_unique<LauParameter>("DeltaX", rhs.deltaX_->value(), minDelta_, maxDelta_, rhs.deltaX_->fixed());
deltaY_ = std::make_unique<LauParameter>("DeltaY", rhs.deltaY_->value(), minDelta_, maxDelta_, rhs.deltaY_->fixed());
deltaXPrime_ = std::make_unique<LauParameter>("DeltaXPrime", rhs.deltaXPrime_->value(), minDelta_, maxDelta_, rhs.deltaXPrime_->fixed());
deltaYPrime_ = std::make_unique<LauParameter>("DeltaYPrime", rhs.deltaYPrime_->value(), minDelta_, maxDelta_, rhs.deltaYPrime_->fixed());
if ( rhs.deltaX_->secondStage() && !rhs.deltaX_->fixed() ) {
deltaX_->secondStage(kTRUE);
deltaX_->initValue(0.0);
}
if ( rhs.deltaY_->secondStage() && !rhs.deltaY_->fixed() ) {
deltaY_->secondStage(kTRUE);
deltaY_->initValue(0.0);
}
if ( rhs.deltaXPrime_->secondStage() && !rhs.deltaXPrime_->fixed() ) {
deltaXPrime_->secondStage(kTRUE);
deltaXPrime_->initValue(0.0);
}
if ( rhs.deltaYPrime_->secondStage() && !rhs.deltaYPrime_->fixed() ) {
deltaYPrime_->secondStage(kTRUE);
deltaYPrime_->initValue(0.0);
}
}
}
std::vector<LauParameter*> LauNSCCartesianCPCoeffSet::getParameters()
{
return { x_.get(), y_.get(), deltaX_.get(), deltaY_.get(), xPrime_.get(), yPrime_.get(), deltaXPrime_.get(), deltaYPrime_.get() };
}
std::vector<const LauParameter*> LauNSCCartesianCPCoeffSet::getParameters() const
{
return { x_.get(), y_.get(), deltaX_.get(), deltaY_.get(), xPrime_.get(), yPrime_.get(), deltaXPrime_.get(), deltaYPrime_.get() };
}
void LauNSCCartesianCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauNSCCartesianCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"x = "<<x_->value()<<",\t";
std::cout<<"y = "<<y_->value()<<",\t";
std::cout<<"Delta x = "<<deltaX_->value()<<",\t";
std::cout<<"Delta y = "<<deltaY_->value()<<",\t";
std::cout<<"xPrime = "<<xPrime_->value()<<",\t";
std::cout<<"yPrime = "<<yPrime_->value()<<",\t";
std::cout<<"Delta xPrime = "<<deltaXPrime_->value()<<",\t";
std::cout<<"Delta yPrime = "<<deltaYPrime_->value()<<"."<<std::endl;
}
void LauNSCCartesianCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $\\Delta$ Real Part & $\\Delta$ Imaginary Part & Real Part Prime & Imaginary Part Prime & $\\Delta$ Real Part Prime & $\\Delta$ Imaginary Part Prime \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauNSCCartesianCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
print.printFormat(stream, deltaX_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaX_->error());
stream<<"$ & $";
print.printFormat(stream, deltaY_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaY_->error());
stream<<"$ & $";
print.printFormat(stream, xPrime_->value());
stream<<" \\pm ";
print.printFormat(stream, xPrime_->error());
stream<<"$ & $";
print.printFormat(stream, yPrime_->value());
stream<<" \\pm ";
print.printFormat(stream, yPrime_->error());
stream<<"$ & $";
print.printFormat(stream, deltaXPrime_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaXPrime_->error());
stream<<"$ & $";
print.printFormat(stream, deltaYPrime_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaYPrime_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauNSCCartesianCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*6.0 - 3.0 };
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*6.0 - 3.0 };
y_->initValue(value); y_->value(value);
}
if (deltaX_->fixed() == kFALSE && deltaX_->secondStage() == kFALSE) {
// Choose a value for "Delta X" between -0.5 and 0.5
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*1.0 - 0.5 };
deltaX_->initValue(value); deltaX_->value(value);
}
if (deltaY_->fixed() == kFALSE && deltaY_->secondStage() == kFALSE) {
// Choose a value for "Delta Y" between -0.5 and 0.5
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*1.0 - 0.5 };
deltaY_->initValue(value); deltaY_->value(value);
}
if (xPrime_->fixed() == kFALSE) {
// Choose a value for "XPrime" between -3.0 and 3.0
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*6.0 - 3.0 };
xPrime_->initValue(value); xPrime_->value(value);
}
if (yPrime_->fixed() == kFALSE) {
// Choose a value for "YPrime" between -3.0 and 3.0
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*6.0 - 3.0 };
yPrime_->initValue(value); yPrime_->value(value);
}
if (deltaXPrime_->fixed() == kFALSE && deltaXPrime_->secondStage() == kFALSE) {
// Choose a value for "Delta XPrime" between -0.5 and 0.5
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*1.0 - 0.5 };
deltaXPrime_->initValue(value); deltaXPrime_->value(value);
}
if (deltaYPrime_->fixed() == kFALSE && deltaYPrime_->secondStage() == kFALSE) {
// Choose a value for "Delta YPrime" between -0.5 and 0.5
const Double_t value { LauRandom::zeroSeedRandom()->Rndm()*1.0 - 0.5 };
deltaYPrime_->initValue(value); deltaYPrime_->value(value);
}
}
void LauNSCCartesianCPCoeffSet::finaliseValues()
{
// update the pulls
x_->updatePull();
y_->updatePull();
deltaX_->updatePull();
deltaY_->updatePull();
xPrime_->updatePull();
yPrime_->updatePull();
deltaXPrime_->updatePull();
deltaYPrime_->updatePull();
}
const LauComplex& LauNSCCartesianCPCoeffSet::particleCoeff()
{
if ( finalStateIsF_ ) {
coeffAf_.setRealImagPart( x_->value() + deltaX_->value(), y_->value() + deltaY_->value() );
return coeffAf_;
} else {
coeffAfbar_.setRealImagPart( xPrime_->value() + deltaXPrime_->value(), yPrime_->value() + deltaYPrime_->value() );
return coeffAfbar_;
}
}
const LauComplex& LauNSCCartesianCPCoeffSet::antiparticleCoeff()
{
if ( finalStateIsF_ ) {
coeffAbarf_.setRealImagPart( xPrime_->value() - deltaXPrime_->value(), yPrime_->value() - deltaYPrime_->value() );
return coeffAbarf_;
} else {
coeffAbarfbar_.setRealImagPart( x_->value() - deltaX_->value(), y_->value() - deltaY_->value() );
return coeffAbarfbar_;
}
}
void LauNSCCartesianCPCoeffSet::setCoeffValues( const LauComplex&, const LauComplex&, Bool_t )
{
std::cerr << "ERROR in LauNSCCartesianCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters" << std::endl;
}
LauParameter LauNSCCartesianCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// a single ACP parameter doesn't really make much sense here
const Double_t value{0.0};
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauNSCCartesianCPCoeffSet* LauNSCCartesianCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauNSCCartesianCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauNSCCartesianCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
-
+//! \cond DOXYGEN_IGNORE
LauNSCCartesianCPCoeffSet nlohmann::adl_serializer<LauNSCCartesianCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::NSCCartesianCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauNSCCartesianCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
const Bool_t finalStateIsF { j.at("finalStateIsF").get<Bool_t>() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t deltaX { j.at("deltaX").get<Double_t>() };
const Double_t deltaY { j.at("deltaY").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t deltaXFixed { j.at("deltaXFixed").get<Bool_t>() };
const Bool_t deltaYFixed { j.at("deltaYFixed").get<Bool_t>() };
const Bool_t deltaXSecondStage { j.at("deltaXSecondStage").get<Bool_t>() };
const Bool_t deltaYSecondStage { j.at("deltaYSecondStage").get<Bool_t>() };
const Double_t xPrime { j.at("xPrime").get<Double_t>() };
const Double_t yPrime { j.at("yPrime").get<Double_t>() };
const Double_t deltaXPrime { j.at("deltaXPrime").get<Double_t>() };
const Double_t deltaYPrime { j.at("deltaYPrime").get<Double_t>() };
const Bool_t xPrimeFixed { j.at("xPrimeFixed").get<Bool_t>() };
const Bool_t yPrimeFixed { j.at("yPrimeFixed").get<Bool_t>() };
const Bool_t deltaXPrimeFixed { j.at("deltaXPrimeFixed").get<Bool_t>() };
const Bool_t deltaYPrimeFixed { j.at("deltaYPrimeFixed").get<Bool_t>() };
const Bool_t deltaXPrimeSecondStage { j.at("deltaXPrimeSecondStage").get<Bool_t>() };
const Bool_t deltaYPrimeSecondStage { j.at("deltaYPrimeSecondStage").get<Bool_t>() };
return { name, finalStateIsF,
x, y, deltaX, deltaY,
xFixed, yFixed, deltaXFixed, deltaYFixed,
deltaXSecondStage, deltaYSecondStage,
xPrime, yPrime, deltaXPrime, deltaYPrime,
xPrimeFixed, yPrimeFixed, deltaXPrimeFixed, deltaYPrimeFixed,
deltaXPrimeSecondStage, deltaYPrimeSecondStage };
}
void nlohmann::adl_serializer<LauNSCCartesianCPCoeffSet>::to_json(json& j, const LauNSCCartesianCPCoeffSet& t)
{
j["type"] = LauCoeffType::NSCCartesianCP;
j["name"] = t.name();
j["finalStateIsF"] = t.finalStateIsF();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["deltaX"] = pars[2]->value();
j["deltaY"] = pars[3]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["deltaXFixed"] = pars[2]->fixed();
j["deltaYFixed"] = pars[3]->fixed();
j["deltaXSecondStage"] = pars[2]->secondStage();
j["deltaYSecondStage"] = pars[3]->secondStage();
j["xPrime"] = pars[4]->value();
j["yPrime"] = pars[5]->value();
j["deltaXPrime"] = pars[6]->value();
j["deltaYPrime"] = pars[7]->value();
j["xPrimeFixed"] = pars[4]->fixed();
j["yPrimeFixed"] = pars[5]->fixed();
j["deltaXPrimeFixed"] = pars[6]->fixed();
j["deltaYPrimeFixed"] = pars[7]->fixed();
j["deltaXPrimeSecondStage"] = pars[6]->secondStage();
j["deltaYPrimeSecondStage"] = pars[7]->secondStage();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauPolarGammaCPCoeffSet.cc b/src/LauPolarGammaCPCoeffSet.cc
index d86e64a..3e80b78 100644
--- a/src/LauPolarGammaCPCoeffSet.cc
+++ b/src/LauPolarGammaCPCoeffSet.cc
@@ -1,890 +1,892 @@
/*
Copyright 2014 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 LauPolarGammaCPCoeffSet.cc
\brief File containing implementation of LauPolarGammaCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauPolarGammaCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
std::shared_ptr<LauParameter> LauPolarGammaCPCoeffSet::gammaGlobal_;
std::shared_ptr<LauParameter> LauPolarGammaCPCoeffSet::rDGlobal_;
std::shared_ptr<LauParameter> LauPolarGammaCPCoeffSet::deltaDGlobal_;
ClassImp(LauPolarGammaCPCoeffSet)
LauPolarGammaCPCoeffSet::LauPolarGammaCPCoeffSet(const TString& compName, const DecayType decayType,
const Double_t x, const Double_t y,
const Double_t rB, const Double_t deltaB, const Double_t gamma,
const Double_t rD, const Double_t deltaD,
const Bool_t xFixed, const Bool_t yFixed,
const Bool_t rBFixed, const Bool_t deltaBFixed, const Bool_t gammaFixed,
const Bool_t rDFixed, const Bool_t deltaDFixed,
const Bool_t rBSecondStage, const Bool_t deltaBSecondStage, const Bool_t gammaSecondStage,
const Bool_t rDSecondStage, const Bool_t deltaDSecondStage,
const Bool_t useGlobalGamma,
const Bool_t useGlobalADSPars) :
LauAbsCoeffSet{ compName },
decayType_{ decayType },
useGlobalGamma_{ useGlobalGamma },
useGlobalADSPars_{ useGlobalADSPars },
acp_{ "ACP", 0.0, -1.0, 1.0 }
{
// All of the possible D decay types need these two parameters
x_ = std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed);
y_ = std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed);
// if we're using a global gamma, create it if it doesn't already exist then set gamma_ to point to it
// otherwise create our individual copy of gamma
if (useGlobalGamma_) {
if (!gammaGlobal_) {
gammaGlobal_ = std::make_shared<LauParameter>("gamma", gamma, minPhase_, maxPhase_, gammaFixed);
gamma_ = gammaGlobal_;
} else {
gamma_.reset( gammaGlobal_->createClone() );
}
} else {
gamma_ = std::make_shared<LauParameter>("gamma", gamma, minPhase_, maxPhase_, gammaFixed);
}
if (gammaSecondStage && !gammaFixed) {
gamma_->secondStage(kTRUE);
gamma_->initValue(0.0);
}
// which of the other parameter we need depends on the D-decay type
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
rB_ = std::make_unique<LauParameter>("rB", rB, minMagnitude_, maxMagnitude_, rBFixed);
deltaB_ = std::make_unique<LauParameter>("deltaB", deltaB, minPhase_, maxPhase_, deltaBFixed);
if (rBSecondStage && !rBFixed) {
rB_->secondStage(kTRUE);
rB_->initValue(0.0);
}
if (deltaBSecondStage && !deltaBFixed) {
deltaB_->secondStage(kTRUE);
deltaB_->initValue(0.0);
}
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
if (useGlobalADSPars_) {
if ( !rDGlobal_ ) {
rDGlobal_ = std::make_shared<LauParameter>("rD", rD, minMagnitude_, maxMagnitude_, rDFixed);
deltaDGlobal_ = std::make_shared<LauParameter>("deltaD", deltaD, minPhase_, maxPhase_, deltaDFixed);
rD_ = rDGlobal_;
deltaD_ = deltaDGlobal_;
} else {
rD_.reset( rDGlobal_->createClone() );
deltaD_.reset( deltaDGlobal_->createClone() );
}
} else {
rD_ = std::make_shared<LauParameter>("rD", rD, minMagnitude_, maxMagnitude_, rDFixed);
deltaD_ = std::make_shared<LauParameter>("deltaD", deltaD, minPhase_, maxPhase_, deltaDFixed);
}
if (rDSecondStage && !rDFixed) {
rD_->secondStage(kTRUE);
rD_->initValue(0.0);
}
if (deltaDSecondStage && !deltaDFixed) {
deltaD_->secondStage(kTRUE);
deltaD_->initValue(0.0);
}
}
}
LauPolarGammaCPCoeffSet::LauPolarGammaCPCoeffSet(const LauPolarGammaCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
decayType_{ rhs.decayType_ },
useGlobalGamma_{ rhs.useGlobalGamma_ },
useGlobalADSPars_{ rhs.useGlobalADSPars_ },
nonCPPart_{ rhs.nonCPPart_ },
cpPart_{ rhs.cpPart_ },
cpAntiPart_{ rhs.cpAntiPart_ },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
gamma_.reset( rhs.gamma_->createClone(constFactor) );
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
rB_.reset( rhs.rB_->createClone(constFactor) );
deltaB_.reset( rhs.deltaB_->createClone(constFactor) );
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
rD_.reset( rhs.rD_->createClone(constFactor) );
deltaD_.reset( rhs.deltaD_->createClone(constFactor) );
}
} else {
if (useGlobalGamma_) {
gamma_.reset( gammaGlobal_->createClone() );
} else {
gamma_ = std::make_shared<LauParameter>("gamma", rhs.gamma_->value(), minPhase_, maxPhase_, rhs.gamma_->fixed());
if ( rhs.gamma_->blind() ) {
const LauBlind* blinder { rhs.gamma_->blinder() };
gamma_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
if ( rhs.gamma_->secondStage() && !rhs.gamma_->fixed() ) {
gamma_->secondStage(kTRUE);
gamma_->initValue(0.0);
}
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
rB_ = std::make_unique<LauParameter>("rB", rhs.rB_->value(), minMagnitude_, maxMagnitude_, rhs.rB_->fixed());
if ( rhs.rB_->blind() ) {
const LauBlind* blinder { rhs.rB_->blinder() };
rB_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
deltaB_ = std::make_unique<LauParameter>("deltaB", rhs.deltaB_->value(), minPhase_, maxPhase_, rhs.deltaB_->fixed());
if ( rhs.deltaB_->blind() ) {
const LauBlind* blinder { rhs.deltaB_->blinder() };
deltaB_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
if ( rhs.rB_->secondStage() && !rhs.rB_->fixed() ) {
rB_->secondStage(kTRUE);
rB_->initValue(0.0);
}
if ( rhs.deltaB_->secondStage() && !rhs.deltaB_->fixed() ) {
deltaB_->secondStage(kTRUE);
deltaB_->initValue(0.0);
}
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
if ( useGlobalADSPars_ ) {
rD_.reset( rDGlobal_->createClone() );
deltaD_.reset( deltaDGlobal_->createClone() );
} else {
rD_ = std::make_shared<LauParameter>("rD", rhs.rD_->value(), minMagnitude_, maxMagnitude_, rhs.rD_->fixed());
if ( rhs.rD_->blind() ) {
const LauBlind* blinder { rhs.rD_->blinder() };
rD_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
deltaD_ = std::make_shared<LauParameter>("deltaD", rhs.deltaD_->value(), minPhase_, maxPhase_, rhs.deltaD_->fixed());
if ( rhs.deltaD_->blind() ) {
const LauBlind* blinder { rhs.deltaD_->blinder() };
deltaD_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
if ( rhs.rD_->secondStage() && !rhs.rD_->fixed() ) {
rD_->secondStage(kTRUE);
rD_->initValue(0.0);
}
if ( rhs.deltaD_->secondStage() && !rhs.deltaD_->fixed() ) {
deltaD_->secondStage(kTRUE);
deltaD_->initValue(0.0);
}
}
}
}
}
void LauPolarGammaCPCoeffSet::adjustName(LauParameter& par, const TString& oldBaseName)
{
if ( ( &par == gamma_.get() && useGlobalGamma_ ) ||
( &par == rD_.get() && useGlobalADSPars_ ) ||
( &par == deltaD_.get() && useGlobalADSPars_ ) ) {
// for global parameters we do not want to adjust their names
return;
} else {
LauAbsCoeffSet::adjustName(par,oldBaseName);
}
}
std::vector<LauParameter*> LauPolarGammaCPCoeffSet::getParameters()
{
std::vector<LauParameter*> pars;
pars.reserve(7);
pars.push_back(x_.get());
pars.push_back(y_.get());
pars.push_back(gamma_.get());
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
pars.push_back(rB_.get());
pars.push_back(deltaB_.get());
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
pars.push_back(rD_.get());
pars.push_back(deltaD_.get());
}
return pars;
}
std::vector<const LauParameter*> LauPolarGammaCPCoeffSet::getParameters() const
{
std::vector<const LauParameter*> pars;
pars.reserve(7);
pars.push_back(x_.get());
pars.push_back(y_.get());
pars.push_back(gamma_.get());
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
pars.push_back(rB_.get());
pars.push_back(deltaB_.get());
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
pars.push_back(rD_.get());
pars.push_back(deltaD_.get());
}
return pars;
}
void LauPolarGammaCPCoeffSet::printParValues() const
{
std::cout<<"INFO in LauPolarGammaCPCoeffSet::printParValues : Component \""<<this->name()<<"\" has ";
std::cout<<"x = "<<x_->value()<<",\t";
std::cout<<"y = "<<y_->value()<<",\t";
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
std::cout<<"rB = "<<rB_->value()<<",\t";
std::cout<<"deltaB = "<<deltaB_->value()<<",\t";
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
std::cout<<"rD = "<<rD_->value()<<",\t";
std::cout<<"deltaD = "<<deltaD_->value()<<",\t";
}
std::cout<<"gamma = "<<gamma_->value()<<"."<<std::endl;
}
void LauPolarGammaCPCoeffSet::printTableHeading(std::ostream& stream) const
{
switch ( decayType_ ) {
case DecayType::GLW_CPOdd :
stream<<"\\begin{tabular}{|l|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $r_B$ & $\\delta_B$ & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::GLW_CPEven :
stream<<"\\begin{tabular}{|l|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $r_B$ & $\\delta_B$ & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::ADS_Favoured :
stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $r_B$ & $\\delta_B$ & $r_D$ & $\\delta_D$ & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::ADS_Suppressed :
stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $r_B$ & $\\delta_B$ & $r_D$ & $\\delta_D$ & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::GLW_CPOdd_btouOnly :
stream<<"\\begin{tabular}{|l|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::GLW_CPEven_btouOnly :
stream<<"\\begin{tabular}{|l|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::ADS_Favoured_btouOnly :
stream<<"\\begin{tabular}{|l|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $r_D$ & $\\delta_D$ & $\\gamma$ \\\\"<<std::endl;
break;
case DecayType::ADS_Suppressed_btouOnly :
stream<<"\\begin{tabular}{|l|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & $\\gamma$ \\\\"<<std::endl;
break;
}
stream<<"\\hline"<<std::endl;
}
void LauPolarGammaCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
print.printFormat(stream, rB_->value());
stream<<" \\pm ";
print.printFormat(stream, rB_->error());
stream<<"$ & $";
print.printFormat(stream, deltaB_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaB_->error());
stream<<"$ & $";
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
print.printFormat(stream, rD_->value());
stream<<" \\pm ";
print.printFormat(stream, rD_->error());
stream<<"$ & $";
print.printFormat(stream, deltaD_->value());
stream<<" \\pm ";
print.printFormat(stream, deltaD_->error());
stream<<"$ & $";
}
print.printFormat(stream, gamma_->value());
stream<<" \\pm ";
print.printFormat(stream, gamma_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauPolarGammaCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
y_->initValue(value); y_->value(value);
}
if (gamma_->fixed() == kFALSE && gamma_->secondStage() == kFALSE) {
// Choose a value for "gamma" between +-pi
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
gamma_->initValue(value); gamma_->value(value);
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
if (rB_->fixed() == kFALSE && rB_->secondStage() == kFALSE) {
// Choose a value for "rB" between 0.0 and 2.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
rB_->initValue(value); rB_->value(value);
}
if (deltaB_->fixed() == kFALSE && deltaB_->secondStage() == kFALSE) {
// Choose a value for "deltaB" between +- pi
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
deltaB_->initValue(value); deltaB_->value(value);
}
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
if (rD_->fixed() == kFALSE && rD_->secondStage() == kFALSE) {
// Choose a value for "rD" between 0.0 and 2.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*2.0 };
rD_->initValue(value); rD_->value(value);
}
if (deltaD_->fixed() == kFALSE && deltaD_->secondStage() == kFALSE) {
// Choose a value for "deltaD" between +- pi
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*LauConstants::twoPi - LauConstants::pi };
deltaD_->initValue(value); deltaD_->value(value);
}
}
}
void LauPolarGammaCPCoeffSet::finaliseValues()
{
// retrieve the current values from the parameters
Double_t gammaVal { gamma_->value() };
Double_t rBVal { 0.0 };
Double_t deltaBVal { 0.0 };
Double_t genDeltaB { 0.0 };
Double_t rDVal { 0.0 };
Double_t deltaDVal { 0.0 };
Double_t genDeltaD { 0.0 };
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
rBVal = rB_->value();
deltaBVal = deltaB_->value();
genDeltaB = deltaB_->genValue();
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
rDVal = rD_->value();
deltaDVal = deltaD_->value();
genDeltaD = deltaD_->genValue();
}
// Check whether we have a negative magnitude.
// If so make it positive and add pi to the phases.
if (rBVal < 0.0) {
rBVal *= -1.0;
deltaBVal += LauConstants::pi;
}
if (rDVal < 0.0) {
rDVal *= -1.0;
deltaDVal += LauConstants::pi;
}
// Check now whether the phases lie in the right range (-pi to pi).
Bool_t deltaBWithinRange{kFALSE};
Bool_t deltaDWithinRange{kFALSE};
Bool_t gammaWithinRange{kFALSE};
while ( deltaBWithinRange == kFALSE ) {
if (deltaBVal > -LauConstants::pi && deltaBVal <= LauConstants::pi) {
deltaBWithinRange = kTRUE;
} else {
// Not within the specified range
if (deltaBVal > LauConstants::pi) {
deltaBVal -= LauConstants::twoPi;
} else if (deltaBVal <= -LauConstants::pi) {
deltaBVal += LauConstants::twoPi;
}
}
}
while ( deltaDWithinRange == kFALSE ) {
if (deltaDVal > -LauConstants::pi && deltaDVal <= LauConstants::pi) {
deltaDWithinRange = kTRUE;
} else {
// Not within the specified range
if (deltaDVal > LauConstants::pi) {
deltaDVal -= LauConstants::twoPi;
} else if (deltaDVal <= -LauConstants::pi) {
deltaDVal += LauConstants::twoPi;
}
}
}
while ( gammaWithinRange == kFALSE ) {
if (gammaVal > -LauConstants::pi && gammaVal <= LauConstants::pi) {
gammaWithinRange = kTRUE;
} else {
// Not within the specified range
if (gammaVal > LauConstants::pi) {
gammaVal -= LauConstants::twoPi;
} else if (gammaVal <= -LauConstants::pi) {
gammaVal += LauConstants::twoPi;
}
}
}
// To resolve the two-fold ambiguity in gamma and deltaB we require gamma to be in the range 0-pi
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
if (gammaVal < 0.0) {
if (deltaBVal <= 0.0) {
gammaVal += LauConstants::pi;
deltaBVal += LauConstants::pi;
} else {
gammaVal += LauConstants::pi;
deltaBVal -= LauConstants::pi;
}
}
}
// A further problem can occur when the generated phase is close to -pi or pi.
// The phase can wrap over to the other end of the scale -
// this leads to artificially large pulls so we wrap it back.
Double_t diff { deltaBVal - genDeltaB };
if (diff > LauConstants::pi) {
deltaBVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
deltaBVal += LauConstants::twoPi;
}
diff = deltaDVal - genDeltaD;
if (diff > LauConstants::pi) {
deltaDVal -= LauConstants::twoPi;
} else if (diff < -LauConstants::pi) {
deltaDVal += LauConstants::twoPi;
}
// finally store the new values in the parameters
// and update the pulls
gamma_->value(gammaVal);
gamma_->updatePull();
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) {
rB_->value(rBVal);
rB_->updatePull();
deltaB_->value(deltaBVal);
deltaB_->updatePull();
}
if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) {
rD_->value(rDVal);
rD_->updatePull();
deltaD_->value(deltaDVal);
deltaD_->updatePull();
}
}
const LauComplex& LauPolarGammaCPCoeffSet::particleCoeff()
{
this->updateAmplitudes();
return particleCoeff_;
}
const LauComplex& LauPolarGammaCPCoeffSet::antiparticleCoeff()
{
this->updateAmplitudes();
return antiparticleCoeff_;
}
void LauPolarGammaCPCoeffSet::updateAmplitudes()
{
nonCPPart_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
const Double_t gammaVal { gamma_->unblindValue() };
switch ( decayType_ ) {
case DecayType::GLW_CPOdd :
{
const Double_t rBVal { rB_->unblindValue() };
const Double_t deltaBVal { deltaB_->unblindValue() };
cpPart_.setRealImagPart( 1.0 - rBVal*TMath::Cos(deltaBVal + gammaVal), -rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart_.setRealImagPart( 1.0 - rBVal*TMath::Cos(deltaBVal - gammaVal), -rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::GLW_CPEven :
{
const Double_t rBVal { rB_->unblindValue() };
const Double_t deltaBVal { deltaB_->unblindValue() };
cpPart_.setRealImagPart( 1.0 + rBVal*TMath::Cos(deltaBVal + gammaVal), rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart_.setRealImagPart( 1.0 + rBVal*TMath::Cos(deltaBVal - gammaVal), rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::ADS_Favoured :
{
const Double_t rBVal { rB_->unblindValue() };
const Double_t deltaBVal { deltaB_->unblindValue() };
const Double_t rDVal { rD_->unblindValue() };
const Double_t deltaDVal { deltaD_->unblindValue() };
cpPart_.setRealImagPart( 1.0 + rBVal*rDVal*TMath::Cos(deltaBVal - deltaDVal + gammaVal), rBVal*rDVal*TMath::Sin(deltaBVal - deltaDVal + gammaVal) );
cpAntiPart_.setRealImagPart( 1.0 + rBVal*rDVal*TMath::Cos(deltaBVal - deltaDVal - gammaVal), rBVal*rDVal*TMath::Sin(deltaBVal - deltaDVal - gammaVal) );
break;
}
case DecayType::ADS_Suppressed :
{
const Double_t rBVal { rB_->unblindValue() };
const Double_t deltaBVal { deltaB_->unblindValue() };
const Double_t rDVal { rD_->unblindValue() };
const Double_t deltaDVal { deltaD_->unblindValue() };
cpPart_.setRealImagPart( rDVal*TMath::Cos(-deltaDVal) + rBVal*TMath::Cos(deltaBVal + gammaVal), rDVal*TMath::Sin(-deltaDVal) + rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart_.setRealImagPart( rDVal*TMath::Cos(-deltaDVal) + rBVal*TMath::Cos(deltaBVal - gammaVal), rDVal*TMath::Sin(-deltaDVal) + rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::GLW_CPOdd_btouOnly :
nonCPPart_.rescale(-1.0);
cpPart_.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart_.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
case DecayType::GLW_CPEven_btouOnly :
cpPart_.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart_.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
case DecayType::ADS_Favoured_btouOnly :
{
const Double_t rDVal { rD_->unblindValue() };
const Double_t deltaDVal { deltaD_->unblindValue() };
cpPart_.setRealImagPart( rDVal * TMath::Cos( -deltaDVal + gammaVal ), rDVal * TMath::Sin( -deltaDVal + gammaVal ) );
cpAntiPart_.setRealImagPart( rDVal * TMath::Cos( -deltaDVal - gammaVal ), rDVal * TMath::Sin( -deltaDVal - gammaVal ) );
break;
}
case DecayType::ADS_Suppressed_btouOnly :
cpPart_.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart_.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
}
particleCoeff_ = nonCPPart_ * cpPart_;
antiparticleCoeff_ = nonCPPart_ * cpAntiPart_;
}
void LauPolarGammaCPCoeffSet::setCoeffValues( const LauComplex&, const LauComplex&, Bool_t )
{
std::cerr << "ERROR in LauPolarGammaCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters" << std::endl;
}
LauParameter LauPolarGammaCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
LauComplex nonCPPart{ x_->value(), y_->value() };
LauComplex cpPart;
LauComplex cpAntiPart;
const Double_t gammaVal { gamma_->value() };
switch ( decayType_ ) {
case DecayType::GLW_CPOdd :
{
const Double_t rBVal { rB_->value() };
const Double_t deltaBVal { deltaB_->value() };
cpPart.setRealImagPart( 1.0 - rBVal*TMath::Cos(deltaBVal + gammaVal), -rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart.setRealImagPart( 1.0 - rBVal*TMath::Cos(deltaBVal - gammaVal), -rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::GLW_CPEven :
{
const Double_t rBVal { rB_->value() };
const Double_t deltaBVal { deltaB_->value() };
cpPart.setRealImagPart( 1.0 + rBVal*TMath::Cos(deltaBVal + gammaVal), rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart.setRealImagPart( 1.0 + rBVal*TMath::Cos(deltaBVal - gammaVal), rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::ADS_Favoured :
{
const Double_t rBVal { rB_->value() };
const Double_t deltaBVal { deltaB_->value() };
const Double_t rDVal { rD_->value() };
const Double_t deltaDVal { deltaD_->value() };
cpPart.setRealImagPart( 1.0 + rBVal*rDVal*TMath::Cos(deltaBVal - deltaDVal + gammaVal), rBVal*rDVal*TMath::Sin(deltaBVal - deltaDVal + gammaVal) );
cpAntiPart.setRealImagPart( 1.0 + rBVal*rDVal*TMath::Cos(deltaBVal - deltaDVal - gammaVal), rBVal*rDVal*TMath::Sin(deltaBVal - deltaDVal - gammaVal) );
break;
}
case DecayType::ADS_Suppressed :
{
const Double_t rBVal { rB_->value() };
const Double_t deltaBVal { deltaB_->value() };
const Double_t rDVal { rD_->value() };
const Double_t deltaDVal { deltaD_->value() };
cpPart.setRealImagPart( rDVal*TMath::Cos(-deltaDVal) + rBVal*TMath::Cos(deltaBVal + gammaVal), rDVal*TMath::Sin(-deltaDVal) + rBVal*TMath::Sin(deltaBVal + gammaVal) );
cpAntiPart.setRealImagPart( rDVal*TMath::Cos(-deltaDVal) + rBVal*TMath::Cos(deltaBVal - gammaVal), rDVal*TMath::Sin(-deltaDVal) + rBVal*TMath::Sin(deltaBVal - gammaVal) );
break;
}
case DecayType::GLW_CPOdd_btouOnly :
nonCPPart.rescale(-1.0);
cpPart.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
case DecayType::GLW_CPEven_btouOnly :
cpPart.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
case DecayType::ADS_Favoured_btouOnly :
{
const Double_t rDVal { rD_->value() };
const Double_t deltaDVal { deltaD_->value() };
cpPart.setRealImagPart( rDVal * TMath::Cos( -deltaDVal + gammaVal ), rDVal * TMath::Sin( -deltaDVal + gammaVal ) );
cpAntiPart.setRealImagPart( rDVal * TMath::Cos( -deltaDVal - gammaVal ), rDVal * TMath::Sin( -deltaDVal - gammaVal ) );
break;
}
case DecayType::ADS_Suppressed_btouOnly :
cpPart.setRealImagPart( 1.0 * TMath::Cos( gammaVal ), 1.0 * TMath::Sin( gammaVal ) );
cpAntiPart.setRealImagPart( 1.0 * TMath::Cos( -gammaVal ), 1.0 * TMath::Sin( -gammaVal ) );
break;
}
const LauComplex partCoeff { nonCPPart * cpPart };
const LauComplex antiCoeff { nonCPPart * cpAntiPart };
const Double_t numer { antiCoeff.abs2() - partCoeff.abs2() };
const Double_t denom { antiCoeff.abs2() + partCoeff.abs2() };
const Double_t value { numer/denom };
// is it fixed?
const Bool_t fixed { gamma_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauPolarGammaCPCoeffSet* LauPolarGammaCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauPolarGammaCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauPolarGammaCPCoeffSet( *this, cloneOption, constFactor );
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauPolarGammaCPCoeffSet nlohmann::adl_serializer<LauPolarGammaCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::PolarGammaCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauPolarGammaCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
const LauPolarGammaCPCoeffSet::DecayType decayType { j.at("decayType").get<LauPolarGammaCPCoeffSet::DecayType>() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t rB { j.at("rB").get<Double_t>() };
const Double_t deltaB { j.at("deltaB").get<Double_t>() };
const Double_t gamma { j.at("gamma").get<Double_t>() };
const Double_t rD { j.at("rD").get<Double_t>() };
const Double_t deltaD { j.at("deltaD").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t rBFixed { j.at("rBFixed").get<Bool_t>() };
const Bool_t deltaBFixed { j.at("deltaBFixed").get<Bool_t>() };
const Bool_t gammaFixed { j.at("gammaFixed").get<Bool_t>() };
const Bool_t rDFixed { j.at("rDFixed").get<Bool_t>() };
const Bool_t deltaDFixed { j.at("deltaDFixed").get<Bool_t>() };
const Bool_t rBSecondStage { j.at("rBSecondStage").get<Bool_t>() };
const Bool_t deltaBSecondStage { j.at("deltaBSecondStage").get<Bool_t>() };
const Bool_t gammaSecondStage { j.at("gammaSecondStage").get<Bool_t>() };
const Bool_t rDSecondStage { j.at("rDSecondStage").get<Bool_t>() };
const Bool_t deltaDSecondStage { j.at("deltaDSecondStage").get<Bool_t>() };
const Bool_t useGlobalGamma { j.at("useGlobalGamma").get<Bool_t>() };
const Bool_t useGlobalADSPars { j.at("useGlobalADSPars").get<Bool_t>() };
return { name, decayType,
x, y, rB, deltaB, gamma, rD, deltaD,
xFixed, yFixed, rBFixed, deltaBFixed, gammaFixed, rDFixed, deltaDFixed,
rBSecondStage, deltaBSecondStage, gammaSecondStage, rDSecondStage, deltaDSecondStage,
useGlobalGamma, useGlobalADSPars };
}
void nlohmann::adl_serializer<LauPolarGammaCPCoeffSet>::to_json(json& j, const LauPolarGammaCPCoeffSet& t)
{
j["type"] = LauCoeffType::PolarGammaCP;
j["name"] = t.name();
const LauPolarGammaCPCoeffSet::DecayType decayType { t.decayType() };
j["decayType"] = decayType;
j["useGlobalGamma"] = t.useGlobalGamma();
j["useGlobalADSPars"] = t.useGlobalADSPars();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["gamma"] = pars[2]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["gammaFixed"] = pars[2]->fixed();
j["gammaSecondStage"] = pars[2]->secondStage();
switch ( decayType ) {
case LauPolarGammaCPCoeffSet::DecayType::GLW_CPOdd :
case LauPolarGammaCPCoeffSet::DecayType::GLW_CPEven :
j["rB"] = pars[3]->value();
j["deltaB"] = pars[4]->value();
j["rBFixed"] = pars[3]->fixed();
j["deltaBFixed"] = pars[4]->fixed();
j["rBSecondStage"] = pars[3]->secondStage();
j["deltaBSecondStage"] = pars[4]->secondStage();
break;
case LauPolarGammaCPCoeffSet::DecayType::ADS_Favoured :
case LauPolarGammaCPCoeffSet::DecayType::ADS_Suppressed :
j["rB"] = pars[3]->value();
j["deltaB"] = pars[4]->value();
j["rD"] = pars[5]->value();
j["deltaD"] = pars[6]->value();
j["rBFixed"] = pars[3]->fixed();
j["deltaBFixed"] = pars[4]->fixed();
j["rDFixed"] = pars[5]->fixed();
j["deltaDFixed"] = pars[6]->fixed();
j["rBSecondStage"] = pars[3]->secondStage();
j["deltaBSecondStage"] = pars[4]->secondStage();
j["rDSecondStage"] = pars[5]->secondStage();
j["deltaDSecondStage"] = pars[6]->secondStage();
break;
case LauPolarGammaCPCoeffSet::DecayType::GLW_CPOdd_btouOnly :
case LauPolarGammaCPCoeffSet::DecayType::GLW_CPEven_btouOnly :
case LauPolarGammaCPCoeffSet::DecayType::ADS_Suppressed_btouOnly :
break;
case LauPolarGammaCPCoeffSet::DecayType::ADS_Favoured_btouOnly :
j["rD"] = pars[3]->value();
j["deltaD"] = pars[4]->value();
j["rDFixed"] = pars[3]->fixed();
j["deltaDFixed"] = pars[4]->fixed();
j["rDSecondStage"] = pars[3]->secondStage();
j["deltaDSecondStage"] = pars[4]->secondStage();
break;
};
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauRealImagCPCoeffSet.cc b/src/LauRealImagCPCoeffSet.cc
index 2962005..9375e91 100644
--- a/src/LauRealImagCPCoeffSet.cc
+++ b/src/LauRealImagCPCoeffSet.cc
@@ -1,299 +1,300 @@
/*
Copyright 2014 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 LauRealImagCPCoeffSet.cc
\brief File containing implementation of LauRealImagCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauRealImagCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauRealImagCPCoeffSet)
LauRealImagCPCoeffSet::LauRealImagCPCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Double_t xbar, const Double_t ybar,
const Bool_t xFixed, const Bool_t yFixed, const Bool_t xbarFixed, const Bool_t ybarFixed) :
LauAbsCoeffSet{ compName },
x_{ std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed) },
y_{ std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) },
xbar_{ std::make_unique<LauParameter>("Xbar", xbar, minRealImagPart_, maxRealImagPart_, xbarFixed) },
ybar_{ std::make_unique<LauParameter>("Ybar", ybar, minRealImagPart_, maxRealImagPart_, ybarFixed) },
particleCoeff_{ x,y },
antiparticleCoeff_{ xbar,ybar },
acp_{ "ACP", 0.0, -1.0, 1.0, kTRUE }
{
}
LauRealImagCPCoeffSet::LauRealImagCPCoeffSet(const LauRealImagCPCoeffSet& rhs, CloneOption cloneOption, Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
xbar_.reset( rhs.xbar_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
xbar_ = std::make_unique<LauParameter>("Xbar", rhs.xbar_->value(), minRealImagPart_, maxRealImagPart_, rhs.xbar_->fixed());
if ( rhs.xbar_->blind() ) {
const LauBlind* blinder { rhs.xbar_->blinder() };
xbar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
ybar_.reset( rhs.ybar_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
ybar_ = std::make_unique<LauParameter>("Ybar", rhs.ybar_->value(), minRealImagPart_, maxRealImagPart_, rhs.ybar_->fixed());
if ( rhs.ybar_->blind() ) {
const LauBlind* blinder { rhs.ybar_->blinder() };
ybar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauRealImagCPCoeffSet::getParameters()
{
return { x_.get(), y_.get(), xbar_.get(), ybar_.get() };
}
std::vector<const LauParameter*> LauRealImagCPCoeffSet::getParameters() const
{
return { x_.get(), y_.get(), xbar_.get(), ybar_.get() };
}
void LauRealImagCPCoeffSet::printParValues() const
{
std::cout << "INFO in LauRealImagCPCoeffSet::printParValues : Component \"" << this->name() << "\" has ";
std::cout << "x = " << x_->value() << ",\t";
std::cout << "y = " << y_->value() << ",\t";
std::cout << "xbar = " << xbar_->value() << ",\t";
std::cout << "ybar = " << ybar_->value() << "." << std::endl;
}
void LauRealImagCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Particle Real Part & Particle Imaginary Part & Antiparticle Real Part & Antiparticle Imaginary Part \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauRealImagCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
print.printFormat(stream, xbar_->value());
stream<<" \\pm ";
print.printFormat(stream, xbar_->error());
stream<<"$ & $";
print.printFormat(stream, ybar_->value());
stream<<" \\pm ";
print.printFormat(stream, ybar_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauRealImagCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
y_->initValue(value); y_->value(value);
}
if (xbar_->fixed() == kFALSE) {
// Choose a value for "Xbar" between -3.0 and 3.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
xbar_->initValue(value); xbar_->value(value);
}
if (ybar_->fixed() == kFALSE) {
// Choose a value for "Ybar" between -3.0 and 3.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
ybar_->initValue(value); ybar_->value(value);
}
}
void LauRealImagCPCoeffSet::finaliseValues()
{
// update the pulls
x_->updatePull();
y_->updatePull();
xbar_->updatePull();
ybar_->updatePull();
}
const LauComplex& LauRealImagCPCoeffSet::particleCoeff()
{
particleCoeff_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
return particleCoeff_;
}
const LauComplex& LauRealImagCPCoeffSet::antiparticleCoeff()
{
antiparticleCoeff_.setRealImagPart( xbar_->unblindValue(), ybar_->unblindValue() );
return antiparticleCoeff_;
}
void LauRealImagCPCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, Bool_t init )
{
const Double_t xVal{ coeff.re() };
const Double_t yVal{ coeff.im() };
const Double_t xBarVal{ coeffBar.re() };
const Double_t yBarVal{ coeffBar.im() };
x_->value( xVal );
y_->value( yVal );
xbar_->value( xBarVal );
ybar_->value( yBarVal );
if ( init ) {
x_->genValue( xVal );
y_->genValue( yVal );
xbar_->genValue( xBarVal );
ybar_->genValue( yBarVal );
x_->initValue( xVal );
y_->initValue( yVal );
xbar_->initValue( xBarVal );
ybar_->initValue( yBarVal );
}
}
LauParameter LauRealImagCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const Double_t csq { x_->value()*x_->value() + y_->value()*y_->value() };
const Double_t cbarsq { xbar_->value()*xbar_->value() + ybar_->value()*ybar_->value() };
const Double_t numer { cbarsq - csq };
const Double_t denom { cbarsq + csq };
const Double_t value { numer/denom };
// is it fixed?
const Bool_t fixed { x_->fixed() && y_->fixed() && xbar_->fixed() && ybar_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauRealImagCPCoeffSet* LauRealImagCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart ) ) {
std::cerr << "ERROR in LauRealImagCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauRealImagCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
-
+//! \cond DOXYGEN_IGNORE
LauRealImagCPCoeffSet nlohmann::adl_serializer<LauRealImagCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::RealImagCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t xbar { j.at("xbar").get<Double_t>() };
const Double_t ybar { j.at("ybar").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t xbarFixed { j.at("xbarFixed").get<Bool_t>() };
const Bool_t ybarFixed { j.at("ybarFixed").get<Bool_t>() };
return { name, x, y, xbar, ybar, xFixed, yFixed, xbarFixed, ybarFixed };
}
void nlohmann::adl_serializer<LauRealImagCPCoeffSet>::to_json(json& j, const LauRealImagCPCoeffSet& t)
{
j["type"] = LauCoeffType::RealImagCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["xbar"] = pars[2]->value();
j["ybar"] = pars[3]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["xbarFixed"] = pars[2]->fixed();
j["ybarFixed"] = pars[3]->fixed();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauRealImagCoeffSet.cc b/src/LauRealImagCoeffSet.cc
index 3784191..f41e462 100644
--- a/src/LauRealImagCoeffSet.cc
+++ b/src/LauRealImagCoeffSet.cc
@@ -1,209 +1,211 @@
/*
Copyright 2006 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 LauRealImagCoeffSet.cc
\brief File containing implementation of LauRealImagCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauRealImagCoeffSet.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauRealImagCoeffSet)
LauRealImagCoeffSet::LauRealImagCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Bool_t xFixed, const Bool_t yFixed) :
LauAbsCoeffSet{ compName },
x_{ std::make_unique<LauParameter>("X",x,minRealImagPart_,maxRealImagPart_,xFixed) },
y_{ std::make_unique<LauParameter>("Y",y,minRealImagPart_,maxRealImagPart_,yFixed) },
coeff_{ x,y }
{
}
LauRealImagCoeffSet::LauRealImagCoeffSet(const LauRealImagCoeffSet& rhs, CloneOption cloneOption, Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
coeff_{ rhs.coeff_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
void LauRealImagCoeffSet::printParValues() const
{
std::cout<<"INFO in LauRealImagCoeffSet::printParValues : Component \""<<this->name()<<"\" has real part = "<<x_->value()<<" and imaginary part = "<<y_->value()<<"."<<std::endl;
}
void LauRealImagCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauRealImagCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauRealImagCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value between -10.0 and 10.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*20.0 - 10.0;
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value between -10.0 and 10.0
Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*20.0 - 10.0;
y_->initValue(value); y_->value(value);
}
}
void LauRealImagCoeffSet::finaliseValues()
{
x_->updatePull();
y_->updatePull();
}
const LauComplex& LauRealImagCoeffSet::particleCoeff()
{
coeff_.setRealImagPart(x_->unblindValue(), y_->unblindValue());
return coeff_;
}
const LauComplex& LauRealImagCoeffSet::antiparticleCoeff()
{
return this->particleCoeff();
}
void LauRealImagCoeffSet::setCoeffValues( const LauComplex& coeff, const LauComplex& coeffBar, Bool_t init )
{
LauComplex average{ coeff };
average += coeffBar;
average.rescale( 0.5 );
const Double_t xVal{ average.re() };
const Double_t yVal{ average.im() };
x_->value( xVal );
y_->value( yVal );
if ( init ) {
x_->genValue( xVal );
y_->genValue( yVal );
x_->initValue( xVal );
y_->initValue( yVal );
}
}
LauParameter LauRealImagCoeffSet::acp()
{
const TString parName{ this->baseName() + "_ACP" };
return LauParameter{parName,0.0};
}
LauRealImagCoeffSet* LauRealImagCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart ) ) {
std::cerr << "ERROR in LauRealImagCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauRealImagCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauRealImagCoeffSet nlohmann::adl_serializer<LauRealImagCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::RealImag ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
return { name, x, y, xFixed, yFixed };
}
void nlohmann::adl_serializer<LauRealImagCoeffSet>::to_json(json& j, const LauRealImagCoeffSet& t)
{
j["type"] = LauCoeffType::RealImag;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
}
+//! \endcond DOXYGEN_IGNORE
diff --git a/src/LauRealImagGammaCPCoeffSet.cc b/src/LauRealImagGammaCPCoeffSet.cc
index 8c6eb5e..c0050b7 100644
--- a/src/LauRealImagGammaCPCoeffSet.cc
+++ b/src/LauRealImagGammaCPCoeffSet.cc
@@ -1,357 +1,359 @@
/*
Copyright 2014 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 LauRealImagGammaCPCoeffSet.cc
\brief File containing implementation of LauRealImagGammaCPCoeffSet class.
*/
#include <iostream>
#include <fstream>
#include <vector>
#include "TMath.h"
#include "TRandom.h"
#include "LauRealImagGammaCPCoeffSet.hh"
#include "LauComplex.hh"
#include "LauConstants.hh"
#include "LauParameter.hh"
#include "LauPrint.hh"
ClassImp(LauRealImagGammaCPCoeffSet)
LauRealImagGammaCPCoeffSet::LauRealImagGammaCPCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Double_t xCP, const Double_t yCP, const Double_t xbarCP, const Double_t ybarCP,
const Bool_t xFixed, const Bool_t yFixed, const Bool_t xCPFixed, const Bool_t yCPFixed, const Bool_t xbarCPFixed, const Bool_t ybarCPFixed) :
LauAbsCoeffSet{ compName },
x_{ std::make_unique<LauParameter>("X", x, minRealImagPart_, maxRealImagPart_, xFixed) },
y_{ std::make_unique<LauParameter>("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) },
xCP_{ std::make_unique<LauParameter>("XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed) },
yCP_{ std::make_unique<LauParameter>("YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed) },
xbarCP_{ std::make_unique<LauParameter>("XbarCP", xbarCP, minRealImagPart_, maxRealImagPart_, xbarCPFixed) },
ybarCP_{ std::make_unique<LauParameter>("YbarCP", ybarCP, minRealImagPart_, maxRealImagPart_, ybarCPFixed) },
nonCPPart_{ x, y },
cpPart_{ 1+xCP, yCP },
cpAntiPart_{ 1+xbarCP, ybarCP },
particleCoeff_{ nonCPPart_ * cpPart_ },
antiparticleCoeff_{ nonCPPart_ * cpAntiPart_ },
acp_{ "ACP", (antiparticleCoeff_.abs2()-particleCoeff_.abs2())/(antiparticleCoeff_.abs2()+particleCoeff_.abs2()), -1.0, 1.0, xCPFixed&&yCPFixed&&xbarCPFixed&&ybarCPFixed }
{
}
LauRealImagGammaCPCoeffSet::LauRealImagGammaCPCoeffSet(const LauRealImagGammaCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor) :
LauAbsCoeffSet{ rhs.name() },
nonCPPart_{ rhs.nonCPPart_ },
cpPart_{ rhs.cpPart_ },
cpAntiPart_{ rhs.cpAntiPart_ },
particleCoeff_{ rhs.particleCoeff_ },
antiparticleCoeff_{ rhs.antiparticleCoeff_ },
acp_{ rhs.acp_ }
{
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) {
x_.reset( rhs.x_->createClone(constFactor) );
} else {
x_ = std::make_unique<LauParameter>("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
if ( rhs.x_->blind() ) {
const LauBlind* blinder { rhs.x_->blinder() };
x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) {
y_.reset( rhs.y_->createClone(constFactor) );
} else {
y_ = std::make_unique<LauParameter>("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
if ( rhs.y_->blind() ) {
const LauBlind* blinder { rhs.y_->blinder() };
y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieCPPars ) {
xCP_.reset( rhs.xCP_->createClone(constFactor) );
yCP_.reset( rhs.yCP_->createClone(constFactor) );
xbarCP_.reset( rhs.xbarCP_->createClone(constFactor) );
ybarCP_.reset( rhs.ybarCP_->createClone(constFactor) );
} else {
xCP_ = std::make_unique<LauParameter>("XCP", rhs.xCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.xCP_->fixed());
if ( rhs.xCP_->blind() ) {
const LauBlind* blinder { rhs.xCP_->blinder() };
xCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
yCP_ = std::make_unique<LauParameter>("YCP", rhs.yCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.yCP_->fixed());
if ( rhs.yCP_->blind() ) {
const LauBlind* blinder { rhs.yCP_->blinder() };
yCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
xbarCP_ = std::make_unique<LauParameter>("XbarCP", rhs.xbarCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.xbarCP_->fixed());
if ( rhs.xbarCP_->blind() ) {
const LauBlind* blinder { rhs.xbarCP_->blinder() };
xbarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
ybarCP_ = std::make_unique<LauParameter>("YbarCP", rhs.ybarCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.ybarCP_->fixed());
if ( rhs.ybarCP_->blind() ) {
const LauBlind* blinder { rhs.ybarCP_->blinder() };
ybarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
}
}
}
std::vector<LauParameter*> LauRealImagGammaCPCoeffSet::getParameters()
{
return {
x_.get(),
y_.get(),
xCP_.get(),
yCP_.get(),
xbarCP_.get(),
ybarCP_.get()
};
}
std::vector<const LauParameter*> LauRealImagGammaCPCoeffSet::getParameters() const
{
return {
x_.get(),
y_.get(),
xCP_.get(),
yCP_.get(),
xbarCP_.get(),
ybarCP_.get()
};
}
void LauRealImagGammaCPCoeffSet::printParValues() const
{
std::cout << "INFO in LauRealImagGammaCPCoeffSet::printParValues : Component \"" << this->name() << "\" has ";
std::cout << "x = " << x_->value() << ",\t";
std::cout << "y = " << y_->value() << ",\t";
std::cout << "xCP = " << xCP_->value() << ",\t";
std::cout << "yCP = " << yCP_->value() << ",\t";
std::cout << "xbarCP = " << xbarCP_->value() << ",\t";
std::cout << "ybarCP = " << ybarCP_->value() << "." << std::endl;
}
void LauRealImagGammaCPCoeffSet::printTableHeading(std::ostream& stream) const
{
stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|}"<<std::endl;
stream<<"\\hline"<<std::endl;
stream<<"Component & Real Part & Imaginary Part & Particle CP Real Part & Particle CP Imaginary Part & Antiparticle CP Real Part & Antiparticle CP Imaginary Part \\\\"<<std::endl;
stream<<"\\hline"<<std::endl;
}
void LauRealImagGammaCPCoeffSet::printTableRow(std::ostream& stream) const
{
LauPrint print;
TString resName { this->name() };
resName = resName.ReplaceAll("_", "\\_");
stream<<resName<<" & $";
print.printFormat(stream, x_->value());
stream<<" \\pm ";
print.printFormat(stream, x_->error());
stream<<"$ & $";
print.printFormat(stream, y_->value());
stream<<" \\pm ";
print.printFormat(stream, y_->error());
stream<<"$ & $";
print.printFormat(stream, xCP_->value());
stream<<" \\pm ";
print.printFormat(stream, xCP_->error());
stream<<"$ & $";
print.printFormat(stream, yCP_->value());
stream<<" \\pm ";
print.printFormat(stream, yCP_->error());
stream<<"$ & $";
print.printFormat(stream, xbarCP_->value());
stream<<" \\pm ";
print.printFormat(stream, xbarCP_->error());
stream<<"$ & $";
print.printFormat(stream, ybarCP_->value());
stream<<" \\pm ";
print.printFormat(stream, ybarCP_->error());
stream<<"$ \\\\"<<std::endl;
}
void LauRealImagGammaCPCoeffSet::randomiseInitValues()
{
if (x_->fixed() == kFALSE) {
// Choose a value for "X" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
x_->initValue(value); x_->value(value);
}
if (y_->fixed() == kFALSE) {
// Choose a value for "Y" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
y_->initValue(value); y_->value(value);
}
if (xCP_->fixed() == kFALSE) {
// Choose a value for "XCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
xCP_->initValue(value); xCP_->value(value);
}
if (yCP_->fixed() == kFALSE) {
// Choose a value for "YCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
yCP_->initValue(value); yCP_->value(value);
}
if (xbarCP_->fixed() == kFALSE) {
// Choose a value for "XbarCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
xbarCP_->initValue(value); xbarCP_->value(value);
}
if (ybarCP_->fixed() == kFALSE) {
// Choose a value for "YbarCP" between -3.0 and 3.0
const Double_t value { LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0 };
ybarCP_->initValue(value); ybarCP_->value(value);
}
}
void LauRealImagGammaCPCoeffSet::finaliseValues()
{
// update the pulls
x_->updatePull();
y_->updatePull();
xCP_->updatePull();
yCP_->updatePull();
xbarCP_->updatePull();
ybarCP_->updatePull();
}
const LauComplex& LauRealImagGammaCPCoeffSet::particleCoeff()
{
nonCPPart_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
cpPart_.setRealImagPart( 1.0+xCP_->unblindValue(), yCP_->unblindValue() );
particleCoeff_ = nonCPPart_ * cpPart_;
return particleCoeff_;
}
const LauComplex& LauRealImagGammaCPCoeffSet::antiparticleCoeff()
{
nonCPPart_.setRealImagPart( x_->unblindValue(), y_->unblindValue() );
cpAntiPart_.setRealImagPart( 1.0+xbarCP_->unblindValue(), ybarCP_->unblindValue() );
antiparticleCoeff_ = nonCPPart_ * cpAntiPart_;
return antiparticleCoeff_;
}
void LauRealImagGammaCPCoeffSet::setCoeffValues( const LauComplex&, const LauComplex&, Bool_t )
{
std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters" << std::endl;
}
LauParameter LauRealImagGammaCPCoeffSet::acp()
{
// set the name
const TString parName{ this->baseName() + "_ACP" };
acp_.name(parName);
// work out the ACP value
const LauComplex nonCPPart{ x_->value(), y_->value() };
const LauComplex cpPart{ 1.0+xCP_->value(), yCP_->value() };
const LauComplex cpAntiPart{ 1.0+xbarCP_->value(), ybarCP_->value() };
const LauComplex partCoeff { nonCPPart * cpPart };
const LauComplex antiCoeff { nonCPPart * cpAntiPart };
const Double_t numer { antiCoeff.abs2() - partCoeff.abs2() };
const Double_t denom { antiCoeff.abs2() + partCoeff.abs2() };
const Double_t value { numer/denom };
// is it fixed?
const Bool_t fixed { xCP_->fixed() && yCP_->fixed() && xbarCP_->fixed() && ybarCP_->fixed() };
acp_.fixed(fixed);
// we can't work out the error without the covariance matrix
const Double_t error{0.0};
// set the value and error
acp_.valueAndErrors(value,error);
return acp_;
}
LauRealImagGammaCPCoeffSet* LauRealImagGammaCPCoeffSet::createClone_impl(const TString& newName, const CloneOption cloneOption, const Double_t constFactor)
{
if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart || cloneOption == CloneOption::TieCPPars ) ) {
std::cerr << "ERROR in LauRealImagGammaCPCoeffSet::createClone : Invalid clone option" << std::endl;
return nullptr;
}
auto clone = new LauRealImagGammaCPCoeffSet{ *this, cloneOption, constFactor };
clone->name( newName );
return clone;
}
+//! \cond DOXYGEN_IGNORE
LauRealImagGammaCPCoeffSet nlohmann::adl_serializer<LauRealImagGammaCPCoeffSet>::from_json(const json& j)
{
const LauCoeffType type { j.at("type").get<LauCoeffType>() };
if ( type != LauCoeffType::RealImagGammaCP ) {
throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagGammaCPCoeffSet");
}
const TString name { j.at("name").get<std::string>().c_str() };
// TODO - handle cloned coeffs
const Double_t x { j.at("x").get<Double_t>() };
const Double_t y { j.at("y").get<Double_t>() };
const Double_t xCP { j.at("xCP").get<Double_t>() };
const Double_t yCP { j.at("yCP").get<Double_t>() };
const Double_t xbarCP { j.at("xbarCP").get<Double_t>() };
const Double_t ybarCP { j.at("ybarCP").get<Double_t>() };
const Bool_t xFixed { j.at("xFixed").get<Bool_t>() };
const Bool_t yFixed { j.at("yFixed").get<Bool_t>() };
const Bool_t xCPFixed { j.at("xCPFixed").get<Bool_t>() };
const Bool_t yCPFixed { j.at("yCPFixed").get<Bool_t>() };
const Bool_t xbarCPFixed { j.at("xbarCPFixed").get<Bool_t>() };
const Bool_t ybarCPFixed { j.at("ybarCPFixed").get<Bool_t>() };
return { name,
x, y, xCP, yCP, xbarCP, ybarCP,
xFixed, yFixed, xCPFixed, yCPFixed, xbarCPFixed, ybarCPFixed
};
}
void nlohmann::adl_serializer<LauRealImagGammaCPCoeffSet>::to_json(json& j, const LauRealImagGammaCPCoeffSet& t)
{
j["type"] = LauCoeffType::RealImagGammaCP;
j["name"] = t.name();
// TODO - handle cloned coeffs
auto pars = t.getParameters();
j["x"] = pars[0]->value();
j["y"] = pars[1]->value();
j["xCP"] = pars[2]->value();
j["yCP"] = pars[3]->value();
j["xbarCP"] = pars[4]->value();
j["ybarCP"] = pars[5]->value();
j["xFixed"] = pars[0]->fixed();
j["yFixed"] = pars[1]->fixed();
j["xCPFixed"] = pars[2]->fixed();
j["yCPFixed"] = pars[3]->fixed();
j["xbarCPFixed"] = pars[4]->fixed();
j["ybarCPFixed"] = pars[5]->fixed();
}
+//! \endcond DOXYGEN_IGNORE

File Metadata

Mime Type
text/x-diff
Expires
Wed, May 14, 11:54 AM (1 h, 18 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
5106146
Default Alt Text
(165 KB)

Event Timeline