diff --git a/examples/Test_Dpipi.cc b/examples/Test_Dpipi.cc index 3b55d45..369e120 100644 --- a/examples/Test_Dpipi.cc +++ b/examples/Test_Dpipi.cc @@ -1,341 +1,382 @@ #include using std::cout; using std::cerr; using std::endl; #include #include +#include #include "TFile.h" #include "TH2.h" #include "TRandom.h" #include "TString.h" #include "TSystem.h" #include "LauDaughters.hh" #include "LauDecayTimePdf.hh" #include "LauEffModel.hh" #include "LauIsobarDynamics.hh" #include "LauMagPhaseCoeffSet.hh" #include "LauRandom.hh" #include "LauRealImagCoeffSet.hh" #include "LauTimeDepFitModel.hh" #include "LauVetoes.hh" #include "LauFlavTag.hh" #include "Lau1DHistPdf.hh" #include "Lau1DCubicSpline.hh" void usage(const TString& progName) { cerr<<"Usage:"< vector_argv{argv, argc+argv}; + + enum efficiencyInputMethod {Spline, Hist}; + efficiencyInputMethod effiMethod = efficiencyInputMethod::Spline; + + auto histTest = std::find(std::begin(vector_argv), std::end(vector_argv), "--hist"); + if(histTest != std::end(vector_argv)) + { + effiMethod = efficiencyInputMethod::Spline; + vector_argv.erase(histTest); + } + + const TString dtype = vector_argv[2]; Int_t port = 0; Int_t iFit = 0; Int_t firstExpt = 0; Int_t firstExptGen = 0; Int_t nExpt = 1; Int_t nExptGen = 1; LauTimeDepFitModel::CPEigenvalue eigenvalue = LauTimeDepFitModel::CPEven; Bool_t fixPhiMix(kFALSE); Bool_t useSinCos(kTRUE); // check the command line arguments if (argc<1) { - usage(argv[0]); + usage(vector_argv[0]); return EXIT_FAILURE; } - TString command = argv[1]; + TString command = vector_argv[1]; if (command != "gen" && command != "fit") { - usage(argv[0]); + usage(vector_argv[0]); return EXIT_FAILURE; } if (command == "fit") { if (argc>3) { - port = atoi(argv[3]); + port = atoi(vector_argv[3]); if (argc>4) { - iFit = atoi(argv[4]); + iFit = atoi(vector_argv[4]); if (argc>5) { - firstExpt = atoi(argv[5]); + firstExpt = atoi(vector_argv[5]); if (argc>6) { - nExpt = atoi(argv[6]); + nExpt = atoi(vector_argv[6]); if (argc>7) { - nExptGen = atoi(argv[7]); + nExptGen = atoi(vector_argv[7]); } } } } } for (firstExptGen = 0; firstExptGen<(firstExpt+nExpt); firstExptGen+=nExptGen) { } firstExptGen -= nExptGen; if ( (nExpt > nExptGen) || (nExptGen%nExpt != 0) ) { cerr<<"Error, nExpt must be a factor of nExptGen."<3) { - firstExptGen = atoi(argv[3]); + firstExptGen = atoi(vector_argv[3]); if (argc>4) { - nExptGen = atoi(argv[4]); + nExptGen = atoi(vector_argv[4]); if (argc>5) { - Int_t eigval = atoi(argv[5]); + Int_t eigval = atoi(vector_argv[5]); if ( eigval == 1 ) { eigenvalue = LauTimeDepFitModel::CPOdd; } else { eigenvalue = LauTimeDepFitModel::CPEven; } } } } } cout<<"dtype "<addMassVeto( 2, 2.00776, 2.01276 ); daughtersB0bar = new LauDaughters("B0_bar", "pi+", "pi-", "D0"); daughtersB0 = new LauDaughters("B0", "pi+", "pi-", "D0_bar"); // efficiency effModelB0bar = new LauEffModel(daughtersB0bar, vetoes); effModelB0 = new LauEffModel(daughtersB0, vetoes); LauIsobarDynamics* sigModelB0bar(0); LauIsobarDynamics* sigModelB0(0); LauTimeDepFitModel* fitModel(0); std::vector params; LauFlavTag* flavTag = new LauFlavTag(); flavTag->setTrueTagVarName("trueTag"); // Use alternative tagging calibration parameters? //flavTag->useAveOmegaDeltaOmega(); TFile* etafile(0); TH1* etahist(0); Lau1DHistPdf* etahistpdf(0); - //if (command == "gen"){ + //if (command == "gen"){ etafile = TFile::Open("histogram.root"); etahist = dynamic_cast(etafile->Get("htemp")); etahistpdf = new Lau1DHistPdf("eta",etahist,0.0,0.54,kTRUE,kFALSE); //} // signal dynamics sigModelB0bar = new LauIsobarDynamics(daughtersB0bar, effModelB0bar); sigModelB0bar->setIntFileName("integ_B0bar.dat"); sigModelB0bar->addResonance("D*+_2", 2, LauAbsResonance::RelBW); sigModelB0bar->addResonance("D*+_0", 2, LauAbsResonance::RelBW); sigModelB0bar->addResonance("rho0(770)", 3, LauAbsResonance::RelBW); sigModelB0bar->addResonance("f_0(980)", 3, LauAbsResonance::RelBW); sigModelB0bar->addResonance("f_2(1270)", 3, LauAbsResonance::RelBW); sigModelB0 = new LauIsobarDynamics(daughtersB0, effModelB0); sigModelB0->setIntFileName("integ_B0.dat"); sigModelB0->addResonance("D*-_2", 1, LauAbsResonance::RelBW); sigModelB0->addResonance("D*-_0", 1, LauAbsResonance::RelBW); sigModelB0->addResonance("rho0(770)", 3, LauAbsResonance::RelBW); sigModelB0->addResonance("f_0(980)", 3, LauAbsResonance::RelBW); sigModelB0->addResonance("f_2(1270)", 3, LauAbsResonance::RelBW); // fit model fitModel = new LauTimeDepFitModel(sigModelB0bar,sigModelB0,flavTag); fitModel->setASqMaxValue(1.45); std::vector coeffset; coeffset.push_back( new LauRealImagCoeffSet("D*-_2", 1.00, 0.00, kTRUE, kTRUE) ); coeffset.push_back( new LauRealImagCoeffSet("D*-_0", 0.53*TMath::Cos( 3.00), 0.53*TMath::Sin( 3.00), kFALSE, kFALSE) ); coeffset.push_back( new LauRealImagCoeffSet("rho0(770)", 1.22*TMath::Cos( 2.25), 1.22*TMath::Sin( 2.25), kFALSE, kFALSE) ); coeffset.push_back( new LauRealImagCoeffSet("f_0(980)", 0.19*TMath::Cos(-2.48), 0.19*TMath::Sin(-2.48), kFALSE, kFALSE) ); coeffset.push_back( new LauRealImagCoeffSet("f_2(1270)", 0.75*TMath::Cos( 2.97), 0.75*TMath::Sin( 2.97), kFALSE, kFALSE) ); for (std::vector::iterator iter=coeffset.begin(); iter!=coeffset.end(); ++iter) { fitModel->setAmpCoeffSet(*iter); } fitModel->setCPEigenvalue( eigenvalue ); fitModel->setPhiMix( 2.0*LauConstants::beta, fixPhiMix, useSinCos ); fitModel->setAsymmetries(0.0,kTRUE); // Tag cat fractions, dilutions and Delta dilutions flavTag->addTagger("OSTagger","tagVal_OS","mistagVal_OS", etahistpdf,0.20,0.5,1.0,0.20,0.5,1.0); flavTag->addTagger("SSTagger","tagVal_SS","mistagVal_SS", etahistpdf,0.30,0.5,1.0,0.30,0.5,1.0); // Delta t PDFs const Double_t minDt(0.0); const Double_t maxDt(20.0); const Double_t minDtErr(0.0); const Double_t maxDtErr(2.5); const Int_t nGauss(3); std::vector scale(nGauss); scale[0] = kTRUE; scale[1] = kTRUE; scale[2] = kFALSE; std::vector dtPars(10); TString mean0Name("dt_"); mean0Name += "_mean_0"; TString mean1Name("dt_"); mean1Name += "_mean_1"; TString mean2Name("dt_"); mean2Name += "_mean_2"; TString sigma0Name("dt_"); sigma0Name += "_sigma_0"; TString sigma1Name("dt_"); sigma1Name += "_sigma_1"; TString sigma2Name("dt_"); sigma2Name += "_sigma_2"; TString frac1Name("dt_"); frac1Name += "_frac_1"; TString frac2Name("dt_"); frac2Name += "_frac_2"; TString tauName("dt_"); tauName += "_tau"; TString freqName("dt_"); freqName += "_deltaM"; LauParameter * mean1 = new LauParameter(mean1Name, -1.27); LauParameter * mean2 = new LauParameter(mean2Name, 0.0); LauParameter * sigma1 = new LauParameter(sigma1Name, 3.0); LauParameter * sigma2 = new LauParameter(sigma2Name, 8.0); LauParameter * frac1 = new LauParameter(frac1Name, 0.0930); LauParameter * frac2 = new LauParameter(frac2Name, 0.0036); LauParameter * tau = new LauParameter(tauName, 1.520); LauParameter * freq = new LauParameter(freqName, 0.5064); TString mean0tagcat63Name("dt_"); mean0tagcat63Name += 63; mean0tagcat63Name += "_mean_0"; TString sigma0tagcat63Name("dt_"); sigma0tagcat63Name += 63; sigma0tagcat63Name += "_sigma_0"; LauParameter * mean0tagcat63 = new LauParameter(mean0tagcat63Name, -0.031); LauParameter * sigma0tagcat63 = new LauParameter(sigma0tagcat63Name, 0.972); //Decay time acceptance spline - same for all tag cats (though doesn't have to be) std::vector dtvals; dtvals.push_back(0.0); dtvals.push_back(1.0); dtvals.push_back(2.0); dtvals.push_back(4.0); dtvals.push_back(6.0); dtvals.push_back(8.0); dtvals.push_back(11.0); dtvals.push_back(14.0); dtvals.push_back(17.0); dtvals.push_back(20.0); std::vector effvals; effvals.push_back(0.0); effvals.push_back(0.15); effvals.push_back(0.25); effvals.push_back(0.33); effvals.push_back(0.38); effvals.push_back(0.4); effvals.push_back(0.43); effvals.push_back(0.45); effvals.push_back(0.47); effvals.push_back(0.50); //Lau1DCubicSpline* dtEffSpline = new Lau1DCubicSpline(dtvals,effvals); - Lau1DCubicSpline* dtEffSpline = new Lau1DCubicSpline(dtvals,effvals,Lau1DCubicSpline::StandardSpline,Lau1DCubicSpline::Natural,Lau1DCubicSpline::Natural); + Lau1DCubicSpline* dtEffSpline = nullptr; + TH1D* dtEffHist = nullptr; + + switch(effiMethod) + { + case efficiencyInputMethod::Spline: + dtEffSpline = new Lau1DCubicSpline(dtvals,effvals,Lau1DCubicSpline::StandardSpline,Lau1DCubicSpline::Natural,Lau1DCubicSpline::Natural); + break; + + case efficiencyInputMethod::Hist: + const Int_t nBins = 9; + const Double_t edges[nBins + 1] = {0,1,2,4,6,8,11,14,17,20}; + const Double_t binFilling[nBins] = {0.0075,0.02,0.029,0.655,0.69,0.715,0.74,0.76,0.785}; + dtEffHist = new TH1D("dtEffHist","Histogram of efficiencies", nBins, edges); + for( Int_t i = 1 ; i <= nBins ; ++i ){ dtEffHist->SetBinContent(i, binFilling[i-1]); } + break; + } dtPars[0] = mean0tagcat63; dtPars[1] = mean1->createClone(); dtPars[2] = mean2->createClone(); dtPars[3] = sigma0tagcat63; dtPars[4] = sigma1->createClone(); dtPars[5] = sigma2->createClone(); dtPars[6] = frac1->createClone(); dtPars[7] = frac2->createClone(); dtPars[8] = tau->createClone(); dtPars[9] = freq->createClone(); if (dtype=="CPEven"){ LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::SimFitSigBd, nGauss, scale, LauDecayTimePdf::DecayTime ); dtPdf->doSmearing(kFALSE); dtPdf->setEffiSpline(dtEffSpline); + dtPdf->setEffiHist(dtEffHist); fitModel->setSignalDtPdf( dtPdf ); } else { LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::SimFitNormBd, nGauss, scale, LauDecayTimePdf::DecayTime ); dtPdf->doSmearing(kFALSE); dtPdf->setEffiSpline(dtEffSpline); + dtPdf->setEffiHist(dtEffHist); fitModel->setSignalDtPdf( dtPdf ); } // set the number of signal events cout<<"nSigEvents = "<setNSigEvents(nSigPar); // set the number of experiments if (command == "fit") { fitModel->setNExpts(nExpt, firstExpt); } else { fitModel->setNExpts(nExptGen, firstExptGen); } fitModel->useAsymmFitErrors(kFALSE); //fitModel->useRandomInitFitPars(kTRUE); fitModel->useRandomInitFitPars(kFALSE); fitModel->doPoissonSmearing(kFALSE); fitModel->doEMLFit(kFALSE); fitModel->writeLatexTable(kFALSE); TString dataFile(""); TString treeName("fitTree"); TString rootFileName(""); TString tableFileName(""); TString fitToyFileName(""); TString splotFileName(""); if (command == "fit") { dataFile = "TEST-Dpipi_"+dtype; dataFile += "_expts"; dataFile += firstExptGen; dataFile += "-"; dataFile += firstExptGen+nExptGen-1; - dataFile += "_CP"; + dataFile += "_CP"; if ( eigenvalue == LauTimeDepFitModel::CPEven ) { dataFile += "even"; } else { dataFile += "odd"; } dataFile += ".root"; rootFileName = "fits/fit"; rootFileName += iFit; rootFileName += "_expts"; rootFileName += firstExpt; rootFileName += "-"; rootFileName += firstExpt+nExpt-1; rootFileName += ".root"; tableFileName = "fitResults_"; tableFileName += iFit; tableFileName += "_expts"; tableFileName += firstExpt; tableFileName += "-"; tableFileName += firstExpt+nExpt-1; fitToyFileName = "fitToyMC_"+dtype; fitToyFileName += iFit; fitToyFileName += "_expts"; fitToyFileName += firstExpt; fitToyFileName += "-"; fitToyFileName += firstExpt+nExpt-1; fitToyFileName += ".root"; splotFileName = "splot_"; splotFileName += iFit; splotFileName += "_expts"; splotFileName += firstExpt; splotFileName += "-"; splotFileName += firstExpt+nExpt-1; splotFileName += ".root"; } else { - dataFile = "TEST-Dpipi_"+dtype; - dataFile += "_expts"; dataFile += firstExptGen; dataFile += "-"; dataFile += firstExptGen+nExptGen-1; dataFile += "_CP"; + dataFile = "TEST-Dpipi"; + switch(effiMethod) + { + case efficiencyInputMethod::Spline: + dataFile += "_Spline"; + break; + case efficiencyInputMethod::Hist: + dataFile += "_Hist"; + break; + } + dataFile += "_"+dtype+"_expts"; dataFile += firstExptGen; dataFile += "-"; dataFile += firstExptGen+nExptGen-1; dataFile += "_CP"; if ( eigenvalue == LauTimeDepFitModel::CPEven ) { dataFile += "even"; } else { dataFile += "odd"; } dataFile += ".root"; rootFileName = "dummy.root"; tableFileName = "genResults"; } // Generate toy from the fitted parameters fitModel->compareFitData(1, fitToyFileName); // Write out per-event likelihoods and sWeights //fitModel->writeSPlotData(splotFileName, "splot", kFALSE); // Execute the generation/fit if ( command == "fit" ){ fitModel->runSlave( dataFile, treeName, rootFileName, tableFileName, "localhost", port ); } else { fitModel->run( command, dataFile, treeName, rootFileName, tableFileName ); } return EXIT_SUCCESS; } diff --git a/examples/Test_Dpipi_efficiencyHist.cc b/examples/Test_Dpipi_efficiencyHist.cc deleted file mode 100755 index 170184f..0000000 --- a/examples/Test_Dpipi_efficiencyHist.cc +++ /dev/null @@ -1,403 +0,0 @@ - -#include -using std::cout; -using std::cerr; -using std::endl; - -#include -#include - -#include "TFile.h" -#include "TH2.h" -#include "TH1D.h" -#include "TRandom.h" -#include "TString.h" -#include "TSystem.h" - -#include "LauDaughters.hh" -#include "LauDecayTimePdf.hh" -#include "LauEffModel.hh" -#include "LauIsobarDynamics.hh" -#include "LauMagPhaseCoeffSet.hh" -#include "LauRandom.hh" -#include "LauRealImagCoeffSet.hh" -#include "LauTimeDepFitModel.hh" -#include "LauVetoes.hh" -#include "LauFlavTag.hh" -#include "Lau1DHistPdf.hh" -#include "Lau1DCubicSpline.hh" - -void usage(const TString& progName) -{ - cerr<<"Usage:"<3) { - port = atoi(argv[3]); - if (argc>4) { - iFit = atoi(argv[4]); - if (argc>5) { - firstExpt = atoi(argv[5]); - if (argc>6) { - nExpt = atoi(argv[6]); - if (argc>7) { - nExptGen = atoi(argv[7]); - } - } - } - } - } - - for (firstExptGen = 0; firstExptGen<(firstExpt+nExpt); firstExptGen+=nExptGen) { - } - firstExptGen -= nExptGen; - - if ( (nExpt > nExptGen) || (nExptGen%nExpt != 0) ) { - cerr<<"Error, nExpt must be a factor of nExptGen."<3) { - firstExptGen = atoi(argv[3]); - if (argc>4) { - nExptGen = atoi(argv[4]); - if (argc>5) { - Int_t eigval = atoi(argv[5]); - if ( eigval == 1 ) { - eigenvalue = LauTimeDepFitModel::CPOdd; - } else { - eigenvalue = LauTimeDepFitModel::CPEven; - } - } - } - } - } - - cout<<"dtype "<addMassVeto( 2, 2.00776, 2.01276 ); - - daughtersB0bar = new LauDaughters("B0_bar", "pi+", "pi-", "D0"); - daughtersB0 = new LauDaughters("B0", "pi+", "pi-", "D0_bar"); - - // efficiency - effModelB0bar = new LauEffModel(daughtersB0bar, vetoes); - effModelB0 = new LauEffModel(daughtersB0, vetoes); - - LauIsobarDynamics* sigModelB0bar(0); - LauIsobarDynamics* sigModelB0(0); - LauTimeDepFitModel* fitModel(0); - - std::vector params; - - //LauFlavTag* flavTag = new LauFlavTag(params, kTRUE, "tagFlv", "tagCat"); - LauFlavTag* flavTag = new LauFlavTag(); - flavTag->setTrueTagVarName("trueTag"); - - // Use alternative tagging calibration parameters? - //flavTag->useAveOmegaDeltaOmega(); - - TFile* etafile(0); - TH1* etahist(0); - Lau1DHistPdf* etahistpdf(0); - //if (command == "gen"){ - etafile = TFile::Open("histogram.root"); - etahist = dynamic_cast(etafile->Get("htemp")); - etahistpdf = new Lau1DHistPdf("eta",etahist,0,0.54,kFALSE,kFALSE); - //} - - // signal dynamics - sigModelB0bar = new LauIsobarDynamics(daughtersB0bar, effModelB0bar); - sigModelB0bar->setIntFileName("integ_B0bar.dat"); - sigModelB0bar->addResonance("D*+_2", 2, LauAbsResonance::RelBW); - sigModelB0bar->addResonance("D*+_0", 2, LauAbsResonance::RelBW); - sigModelB0bar->addResonance("rho0(770)", 3, LauAbsResonance::RelBW); - sigModelB0bar->addResonance("f_0(980)", 3, LauAbsResonance::RelBW); - sigModelB0bar->addResonance("f_2(1270)", 3, LauAbsResonance::RelBW); - - sigModelB0 = new LauIsobarDynamics(daughtersB0, effModelB0); - sigModelB0->setIntFileName("integ_B0.dat"); - sigModelB0->addResonance("D*-_2", 1, LauAbsResonance::RelBW); - sigModelB0->addResonance("D*-_0", 1, LauAbsResonance::RelBW); - sigModelB0->addResonance("rho0(770)", 3, LauAbsResonance::RelBW); - sigModelB0->addResonance("f_0(980)", 3, LauAbsResonance::RelBW); - sigModelB0->addResonance("f_2(1270)", 3, LauAbsResonance::RelBW); - - // fit model - fitModel = new LauTimeDepFitModel(sigModelB0bar,sigModelB0,flavTag); - fitModel->setASqMaxValue(1.45); - - std::vector coeffset; - coeffset.push_back( new LauRealImagCoeffSet("D*-_2", 1.00, 0.00, kTRUE, kTRUE) ); - coeffset.push_back( new LauRealImagCoeffSet("D*-_0", 0.53*TMath::Cos( 3.00), 0.53*TMath::Sin( 3.00), kFALSE, kFALSE) ); - coeffset.push_back( new LauRealImagCoeffSet("rho0(770)", 1.22*TMath::Cos( 2.25), 1.22*TMath::Sin( 2.25), kFALSE, kFALSE) ); - coeffset.push_back( new LauRealImagCoeffSet("f_0(980)", 0.19*TMath::Cos(-2.48), 0.19*TMath::Sin(-2.48), kFALSE, kFALSE) ); - coeffset.push_back( new LauRealImagCoeffSet("f_2(1270)", 0.75*TMath::Cos( 2.97), 0.75*TMath::Sin( 2.97), kFALSE, kFALSE) ); - - for (std::vector::iterator iter=coeffset.begin(); iter!=coeffset.end(); ++iter) { - fitModel->setAmpCoeffSet(*iter); - } - - fitModel->setCPEigenvalue( eigenvalue ); - fitModel->setPhiMix( 2.0*LauConstants::beta, fixPhiMix, useSinCos ); - - fitModel->setAsymmetries(0.0,kFALSE); - - int tagCat(-1); - // Tag cat fractions, dilutions and Delta dilutions -// if (dtype=="CPEven"){ -// flavTag->addValidTagCat(63); -// flavTag->setSignalTagCatPars(63, 0.6, 1.0, 0.0, kTRUE, kTRUE,"PerEvtMistag",0.5,0.5,1.0,1.0,0.0); -// tagCat=63; -// } else { -// // Still need tagCat 63 in QFS channels for calibration etc -// flavTag->addValidTagCat(63); -// flavTag->setSignalTagCatPars(63, 0.6, 1.0, 0.0, kTRUE, kTRUE,"PerEvtMistag",0.5,0.5,1.0,1.0,0.0); -// tagCat=63; -// } - - // Tag cat fractions, dilutions and Delta dilutions - flavTag->addTagger("OSTagger","tagVal_OS","mistagVal_OS", etahistpdf,0.20,0.5,1.0,0.20,0.5,1.0); - flavTag->addTagger("SSTagger","tagVal_SS","mistagVal_SS", etahistpdf,0.30,0.5,1.0,0.30,0.5,1.0); - - //if (command == "gen"){ - // fitModel->setSignalFlavTagPdfs(63,etahistpdf); - //} - - - // Delta t PDFs - - const Double_t minDt(0.0); - const Double_t maxDt(20.0); - const Double_t minDtErr(0.0); - const Double_t maxDtErr(2.5); - const Int_t nGauss(3); - std::vector scale(nGauss); - scale[0] = kTRUE; - scale[1] = kTRUE; - scale[2] = kFALSE; - - std::vector dtPars(10); - std::vector dtParsExp(9); - TString mean0Name("dt_"); mean0Name += "_mean_0"; - TString mean1Name("dt_"); mean1Name += "_mean_1"; - TString mean2Name("dt_"); mean2Name += "_mean_2"; - TString sigma0Name("dt_"); sigma0Name += "_sigma_0"; - TString sigma1Name("dt_"); sigma1Name += "_sigma_1"; - TString sigma2Name("dt_"); sigma2Name += "_sigma_2"; - TString frac1Name("dt_"); frac1Name += "_frac_1"; - TString frac2Name("dt_"); frac2Name += "_frac_2"; - TString tauName("dt_"); tauName += "_tau"; - TString freqName("dt_"); freqName += "_deltaM"; - - LauParameter * mean0 = new LauParameter(mean0Name, -0.181); - LauParameter * mean1 = new LauParameter(mean1Name, -1.27); - LauParameter * mean2 = new LauParameter(mean2Name, 0.0); - LauParameter * sigma0 = new LauParameter(sigma0Name, 1.067); - LauParameter * sigma1 = new LauParameter(sigma1Name, 3.0); - LauParameter * sigma2 = new LauParameter(sigma2Name, 8.0); - LauParameter * frac1 = new LauParameter(frac1Name, 0.0930); - LauParameter * frac2 = new LauParameter(frac2Name, 0.0036); - LauParameter * tau = new LauParameter(tauName, 1.520); - LauParameter * freq = new LauParameter(freqName, 0.5064); - - TString mean0tagcat63Name("dt_"); mean0tagcat63Name += 63; mean0tagcat63Name += "_mean_0"; - TString sigma0tagcat63Name("dt_"); sigma0tagcat63Name += 63; sigma0tagcat63Name += "_sigma_0"; - - LauParameter * mean0tagcat63 = new LauParameter(mean0tagcat63Name, -0.031); - LauParameter * sigma0tagcat63 = new LauParameter(sigma0tagcat63Name, 0.972); - - dtPars[0] = mean0; dtParsExp[0] = mean0; - dtPars[1] = mean1; dtParsExp[1] = mean1; - dtPars[2] = mean2; dtParsExp[2] = mean2; - dtPars[3] = sigma0; dtParsExp[3] = sigma0; - dtPars[4] = sigma1; dtParsExp[4] = sigma1; - dtPars[5] = sigma2; dtParsExp[5] = sigma2; - dtPars[6] = frac1; dtParsExp[6] = frac1; - dtPars[7] = frac2; dtParsExp[7] = frac2; - dtPars[8] = tau; dtParsExp[8] = tau; - dtPars[9] = freq; - - //Decay time acceptance spline - same for all tag cats (though doesn't have to be) - std::vector dtvals; - dtvals.push_back(0.0); dtvals.push_back(1.0); dtvals.push_back(2.0); dtvals.push_back(4.0); dtvals.push_back(6.0); - dtvals.push_back(8.0); dtvals.push_back(11.0); dtvals.push_back(14.0); dtvals.push_back(17.0); dtvals.push_back(20.0); - std::vector effvals; - effvals.push_back(0.0); effvals.push_back(0.15); effvals.push_back(0.25); effvals.push_back(0.33); effvals.push_back(0.38); - effvals.push_back(0.4); effvals.push_back(0.43); effvals.push_back(0.45); effvals.push_back(0.47); effvals.push_back(0.50); - - //Lau1DCubicSpline* dtEffSpline = new Lau1DCubicSpline(dtvals,effvals); -// Lau1DCubicSpline* dtEffSpline = new Lau1DCubicSpline(dtvals,effvals,Lau1DCubicSpline::StandardSpline,Lau1DCubicSpline::Natural,Lau1DCubicSpline::Natural); - - const Int_t nBins = 9; - const Double_t edges[nBins + 1] = {0,1,2,4,6,8,11,14,17,20}; - const Double_t binFilling[nBins] = {0.0075,0.02,0.029,0.655,0.69,0.715,0.74,0.76,0.785}; - TH1D effHist("effHist","Histogram of efficiencies", nBins, edges); - - for( Int_t i = 1 ; i <= nBins ; ++i ){ effHist.SetBinContent(i, binFilling[i-1]); } - - if (dtype=="CPEven"){ - LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::ExpTrig, nGauss, scale, LauDecayTimePdf::DecayTime, LauDecayTimePdf::Flat ); - dtPdf->doSmearing(kFALSE); - dtPdf->setEffiHist(&effHist); -// dtPdf->setEffiSpline(dtEffSpline, effPars0); - fitModel->setSignalDtPdf( dtPdf ); - } else { - LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::ExpTrig, nGauss, scale, LauDecayTimePdf::DecayTime, LauDecayTimePdf::Flat ); - dtPdf->doSmearing(kFALSE); - dtPdf->setEffiHist(&effHist); -// dtPdf->setEffiSpline(dtEffSpline, effPars0); - fitModel->setSignalDtPdf( dtPdf ); - } - - //if (tagCat==63){ - dtPars[0] = mean0tagcat63; - dtPars[1] = mean1->createClone(); - dtPars[2] = mean2->createClone(); - dtPars[3] = sigma0tagcat63; - dtPars[4] = sigma1->createClone(); - dtPars[5] = sigma2->createClone(); - dtPars[6] = frac1->createClone(); - dtPars[7] = frac2->createClone(); - dtPars[8] = tau->createClone(); - dtPars[9] = freq->createClone(); - - //LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::ExpTrig, nGauss, scale, LauDecayTimePdf::DecayTime ); - if (dtype=="CPEven"){ - LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::SimFitSigBd, nGauss, scale, LauDecayTimePdf::DecayTime , LauDecayTimePdf::Flat ); - dtPdf->doSmearing(kFALSE); -// dtPdf->setEffiSpline(dtEffSpline, effPars63); - dtPdf->setEffiHist(&effHist); - fitModel->setSignalDtPdf( dtPdf ); - } else { - LauDecayTimePdf * dtPdf = new LauDecayTimePdf( "deltaTAvg", "deltaTAvgErr", dtPars, minDt, maxDt, minDtErr, maxDtErr, LauDecayTimePdf::SimFitNormBd, nGauss, scale, LauDecayTimePdf::DecayTime , LauDecayTimePdf::Flat); - dtPdf->doSmearing(kFALSE); -// dtPdf->setEffiSpline(dtEffSpline, effPars63); - dtPdf->setEffiHist(&effHist); - fitModel->setSignalDtPdf( dtPdf ); - } - - //} - - // set the number of signal events - cout<<"nSigEvents = "<setNSigEvents(nSigPar); - - // set the number of experiments - if (command == "fit") { - fitModel->setNExpts(nExpt, firstExpt); - } else { - fitModel->setNExpts(nExptGen, firstExptGen); - } - - fitModel->useAsymmFitErrors(kFALSE); - //fitModel->useRandomInitFitPars(kTRUE); - fitModel->useRandomInitFitPars(kFALSE); - fitModel->doPoissonSmearing(kFALSE); - fitModel->doEMLFit(kFALSE); - fitModel->writeLatexTable(kFALSE); - - - TString dataFile(""); - TString treeName("fitTree"); - TString rootFileName(""); - TString tableFileName(""); - TString fitToyFileName(""); - TString splotFileName(""); - if (command == "fit") { - dataFile = "TEST-Dpipi_"+dtype; - dataFile += "_expts"; dataFile += firstExptGen; dataFile += "-"; dataFile += firstExptGen+nExptGen-1; - dataFile += "_CP"; - if ( eigenvalue == LauTimeDepFitModel::CPEven ) { - dataFile += "even"; - } else { - dataFile += "odd"; - } - dataFile += ".root"; - - rootFileName = "fits/fit"; rootFileName += iFit; - rootFileName += "_expts"; rootFileName += firstExpt; rootFileName += "-"; rootFileName += firstExpt+nExpt-1; - rootFileName += ".root"; - - tableFileName = "fitResults_"; tableFileName += iFit; - tableFileName += "_expts"; tableFileName += firstExpt; tableFileName += "-"; tableFileName += firstExpt+nExpt-1; - - fitToyFileName = "fitToyMC_"+dtype; fitToyFileName += iFit; - fitToyFileName += "_expts"; fitToyFileName += firstExpt; fitToyFileName += "-"; fitToyFileName += firstExpt+nExpt-1; - fitToyFileName += ".root"; - - splotFileName = "splot_"; splotFileName += iFit; - splotFileName += "_expts"; splotFileName += firstExpt; splotFileName += "-"; splotFileName += firstExpt+nExpt-1; - splotFileName += ".root"; - } else { - dataFile = "TEST-Dpipi_"+dtype; - dataFile += "_expts"; dataFile += firstExptGen; dataFile += "-"; dataFile += firstExptGen+nExptGen-1; dataFile += "_CP"; - if ( eigenvalue == LauTimeDepFitModel::CPEven ) { - dataFile += "even"; - } else { - dataFile += "odd"; - } - dataFile += ".root"; - - rootFileName = "dummy.root"; - tableFileName = "genResults"; - } - - // Generate toy from the fitted parameters - fitModel->compareFitData(1, fitToyFileName); - - // Write out per-event likelihoods and sWeights - //fitModel->writeSPlotData(splotFileName, "splot", kFALSE); - - // Execute the generation/fit - if ( command == "fit" ){ - fitModel->runSlave( dataFile, treeName, rootFileName, tableFileName, "localhost", port ); - } else { - fitModel->run( command, dataFile, treeName, rootFileName, tableFileName ); - } - - return EXIT_SUCCESS; -} diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 237c787..86cd4e6 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -1,54 +1,54 @@ # Set the include directory (unfortunately this old-style stuff is necessary for the ROOT_GENERATE_DICTIONARY macro) include_directories(${PROJECT_SOURCE_DIR}/inc) # Use glob to find the headers and sources file(GLOB LAURA_HEADERS ${PROJECT_SOURCE_DIR}/inc/*.hh) file(GLOB LAURA_SOURCES ${CMAKE_CURRENT_SOURCE_DIR}/*.cc) if (NOT LAURA_BUILD_ROOFIT_SLAVE) list(REMOVE_ITEM LAURA_HEADERS ${PROJECT_SOURCE_DIR}/inc/LauRooFitSlave.hh) list(REMOVE_ITEM LAURA_SOURCES ${CMAKE_CURRENT_SOURCE_DIR}/LauRooFitSlave.cc) endif() # Generate the rootcint file set(LAURA_LINKDEF ${PROJECT_SOURCE_DIR}/inc/Laura++_LinkDef.h) set(LAURA_DICTIONARY_ROOT G__Laura++) set(LAURA_DICTIONARY ${LAURA_DICTIONARY_ROOT}.cxx) if (LAURA_BUILD_ROOFIT_SLAVE) ROOT_GENERATE_DICTIONARY( ${LAURA_DICTIONARY_ROOT} ${LAURA_HEADERS} LINKDEF ${LAURA_LINKDEF} OPTIONS -DDOLAUROOFITSLAVE ) else() ROOT_GENERATE_DICTIONARY( ${LAURA_DICTIONARY_ROOT} ${LAURA_HEADERS} LINKDEF ${LAURA_LINKDEF} ) endif() # Build the shared library add_library(Laura++ SHARED ${LAURA_SOURCES} ${LAURA_DICTIONARY}) set_target_properties(Laura++ PROPERTIES OUTPUT_NAME Laura++) set_target_properties(Laura++ PROPERTIES VERSION ${CMAKE_PROJECT_VERSION} SOVERSION ${CMAKE_PROJECT_VERSION_MAJOR}) set_target_properties(Laura++ PROPERTIES LIBRARY_OUTPUT_DIRECTORY ${PROJECT_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}) target_include_directories(Laura++ PUBLIC $ $) -target_link_libraries(Laura++ ROOT::Core ROOT::Hist ROOT::Matrix ROOT::Physics ROOT::Minuit ROOT::EG ROOT::Tree) +target_link_libraries(Laura++ ROOT::Core ROOT::Hist ROOT::Matrix ROOT::Physics ROOT::Minuit ROOT::EG ROOT::Tree ROOT::RooFitCore) if (LAURA_BUILD_ROOFIT_SLAVE) - target_link_libraries(Laura++ ROOT::RooFit ROOT::RooFitCore) + target_link_libraries(Laura++ ROOT::RooFit) endif() # Install the libraries install( TARGETS Laura++ EXPORT "LauraTargets" LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} ) # Install the pcm and rootmap files generated by ROOT_GENERATE_DICTIONARY install( FILES ${CMAKE_CURRENT_BINARY_DIR}/libLaura++.rootmap ${CMAKE_CURRENT_BINARY_DIR}/libLaura++_rdict.pcm DESTINATION ${CMAKE_INSTALL_LIBDIR} )