Changeset View
Changeset View
Standalone View
Standalone View
src/LauCPFitModel.cc
Show First 20 Lines • Show All 1,417 Lines • ▼ Show 20 Lines | std::pair<LauCPFitModel::LauGenInfo,Bool_t> LauCPFitModel::eventsToGenerate() | ||||
// If we're smearing then smear each one individually | // If we're smearing then smear each one individually | ||||
LauGenInfo nEvtsGen; | LauGenInfo nEvtsGen; | ||||
// Keep track of whether any yield or asymmetry parameters are blinded | // Keep track of whether any yield or asymmetry parameters are blinded | ||||
Bool_t blind = kFALSE; | Bool_t blind = kFALSE; | ||||
// Signal | // Signal | ||||
if ( signalEvents_->blind() ) { | |||||
blind = kTRUE; | |||||
} | |||||
Double_t evtWeight(1.0); | Double_t evtWeight(1.0); | ||||
Double_t nEvts = signalEvents_->genValue(); | Double_t nEvts = signalEvents_->genValue(); | ||||
if ( nEvts < 0.0 ) { | if ( nEvts < 0.0 ) { | ||||
evtWeight = -1.0; | evtWeight = -1.0; | ||||
nEvts = TMath::Abs( nEvts ); | nEvts = TMath::Abs( nEvts ); | ||||
} | } | ||||
if ( signalEvents_->blind() ) { | |||||
blind = kTRUE; | |||||
} | |||||
Double_t asym(0.0); | |||||
Double_t sigAsym(0.0); | Double_t sigAsym(0.0); | ||||
// need to include this as an alternative in case the DP isn't in the model | // need to include this as an alternative in case the DP isn't in the model | ||||
if ( !this->useDP() || forceAsym_ ) { | if ( !this->useDP() || forceAsym_ ) { | ||||
sigAsym = signalAsym_->genValue(); | sigAsym = signalAsym_->genValue(); | ||||
if ( signalAsym_->blind() ) { | if ( signalAsym_->blind() ) { | ||||
blind = kTRUE; | blind = kTRUE; | ||||
} | } | ||||
} else { | } else { | ||||
Double_t negRate = negSigModel_->getDPNorm(); | Double_t negRate = negSigModel_->getDPNorm(); | ||||
Double_t posRate = posSigModel_->getDPNorm(); | Double_t posRate = posSigModel_->getDPNorm(); | ||||
if (negRate+posRate>1e-30) { | if (negRate+posRate>1e-30) { | ||||
sigAsym = (negRate-posRate)/(negRate+posRate); | sigAsym = (negRate-posRate)/(negRate+posRate); | ||||
} | } | ||||
} | } | ||||
asym = sigAsym; | Double_t nPosEvts = (nEvts/2.0 * (1.0 - sigAsym)); | ||||
Int_t nPosEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 - asym)) + 0.5); | Double_t nNegEvts = (nEvts/2.0 * (1.0 + sigAsym)); | ||||
Int_t nNegEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 + asym)) + 0.5); | |||||
Int_t nPosEvtsToGen { static_cast<Int_t>(nPosEvts) }; | |||||
Int_t nNegEvtsToGen { static_cast<Int_t>(nNegEvts) }; | |||||
if (this->doPoissonSmearing()) { | if (this->doPoissonSmearing()) { | ||||
nNegEvts = LauRandom::randomFun()->Poisson(nNegEvts); | nPosEvtsToGen = LauRandom::randomFun()->Poisson(nPosEvts); | ||||
nPosEvts = LauRandom::randomFun()->Poisson(nPosEvts); | nNegEvtsToGen = LauRandom::randomFun()->Poisson(nNegEvts); | ||||
} | } | ||||
nEvtsGen[std::make_pair("signal",-1)] = std::make_pair(nNegEvts,evtWeight); | |||||
nEvtsGen[std::make_pair("signal",+1)] = std::make_pair(nPosEvts,evtWeight); | nEvtsGen[std::make_pair("signal",+1)] = std::make_pair(nPosEvtsToGen,evtWeight); | ||||
nEvtsGen[std::make_pair("signal",-1)] = std::make_pair(nNegEvtsToGen,evtWeight); | |||||
// backgrounds | // backgrounds | ||||
const UInt_t nBkgnds = this->nBkgndClasses(); | const UInt_t nBkgnds = this->nBkgndClasses(); | ||||
for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) { | for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) { | ||||
const TString& bkgndClass = this->bkgndClassName(bkgndID); | |||||
const LauAbsRValue* evtsPar = bkgndEvents_[bkgndID]; | const LauAbsRValue* evtsPar = bkgndEvents_[bkgndID]; | ||||
const LauAbsRValue* asymPar = bkgndAsym_[bkgndID]; | const LauAbsRValue* asymPar = bkgndAsym_[bkgndID]; | ||||
if ( evtsPar->blind() || asymPar->blind() ) { | if ( evtsPar->blind() || asymPar->blind() ) { | ||||
blind = kTRUE; | blind = kTRUE; | ||||
} | } | ||||
evtWeight = 1.0; | evtWeight = 1.0; | ||||
nEvts = TMath::FloorNint( evtsPar->genValue() ); | nEvts = evtsPar->genValue(); | ||||
if ( nEvts < 0 ) { | if ( nEvts < 0 ) { | ||||
evtWeight = -1.0; | evtWeight = -1.0; | ||||
nEvts = TMath::Abs( nEvts ); | nEvts = TMath::Abs( nEvts ); | ||||
} | } | ||||
asym = asymPar->genValue(); | |||||
nPosEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 - asym)) + 0.5); | const Double_t asym = asymPar->genValue(); | ||||
nNegEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 + asym)) + 0.5); | nPosEvts = (nEvts/2.0 * (1.0 - asym)); | ||||
nNegEvts = (nEvts/2.0 * (1.0 + asym)); | |||||
nPosEvtsToGen = static_cast<Int_t>(nPosEvts); | |||||
nNegEvtsToGen = static_cast<Int_t>(nNegEvts); | |||||
if (this->doPoissonSmearing()) { | if (this->doPoissonSmearing()) { | ||||
nNegEvts = LauRandom::randomFun()->Poisson(nNegEvts); | nPosEvtsToGen = LauRandom::randomFun()->Poisson(nPosEvts); | ||||
nPosEvts = LauRandom::randomFun()->Poisson(nPosEvts); | nNegEvtsToGen = LauRandom::randomFun()->Poisson(nNegEvts); | ||||
} | } | ||||
nEvtsGen[std::make_pair(bkgndClass,-1)] = std::make_pair(nNegEvts,evtWeight); | |||||
nEvtsGen[std::make_pair(bkgndClass,+1)] = std::make_pair(nPosEvts,evtWeight); | const TString& bkgndClass = this->bkgndClassName(bkgndID); | ||||
nEvtsGen[std::make_pair(bkgndClass,+1)] = std::make_pair(nPosEvtsToGen,evtWeight); | |||||
nEvtsGen[std::make_pair(bkgndClass,-1)] = std::make_pair(nNegEvtsToGen,evtWeight); | |||||
} | } | ||||
// Print out the information on what we're generating, but only if none of the parameters are blind (otherwise we risk unblinding them!) | // Print out the information on what we're generating, but only if none of the parameters are blind (otherwise we risk unblinding them!) | ||||
if ( !blind ) { | if ( !blind ) { | ||||
std::cout << "INFO in LauCPFitModel::eventsToGenerate : Generating toy MC with:" << std::endl; | std::cout << "INFO in LauCPFitModel::eventsToGenerate : Generating toy MC with:" << std::endl; | ||||
std::cout << " : Signal asymmetry = " << sigAsym << " and number of signal events = " << signalEvents_->genValue() << std::endl; | std::cout << " : Signal asymmetry = " << sigAsym << " and number of signal events = " << signalEvents_->genValue() << std::endl; | ||||
for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) { | for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) { | ||||
const TString& bkgndClass = this->bkgndClassName(bkgndID); | const TString& bkgndClass = this->bkgndClassName(bkgndID); | ||||
▲ Show 20 Lines • Show All 1,949 Lines • Show Last 20 Lines |