diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index b0646ef..2456f73 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -1,44 +1,48 @@ list(APPEND EXAMPLE_SOURCES B3piKMatrixMassProj B3piKMatrixPlots CalcChiSq + FlatPhaseSpace + FlatPhaseSpace-CustomMasses + FlatSqDalitz + FlatSqDalitz-CustomMasses GenFit3K GenFit3KS GenFit3pi GenFitBelleCPKpipi GenFitDpipi GenFitDs2KKpi GenFitEFKLLM GenFitKpipi GenFitNoDP GenFitNoDPMultiDim KMatrixDto3pi KMatrixExample MergeDataFiles mixedSampleTest PlotKMatrixTAmp PlotResults point2PointTestSample QuasiFlatSqDalitz QuasiFlatSqDalitz-CustomMasses ResultsExtractor SimFitCoordinator SimFitTask SimFitTaskRooFit ) if(NOT LAURA_BUILD_ROOFIT_TASK) list(REMOVE_ITEM EXAMPLE_SOURCES SimFitTaskRooFit) endif() foreach( _example ${EXAMPLE_SOURCES}) add_executable(${_example} ${_example}.cc) target_link_libraries(${_example} PRIVATE Laura++) install(TARGETS ${_example} RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) endforeach() # Also install the python script version of GenFit3pi configure_file(GenFit3pi.py.in GenFit3pi.py @ONLY) install(PROGRAMS ${CMAKE_CURRENT_BINARY_DIR}/GenFit3pi.py DESTINATION ${CMAKE_INSTALL_BINDIR}) diff --git a/examples/QuasiFlatSqDalitz-CustomMasses.cc b/examples/FlatPhaseSpace-CustomMasses.cc similarity index 64% copy from examples/QuasiFlatSqDalitz-CustomMasses.cc copy to examples/FlatPhaseSpace-CustomMasses.cc index 017678f..a771b8c 100644 --- a/examples/QuasiFlatSqDalitz-CustomMasses.cc +++ b/examples/FlatPhaseSpace-CustomMasses.cc @@ -1,114 +1,100 @@ /* -Copyright 2018 University of Warwick +Copyright 2022 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauKinematics.hh" #include "LauRandom.hh" /* - * Generates toy MC according to EvtGen's FLATSQDALITZ model + * Generates toy MC uniformly in the Dalitz plot * - * Usage: QuasiFlatSqDalitz + * Usage: FlatPhaseSpace-CustomMasses */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 7 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const Double_t parent_mass { std::stod( cmdLineArgs[1] ) }; const Double_t child_1_mass { std::stod( cmdLineArgs[2] ) }; const Double_t child_2_mass { std::stod( cmdLineArgs[3] ) }; const Double_t child_3_mass { std::stod( cmdLineArgs[4] ) }; const ULong_t nEvents { std::stoul( cmdLineArgs[5] ) }; const TString outputFileName { cmdLineArgs[6] }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE const Bool_t squareDP { kTRUE }; LauKinematics kinematics { child_1_mass, child_2_mass, child_3_mass, parent_mass, squareDP }; - Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; + Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}; TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); - tree->Branch("mPrime",&mPrime,"mPrime/D"); - tree->Branch("thPrime",&thPrime,"thPrime/D"); - tree->Branch("jacobian",&jacobian,"jacobian/D"); - tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); - tree->Branch("randomNum",&randomNum,"randomNum/D"); + if ( squareDP ) { + tree->Branch("mPrime",&mPrime,"mPrime/D"); + tree->Branch("thPrime",&thPrime,"thPrime/D"); + } - Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } - evtAccepted = kFALSE; - - while ( !evtAccepted ) { - kinematics.genFlatPhaseSpace( m13Sq, m23Sq ); - kinematics.updateKinematics( m13Sq, m23Sq ); - jacobian = kinematics.calcSqDPJacobian(); - invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; - - randomNum = LauRandom::randomFun()->Rndm(); - if ( randomNum < invJacobian ) { - evtAccepted = kTRUE; - m12Sq = kinematics.getm12Sq(); - mPrime = kinematics.getmPrime(); - thPrime = kinematics.getThetaPrime(); - - tree->Fill(); - } + kinematics.genFlatPhaseSpace( m13Sq, m23Sq ); + kinematics.updateKinematics( m13Sq, m23Sq ); + m12Sq = kinematics.getm12Sq(); + if ( squareDP ) { + mPrime = kinematics.getmPrime(); + thPrime = kinematics.getThetaPrime(); } + + tree->Fill(); } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; } diff --git a/examples/QuasiFlatSqDalitz.cc b/examples/FlatPhaseSpace.cc similarity index 54% copy from examples/QuasiFlatSqDalitz.cc copy to examples/FlatPhaseSpace.cc index 272039f..c505a9b 100644 --- a/examples/QuasiFlatSqDalitz.cc +++ b/examples/FlatPhaseSpace.cc @@ -1,122 +1,102 @@ /* -Copyright 2018 University of Warwick +Copyright 2022 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauDaughters.hh" #include "LauKinematics.hh" #include "LauRandom.hh" /* - * Generates toy MC according to EvtGen's FLATSQDALITZ model + * Generates toy MC uniformly in the Dalitz plot * - * Usage: QuasiFlatSqDalitz + * Usage: FlatPhaseSpace */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 6 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const TString& parent { cmdLineArgs[1] }; const TString& child_1 { cmdLineArgs[2] }; const TString& child_2 { cmdLineArgs[3] }; const TString& child_3 { cmdLineArgs[4] }; const ULong_t nEvents { std::stoul( cmdLineArgs[5].Data() ) }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE - Bool_t squareDP { kTRUE }; + const Bool_t squareDP { kTRUE }; - // This defines the DP => decay is B+ -> pi+ pi+ pi- - // Particle 1 = pi+ - // Particle 2 = pi+ - // Particle 3 = pi- - // The DP is defined in terms of m13Sq and m23Sq LauDaughters* daughters { new LauDaughters(parent, child_1, child_2, child_3, squareDP) }; - LauKinematics* kinematics { daughters->getKinematics() }; - Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; + Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}; - TString fileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; - TFile * file { TFile::Open(fileName,"recreate") }; + const TString outputFileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; + TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); - tree->Branch("mPrime",&mPrime,"mPrime/D"); - tree->Branch("thPrime",&thPrime,"thPrime/D"); - tree->Branch("jacobian",&jacobian,"jacobian/D"); - tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); - tree->Branch("randomNum",&randomNum,"randomNum/D"); + if ( squareDP ) { + tree->Branch("mPrime",&mPrime,"mPrime/D"); + tree->Branch("thPrime",&thPrime,"thPrime/D"); + } - Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } - evtAccepted = kFALSE; - - while ( !evtAccepted ) { - kinematics->genFlatPhaseSpace( m13Sq, m23Sq ); - kinematics->updateKinematics( m13Sq, m23Sq ); - jacobian = kinematics->calcSqDPJacobian(); - invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; - - randomNum = LauRandom::randomFun()->Rndm(); - if ( randomNum < invJacobian ) { - evtAccepted = kTRUE; - m12Sq = kinematics->getm12Sq(); - mPrime = kinematics->getmPrime(); - thPrime = kinematics->getThetaPrime(); - - tree->Fill(); - } + kinematics->genFlatPhaseSpace( m13Sq, m23Sq ); + kinematics->updateKinematics( m13Sq, m23Sq ); + m12Sq = kinematics->getm12Sq(); + if ( squareDP ) { + mPrime = kinematics->getmPrime(); + thPrime = kinematics->getThetaPrime(); } + + tree->Fill(); } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; } diff --git a/examples/QuasiFlatSqDalitz-CustomMasses.cc b/examples/FlatSqDalitz-CustomMasses.cc similarity index 66% copy from examples/QuasiFlatSqDalitz-CustomMasses.cc copy to examples/FlatSqDalitz-CustomMasses.cc index 017678f..a63a2f8 100644 --- a/examples/QuasiFlatSqDalitz-CustomMasses.cc +++ b/examples/FlatSqDalitz-CustomMasses.cc @@ -1,114 +1,98 @@ /* -Copyright 2018 University of Warwick +Copyright 2022 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauKinematics.hh" #include "LauRandom.hh" /* - * Generates toy MC according to EvtGen's FLATSQDALITZ model + * Generates toy MC uniformly in the square Dalitz plot * - * Usage: QuasiFlatSqDalitz + * Usage: FlatSqDalitz-CustomMasses */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 7 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const Double_t parent_mass { std::stod( cmdLineArgs[1] ) }; const Double_t child_1_mass { std::stod( cmdLineArgs[2] ) }; const Double_t child_2_mass { std::stod( cmdLineArgs[3] ) }; const Double_t child_3_mass { std::stod( cmdLineArgs[4] ) }; const ULong_t nEvents { std::stoul( cmdLineArgs[5] ) }; const TString outputFileName { cmdLineArgs[6] }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE + // Do not change this! const Bool_t squareDP { kTRUE }; LauKinematics kinematics { child_1_mass, child_2_mass, child_3_mass, parent_mass, squareDP }; - Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; + Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}; TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); tree->Branch("mPrime",&mPrime,"mPrime/D"); tree->Branch("thPrime",&thPrime,"thPrime/D"); - tree->Branch("jacobian",&jacobian,"jacobian/D"); - tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); - tree->Branch("randomNum",&randomNum,"randomNum/D"); - Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } - evtAccepted = kFALSE; + kinematics.genFlatSqDP( mPrime, thPrime ); + kinematics.updateSqDPKinematics( mPrime, thPrime ); - while ( !evtAccepted ) { - kinematics.genFlatPhaseSpace( m13Sq, m23Sq ); - kinematics.updateKinematics( m13Sq, m23Sq ); - jacobian = kinematics.calcSqDPJacobian(); - invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; + m13Sq = kinematics.getm13Sq(); + m23Sq = kinematics.getm23Sq(); + m12Sq = kinematics.getm12Sq(); - randomNum = LauRandom::randomFun()->Rndm(); - if ( randomNum < invJacobian ) { - evtAccepted = kTRUE; - m12Sq = kinematics.getm12Sq(); - mPrime = kinematics.getmPrime(); - thPrime = kinematics.getThetaPrime(); - - tree->Fill(); - } - } + tree->Fill(); } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; } diff --git a/examples/QuasiFlatSqDalitz.cc b/examples/FlatSqDalitz.cc similarity index 57% copy from examples/QuasiFlatSqDalitz.cc copy to examples/FlatSqDalitz.cc index 272039f..e43dcf8 100644 --- a/examples/QuasiFlatSqDalitz.cc +++ b/examples/FlatSqDalitz.cc @@ -1,122 +1,100 @@ /* -Copyright 2018 University of Warwick +Copyright 2022 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauDaughters.hh" #include "LauKinematics.hh" #include "LauRandom.hh" /* - * Generates toy MC according to EvtGen's FLATSQDALITZ model + * Generates toy MC uniformly in the square Dalitz plot * - * Usage: QuasiFlatSqDalitz + * Usage: FlatSqDalitz */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 6 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const TString& parent { cmdLineArgs[1] }; const TString& child_1 { cmdLineArgs[2] }; const TString& child_2 { cmdLineArgs[3] }; const TString& child_3 { cmdLineArgs[4] }; const ULong_t nEvents { std::stoul( cmdLineArgs[5].Data() ) }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE - Bool_t squareDP { kTRUE }; + // Do not change this! + const Bool_t squareDP { kTRUE }; - // This defines the DP => decay is B+ -> pi+ pi+ pi- - // Particle 1 = pi+ - // Particle 2 = pi+ - // Particle 3 = pi- - // The DP is defined in terms of m13Sq and m23Sq LauDaughters* daughters { new LauDaughters(parent, child_1, child_2, child_3, squareDP) }; - LauKinematics* kinematics { daughters->getKinematics() }; - Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; + Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}; - TString fileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; - TFile * file { TFile::Open(fileName,"recreate") }; + const TString outputFileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; + TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); tree->Branch("mPrime",&mPrime,"mPrime/D"); tree->Branch("thPrime",&thPrime,"thPrime/D"); - tree->Branch("jacobian",&jacobian,"jacobian/D"); - tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); - tree->Branch("randomNum",&randomNum,"randomNum/D"); - Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } - evtAccepted = kFALSE; - - while ( !evtAccepted ) { - kinematics->genFlatPhaseSpace( m13Sq, m23Sq ); - kinematics->updateKinematics( m13Sq, m23Sq ); - jacobian = kinematics->calcSqDPJacobian(); - invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; + kinematics->genFlatSqDP( mPrime, thPrime ); + kinematics->updateSqDPKinematics( mPrime, thPrime ); - randomNum = LauRandom::randomFun()->Rndm(); - if ( randomNum < invJacobian ) { - evtAccepted = kTRUE; - m12Sq = kinematics->getm12Sq(); - mPrime = kinematics->getmPrime(); - thPrime = kinematics->getThetaPrime(); + m13Sq = kinematics->getm13Sq(); + m23Sq = kinematics->getm23Sq(); + m12Sq = kinematics->getm12Sq(); - tree->Fill(); - } - } + tree->Fill(); } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; } diff --git a/examples/QuasiFlatSqDalitz-CustomMasses.cc b/examples/QuasiFlatSqDalitz-CustomMasses.cc index 017678f..f5d2667 100644 --- a/examples/QuasiFlatSqDalitz-CustomMasses.cc +++ b/examples/QuasiFlatSqDalitz-CustomMasses.cc @@ -1,114 +1,119 @@ /* Copyright 2018 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauKinematics.hh" #include "LauRandom.hh" /* * Generates toy MC according to EvtGen's FLATSQDALITZ model * - * Usage: QuasiFlatSqDalitz + * In LHCb simulation versions Sim08 and Sim09, and also in the global EvtGen + * version 2.0.0 the FLATSQDALITZ model does not generate perfectly uniformly + * in the square Dalitz plot. + * + * This issue was fixed in EvtGen 2.1.0 onwards, and hence in LHCb Sim10 since + * that uses a later tag of the global EvtGen repository. + * + * Usage: QuasiFlatSqDalitz-CustomMasses */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 7 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const Double_t parent_mass { std::stod( cmdLineArgs[1] ) }; const Double_t child_1_mass { std::stod( cmdLineArgs[2] ) }; const Double_t child_2_mass { std::stod( cmdLineArgs[3] ) }; const Double_t child_3_mass { std::stod( cmdLineArgs[4] ) }; const ULong_t nEvents { std::stoul( cmdLineArgs[5] ) }; const TString outputFileName { cmdLineArgs[6] }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE + // Do not change this! const Bool_t squareDP { kTRUE }; LauKinematics kinematics { child_1_mass, child_2_mass, child_3_mass, parent_mass, squareDP }; Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); tree->Branch("mPrime",&mPrime,"mPrime/D"); tree->Branch("thPrime",&thPrime,"thPrime/D"); tree->Branch("jacobian",&jacobian,"jacobian/D"); tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); tree->Branch("randomNum",&randomNum,"randomNum/D"); Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } evtAccepted = kFALSE; while ( !evtAccepted ) { kinematics.genFlatPhaseSpace( m13Sq, m23Sq ); kinematics.updateKinematics( m13Sq, m23Sq ); jacobian = kinematics.calcSqDPJacobian(); invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; randomNum = LauRandom::randomFun()->Rndm(); if ( randomNum < invJacobian ) { evtAccepted = kTRUE; m12Sq = kinematics.getm12Sq(); mPrime = kinematics.getmPrime(); thPrime = kinematics.getThetaPrime(); tree->Fill(); } } } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; } diff --git a/examples/QuasiFlatSqDalitz.cc b/examples/QuasiFlatSqDalitz.cc index 272039f..fe50222 100644 --- a/examples/QuasiFlatSqDalitz.cc +++ b/examples/QuasiFlatSqDalitz.cc @@ -1,122 +1,121 @@ /* Copyright 2018 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 */ #include #include #include #include #include "TFile.h" #include "TRandom.h" #include "TTree.h" #include "LauDaughters.hh" #include "LauKinematics.hh" #include "LauRandom.hh" /* * Generates toy MC according to EvtGen's FLATSQDALITZ model * + * In LHCb simulation versions Sim08 and Sim09, and also in the global EvtGen + * version 2.0.0 the FLATSQDALITZ model does not generate perfectly uniformly + * in the square Dalitz plot. + * + * This issue was fixed in EvtGen 2.1.0 onwards, and hence in LHCb Sim10 since + * that uses a later tag of the global EvtGen repository. + * * Usage: QuasiFlatSqDalitz */ int main( int argc, char** argv ) { // Convert the command-line arguments into a useful format const std::vector cmdLineArgs{ argv, argv+argc }; if ( cmdLineArgs.size() != 6 ) { std::cout << "Usage: " << cmdLineArgs[0] << " " << std::endl; return 0; } const TString& parent { cmdLineArgs[1] }; const TString& child_1 { cmdLineArgs[2] }; const TString& child_2 { cmdLineArgs[3] }; const TString& child_3 { cmdLineArgs[4] }; const ULong_t nEvents { std::stoul( cmdLineArgs[5].Data() ) }; - // If you want to use square DP histograms for efficiency, - // backgrounds or you just want the square DP co-ordinates - // stored in the toy MC ntuple then set this to kTRUE - Bool_t squareDP { kTRUE }; + // Do not change this! + const Bool_t squareDP { kTRUE }; - // This defines the DP => decay is B+ -> pi+ pi+ pi- - // Particle 1 = pi+ - // Particle 2 = pi+ - // Particle 3 = pi- - // The DP is defined in terms of m13Sq and m23Sq LauDaughters* daughters { new LauDaughters(parent, child_1, child_2, child_3, squareDP) }; - LauKinematics* kinematics { daughters->getKinematics() }; Double_t m12Sq{0.0}, m13Sq{0.0}, m23Sq{0.0}, mPrime{0.0}, thPrime{0.0}, jacobian{0.0}, invJacobian{0.0}, randomNum{0.0}; - TString fileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; - TFile * file { TFile::Open(fileName,"recreate") }; + const TString outputFileName { TString::Format( "%s_%s_%s_%s.root", daughters->getSanitisedNameParent().Data(), daughters->getSanitisedNameDaug1().Data(), daughters->getSanitisedNameDaug2().Data(), daughters->getSanitisedNameDaug3().Data() ) }; + TFile * file { TFile::Open(outputFileName,"recreate") }; TTree * tree { new TTree("genResults","genResults") }; tree->Branch("m12Sq",&m12Sq,"m12Sq/D"); tree->Branch("m13Sq",&m13Sq,"m13Sq/D"); tree->Branch("m23Sq",&m23Sq,"m23Sq/D"); tree->Branch("mPrime",&mPrime,"mPrime/D"); tree->Branch("thPrime",&thPrime,"thPrime/D"); tree->Branch("jacobian",&jacobian,"jacobian/D"); tree->Branch("invJacobian",&invJacobian,"invJacobian/D"); tree->Branch("randomNum",&randomNum,"randomNum/D"); Bool_t evtAccepted{kFALSE}; for ( ULong_t i{0}; i < nEvents; ++i ) { if ( nEvents > 1000 && i%(nEvents/1000)==0 ) { std::cout << "Generating event " << i << std::endl; } evtAccepted = kFALSE; while ( !evtAccepted ) { kinematics->genFlatPhaseSpace( m13Sq, m23Sq ); kinematics->updateKinematics( m13Sq, m23Sq ); jacobian = kinematics->calcSqDPJacobian(); invJacobian = (jacobian<1.0) ? 1.0 : 1.0/jacobian; randomNum = LauRandom::randomFun()->Rndm(); if ( randomNum < invJacobian ) { evtAccepted = kTRUE; m12Sq = kinematics->getm12Sq(); mPrime = kinematics->getmPrime(); thPrime = kinematics->getThetaPrime(); tree->Fill(); } } } file->Write(); file->Close(); delete file; return EXIT_SUCCESS; }