diff --git a/src/LauBelleCPCoeffSet.cc b/src/LauBelleCPCoeffSet.cc index b545052..b4ef146 100644 --- a/src/LauBelleCPCoeffSet.cc +++ b/src/LauBelleCPCoeffSet.cc @@ -1,444 +1,439 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauBelleCPCoeffSet.hh" #include "LauJsonTools.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("A", a, minMagnitude_, maxMagnitude_, aFixed) }, b_{ std::make_unique("B", b, minMagnitude_, maxMagnitude_, bFixed) }, delta_{ std::make_unique("Delta", delta, minPhase_, maxPhase_, deltaFixed) }, phi_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) { a_.reset( rhs.a_->createClone(constFactor) ); } else { const auto aVal { getOptionalValue( coeffInfo, "A", JsonType::Number ).value_or( rhs.a_->value() ) }; const auto aFixed { getOptionalValue( coeffInfo, "AFixed", JsonType::Boolean ).value_or( rhs.a_->fixed() ) }; a_ = std::make_unique("A", aVal, minMagnitude_, maxMagnitude_, aFixed); if ( getOptionalValue( coeffInfo, "ABlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "ABlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "ABlindingWidth" ) }; a_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto bVal { getOptionalValue( coeffInfo, "B", JsonType::Number ).value_or( rhs.b_->value() ) }; const auto bFixed { getOptionalValue( coeffInfo, "BFixed", JsonType::Boolean ).value_or( rhs.b_->fixed() ) }; const auto bSecondStage { getOptionalValue( coeffInfo, "BSecondStage", JsonType::Boolean ).value_or( rhs.b_->secondStage() ) }; b_ = std::make_unique("B", bVal, minMagnitude_, maxMagnitude_, bFixed); if (bSecondStage && !bFixed) { b_->secondStage(kTRUE); b_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "BBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "BBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "BBlindingWidth" ) }; b_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto deltaVal { getOptionalValue( coeffInfo, "Delta", JsonType::Number ).value_or( rhs.delta_->value() ) }; const auto deltaFixed { getOptionalValue( coeffInfo, "DeltaFixed", JsonType::Boolean ).value_or( rhs.delta_->fixed() ) }; delta_ = std::make_unique("Delta", deltaVal, minPhase_, maxPhase_, deltaFixed); if ( getOptionalValue( coeffInfo, "DeltaBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaBlindingWidth" ) }; delta_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto phiVal { getOptionalValue( coeffInfo, "Phi", JsonType::Number ).value_or( rhs.phi_->value() ) }; const auto phiFixed { getOptionalValue( coeffInfo, "PhiFixed", JsonType::Boolean ).value_or( rhs.phi_->fixed() ) }; const auto phiSecondStage { getOptionalValue( coeffInfo, "PhiSecondStage", JsonType::Boolean ).value_or( rhs.phi_->secondStage() ) }; phi_ = std::make_unique("Phi", phiVal, minPhase_, maxPhase_, phiFixed); if (phiSecondStage && !phiFixed) { phi_->secondStage(kTRUE); phi_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "PhiBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "PhiBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "PhiBlindingWidth" ) }; phi_->blindParameter( blindingString, blindingWidth ); } else 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 \""<name()<<"\" has "; std::cout<<"A-magnitude = "<value()<<",\t"; std::cout<<"Delta = "<value()<<",\t"; std::cout<<"B-magnitude = "<value()<<",\t"; std::cout<<"Phi = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauBelleCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauBelleCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauBelleCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::BelleCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauBelleCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauBelleCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("A", JsonType::Number), std::make_pair("Delta", JsonType::Number), std::make_pair("B", JsonType::Number), std::make_pair("Phi", JsonType::Number), std::make_pair("AFixed", JsonType::Boolean), std::make_pair("DeltaFixed", JsonType::Boolean), std::make_pair("BFixed", JsonType::Boolean), std::make_pair("PhiFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauBelleCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto a { getValue( j, "A" ) }; const auto delta { getValue( j, "Delta" ) }; const auto b { getValue( j, "B" ) }; const auto phi { getValue( j, "Phi" ) }; const auto aFixed { getValue( j, "AFixed" ) }; const auto deltaFixed { getValue( j, "DeltaFixed" ) }; const auto bFixed { getValue( j, "BFixed" ) }; const auto phiFixed { getValue( j, "PhiFixed" ) }; const auto bSecondStage { getOptionalValue( j, "BSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto phiSecondStage { getOptionalValue( j, "PhiSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; LauBelleCPCoeffSet coeff{ name, a, delta, b, phi, aFixed, deltaFixed, bFixed, phiFixed, bSecondStage, phiSecondStage }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauCartesianCPCoeffSet.cc b/src/LauCartesianCPCoeffSet.cc index 49285a7..5591328 100644 --- a/src/LauCartesianCPCoeffSet.cc +++ b/src/LauCartesianCPCoeffSet.cc @@ -1,366 +1,361 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauCartesianCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("X", x, minRealImagPart_, maxRealImagPart_, xFixed) }, y_{ std::make_unique("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) }, deltaX_{ std::make_unique("DeltaX", deltaX, minDelta_, maxDelta_, deltaXFixed) }, deltaY_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto deltaXVal { getOptionalValue( coeffInfo, "DeltaX", JsonType::Number ).value_or( rhs.deltaX_->value() ) }; const auto deltaXFixed { getOptionalValue( coeffInfo, "DeltaXFixed", JsonType::Boolean ).value_or( rhs.deltaX_->fixed() ) }; const auto deltaXSecondStage { getOptionalValue( coeffInfo, "DeltaXSecondStage", JsonType::Boolean ).value_or( rhs.deltaX_->secondStage() ) }; const auto deltaYVal { getOptionalValue( coeffInfo, "DeltaY", JsonType::Number ).value_or( rhs.deltaY_->value() ) }; const auto deltaYFixed { getOptionalValue( coeffInfo, "DeltaYFixed", JsonType::Boolean ).value_or( rhs.deltaY_->fixed() ) }; const auto deltaYSecondStage { getOptionalValue( coeffInfo, "DeltaYSecondStage", JsonType::Boolean ).value_or( rhs.deltaY_->secondStage() ) }; deltaX_ = std::make_unique("DeltaX", deltaXVal, minDelta_, maxDelta_, deltaXFixed); deltaY_ = std::make_unique("DeltaY", deltaYVal, minDelta_, maxDelta_, deltaYFixed); if ( deltaXSecondStage && !deltaXFixed ) { deltaX_->secondStage(kTRUE); deltaX_->initValue(0.0); } if ( deltaYSecondStage && !deltaYFixed ) { deltaY_->secondStage(kTRUE); deltaY_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "DeltaXBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaXBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaXBlindingWidth" ) }; deltaX_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaX_->blind() ) { const LauBlind* blinder { rhs.deltaX_->blinder() }; deltaX_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } if ( getOptionalValue( coeffInfo, "DeltaYBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaYBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaYBlindingWidth" ) }; deltaY_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaY_->blind() ) { const LauBlind* blinder { rhs.deltaY_->blinder() }; deltaY_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } void LauCartesianCPCoeffSet::printParValues() const { std::cout<<"INFO in LauCartesianCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"X = "<value()<<",\t"; std::cout<<"Y = "<value()<<",\t"; std::cout<<"Delta x = "<value()<<",\t"; std::cout<<"Delta y = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauCartesianCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauCartesianCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauCartesianCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::CartesianCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauCartesianCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauCartesianCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("DeltaX", JsonType::Number), std::make_pair("DeltaY", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("DeltaXFixed", JsonType::Boolean), std::make_pair("DeltaYFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauCartesianCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto deltaX { getValue( j, "DeltaX" ) }; const auto deltaY { getValue( j, "DeltaY" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto deltaXFixed { getValue( j, "DeltaXFixed" ) }; const auto deltaYFixed { getValue( j, "DeltaYFixed" ) }; const auto deltaXSecondStage { getOptionalValue( j, "DeltaXSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto deltaYSecondStage { getOptionalValue( j, "DeltaYSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; LauCartesianCPCoeffSet coeff{ name, x, y, deltaX, deltaY, xFixed, yFixed, deltaXFixed, deltaYFixed, deltaXSecondStage, deltaYSecondStage }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauCartesianGammaCPCoeffSet.cc b/src/LauCartesianGammaCPCoeffSet.cc index 1da3124..802f653 100644 --- a/src/LauCartesianGammaCPCoeffSet.cc +++ b/src/LauCartesianGammaCPCoeffSet.cc @@ -1,422 +1,417 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauCartesianGammaCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("X", x, minRealImagPart_, maxRealImagPart_, xFixed) }, y_{ std::make_unique("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) }, xCP_{ std::make_unique("XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed) }, yCP_{ std::make_unique("YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed) }, deltaXCP_{ std::make_unique("DeltaXCP", deltaXCP, minDelta_, maxDelta_, deltaXCPFixed) }, deltaYCP_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, nonCPPart_{ rhs.nonCPPart_ }, cpPart_{ rhs.cpPart_ }, cpAntiPart_{ rhs.cpAntiPart_ }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto xCPVal { getOptionalValue( coeffInfo, "XCP", JsonType::Number ).value_or( rhs.xCP_->value() ) }; const auto xCPFixed { getOptionalValue( coeffInfo, "XCPFixed", JsonType::Boolean ).value_or( rhs.xCP_->fixed() ) }; xCP_ = std::make_unique("XCP", xCPVal, minRealImagPart_, maxRealImagPart_, xCPFixed); if ( getOptionalValue( coeffInfo, "XCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XCPBlindingWidth" ) }; xCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.xCP_->blind() ) { const LauBlind* blinder { rhs.xCP_->blinder() }; xCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto yCPVal { getOptionalValue( coeffInfo, "YCP", JsonType::Number ).value_or( rhs.yCP_->value() ) }; const auto yCPFixed { getOptionalValue( coeffInfo, "YCPFixed", JsonType::Boolean ).value_or( rhs.yCP_->fixed() ) }; yCP_ = std::make_unique("YCP", yCPVal, minRealImagPart_, maxRealImagPart_, yCPFixed); if ( getOptionalValue( coeffInfo, "YCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YCPBlindingWidth" ) }; yCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.yCP_->blind() ) { const LauBlind* blinder { rhs.yCP_->blinder() }; yCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto deltaXCPVal { getOptionalValue( coeffInfo, "DeltaXCP", JsonType::Number ).value_or( rhs.deltaXCP_->value() ) }; const auto deltaXCPFixed { getOptionalValue( coeffInfo, "DeltaXCPFixed", JsonType::Boolean ).value_or( rhs.deltaXCP_->fixed() ) }; const auto deltaXCPSecondStage { getOptionalValue( coeffInfo, "DeltaXCPSecondStage", JsonType::Boolean ).value_or( rhs.deltaXCP_->secondStage() ) }; const auto deltaYCPVal { getOptionalValue( coeffInfo, "DeltaYCP", JsonType::Number ).value_or( rhs.deltaYCP_->value() ) }; const auto deltaYCPFixed { getOptionalValue( coeffInfo, "DeltaYCPFixed", JsonType::Boolean ).value_or( rhs.deltaYCP_->fixed() ) }; const auto deltaYCPSecondStage { getOptionalValue( coeffInfo, "DeltaYCPSecondStage", JsonType::Boolean ).value_or( rhs.deltaYCP_->secondStage() ) }; deltaXCP_ = std::make_unique("DeltaXCP", deltaXCPVal, minDelta_, maxDelta_, deltaXCPFixed); deltaYCP_ = std::make_unique("DeltaYCP", deltaYCPVal, minDelta_, maxDelta_, deltaYCPFixed); if ( deltaXCPSecondStage && !deltaXCPFixed ) { deltaXCP_->secondStage(kTRUE); deltaXCP_->initValue(0.0); } if ( deltaYCPSecondStage && !deltaYCPFixed ) { deltaYCP_->secondStage(kTRUE); deltaYCP_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "DeltaXCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaXCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaXCPBlindingWidth" ) }; deltaXCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaXCP_->blind() ) { const LauBlind* blinder { rhs.deltaXCP_->blinder() }; deltaXCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } if ( getOptionalValue( coeffInfo, "DeltaYCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaYCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaYCPBlindingWidth" ) }; deltaYCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaYCP_->blind() ) { const LauBlind* blinder { rhs.deltaYCP_->blinder() }; deltaYCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } void LauCartesianGammaCPCoeffSet::printParValues() const { std::cout<<"INFO in LauCartesianGammaCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"X = "<value()<<",\t"; std::cout<<"Y = "<value()<<",\t"; std::cout<<"XCP = "<value()<<",\t"; std::cout<<"YCP = "<value()<<",\t"; std::cout<<"Delta xCP = "<value()<<",\t"; std::cout<<"Delta yCP = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauCartesianGammaCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauCartesianGammaCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauCartesianGammaCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::CartesianGammaCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauCartesianGammaCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauCartesianGammaCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("XCP", JsonType::Number), std::make_pair("YCP", JsonType::Number), std::make_pair("DeltaXCP", JsonType::Number), std::make_pair("DeltaYCP", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("XCPFixed", JsonType::Boolean), std::make_pair("YCPFixed", JsonType::Boolean), std::make_pair("DeltaXCPFixed", JsonType::Boolean), std::make_pair("DeltaYCPFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauCartesianGammaCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto xCP { getValue( j, "XCP" ) }; const auto yCP { getValue( j, "YCP" ) }; const auto deltaXCP { getValue( j, "DeltaXCP" ) }; const auto deltaYCP { getValue( j, "DeltaYCP" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto xCPFixed { getValue( j, "XCPFixed" ) }; const auto yCPFixed { getValue( j, "YCPFixed" ) }; const auto deltaXCPFixed { getValue( j, "DeltaXCPFixed" ) }; const auto deltaYCPFixed { getValue( j, "DeltaYCPFixed" ) }; const auto deltaXCPSecondStage { getOptionalValue( j, "DeltaXCPSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto deltaYCPSecondStage { getOptionalValue( j, "DeltaYCPSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; LauCartesianGammaCPCoeffSet coeff{ name, x, y, xCP, yCP, deltaXCP, deltaYCP, xFixed, yFixed, xCPFixed, yCPFixed, deltaXCPFixed, deltaYCPFixed, deltaXCPSecondStage, deltaYCPSecondStage }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauCleoCPCoeffSet.cc b/src/LauCleoCPCoeffSet.cc index e0c843e..62b46f1 100644 --- a/src/LauCleoCPCoeffSet.cc +++ b/src/LauCleoCPCoeffSet.cc @@ -1,437 +1,432 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauCleoCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("A", a, minMagnitude_, maxMagnitude_, aFixed) }, b_{ std::make_unique("B", b, minMagnitude_, maxMagnitude_, bFixed) }, delta_{ std::make_unique("Delta", delta, minPhase_, maxPhase_, deltaFixed) }, phi_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) { a_.reset( rhs.a_->createClone(constFactor) ); } else { const auto aVal { getOptionalValue( coeffInfo, "A", JsonType::Number ).value_or( rhs.a_->value() ) }; const auto aFixed { getOptionalValue( coeffInfo, "AFixed", JsonType::Boolean ).value_or( rhs.a_->fixed() ) }; a_ = std::make_unique("A", aVal, minMagnitude_, maxMagnitude_, aFixed); if ( getOptionalValue( coeffInfo, "ABlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "ABlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "ABlindingWidth" ) }; a_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto bVal { getOptionalValue( coeffInfo, "B", JsonType::Number ).value_or( rhs.b_->value() ) }; const auto bFixed { getOptionalValue( coeffInfo, "BFixed", JsonType::Boolean ).value_or( rhs.b_->fixed() ) }; const auto bSecondStage { getOptionalValue( coeffInfo, "BSecondStage", JsonType::Boolean ).value_or( rhs.b_->secondStage() ) }; b_ = std::make_unique("B", bVal, minMagnitude_, maxMagnitude_, bFixed); if ( bSecondStage && !bFixed ) { b_->secondStage(kTRUE); b_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "BBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "BBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "BBlindingWidth" ) }; b_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto deltaVal { getOptionalValue( coeffInfo, "Delta", JsonType::Number ).value_or( rhs.delta_->value() ) }; const auto deltaFixed { getOptionalValue( coeffInfo, "DeltaFixed", JsonType::Boolean ).value_or( rhs.delta_->fixed() ) }; delta_ = std::make_unique("Delta", deltaVal, minPhase_, maxPhase_, deltaFixed); if ( getOptionalValue( coeffInfo, "DeltaBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaBlindingWidth" ) }; delta_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto phiVal { getOptionalValue( coeffInfo, "Phi", JsonType::Number ).value_or( rhs.phi_->value() ) }; const auto phiFixed { getOptionalValue( coeffInfo, "PhiFixed", JsonType::Boolean ).value_or( rhs.phi_->fixed() ) }; const auto phiSecondStage { getOptionalValue( coeffInfo, "PhiSecondStage", JsonType::Boolean ).value_or( rhs.phi_->secondStage() ) }; phi_ = std::make_unique("Phi", phiVal, minPhase_, maxPhase_, phiFixed); if ( phiSecondStage && !phiFixed ) { phi_->secondStage(kTRUE); phi_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "PhiBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "PhiBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "PhiBlindingWidth" ) }; phi_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.phi_->blind() ) { const LauBlind* blinder { rhs.phi_->blinder() }; phi_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } void LauCleoCPCoeffSet::printParValues() const { std::cout<<"INFO in LauCleoCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"A-magnitude = "<value()<<",\t"; std::cout<<"Delta = "<value()<<",\t"; std::cout<<"B-magnitude = "<value()<<",\t"; std::cout<<"Phi = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauCleoCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauCleoCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauCleoCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::CleoCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauCleoCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauCleoCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("A", JsonType::Number), std::make_pair("Delta", JsonType::Number), std::make_pair("B", JsonType::Number), std::make_pair("Phi", JsonType::Number), std::make_pair("AFixed", JsonType::Boolean), std::make_pair("DeltaFixed", JsonType::Boolean), std::make_pair("BFixed", JsonType::Boolean), std::make_pair("PhiFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauCleoCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto a { getValue( j, "A" ) }; const auto delta { getValue( j, "Delta" ) }; const auto b { getValue( j, "B" ) }; const auto phi { getValue( j, "Phi" ) }; const auto aFixed { getValue( j, "AFixed" ) }; const auto deltaFixed { getValue( j, "DeltaFixed" ) }; const auto bFixed { getValue( j, "BFixed" ) }; const auto phiFixed { getValue( j, "PhiFixed" ) }; const auto bSecondStage { getOptionalValue( j, "BSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto phiSecondStage { getOptionalValue( j, "PhiSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; LauCleoCPCoeffSet coeff{ name, a, delta, b, phi, aFixed, deltaFixed, bFixed, phiFixed, bSecondStage, phiSecondStage }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauMagPhaseCPCoeffSet.cc b/src/LauMagPhaseCPCoeffSet.cc index bd5ca55..6a768f3 100644 --- a/src/LauMagPhaseCPCoeffSet.cc +++ b/src/LauMagPhaseCPCoeffSet.cc @@ -1,400 +1,395 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauMagPhaseCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("Mag", mag, minMagnitude_, maxMagnitude_, magFixed) }, phase_{ std::make_unique("Phase", phase, minPhase_, maxPhase_, phaseFixed) }, magBar_{ std::make_unique("MagBar", magBar, minMagnitude_, maxMagnitude_, magBarFixed) }, phaseBar_{ std::make_unique("PhaseBar", phaseBar, minPhase_, maxPhase_, phaseBarFixed) }, acp_{ "ACP", (magBar*magBar - mag*mag)/(magBar*magBar + mag*mag), -1.0, 1.0 } { } LauMagPhaseCPCoeffSet::LauMagPhaseCPCoeffSet(const LauMagPhaseCPCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) { mag_.reset( rhs.mag_->createClone(constFactor) ); magBar_.reset( rhs.magBar_->createClone(constFactor) ); } else { const auto magVal { getOptionalValue( coeffInfo, "Mag", JsonType::Number ).value_or( rhs.mag_->value() ) }; const auto magFixed { getOptionalValue( coeffInfo, "MagFixed", JsonType::Boolean ).value_or( rhs.mag_->fixed() ) }; mag_ = std::make_unique("Mag", magVal, minMagnitude_, maxMagnitude_, magFixed); if ( getOptionalValue( coeffInfo, "MagBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "MagBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "MagBlindingWidth" ) }; mag_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.mag_->blind() ) { const LauBlind* blinder { rhs.mag_->blinder() }; mag_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto magBarVal { getOptionalValue( coeffInfo, "MagBar", JsonType::Number ).value_or( rhs.magBar_->value() ) }; const auto magBarFixed { getOptionalValue( coeffInfo, "MagBarFixed", JsonType::Boolean ).value_or( rhs.magBar_->fixed() ) }; magBar_ = std::make_unique("MagBar", magBarVal, minMagnitude_, maxMagnitude_, magBarFixed); if ( getOptionalValue( coeffInfo, "MagBarBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "MagBarBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "MagBarBlindingWidth" ) }; magBar_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto phaseVal { getOptionalValue( coeffInfo, "Phase", JsonType::Number ).value_or( rhs.phase_->value() ) }; const auto phaseFixed { getOptionalValue( coeffInfo, "PhaseFixed", JsonType::Boolean ).value_or( rhs.phase_->fixed() ) }; phase_ = std::make_unique("Phase", phaseVal, minPhase_, maxPhase_, phaseFixed); if ( getOptionalValue( coeffInfo, "PhaseBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "PhaseBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "PhaseBlindingWidth" ) }; phase_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.phase_->blind() ) { const LauBlind* blinder { rhs.phase_->blinder() }; phase_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto phaseBarVal { getOptionalValue( coeffInfo, "PhaseBar", JsonType::Number ).value_or( rhs.phaseBar_->value() ) }; const auto phaseBarFixed { getOptionalValue( coeffInfo, "PhaseBarFixed", JsonType::Boolean ).value_or( rhs.phaseBar_->fixed() ) }; phaseBar_ = std::make_unique("PhaseBar", phaseBarVal, minPhase_, maxPhase_, phaseBarFixed); if ( getOptionalValue( coeffInfo, "PhaseBarBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "PhaseBarBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "PhaseBarBlindingWidth" ) }; phaseBar_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.phaseBar_->blind() ) { const LauBlind* blinder { rhs.phaseBar_->blinder() }; phaseBar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } void LauMagPhaseCPCoeffSet::printParValues() const { std::cout<<"INFO in LauMagPhaseCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"Mag = "<value()<<",\t"; std::cout<<"Phase = "<value()<<",\t"; std::cout<<"MagBar = "<value()<<",\t"; std::cout<<"PhaseBar = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude ) ) { std::cerr << "ERROR in LauMagPhaseCPCoeffSet::createClone : Invalid clone option" << std::endl; return nullptr; } - if ( this->clone() ) { - const LauMagPhaseCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauMagPhaseCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauMagPhaseCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::MagPhaseCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauMagPhaseCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauMagPhaseCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("Mag", JsonType::Number), std::make_pair("Phase", JsonType::Number), std::make_pair("MagBar", JsonType::Number), std::make_pair("PhaseBar", JsonType::Number), std::make_pair("MagFixed", JsonType::Boolean), std::make_pair("PhaseFixed", JsonType::Boolean), std::make_pair("MagBarFixed", JsonType::Boolean), std::make_pair("PhaseBarFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauMagPhaseCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto mag { getValue( j, "Mag" ) }; const auto phase { getValue( j, "Phase" ) }; const auto magBar { getValue( j, "MagBar" ) }; const auto phaseBar { getValue( j, "PhaseBar" ) }; const auto magFixed { getValue( j, "MagFixed" ) }; const auto phaseFixed { getValue( j, "PhaseFixed" ) }; const auto magBarFixed { getValue( j, "MagBarFixed" ) }; const auto phaseBarFixed { getValue( j, "PhaseBarFixed" ) }; LauMagPhaseCPCoeffSet coeff{ name, mag, phase, magBar, phaseBar, magFixed, phaseFixed, magBarFixed, phaseBarFixed }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauMagPhaseCoeffSet.cc b/src/LauMagPhaseCoeffSet.cc index b0ae7b8..5898f63 100644 --- a/src/LauMagPhaseCoeffSet.cc +++ b/src/LauMagPhaseCoeffSet.cc @@ -1,283 +1,278 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("Mag", magnitude, minMagnitude_, maxMagnitude_, magFixed) }, phase_{ std::make_unique("Phase", phase, minPhase_, maxPhase_, phaseFixed) } { } LauMagPhaseCoeffSet::LauMagPhaseCoeffSet(const LauMagPhaseCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, coeff_{ rhs.coeff_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieMagnitude ) { magnitude_.reset( rhs.magnitude_->createClone(constFactor) ); } else { const auto magVal { getOptionalValue( coeffInfo, "Mag", JsonType::Number ).value_or( rhs.magnitude_->value() ) }; const auto magFixed { getOptionalValue( coeffInfo, "MagFixed", JsonType::Boolean ).value_or( rhs.magnitude_->fixed() ) }; magnitude_ = std::make_unique("Mag", magVal, minMagnitude_, maxMagnitude_, magFixed); if ( getOptionalValue( coeffInfo, "MagBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "MagBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "MagBlindingWidth" ) }; magnitude_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto phaseVal { getOptionalValue( coeffInfo, "Phase", JsonType::Number ).value_or( rhs.phase_->value() ) }; const auto phaseFixed { getOptionalValue( coeffInfo, "PhaseFixed", JsonType::Boolean ).value_or( rhs.phase_->fixed() ) }; phase_ = std::make_unique("Phase", phaseVal, minPhase_, maxPhase_, phaseFixed); if ( getOptionalValue( coeffInfo, "PhaseBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "PhaseBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "PhaseBlindingWidth" ) }; phase_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.phase_->blind() ) { const LauBlind* blinder { rhs.phase_->blinder() }; phase_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } void LauMagPhaseCoeffSet::printParValues() const { std::cout<<"INFO in LauMagPhaseCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"Mag = "<value()<<",\t"; std::cout<<"Phase = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<value()); stream<<" \\pm "; print.printFormat(stream, magnitude_->error()); stream<<"$ & $"; print.printFormat(stream, phase_->value()); stream<<" \\pm "; print.printFormat(stream, phase_->error()); stream<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TiePhase || cloneOption == CloneOption::TieMagnitude ) ) { std::cerr << "ERROR in LauMagPhaseCoeffSet::createClone : Invalid clone option" << std::endl; return nullptr; } - if ( this->clone() ) { - const LauMagPhaseCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauMagPhaseCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauMagPhaseCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::MagPhase ) { throw LauWrongCoeffType{"Wrong coefficient type given to construct LauMagPhaseCoeffSet"}; } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauMagPhaseCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("Mag", JsonType::Number), std::make_pair("Phase", JsonType::Number), std::make_pair("MagFixed", JsonType::Boolean), std::make_pair("PhaseFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauMagPhaseCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto mag { getValue( j, "Mag" ) }; const auto phase { getValue( j, "Phase" ) }; const auto magFixed { getValue( j, "MagFixed" ) }; const auto phaseFixed { getValue( j, "PhaseFixed" ) }; LauMagPhaseCoeffSet coeff{ name, mag, phase, magFixed, phaseFixed }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauNSCCartesianCPCoeffSet.cc b/src/LauNSCCartesianCPCoeffSet.cc index db94f5b..766f788 100644 --- a/src/LauNSCCartesianCPCoeffSet.cc +++ b/src/LauNSCCartesianCPCoeffSet.cc @@ -1,516 +1,511 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauNSCCartesianCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("X", x, minRealImagPart_, maxRealImagPart_, xFixed) }, y_{ std::make_unique("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) }, deltaX_{ std::make_unique("DeltaX", deltaX, minDelta_, maxDelta_, deltaXFixed) }, deltaY_{ std::make_unique("DeltaY", deltaY, minDelta_, maxDelta_, deltaYFixed) }, xPrime_{ std::make_unique("XPrime", xPrime, minRealImagPart_, maxRealImagPart_, xPrimeFixed) }, yPrime_{ std::make_unique("YPrime", yPrime, minRealImagPart_, maxRealImagPart_, yPrimeFixed) }, deltaXPrime_{ std::make_unique("DeltaXPrime", deltaXPrime, minDelta_, maxDelta_, deltaXPrimeFixed) }, deltaYPrime_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, finalStateIsF_{ rhs.finalStateIsF_ }, coeffAf_{ rhs.coeffAf_ }, coeffAfbar_{ rhs.coeffAfbar_ }, coeffAbarf_{ rhs.coeffAbarf_ }, coeffAbarfbar_{ rhs.coeffAbarfbar_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); xPrime_.reset( rhs.xPrime_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.x_->blind() ) { const LauBlind* blinder { rhs.x_->blinder() }; x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto xPrimeVal { getOptionalValue( coeffInfo, "XPrime", JsonType::Number ).value_or( rhs.xPrime_->value() ) }; const auto xPrimeFixed { getOptionalValue( coeffInfo, "XPrimeFixed", JsonType::Boolean ).value_or( rhs.xPrime_->fixed() ) }; xPrime_ = std::make_unique("XPrime", xPrimeVal, minRealImagPart_, maxRealImagPart_, xPrimeFixed); if ( getOptionalValue( coeffInfo, "XPrimeBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XPrimeBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XPrimeBlindingWidth" ) }; xPrime_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.xPrime_->blind() ) { const LauBlind* blinder { rhs.xPrime_->blinder() }; xPrime_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieImagPart ) { y_.reset( rhs.y_->createClone(constFactor) ); yPrime_.reset( rhs.yPrime_->createClone(constFactor) ); } else { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.y_->blind() ) { const LauBlind* blinder { rhs.y_->blinder() }; y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto yPrimeVal { getOptionalValue( coeffInfo, "YPrime", JsonType::Number ).value_or( rhs.yPrime_->value() ) }; const auto yPrimeFixed { getOptionalValue( coeffInfo, "YPrimeFixed", JsonType::Boolean ).value_or( rhs.yPrime_->fixed() ) }; yPrime_ = std::make_unique("YPrime", yPrimeVal, minRealImagPart_, maxRealImagPart_, yPrimeFixed); if ( getOptionalValue( coeffInfo, "YPrimeBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YPrimeBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YPrimeBlindingWidth" ) }; yPrime_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.yPrime_->blind() ) { const LauBlind* blinder { rhs.yPrime_->blinder() }; yPrime_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } 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 { const auto deltaXVal { getOptionalValue( coeffInfo, "DeltaX", JsonType::Number ).value_or( rhs.deltaX_->value() ) }; const auto deltaXFixed { getOptionalValue( coeffInfo, "DeltaXFixed", JsonType::Boolean ).value_or( rhs.deltaX_->fixed() ) }; const auto deltaXSecondStage { getOptionalValue( coeffInfo, "DeltaXSecondStage", JsonType::Boolean ).value_or( rhs.deltaX_->secondStage() ) }; const auto deltaYVal { getOptionalValue( coeffInfo, "DeltaY", JsonType::Number ).value_or( rhs.deltaY_->value() ) }; const auto deltaYFixed { getOptionalValue( coeffInfo, "DeltaYFixed", JsonType::Boolean ).value_or( rhs.deltaY_->fixed() ) }; const auto deltaYSecondStage { getOptionalValue( coeffInfo, "DeltaYSecondStage", JsonType::Boolean ).value_or( rhs.deltaY_->secondStage() ) }; deltaX_ = std::make_unique("DeltaX", deltaXVal, minDelta_, maxDelta_, deltaXFixed); deltaY_ = std::make_unique("DeltaY", deltaYVal, minDelta_, maxDelta_, deltaYFixed); if ( getOptionalValue( coeffInfo, "DeltaXBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaXBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaXBlindingWidth" ) }; deltaX_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaX_->blind() ) { const LauBlind* blinder { rhs.deltaX_->blinder() }; deltaX_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } if ( getOptionalValue( coeffInfo, "DeltaYBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaYBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaYBlindingWidth" ) }; deltaY_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaY_->blind() ) { const LauBlind* blinder { rhs.deltaY_->blinder() }; deltaY_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto deltaXPrimeVal { getOptionalValue( coeffInfo, "DeltaXPrime", JsonType::Number ).value_or( rhs.deltaXPrime_->value() ) }; const auto deltaXPrimeFixed { getOptionalValue( coeffInfo, "DeltaXPrimeFixed", JsonType::Boolean ).value_or( rhs.deltaXPrime_->fixed() ) }; const auto deltaXPrimeSecondStage { getOptionalValue( coeffInfo, "DeltaXPrimeSecondStage", JsonType::Boolean ).value_or( rhs.deltaXPrime_->secondStage() ) }; const auto deltaYPrimeVal { getOptionalValue( coeffInfo, "DeltaYPrime", JsonType::Number ).value_or( rhs.deltaYPrime_->value() ) }; const auto deltaYPrimeFixed { getOptionalValue( coeffInfo, "DeltaYPrimeFixed", JsonType::Boolean ).value_or( rhs.deltaYPrime_->fixed() ) }; const auto deltaYPrimeSecondStage { getOptionalValue( coeffInfo, "DeltaYPrimeSecondStage", JsonType::Boolean ).value_or( rhs.deltaYPrime_->secondStage() ) }; deltaXPrime_ = std::make_unique("DeltaXPrime", deltaXPrimeVal, minDelta_, maxDelta_, deltaXPrimeFixed); deltaYPrime_ = std::make_unique("DeltaYPrime", deltaYPrimeVal, minDelta_, maxDelta_, deltaYPrimeFixed); if ( getOptionalValue( coeffInfo, "DeltaXPrimeBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaXPrimeBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaXPrimeBlindingWidth" ) }; deltaXPrime_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaXPrime_->blind() ) { const LauBlind* blinder { rhs.deltaXPrime_->blinder() }; deltaXPrime_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } if ( getOptionalValue( coeffInfo, "DeltaYPrimeBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "DeltaYPrimeBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "DeltaYPrimeBlindingWidth" ) }; deltaYPrime_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaYPrime_->blind() ) { const LauBlind* blinder { rhs.deltaYPrime_->blinder() }; deltaYPrime_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } 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); } } } void LauNSCCartesianCPCoeffSet::printParValues() const { std::cout<<"INFO in LauNSCCartesianCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"X = "<value()<<",\t"; std::cout<<"Y = "<value()<<",\t"; std::cout<<"Delta x = "<value()<<",\t"; std::cout<<"Delta y = "<value()<<",\t"; std::cout<<"XPrime = "<value()<<",\t"; std::cout<<"YPrime = "<value()<<",\t"; std::cout<<"Delta xPrime = "<value()<<",\t"; std::cout<<"Delta yPrime = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauNSCCartesianCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauNSCCartesianCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } void LauNSCCartesianCPCoeffSet::serialiseToJson( nlohmann::json& j ) const { // Call the base class method to do most of the work LauAbsCoeffSet::serialiseToJson(j); j["finalStateIsF"] = this->finalStateIsF(); } //! \cond DOXYGEN_IGNORE LauNSCCartesianCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::NSCCartesianCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauNSCCartesianCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauNSCCartesianCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("DeltaX", JsonType::Number), std::make_pair("DeltaY", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("DeltaXFixed", JsonType::Boolean), std::make_pair("DeltaYFixed", JsonType::Boolean), std::make_pair("XPrime", JsonType::Number), std::make_pair("YPrime", JsonType::Number), std::make_pair("DeltaXPrime", JsonType::Number), std::make_pair("DeltaYPrime", JsonType::Number), std::make_pair("XPrimeFixed", JsonType::Boolean), std::make_pair("YPrimeFixed", JsonType::Boolean), std::make_pair("DeltaXPrimeFixed", JsonType::Boolean), std::make_pair("DeltaYPrimeFixed", JsonType::Boolean), }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauNSCCartesianCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto finalStateIsF { getValue( j, "finalStateIsF" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto deltaX { getValue( j, "DeltaX" ) }; const auto deltaY { getValue( j, "DeltaY" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto deltaXFixed { getValue( j, "DeltaXFixed" ) }; const auto deltaYFixed { getValue( j, "DeltaYFixed" ) }; const auto deltaXSecondStage { getOptionalValue( j, "DeltaXSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto deltaYSecondStage { getOptionalValue( j, "DeltaYSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto xPrime { getValue( j, "XPrime" ) }; const auto yPrime { getValue( j, "YPrime" ) }; const auto deltaXPrime { getValue( j, "DeltaXPrime" ) }; const auto deltaYPrime { getValue( j, "DeltaYPrime" ) }; const auto xPrimeFixed { getValue( j, "XPrimeFixed" ) }; const auto yPrimeFixed { getValue( j, "YPrimeFixed" ) }; const auto deltaXPrimeFixed { getValue( j, "DeltaXPrimeFixed" ) }; const auto deltaYPrimeFixed { getValue( j, "DeltaYPrimeFixed" ) }; const auto deltaXPrimeSecondStage { getOptionalValue( j, "DeltaXPrimeSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; const auto deltaYPrimeSecondStage { getOptionalValue( j, "DeltaYPrimeSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; LauNSCCartesianCPCoeffSet coeff{ name, finalStateIsF, x, y, deltaX, deltaY, xFixed, yFixed, deltaXFixed, deltaYFixed, deltaXSecondStage, deltaYSecondStage, xPrime, yPrime, deltaXPrime, deltaYPrime, xPrimeFixed, yPrimeFixed, deltaXPrimeFixed, deltaYPrimeFixed, deltaXPrimeSecondStage, deltaYPrimeSecondStage }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauPolarGammaCPCoeffSet.cc b/src/LauPolarGammaCPCoeffSet.cc index 0c009f5..6c44d06 100644 --- a/src/LauPolarGammaCPCoeffSet.cc +++ b/src/LauPolarGammaCPCoeffSet.cc @@ -1,1047 +1,1042 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauPolarGammaCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.hh" #include "LauParameter.hh" #include "LauPrint.hh" std::shared_ptr LauPolarGammaCPCoeffSet::gammaGlobal_; std::shared_ptr LauPolarGammaCPCoeffSet::rDGlobal_; std::shared_ptr 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("X", x, minRealImagPart_, maxRealImagPart_, xFixed); y_ = std::make_unique("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("gamma", gamma, minPhase_, maxPhase_, gammaFixed); gamma_ = gammaGlobal_; } else { gamma_.reset( gammaGlobal_->createClone() ); } } else { gamma_ = std::make_shared("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("rB", rB, minMagnitude_, maxMagnitude_, rBFixed); deltaB_ = std::make_unique("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("rD", rD, minMagnitude_, maxMagnitude_, rDFixed); deltaDGlobal_ = std::make_shared("deltaD", deltaD, minPhase_, maxPhase_, deltaDFixed); rD_ = rDGlobal_; deltaD_ = deltaDGlobal_; } else { rD_.reset( rDGlobal_->createClone() ); deltaD_.reset( deltaDGlobal_->createClone() ); } } else { rD_ = std::make_shared("rD", rD, minMagnitude_, maxMagnitude_, rDFixed); deltaD_ = std::make_shared("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, 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_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto gammaVal { getOptionalValue( coeffInfo, "gamma", JsonType::Number ).value_or( rhs.gamma_->value() ) }; const auto gammaFixed { getOptionalValue( coeffInfo, "gammaFixed", JsonType::Boolean ).value_or( rhs.gamma_->fixed() ) }; const auto gammaSecondStage { getOptionalValue( coeffInfo, "gammaSecondStage", JsonType::Boolean ).value_or( rhs.gamma_->secondStage() ) }; gamma_ = std::make_shared("gamma", gammaVal, minPhase_, maxPhase_, gammaFixed); if ( gammaSecondStage && !gammaFixed ) { gamma_->secondStage(kTRUE); gamma_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "gammaBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "gammaBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "gammaBlindingWidth" ) }; gamma_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.gamma_->blind() ) { const LauBlind* blinder { rhs.gamma_->blinder() }; gamma_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) { const auto rBVal { getOptionalValue( coeffInfo, "rB", JsonType::Number ).value_or( rhs.rB_->value() ) }; const auto rBFixed { getOptionalValue( coeffInfo, "rBFixed", JsonType::Boolean ).value_or( rhs.rB_->fixed() ) }; const auto rBSecondStage { getOptionalValue( coeffInfo, "rBSecondStage", JsonType::Boolean ).value_or( rhs.rB_->secondStage() ) }; rB_ = std::make_unique("rB", rBVal, minMagnitude_, maxMagnitude_, rBFixed); if ( rBSecondStage && !rBFixed ) { rB_->secondStage(kTRUE); rB_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "rBBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "rBBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "rBBlindingWidth" ) }; rB_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.rB_->blind() ) { const LauBlind* blinder { rhs.rB_->blinder() }; rB_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto deltaBVal { getOptionalValue( coeffInfo, "deltaB", JsonType::Number ).value_or( rhs.deltaB_->value() ) }; const auto deltaBFixed { getOptionalValue( coeffInfo, "deltaBFixed", JsonType::Boolean ).value_or( rhs.deltaB_->fixed() ) }; const auto deltaBSecondStage { getOptionalValue( coeffInfo, "deltaBSecondStage", JsonType::Boolean ).value_or( rhs.deltaB_->secondStage() ) }; deltaB_ = std::make_unique("deltaB", deltaBVal, minPhase_, maxPhase_, deltaBFixed); if ( deltaBSecondStage && !deltaBFixed ) { deltaB_->secondStage(kTRUE); deltaB_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "deltaBBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "deltaBBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "deltaBBlindingWidth" ) }; deltaB_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaB_->blind() ) { const LauBlind* blinder { rhs.deltaB_->blinder() }; deltaB_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } 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 { const auto rDVal { getOptionalValue( coeffInfo, "rD", JsonType::Number ).value_or( rhs.rD_->value() ) }; const auto rDFixed { getOptionalValue( coeffInfo, "rDFixed", JsonType::Boolean ).value_or( rhs.rD_->fixed() ) }; const auto rDSecondStage { getOptionalValue( coeffInfo, "rDSecondStage", JsonType::Boolean ).value_or( rhs.rD_->secondStage() ) }; rD_ = std::make_unique("rD", rDVal, minMagnitude_, maxMagnitude_, rDFixed); if ( rDSecondStage && !rDFixed ) { rD_->secondStage(kTRUE); rD_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "rDBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "rDBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "rDBlindingWidth" ) }; rD_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.rD_->blind() ) { const LauBlind* blinder { rhs.rD_->blinder() }; rD_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto deltaDVal { getOptionalValue( coeffInfo, "deltaD", JsonType::Number ).value_or( rhs.deltaD_->value() ) }; const auto deltaDFixed { getOptionalValue( coeffInfo, "deltaDFixed", JsonType::Boolean ).value_or( rhs.deltaD_->fixed() ) }; const auto deltaDSecondStage { getOptionalValue( coeffInfo, "deltaDSecondStage", JsonType::Boolean ).value_or( rhs.deltaD_->secondStage() ) }; deltaD_ = std::make_unique("deltaD", deltaDVal, minPhase_, maxPhase_, deltaDFixed); if ( deltaDSecondStage && !deltaDFixed ) { deltaD_->secondStage(kTRUE); deltaD_->initValue(0.0); } if ( getOptionalValue( coeffInfo, "deltaDBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "deltaDBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "deltaDBlindingWidth" ) }; deltaD_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.deltaD_->blind() ) { const LauBlind* blinder { rhs.deltaD_->blinder() }; deltaD_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } } } 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 LauPolarGammaCPCoeffSet::getParameters() { std::vector 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 LauPolarGammaCPCoeffSet::getParameters() const { std::vector 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 LauPolarGammaCPCoeffSet::getParNames() const { std::vector pars; pars.reserve(7); pars.push_back("X"); pars.push_back("Y"); pars.push_back("gamma"); if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) { pars.push_back("rB"); pars.push_back("deltaB"); } if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) { pars.push_back("rD"); pars.push_back("deltaD"); } return pars; } void LauPolarGammaCPCoeffSet::printParValues() const { std::cout<<"INFO in LauPolarGammaCPCoeffSet::printParValues : Component \""<name()<<"\" has "; std::cout<<"X = "<value()<<",\t"; std::cout<<"Y = "<value()<<",\t"; if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::GLW_CPOdd || decayType_ == DecayType::GLW_CPEven ) { std::cout<<"rB = "<value()<<",\t"; std::cout<<"deltaB = "<value()<<",\t"; } if ( decayType_ == DecayType::ADS_Favoured || decayType_ == DecayType::ADS_Suppressed || decayType_ == DecayType::ADS_Favoured_btouOnly ) { std::cout<<"rD = "<value()<<",\t"; std::cout<<"deltaD = "<value()<<",\t"; } std::cout<<"gamma = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauPolarGammaCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauPolarGammaCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } std::ostream& operator<<( std::ostream& os, const LauPolarGammaCPCoeffSet::DecayType type ) { using DecayType = LauPolarGammaCPCoeffSet::DecayType; switch ( type ) { case DecayType::GLW_CPOdd : os << "GLW_CPOdd"; break; case DecayType::GLW_CPEven : os << "GLW_CPEven"; break; case DecayType::ADS_Favoured : os << "ADS_Favoured"; break; case DecayType::ADS_Suppressed : os << "ADS_Suppressed"; break; case DecayType::GLW_CPOdd_btouOnly : os << "GLW_CPOdd_btouOnly"; break; case DecayType::GLW_CPEven_btouOnly : os << "GLW_CPEven_btouOnly"; break; case DecayType::ADS_Favoured_btouOnly : os << "ADS_Favoured_btouOnly"; break; case DecayType::ADS_Suppressed_btouOnly : os << "ADS_Suppressed_btouOnly"; break; } return os; } void LauPolarGammaCPCoeffSet::serialiseToJson( nlohmann::json& j ) const { // Call the base class method to do most of the work LauAbsCoeffSet::serialiseToJson(j); const LauPolarGammaCPCoeffSet::DecayType decayType { this->decayType() }; j["decayType"] = decayType; j["useGlobalGamma"] = this->useGlobalGamma(); j["useGlobalADSPars"] = this->useGlobalADSPars(); if ( this->useGlobalGamma() && gamma_->clone() ) { j["gamma"] = gamma_->value(); j["gammaFixed"] = gamma_->fixed(); j["gammaSecondStage"] = gamma_->secondStage(); } if ( this->useGlobalADSPars() && rD_ && rD_->clone() ) { j["rD"] = rD_->value(); j["rDFixed"] = rD_->fixed(); j["rDSecondStage"] = rD_->secondStage(); } if ( this->useGlobalADSPars() && deltaD_ && deltaD_->clone() ) { j["deltaD"] = deltaD_->value(); j["deltaDFixed"] = deltaD_->fixed(); j["deltaDSecondStage"] = deltaD_->secondStage(); } } //! \cond DOXYGEN_IGNORE LauPolarGammaCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; using DecayType = LauPolarGammaCPCoeffSet::DecayType; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::PolarGammaCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauPolarGammaCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauPolarGammaCPCoeffSet standalone"}; } const auto decayType { getValue( j, "decayType" ) }; std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("Y", JsonType::Number), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("gamma", JsonType::Number), std::make_pair("gammaFixed", JsonType::Boolean) }; if ( decayType == DecayType::ADS_Favoured || decayType == DecayType::ADS_Suppressed || decayType == DecayType::GLW_CPOdd || decayType == DecayType::GLW_CPEven ) { mandatoryElements.push_back( std::make_pair( "rB", JsonType::Number ) ); mandatoryElements.push_back( std::make_pair( "rBFixed", JsonType::Boolean ) ); mandatoryElements.push_back( std::make_pair( "deltaB", JsonType::Number ) ); mandatoryElements.push_back( std::make_pair( "deltaBFixed", JsonType::Boolean ) ); } if ( decayType == DecayType::ADS_Favoured || decayType == DecayType::ADS_Suppressed || decayType == DecayType::ADS_Favoured_btouOnly ) { mandatoryElements.push_back( std::make_pair( "rD", JsonType::Number ) ); mandatoryElements.push_back( std::make_pair( "rDFixed", JsonType::Boolean ) ); mandatoryElements.push_back( std::make_pair( "deltaD", JsonType::Number ) ); mandatoryElements.push_back( std::make_pair( "deltaDFixed", JsonType::Boolean ) ); } if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauPolarGammaCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto useGlobalGamma { getOptionalValue( j, "useGlobalGamma", JsonType::Boolean ).value_or( kFALSE ) }; const auto useGlobalADSPars { getOptionalValue( j, "useGlobalADSPars", JsonType::Boolean ).value_or( kFALSE ) }; const auto x { getValue( j, "X" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto y { getValue( j, "Y" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto gamma { getValue( j, "gamma" ) }; const auto gammaFixed { getValue( j, "gammaFixed" ) }; const auto gammaSecondStage { getOptionalValue( j, "gammaSecondStage", JsonType::Boolean ).value_or( kFALSE ) }; Double_t rB{0.0}; Bool_t rBFixed{kTRUE}; Bool_t rBSecondStage{kFALSE}; Double_t deltaB{0.0}; Bool_t deltaBFixed{kTRUE}; Bool_t deltaBSecondStage{kFALSE}; Double_t rD{0.0}; Bool_t rDFixed{kTRUE}; Bool_t rDSecondStage{kFALSE}; Double_t deltaD{0.0}; Bool_t deltaDFixed{kTRUE}; Bool_t deltaDSecondStage{kFALSE}; if ( decayType == DecayType::ADS_Favoured || decayType == DecayType::ADS_Suppressed || decayType == DecayType::GLW_CPOdd || decayType == DecayType::GLW_CPEven ) { rB = getValue( j, "rB" ); rBFixed = getValue( j, "rBFixed" ); rBSecondStage = getOptionalValue( j, "rBSecondStage", JsonType::Boolean ).value_or(kFALSE); deltaB = getValue( j, "deltaB" ); deltaBFixed = getValue( j, "deltaBFixed" ); deltaBSecondStage = getOptionalValue( j, "deltaBSecondStage", JsonType::Boolean ).value_or(kFALSE); } if ( decayType == DecayType::ADS_Favoured || decayType == DecayType::ADS_Suppressed || decayType == DecayType::ADS_Favoured_btouOnly ) { rD = getValue( j, "rD" ); rDFixed = getValue( j, "rDFixed" ); rDSecondStage = getOptionalValue( j, "rDSecondStage", JsonType::Boolean ).value_or(kFALSE); deltaD = getValue( j, "deltaD" ); deltaDFixed = getValue( j, "deltaDFixed" ); deltaDSecondStage = getOptionalValue( j, "deltaDSecondStage", JsonType::Boolean ).value_or(kFALSE); } LauPolarGammaCPCoeffSet coeff{ name, decayType, x, y, rB, deltaB, gamma, rD, deltaD, xFixed, yFixed, rBFixed, deltaBFixed, gammaFixed, rDFixed, deltaDFixed, rBSecondStage, deltaBSecondStage, gammaSecondStage, rDSecondStage, deltaDSecondStage, useGlobalGamma, useGlobalADSPars }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauRealImagCPCoeffSet.cc b/src/LauRealImagCPCoeffSet.cc index 348b789..cef9c28 100644 --- a/src/LauRealImagCPCoeffSet.cc +++ b/src/LauRealImagCPCoeffSet.cc @@ -1,340 +1,335 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauRealImagCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("X", x, minRealImagPart_, maxRealImagPart_, xFixed) }, y_{ std::make_unique("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) }, xbar_{ std::make_unique("XBar", xbar, minRealImagPart_, maxRealImagPart_, xbarFixed) }, ybar_{ std::make_unique("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, const CloneOption cloneOption, const Double_t constFactor, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); xbar_.reset( rhs.xbar_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.x_->blind() ) { const LauBlind* blinder { rhs.x_->blinder() }; x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto xbarVal { getOptionalValue( coeffInfo, "XBar", JsonType::Number ).value_or( rhs.xbar_->value() ) }; const auto xbarFixed { getOptionalValue( coeffInfo, "XBarFixed", JsonType::Boolean ).value_or( rhs.xbar_->fixed() ) }; xbar_ = std::make_unique("XBar", xbarVal, minRealImagPart_, maxRealImagPart_, xbarFixed); if ( getOptionalValue( coeffInfo, "XBarBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBarBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBarBlindingWidth" ) }; xbar_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.y_->blind() ) { const LauBlind* blinder { rhs.y_->blinder() }; y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto ybarVal { getOptionalValue( coeffInfo, "YBar", JsonType::Number ).value_or( rhs.ybar_->value() ) }; const auto ybarFixed { getOptionalValue( coeffInfo, "YBarFixed", JsonType::Boolean ).value_or( rhs.ybar_->fixed() ) }; ybar_ = std::make_unique("YBar", ybarVal, minRealImagPart_, maxRealImagPart_, ybarFixed); if ( getOptionalValue( coeffInfo, "YBarBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBarBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBarBlindingWidth" ) }; ybar_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.ybar_->blind() ) { const LauBlind* blinder { rhs.ybar_->blinder() }; ybar_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } 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|}"<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart ) ) { std::cerr << "ERROR in LauRealImagCPCoeffSet::createClone : Invalid clone option" << std::endl; return nullptr; } - if ( this->clone() ) { - const LauRealImagCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauRealImagCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauRealImagCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::RealImagCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauRealImagCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("XBar", JsonType::Number), std::make_pair("YBar", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("XBarFixed", JsonType::Boolean), std::make_pair("YBarFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauRealImagCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto xbar { getValue( j, "XBar" ) }; const auto ybar { getValue( j, "YBar" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto xbarFixed { getValue( j, "XBarFixed" ) }; const auto ybarFixed { getValue( j, "YBarFixed" ) }; LauRealImagCPCoeffSet coeff{ name, x, y, xbar, ybar, xFixed, yFixed, xbarFixed, ybarFixed }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauRealImagCoeffSet.cc b/src/LauRealImagCoeffSet.cc index 73002c6..efcad1e 100644 --- a/src/LauRealImagCoeffSet.cc +++ b/src/LauRealImagCoeffSet.cc @@ -1,244 +1,239 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauRealImagCoeffSet.hh" #include "LauJsonTools.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("X",x,minRealImagPart_,maxRealImagPart_,xFixed) }, y_{ std::make_unique("Y",y,minRealImagPart_,maxRealImagPart_,yFixed) }, coeff_{ x,y } { } LauRealImagCoeffSet::LauRealImagCoeffSet(const LauRealImagCoeffSet& rhs, const CloneOption cloneOption, const Double_t constFactor, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, coeff_{ rhs.coeff_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else 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 \""<name()<<"\" has real part = "<value()<<" and imaginary part = "<value()<<"."<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<value()); stream<<" \\pm "; print.printFormat(stream, x_->error()); stream<<"$ & $"; print.printFormat(stream, y_->value()); stream<<" \\pm "; print.printFormat(stream, y_->error()); stream<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { if ( ! ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart || cloneOption == CloneOption::TieImagPart ) ) { std::cerr << "ERROR in LauRealImagCoeffSet::createClone : Invalid clone option" << std::endl; return nullptr; } - if ( this->clone() ) { - const LauRealImagCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauRealImagCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauRealImagCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::RealImag ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauRealImagCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauRealImagCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; LauRealImagCoeffSet coeff{ name, x, y, xFixed, yFixed }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE diff --git a/src/LauRealImagGammaCPCoeffSet.cc b/src/LauRealImagGammaCPCoeffSet.cc index ec9ebb7..10238a9 100644 --- a/src/LauRealImagGammaCPCoeffSet.cc +++ b/src/LauRealImagGammaCPCoeffSet.cc @@ -1,403 +1,398 @@ /* 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 #include #include #include "TMath.h" #include "TRandom.h" #include "LauRealImagGammaCPCoeffSet.hh" #include "LauComplex.hh" #include "LauConstants.hh" #include "LauJsonTools.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("X", x, minRealImagPart_, maxRealImagPart_, xFixed) }, y_{ std::make_unique("Y", y, minRealImagPart_, maxRealImagPart_, yFixed) }, xCP_{ std::make_unique("XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed) }, yCP_{ std::make_unique("YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed) }, xbarCP_{ std::make_unique("XBarCP", xbarCP, minRealImagPart_, maxRealImagPart_, xbarCPFixed) }, ybarCP_{ std::make_unique("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, const nlohmann::json& coeffInfo) : LauAbsCoeffSet{ rhs.name(), rhs.baseName(), &rhs, cloneOption, constFactor }, nonCPPart_{ rhs.nonCPPart_ }, cpPart_{ rhs.cpPart_ }, cpAntiPart_{ rhs.cpAntiPart_ }, particleCoeff_{ rhs.particleCoeff_ }, antiparticleCoeff_{ rhs.antiparticleCoeff_ }, acp_{ rhs.acp_ } { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; if ( cloneOption == CloneOption::All || cloneOption == CloneOption::TieRealPart ) { x_.reset( rhs.x_->createClone(constFactor) ); } else { const auto xVal { getOptionalValue( coeffInfo, "X", JsonType::Number ).value_or( rhs.x_->value() ) }; const auto xFixed { getOptionalValue( coeffInfo, "XFixed", JsonType::Boolean ).value_or( rhs.x_->fixed() ) }; x_ = std::make_unique("X", xVal, minRealImagPart_, maxRealImagPart_, xFixed); if ( getOptionalValue( coeffInfo, "XBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBlindingWidth" ) }; x_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto yVal { getOptionalValue( coeffInfo, "Y", JsonType::Number ).value_or( rhs.y_->value() ) }; const auto yFixed { getOptionalValue( coeffInfo, "YFixed", JsonType::Boolean ).value_or( rhs.y_->fixed() ) }; y_ = std::make_unique("Y", yVal, minRealImagPart_, maxRealImagPart_, yFixed); if ( getOptionalValue( coeffInfo, "YBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBlindingWidth" ) }; y_->blindParameter( blindingString, blindingWidth ); } else 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 { const auto xCPVal { getOptionalValue( coeffInfo, "XCP", JsonType::Number ).value_or( rhs.xCP_->value() ) }; const auto xCPFixed { getOptionalValue( coeffInfo, "XCPFixed", JsonType::Boolean ).value_or( rhs.xCP_->fixed() ) }; xCP_ = std::make_unique("XCP", xCPVal, minRealImagPart_, maxRealImagPart_, xCPFixed); if ( getOptionalValue( coeffInfo, "XCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XCPBlindingWidth" ) }; xCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.xCP_->blind() ) { const LauBlind* blinder { rhs.xCP_->blinder() }; xCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto yCPVal { getOptionalValue( coeffInfo, "YCP", JsonType::Number ).value_or( rhs.yCP_->value() ) }; const auto yCPFixed { getOptionalValue( coeffInfo, "YCPFixed", JsonType::Boolean ).value_or( rhs.yCP_->fixed() ) }; yCP_ = std::make_unique("YCP", yCPVal, minRealImagPart_, maxRealImagPart_, yCPFixed); if ( getOptionalValue( coeffInfo, "YCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YCPBlindingWidth" ) }; yCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.yCP_->blind() ) { const LauBlind* blinder { rhs.yCP_->blinder() }; yCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto xbarCPVal { getOptionalValue( coeffInfo, "XBarCP", JsonType::Number ).value_or( rhs.xbarCP_->value() ) }; const auto xbarCPFixed { getOptionalValue( coeffInfo, "XBarCPFixed", JsonType::Boolean ).value_or( rhs.xbarCP_->fixed() ) }; xbarCP_ = std::make_unique("XBarCP", xbarCPVal, minRealImagPart_, maxRealImagPart_, xbarCPFixed); if ( getOptionalValue( coeffInfo, "XBarCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "XBarCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "XBarCPBlindingWidth" ) }; xbarCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.xbarCP_->blind() ) { const LauBlind* blinder { rhs.xbarCP_->blinder() }; xbarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } const auto ybarCPVal { getOptionalValue( coeffInfo, "YBarCP", JsonType::Number ).value_or( rhs.ybarCP_->value() ) }; const auto ybarCPFixed { getOptionalValue( coeffInfo, "YBarCPFixed", JsonType::Boolean ).value_or( rhs.ybarCP_->fixed() ) }; ybarCP_ = std::make_unique("YBarCP", ybarCPVal, minRealImagPart_, maxRealImagPart_, ybarCPFixed); if ( getOptionalValue( coeffInfo, "YBarCPBlind", JsonType::Boolean ).value_or( kFALSE ) ) { const auto blindingString { getValue( coeffInfo, "YBarCPBlindingString" ) }; const auto blindingWidth { getValue( coeffInfo, "YBarCPBlindingWidth" ) }; ybarCP_->blindParameter( blindingString, blindingWidth ); } else if ( rhs.ybarCP_->blind() ) { const LauBlind* blinder { rhs.ybarCP_->blinder() }; ybarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() ); } } } 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|}"<name() }; resName = resName.ReplaceAll("_", "\\_"); stream<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<<"$ \\\\"<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, const nlohmann::json& coeffInfo) const { 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; } - if ( this->clone() ) { - const LauRealImagGammaCPCoeffSet* parent { static_cast(this->parent()) }; - return parent->createClone_impl( newName, cloneOption, constFactor, coeffInfo ); - } - auto clone = new LauRealImagGammaCPCoeffSet{ *this, cloneOption, constFactor, coeffInfo }; clone->name( newName ); return clone; } //! \cond DOXYGEN_IGNORE LauRealImagGammaCPCoeffSet nlohmann::adl_serializer::from_json(const json& j) { using LauJsonTools::JsonType; using LauJsonTools::getValue; using LauJsonTools::getOptionalValue; const auto type { getValue( j, "type" ) }; if ( type != LauCoeffType::RealImagGammaCP ) { throw LauWrongCoeffType("Wrong coefficient type given to construct LauRealImagGammaCPCoeffSet"); } const auto clone { getValue( j, "clone" ) }; if ( clone ) { throw LauClonedCoeff{"Cannot build a cloned LauRealImagGammaCPCoeffSet standalone"}; } const std::vector mandatoryElements { std::make_pair("X", JsonType::Number), std::make_pair("Y", JsonType::Number), std::make_pair("XCP", JsonType::Number), std::make_pair("YCP", JsonType::Number), std::make_pair("XBarCP", JsonType::Number), std::make_pair("YBarCP", JsonType::Number), std::make_pair("XFixed", JsonType::Boolean), std::make_pair("YFixed", JsonType::Boolean), std::make_pair("XCPFixed", JsonType::Boolean), std::make_pair("YCPFixed", JsonType::Boolean), std::make_pair("XBarCPFixed", JsonType::Boolean), std::make_pair("YBarCPFixed", JsonType::Boolean) }; if ( ! LauJsonTools::checkObjectElements( j, mandatoryElements ) ) { throw LauJsonTools::MissingJsonElement{"Missing elements needed to construct LauRealImagGammaCPCoeffSet"}; } const auto name { getValue( j, "name" ) }; const auto x { getValue( j, "X" ) }; const auto y { getValue( j, "Y" ) }; const auto xCP { getValue( j, "XCP" ) }; const auto yCP { getValue( j, "YCP" ) }; const auto xbarCP { getValue( j, "XBarCP" ) }; const auto ybarCP { getValue( j, "YBarCP" ) }; const auto xFixed { getValue( j, "XFixed" ) }; const auto yFixed { getValue( j, "YFixed" ) }; const auto xCPFixed { getValue( j, "XCPFixed" ) }; const auto yCPFixed { getValue( j, "YCPFixed" ) }; const auto xbarCPFixed { getValue( j, "XBarCPFixed" ) }; const auto ybarCPFixed { getValue( j, "YBarCPFixed" ) }; LauRealImagGammaCPCoeffSet coeff{ name, x, y, xCP, yCP, xbarCP, ybarCP, xFixed, yFixed, xCPFixed, yCPFixed, xbarCPFixed, ybarCPFixed }; coeff.applyBlinding( j ); return coeff; } //! \endcond DOXYGEN_IGNORE