82#define NoHistConst_Low "0"
83#define NoHistConst_High "2000"
100 fSystToFix( measurement.GetConstantParams() ),
101 fParamValues( measurement.GetParamValues() ),
102 fNomLumi( measurement.GetLumi() ),
103 fLumiError( measurement.GetLumi()*measurement.GetLumiRelErr() ),
104 fLowBin( measurement.GetBinLow() ),
105 fHighBin( measurement.GetBinHigh() ),
121 if( proto_config == NULL ) {
122 std::cout <<
"Error: Did not find 'ModelConfig' object in file: " << ws_single->
GetName()
127 std::vector<std::string> poi_list = measurement.
GetPOIList();
128 if( poi_list.size()==0 ) {
129 cxcoutWHF <<
"No Parametetrs of interest are set" << std::endl;
133 std::stringstream sstream;
134 sstream <<
"Setting Parameter(s) of Interest as: ";
135 for(
unsigned int i = 0; i < poi_list.size(); ++i) {
136 sstream << poi_list.at(i) <<
" ";
141 for(
unsigned int i = 0; i < poi_list.size(); ++i ) {
142 std::string poi_name = poi_list.at(i);
148 std::cout <<
"WARNING: Can't find parameter of interest: " << poi_name
149 <<
" in Workspace. Not setting in ModelConfig." << std::endl;
156 std::string NewModelName =
"newSimPdf";
161 std::cout <<
"Error: Failed to find dataset: " << expData
162 <<
" in workspace" << std::endl;
165 if(poi_list.size()!=0){
177 if( !pdf ) pdf = ws_single->
pdf( ModelName );
178 const RooArgSet* observables = ws_single->
set(
"observables");
181 std::string SnapShotName =
"NominalParamValues";
188 std::string AsimovName = asimov.
GetName();
190 cxcoutPHF <<
"Generating additional Asimov Dataset: " << AsimovName << std::endl;
195 cxcoutPHF <<
"Importing Asimov dataset" << std::endl;
196 bool failure = ws_single->
import(*asimov_dataset,
Rename(AsimovName.c_str()));
198 std::cout <<
"Error: Failed to import Asimov dataset: " << AsimovName
200 delete asimov_dataset;
209 delete asimov_dataset;
231 string ch_name = channel.
GetName();
235 if( ws_single == NULL ) {
236 cxcoutF(
HistFactory) <<
"Error: Failed to make Single-Channel workspace for channel: " << ch_name
237 <<
" and measurement: " << measurement.
GetName() << std::endl;
267 vector<std::unique_ptr<RooWorkspace>> channel_workspaces;
268 vector<string> channel_names;
270 for(
unsigned int chanItr = 0; chanItr < measurement.
GetChannels().
size(); ++chanItr ) {
276 <<
" has uninitialized histogram pointers" << std::endl;
280 string ch_name = channel.
GetName();
281 channel_names.push_back(ch_name);
286 channel_workspaces.emplace_back(ws_single);
308 for (
unsigned int idx=0; idx <
fObsNameVec.size(); ++idx) {
320 obs->setBinning(binning);
342 unsigned int histndim(1);
343 std::string classname = hist->
ClassName();
344 if (classname.find(
"TH1")==0) { histndim=1; }
345 else if (classname.find(
"TH2")==0) { histndim=2; }
346 else if (classname.find(
"TH3")==0) { histndim=3; }
349 prefix +=
"_Hist_alphanominal";
351 RooDataHist histDHist((prefix +
"DHist").c_str(),
"",observables,hist);
352 RooHistFunc histFunc(prefix.c_str(),
"",observables,histDHist,0);
363 std::vector<std::string> & constraintTermNames) {
364 std::string paramName = param.
GetName();
365 std::string constraintName = paramName +
"Constraint";
368 if(
proto.pdf(constraintName))
return;
372 const double gaussSigma = isUniform ? 100. : 1.0;
374 cxcoutIHF <<
"Added a uniform constraint for " << paramName <<
" as a Gaussian constraint with a very large sigma " << std::endl;
377 std::stringstream command;
378 command <<
"Gaussian::" << constraintName <<
"(" << paramName <<
",nom_" << paramName <<
"[0.,-10,10],"
379 << gaussSigma <<
")";
380 constraintTermNames.emplace_back(
proto.factory( command.str().c_str() )->GetName());
381 auto * normParam =
proto.var(std::string(
"nom_") + paramName);
382 normParam->setConstant();
383 const_cast<RooArgSet*
>(
proto.set(
"globalObservables"))->add(*normParam);
393 for(
auto const& histoSys : histoSysList) {
394 const std::string histoSysName = histoSys.GetName();
397 temp = (
RooRealVar*)
proto.factory((
"alpha_" + histoSysName + range).c_str());
410 const string& prefix,
415 vector<double> low, high;
418 for(
unsigned int j=0; j<histoSysList.size(); ++j){
419 std::stringstream str;
422 const HistoSys& histoSys = histoSysList.at(j);
423 RooDataHist* lowDHist =
new RooDataHist((prefix+str.str()+
"lowDHist").c_str(),
"",observables, histoSys.GetHistoLow());
424 RooDataHist* highDHist =
new RooDataHist((prefix+str.str()+
"highDHist").c_str(),
"",observables, histoSys.GetHistoHigh());
427 lowSet.
add(*lowFunc);
428 highSet.
add(*highFunc);
433 interp.setPositiveDefinite();
434 interp.setAllInterpCodes(4);
437 interp.setBinIntegrator(observableSet);
438 interp.forceNumInt();
442 auto interpInWS =
proto->arg(prefix.c_str());
453 std::vector<string> prodNames;
456 vector<string> normFactorNames, rangeNames;
459 string overallNorm_times_sigmaEpsilon = sample.
GetName() +
"_" + channel +
"_scaleFactors";
460 auto sigEps =
proto->arg(sigmaEpsilon.c_str());
462 auto normFactor = std::make_unique<RooProduct>(overallNorm_times_sigmaEpsilon.c_str(), overallNorm_times_sigmaEpsilon.c_str(),
RooArgList(*sigEps));
464 if(normList.size() > 0){
466 for(vector<NormFactor>::iterator itr = normList.begin(); itr != normList.end(); ++itr){
472 varname = norm.
GetName() +
"_" + channel;
481 std::stringstream range;
484 if(
proto->obj(varname) == NULL) {
487 proto->factory((varname + range.str()).c_str());
493 cxcoutW(
HistFactory) <<
"Const attribute to <NormFactor> tag is deprecated, will ignore." <<
494 " Instead, add \n\t<ParamSetting Const=\"True\"> " << varname <<
" </ParamSetting>\n" <<
495 " to your top-level XML's <Measurement> entry" << endl;
497 prodNames.push_back(varname);
498 rangeNames.push_back(range.str());
499 normFactorNames.push_back(varname);
503 for (
const auto&
name : prodNames) {
506 normFactor->addTerm(arg);
511 unsigned int rangeIndex=0;
512 for( vector<string>::iterator nit = normFactorNames.begin(); nit!=normFactorNames.end(); ++nit){
513 if( count (normFactorNames.begin(), normFactorNames.end(), *nit) > 1 ){
514 cxcoutI(
HistFactory) <<
"<NormFactor Name =\""<<*nit<<
"\"> is duplicated for <Sample Name=\""
515 << sample.
GetName() <<
"\">, but only one factor will be included. \n Instead, define something like"
516 <<
"\n\t<Function Name=\""<<*nit<<
"Squared\" Expression=\""<<*nit<<
"*"<<*nit<<
"\" Var=\""<<*nit<<rangeNames.at(rangeIndex)
517 <<
"\"> \nin your top-level XML's <Measurment> entry and use <NormFactor Name=\""<<*nit<<
"Squared\" in your channel XML file."<< endl;
527 std::vector<OverallSys>& systList,
528 vector<string>& constraintTermNames,
529 vector<string>& totSystTermNames) {
535 totSystTermNames.push_back(prefix);
538 vector<double> lowVec, highVec;
540 std::map<std::string, double>::iterator itconstr;
541 for(
unsigned int i = 0; i < systList.size(); ++i) {
544 std::string strname = sys.
GetName();
545 const char *
name = strname.c_str();
556 cxcoutI(
HistFactory) <<
"HistoToWorkspaceFast::AddConstraintTerm - zero uncertainty assigned - skip systematic " << sys.
GetName() << std::endl;
559 double tauVal = 1./(relerr*relerr);
560 double sqtau = 1./relerr;
574 alphaOfBeta->
Print(
"t");
577 constraintTermNames.push_back(
gamma->GetName());
581 const_cast<RooArgSet*
>(
proto->set(
"globalObservables"))->add(*yvar);
584 params.
add(*alphaOfBeta);
595 makeGaussianConstraint(*alpha, *
proto, isUniform, constraintTermNames);
606 double tauVal = 1./relerr;
607 std::string tauName =
"tau_" + sys.
GetName();
609 double kappaVal = 1. + relerr;
610 std::string kappaName =
"kappa_" + sys.
GetName();
612 const char * alphaName = alpha->
GetName();
614 std::string alphaOfBetaName =
"alphaOfBeta_" + sys.
GetName();
616 tauName.c_str(),kappaName.c_str(),alphaName,
617 tauName.c_str(),kappaName.c_str(),alphaName ) );
621 alphaOfBeta->
Print(
"t");
622 params.
add(*alphaOfBeta);
627 double low = sys.
GetLow();
629 lowVec.push_back(low);
630 highVec.push_back(high);
634 if(systList.size() > 0){
639 assert(
int(lowVec.size()) == params.
getSize() );
641 FlexibleInterpVar interp( (interpName).c_str(),
"", params, 1., lowVec, highVec);
644 proto->import(interp);
649 proto->import(interp);
659 const vector<RooProduct*>& sampleScaleFactors, std::vector<vector<RooAbsArg*>>& sampleHistFuncs)
const {
660 assert(sampleScaleFactors.size() == sampleHistFuncs.size());
665 throw std::logic_error(
"HistFactory didn't process the observables correctly. Please file a bug report.");
667 auto firstHistFunc =
dynamic_cast<const RooHistFunc*
>(sampleHistFuncs.front().front());
668 if (!firstHistFunc) {
670 firstHistFunc =
dynamic_cast<const RooHistFunc*
>(piecewiseInt->nominalHist());
672 assert(firstHistFunc);
675 const std::string binWidthFunctionName = totName +
"_binWidth";
676 RooBinWidthFunction binWidth(binWidthFunctionName.c_str(),
"Divide by bin width to obtain probability density", *firstHistFunc,
true);
677 proto->import(binWidth);
678 auto binWidthWS =
proto->function(binWidthFunctionName.c_str());
684 for (
unsigned int i=0; i < sampleHistFuncs.size(); ++i) {
685 assert(!sampleHistFuncs[i].empty());
686 coefList.
add(*sampleScaleFactors[i]);
688 std::vector<RooAbsArg*>& thisSampleHistFuncs = sampleHistFuncs[i];
689 thisSampleHistFuncs.push_back(binWidthWS);
691 if (thisSampleHistFuncs.size() == 1) {
693 shapeList.
add(*thisSampleHistFuncs.front());
696 std::string
name = thisSampleHistFuncs.front()->GetName();
697 auto pos =
name.find(
"Hist_alpha");
698 if (pos != std::string::npos) {
699 name =
name.substr(0, pos) +
"shapes";
700 }
else if ( (pos =
name.find(
"nominal")) != std::string::npos) {
701 name =
name.substr(0, pos) +
"shapes";
711 RooRealSumPdf tot(totName.c_str(), totName.c_str(), shapeList, coefList,
true);
729 vector<string>& likelihoodTermNames){
734 for(
Int_t i=lowBin; i<highBin; ++i){
735 std::stringstream str;
738 string command(
"Poisson::"+prefix+str.str()+
"("+obsPrefix+str.str()+
","+expPrefix+str.str()+
",1)");
742 cout <<
"Poisson Term " << command << endl;
746 likelihoodTermNames.push_back( temp->
GetName() );
749 proto->defineSet(prefix.c_str(),Pois);
755 map<string,double> gammaSyst,
756 map<string,double> uniformSyst,
757 map<string,double> logNormSyst,
758 map<string,double> noSyst) {
759 string pdfName(pdfNameChar);
762 if( combined_config==NULL ) {
763 std::cout <<
"Error: Failed to find object 'ModelConfig' in workspace: "
764 <<
proto->GetName() << std::endl;
769 string edit=
"EDIT::newSimPdf("+pdfName+
",";
771 string lastPdf=pdfName;
773 unsigned int numReplacements = 0;
774 unsigned int nskipped = 0;
775 map<string,double>::iterator it;
779 for(it=gammaSyst.begin(); it!=gammaSyst.end(); ++it) {
781 if(!
proto->var(
"alpha_"+it->first)){
788 double relativeUncertainty = it->second;
789 double scale = 1/
sqrt((1+1/
pow(relativeUncertainty,2)));
792 proto->factory(
Form(
"beta_%s[1,0,10]",it->first.c_str()));
793 proto->factory(
Form(
"y_%s[%f]",it->first.c_str(),1./
pow(relativeUncertainty,2))) ;
794 proto->factory(
Form(
"theta_%s[%f]",it->first.c_str(),
pow(relativeUncertainty,2))) ;
795 proto->factory(
Form(
"Gamma::beta_%sConstraint(beta_%s,sum::k_%s(y_%s,one[1]),theta_%s,zero[0])",
800 it->first.c_str())) ;
816 proto->factory(
Form(
"PolyVar::alphaOfBeta_%s(beta_%s,{%f,%f})",it->first.c_str(),it->first.c_str(),-1./scale,1./scale));
819 if(
proto->var(
"alpha_" + it->first)->isConstant()) {
820 proto->var(
"beta_" + it->first)->setConstant(
true);
823 proto->var(
"beta_" + it->first)->setConstant(
false);
829 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
831 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
841 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
842 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
846 cout <<
"Going to issue this edit command\n" << edit<< endl;
847 proto->factory( edit.c_str() );
850 cxcoutWHF <<
"---------------------\n WARNING: failed to make EDIT\n\n" << endl;
856 for(it=uniformSyst.begin(); it!=uniformSyst.end(); ++it) {
857 cout <<
"edit for " << it->first <<
"with rel uncert = " << it->second << endl;
858 if(!
proto->var(
"alpha_"+it->first)){
859 cout <<
"systematic not there" << endl;
866 proto->factory(
Form(
"beta_%s[1,0,10]",it->first.c_str()));
867 proto->factory(
Form(
"Uniform::beta_%sConstraint(beta_%s)",it->first.c_str(),it->first.c_str()));
868 proto->factory(
Form(
"PolyVar::alphaOfBeta_%s(beta_%s,{-1,1})",it->first.c_str(),it->first.c_str()));
871 if(
proto->var(
"alpha_" + it->first)->isConstant())
872 proto->var(
"beta_" + it->first)->setConstant(
true);
874 proto->var(
"beta_" + it->first)->setConstant(
false);
879 cout <<
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint" << endl;
880 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
882 cout <<
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first << endl;
883 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
885 if(
proto->pdf(
"alpha_"+it->first+
"Constraint") &&
proto->var(
"alpha_"+it->first))
886 cout <<
" checked they are there" <<
proto->pdf(
"alpha_"+it->first+
"Constraint") <<
" " <<
proto->var(
"alpha_"+it->first) << endl;
888 cout <<
"NOT THERE" << endl;
891 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
892 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
897 proto->factory( edit.c_str() );
900 cxcoutWHF <<
"---------------------\n WARNING: failed to make EDIT\n\n" << endl;
910 for(it=logNormSyst.begin(); it!=logNormSyst.end(); ++it) {
911 cout <<
"edit for " << it->first <<
"with rel uncert = " << it->second << endl;
912 if(!
proto->var(
"alpha_"+it->first)){
913 cout <<
"systematic not there" << endl;
919 double relativeUncertainty = it->second;
920 double kappa = 1+relativeUncertainty;
924 double scale = relativeUncertainty;
927 const char *
cname = it->first.c_str();
947 cout <<
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint" << endl;
948 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
950 cout <<
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first << endl;
951 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
953 if(
proto->pdf(
"alpha_"+it->first+
"Constraint") &&
proto->var(
"alpha_"+it->first) )
954 cout <<
" checked they are there" <<
proto->pdf(
"alpha_"+it->first+
"Constraint") <<
" " <<
proto->var(
"alpha_"+it->first) << endl;
956 cout <<
"NOT THERE" << endl;
959 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
960 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
965 proto->factory( edit.c_str() );
968 cxcoutWHF <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
975 proto->removeSet(
"globalObservables");
976 proto->defineSet(
"globalObservables",gobsNew);
984 for(it=noSyst.begin(); it!=noSyst.end(); ++it) {
986 cout <<
"remove constraint for parameter" << it->first << endl;
987 if(!
proto->var(
"alpha_"+it->first) || !
proto->pdf(
"alpha_"+it->first+
"Constraint") ) {
988 cout <<
"systematic not there" << endl;
995 if ( !
proto->var(
"one") ) {
proto->factory(
"one[1.0]"); }
996 proto->var(
"one")->setConstant();
999 cout <<
"alpha_"+it->first+
"Constraint=one" << endl;
1000 editList+=precede +
"alpha_"+it->first+
"Constraint=one";
1004 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
1005 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
1009 cout << edit << endl;
1010 proto->factory( edit.c_str() );
1013 cxcoutWHF <<
"---------------------\n WARNING: failed to make EDIT\n\n" << endl;
1021 edit=
"EDIT::newSimPdf("+lastPdf+
","+editList+
")";
1022 cout << edit<< endl;
1023 proto->factory( edit.c_str() );
1028 combined_config->
SetPdf(*newOne);
1031 cxcoutWHF <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
1037 FILE* covFile = fopen ((
filename).c_str(),
"w");
1038 fprintf(covFile,
" ") ;
1039 for (
auto const *myargi : static_range_cast<RooRealVar *>(*params)) {
1040 if(myargi->isConstant())
continue;
1041 fprintf(covFile,
" & %s", myargi->GetName());
1043 fprintf(covFile,
"\\\\ \\hline \n" );
1044 for (
auto const *myargi : static_range_cast<RooRealVar *>(*params)) {
1045 if(myargi->isConstant())
continue;
1046 fprintf(covFile,
"%s", myargi->GetName());
1047 for (
auto const *myargj : static_range_cast<RooRealVar *>(*params)) {
1048 if(myargj->isConstant())
continue;
1049 cout << myargi->GetName() <<
"," << myargj->GetName();
1050 fprintf(covFile,
" & %.2f",
result->correlation(*myargi, *myargj));
1053 fprintf(covFile,
" \\\\\n");
1066 Error(
"MakeSingleChannelWorkspace",
1067 "The input Channel does not contain any sample - return a nullptr");
1071 const TH1* channel_hist_template = channel.
GetSamples().front().GetHisto();
1072 if (channel_hist_template ==
nullptr) {
1074 channel_hist_template = channel.
GetSamples().front().GetHisto();
1076 if (channel_hist_template ==
nullptr) {
1077 std::ostringstream stream;
1078 stream <<
"The sample " << channel.
GetSamples().front().GetName()
1079 <<
" in channel " << channel.
GetName() <<
" does not contain a histogram. This is the channel:\n";
1080 channel.
Print(stream);
1081 Error(
"MakeSingleChannelWorkspace",
"%s", stream.str().c_str());
1086 std::cout <<
"MakeSingleChannelWorkspace: Channel: " << channel.
GetName()
1087 <<
" has uninitialized histogram pointers" << std::endl;
1101 string channel_name = channel.
GetName();
1111 for (
unsigned int idx=0; idx<
fObsNameVec.size(); ++idx ) {
1121 throw hf_exc(
"HistFactory is limited to 1- to 3-dimensional histograms.");
1125 <<
"\tStarting to process '"
1126 << channel_name <<
"' channel with " <<
fObsNameVec.size() <<
" observables"
1127 <<
"\n-----------------------------------------\n" << endl;
1133 auto proto_config = make_unique<ModelConfig>(
"ModelConfig",
proto);
1134 proto_config->SetWorkspace(*
proto);
1140 proto->factory(funcIter->c_str());
1144 RooArgSet likelihoodTerms(
"likelihoodTerms"), constraintTerms(
"constraintTerms");
1145 vector<string> likelihoodTermNames, constraintTermNames, totSystTermNames;
1147 std::vector<std::vector<RooAbsArg*>> allSampleHistFuncs;
1148 std::vector<RooProduct*> sampleScaleFactors;
1150 vector< pair<string,string> > statNamePairs;
1151 vector< pair<const TH1*, const TH1*> > statHistPairs;
1152 const std::string statFuncName =
"mc_stat_" + channel_name;
1154 string prefix, range;
1159 std::stringstream lumiStr;
1162 proto->factory(lumiStr.str().c_str());
1165 std::stringstream lumiErrorStr;
1167 proto->factory((
"Gaussian::lumiConstraint(Lumi,"+lumiErrorStr.str()+
")").c_str());
1168 proto->var(
"nominalLumi")->setConstant();
1169 proto->defineSet(
"globalObservables",
"nominalLumi");
1171 constraintTermNames.push_back(
"lumiConstraint");
1180 string overallSystName = sample.GetName() +
"_" + channel_name +
"_epsilon";
1182 string systSourcePrefix =
"alpha_";
1187 sample.GetOverallSysList(), constraintTermNames , totSystTermNames);
1189 allSampleHistFuncs.emplace_back();
1190 std::vector<RooAbsArg*>& sampleHistFuncs = allSampleHistFuncs.back();
1194 assert(normFactors);
1203 const TH1* nominal = sample.GetHisto();
1215 string expPrefix = sample.
GetName() +
"_" + channel_name;
1219 assert(nominalHistFunc);
1221 if(sample.GetHistoSysList().size() == 0) {
1223 cxcoutI(
HistFactory) << sample.GetName() +
"_" + channel_name +
" has no variation histograms " << endl;
1225 sampleHistFuncs.push_back(nominalHistFunc);
1229 string constraintPrefix = sample.GetName() +
"_" + channel_name +
"_Hist_alpha";
1232 RooArgList interpParams = makeInterpolationParameters(sample.GetHistoSysList(), *
proto);
1235 for(std::size_t i = 0; i < interpParams.
size(); ++i) {
1236 bool isUniform = measurement.
GetUniformSyst().count(sample.GetHistoSysList()[i].GetName()) > 0;
1237 makeGaussianConstraint(interpParams[i], *
proto, isUniform, constraintTermNames);
1241 sampleHistFuncs.push_back( makeLinInterp(interpParams, nominalHistFunc,
proto,
1242 sample.GetHistoSysList(), constraintPrefix, observables) );
1249 if( sample.GetStatError().GetActivate() ) {
1252 cxcoutF(
HistFactory) <<
"Cannot include Stat Error for histograms of more than 3 dimensions."
1261 cxcoutI(
HistFactory) <<
"Sample: " << sample.GetName() <<
" to be included in Stat Error "
1262 <<
"for channel " << channel_name
1265 string UncertName = sample.GetName() +
"_" + channel_name +
"_StatAbsolUncert";
1266 TH1* statErrorHist = NULL;
1268 if( sample.GetStatError().GetErrorHist() == NULL ) {
1271 <<
" Channel: " << channel_name
1272 <<
" Sample: " << sample.GetName()
1279 statErrorHist = (
TH1*) sample.GetStatError().GetErrorHist()->
Clone();
1284 <<
"\tChannel: " << channel_name
1285 <<
"\tSample: " << sample.GetName()
1286 <<
"\tError Histogram: " << statErrorHist->
GetName() << std::endl;
1289 statErrorHist->
Multiply( nominal );
1290 statErrorHist->
SetName( UncertName.c_str() );
1295 statHistPairs.push_back( std::make_pair(nominal, statErrorHist) );
1313 if( paramHist ==
nullptr ) {
1318 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1319 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
1320 theObservables.
add( *
proto->var(*itr) );
1325 std::string ParamSetPrefix =
"gamma_stat_" + channel_name;
1326 double gammaMin = 0.0;
1327 double gammaMax = 10.0;
1329 ParamSetPrefix.c_str(),
1331 gammaMin, gammaMax);
1333 ParamHistFunc statUncertFunc(statFuncName.c_str(), statFuncName.c_str(),
1334 theObservables, statFactorParams );
1342 sampleHistFuncs.push_back(paramHist);
1351 if( sample.GetShapeFactorList().size() > 0 ) {
1354 cxcoutF(
HistFactory) <<
"Cannot include Stat Error for histograms of more than 3 dimensions."
1359 cxcoutI(
HistFactory) <<
"Sample: " << sample.GetName() <<
" in channel: " << channel_name
1360 <<
" to be include a ShapeFactor."
1363 std::vector<ParamHistFunc*> paramHistFuncList;
1365 for(
unsigned int i=0; i < sample.GetShapeFactorList().
size(); ++i) {
1367 ShapeFactor& shapeFactor = sample.GetShapeFactorList().at(i);
1369 std::string funcName = channel_name +
"_" + shapeFactor.
GetName() +
"_shapeFactor";
1371 if( paramHist ==
nullptr ) {
1374 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1375 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
1376 theObservables.
add( *
proto->var(*itr) );
1380 std::string funcParams =
"gamma_" + shapeFactor.
GetName();
1386 theObservables, 0, 1000);
1389 ParamHistFunc shapeFactorFunc( funcName.c_str(), funcName.c_str(),
1390 theObservables, shapeFactorParams );
1396 <<
" to have initial shape from hist: "
1399 shapeFactorFunc.
setShape( initialShape );
1405 <<
" to be constant" << std::endl;
1414 paramHistFuncList.push_back(paramHist);
1421 for (
auto paramHistFunc : paramHistFuncList) {
1423 auto paramHFinWS =
proto->arg(paramHistFunc->GetName());
1424 assert(paramHFinWS);
1425 delete paramHistFunc;
1427 sampleHistFuncs.push_back(paramHFinWS);
1438 if( sample.GetShapeSysList().size() != 0 ) {
1441 cxcoutF(
HistFactory) <<
"Cannot include Stat Error for histograms of more than 3 dimensions."
1447 std::vector<string> ShapeSysNames;
1449 for(
unsigned int i = 0; i < sample.GetShapeSysList().
size(); ++i) {
1462 cxcoutI(
HistFactory) <<
"Sample: " << sample.GetName() <<
" in channel: " << channel_name
1463 <<
" to include a ShapeSys." << std::endl;
1465 std::string funcName = channel_name +
"_" + shapeSys.
GetName() +
"_ShapeSys";
1466 ShapeSysNames.push_back( funcName );
1468 if( paramHist == NULL ) {
1474 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1476 theObservables.
add( *
proto->var(*itr) );
1480 std::string funcParams =
"gamma_" + shapeSys.
GetName();
1483 theObservables, 0, 10);
1486 ParamHistFunc shapeFactorFunc( funcName.c_str(), funcName.c_str(),
1487 theObservables, shapeFactorParams );
1510 double minShapeUncertainty = 0.0;
1512 *paramHist, shapeErrorHist,
1514 minShapeUncertainty);
1523 for(
unsigned int i = 0; i < ShapeSysNames.size(); ++i ) {
1524 auto func =
proto->function(ShapeSysNames.at(i).c_str());
1526 sampleHistFuncs.push_back(func);
1537 auto lumi =
proto->arg(
"Lumi");
1538 if( !sample.GetNormalizeByTheory() ) {
1541 lumiParamString += measurement.
GetLumi();
1543 lumi =
proto->factory((
"Lumi[" + lumiParamString +
"]").
Data());
1549 normFactors->addTerm(lumi);
1555 auto normFactorsInWS =
dynamic_cast<RooProduct*
>(
proto->arg(normFactors->GetName()));
1556 assert(normFactorsInWS);
1558 sampleScaleFactors.push_back(normFactorsInWS);
1563 if( statHistPairs.size() > 0 ) {
1568 if( fracStatError == NULL ) {
1570 << channel_name +
"_StatUncert" +
"_RelErr" << std::endl;
1578 << chanStatUncertFunc->
GetName()
1579 <<
" params: " << chanStatUncertFunc->
paramList()
1598 *chanStatUncertFunc, fracStatError.
get(),
1600 statRelErrorThreshold);
1604 for (
unsigned int i = 0; i < statHistPairs.size() ; ++i )
1605 delete statHistPairs[i].
second;
1607 statHistPairs.clear();
1615 sampleScaleFactors, allSampleHistFuncs);
1616 likelihoodTermNames.push_back(channel_name+
"_model");
1620 for(
unsigned int i=0; i<systToFix.size(); ++i){
1628 if(systToFix.at(i)==
"Lumi"){
1629 auxMeas =
proto->var(
"nominalLumi");
1631 auxMeas =
proto->var(std::string(
"nom_") + temp->
GetName());
1635 const_cast<RooArgSet*
>(
proto->set(
"globalObservables"))->remove(*auxMeas);
1642 <<
" could not set it to constant" << endl;
1648 for(
unsigned int i=0; i<constraintTermNames.size(); ++i){
1649 RooAbsArg* proto_arg = (
proto->arg(constraintTermNames[i].c_str()));
1650 if( proto_arg==NULL ) {
1652 <<
" in workspace: " <<
proto->GetName() << std::endl;
1655 constraintTerms.
add( *proto_arg );
1658 for(
unsigned int i=0; i<likelihoodTermNames.size(); ++i){
1659 RooAbsArg* proto_arg = (
proto->arg(likelihoodTermNames[i].c_str()));
1660 if( proto_arg==NULL ) {
1662 <<
" in workspace: " <<
proto->GetName() << std::endl;
1665 likelihoodTerms.add( *proto_arg );
1667 proto->defineSet(
"constraintTerms",constraintTerms);
1668 proto->defineSet(
"likelihoodTerms",likelihoodTerms);
1672 std::string observablesStr;
1674 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1676 observables.
add( *
proto->var(*itr) );
1677 if (!observablesStr.empty()) { observablesStr +=
","; }
1678 observablesStr += *itr;
1690 <<
"\timport model into workspace"
1691 <<
"\n-----------------------------------------\n" << endl;
1693 auto model = make_unique<RooProdPdf>(
1694 (
"model_"+channel_name).c_str(),
1695 "product of Poissons accross bins for a single channel",
1696 constraintTerms,
Conditional(likelihoodTerms,observables));
1699 proto_config->SetPdf(*model);
1700 proto_config->SetObservables(observables);
1701 proto_config->SetGlobalObservables(*
proto->set(
"globalObservables"));
1705 proto->import(*proto_config,proto_config->GetName());
1706 proto->importClassCode();
1711 const char* weightName=
"weightVar";
1717 int asymcalcPrintLevel = 0;
1728 RooDataSet dataset{
"obsData",
"",*
proto->set(
"obsAndWeight"),weightName};
1730 proto->import(dataset);
1735 if(
data.GetName().empty()) {
1737 <<
" has no name! The name always needs to be set for additional datasets, "
1738 <<
"either via the \"Name\" tag in the XML or via RooStats::HistFactory::Data::SetName()." << std::endl;
1741 std::string
const& dataName =
data.GetName();
1742 TH1 const* mnominal =
data.GetHisto();
1745 <<
" with name: " << dataName <<
" is NULL" << std::endl;
1750 RooDataSet dataset{dataName.c_str(),
"", *
proto->set(
"obsAndWeight"), weightName};
1752 proto->import(dataset);
1764 TH1 const& mnominal,
1766 std::vector<std::string>
const& obsNameVec) {
1772 if (obsNameVec.empty() ) {
1773 Error(
"ConfigureHistFactoryDataset",
"Invalid input - return");
1781 for (
int i=1; i<=ax->
GetNbins(); ++i) {
1784 proto.var( obsNameVec[0] )->setVal( xval );
1786 if(obsNameVec.size()==1) {
1788 obsDataUnbinned.
add( *
proto.set(
"obsAndWeight"), fval );
1791 for(
int j=1; j<=ay->
GetNbins(); ++j) {
1793 proto.var( obsNameVec[1] )->setVal( yval );
1795 if(obsNameVec.size()==2) {
1797 obsDataUnbinned.
add( *
proto.set(
"obsAndWeight"), fval );
1800 for(
int k=1; k<=az->
GetNbins(); ++k) {
1802 proto.var( obsNameVec[2] )->setVal( zval );
1804 obsDataUnbinned.
add( *
proto.set(
"obsAndWeight"), fval );
1817 unsigned int histndim(1);
1818 std::string classname = hist->
ClassName();
1819 if (classname.find(
"TH1")==0) { histndim=1; }
1820 else if (classname.find(
"TH2")==0) { histndim=2; }
1821 else if (classname.find(
"TH3")==0) { histndim=3; }
1823 for (
unsigned int idx=0; idx<histndim; ++idx ) {
1836 if (ch_names.empty() || chs.empty() ) {
1837 Error(
"MakeCombinedModel",
"Input vectors are empty - return a nullptr");
1840 if (chs.size() < ch_names.size() ) {
1841 Error(
"MakeCombinedModel",
"Input vector of workspace has an invalid size - return a nullptr");
1849 map<string, RooAbsPdf*> pdfMap;
1850 vector<RooAbsPdf*> models;
1853 for(
unsigned int i = 0; i< ch_names.size(); ++i){
1860 stringstream channelString;
1861 channelString <<
"channelCat[";
1862 for(
unsigned int i = 0; i< ch_names.size(); ++i){
1863 string channel_name=ch_names[i];
1864 if (i == 0 && isdigit(channel_name[0])) {
1865 throw std::invalid_argument(
"The first channel name for HistFactory cannot start with a digit. Got " + channel_name);
1867 if (channel_name.find(
',') != std::string::npos) {
1868 throw std::invalid_argument(
"Channel names for HistFactory cannot contain ','. Got " + channel_name);
1871 if (i == 0) channelString << channel_name ;
1872 else channelString <<
',' << channel_name ;
1876 if(!model) cout <<
"failed to find model for channel"<<endl;
1878 models.push_back(model);
1879 globalObs.
add(*ch->
set(
"globalObservables"),
true);
1882 pdfMap[channel_name]=model;
1884 channelString <<
"]";
1887 <<
"\tEntering combination"
1888 <<
"\n-----------------------------------------\n" << endl;
1894 if (!channelCat)
throw std::runtime_error(
"Unable to construct a category from string " + channelString.str());
1901 combined->
import(globalObs);
1902 combined->
defineSet(
"globalObservables",globalObs);
1909 <<
"\tcreate toy data for " << channelString.str()
1910 <<
"\n-----------------------------------------\n" << endl;
1916 combined->
factory(
"weightVar[0,-1e10,1e10]");
1917 obsList.
add(*combined->
var(
"weightVar"));
1922 if( asimov_combined ) {
1923 combined->
import( *asimov_combined,
Rename(
"asimovData"));
1926 std::cout <<
"Error: Failed to create combined asimov dataset" << std::endl;
1929 delete asimov_combined;
1941 obsList.
add(*channelCat);
1942 combined->
defineSet(
"observables",obsList);
1949 <<
"\tImporting combined model"
1950 <<
"\n-----------------------------------------\n" << endl;
1956 std::map< std::string, double>::iterator param_itr =
fParamValues.begin();
1959 std::string paramName = param_itr->first;
1960 double paramVal = param_itr->second;
1964 temp->
setVal( paramVal );
1965 cxcoutI(
HistFactory) <<
"setting " << paramName <<
" to the value: " << paramVal << endl;
1967 cxcoutE(
HistFactory) <<
"could not find variable " << paramName <<
" could not set its value" << endl;
1971 for(
unsigned int i=0; i<
fSystToFix.size(); ++i){
1986 combined_config->
SetPdf(*simPdf);
1989 combined->
import(*combined_config,combined_config->
GetName());
1994 delete combined_config;
2002 std::vector<std::unique_ptr<RooWorkspace>>& wspace_vec,
2003 std::vector<std::string>
const& channel_names,
2004 std::string
const& dataSetName,
2013 for(
unsigned int i = 0; i< channel_names.size(); ++i){
2016 cxcoutPHF <<
"Merging data for channel " << channel_names[i].c_str() << std::endl;
2018 if( !obsDataInChannel ) {
2019 std::cout <<
"Error: Can't find DataSet: " << dataSetName
2020 <<
" in channel: " << channel_names.at(i)
2027 obsList,
Index(*channelCat),
2029 Import(channel_names[i].c_str(),*obsDataInChannel));
2031 simData->
append(*tempData);
2042 combined->
import(*simData,
Rename(dataSetName.c_str()));
2047 std::cout <<
"Error: Unable to merge observable datasets" << std::endl;
2066 Int_t binNumber = 0;
2069 for(
Int_t i_bin = 0; i_bin < numBins; ++i_bin) {
2077 double histError = Nominal->
GetBinError( binNumber );
2080 if( histError != histError ) {
2081 std::cout <<
"Warning: In histogram " << Nominal->
GetName()
2082 <<
" bin error for bin " << i_bin
2083 <<
" is NAN. Not using Error!!!"
2091 if( histError < 0 ) {
2092 std::cout <<
"Warning: In histogram " << Nominal->
GetName()
2093 <<
" bin error for bin " << binNumber
2094 <<
" is < 0. Setting Error to 0"
2119 unsigned int numHists = HistVec.size();
2121 if( numHists == 0 ) {
2122 cxcoutE(
HistFactory) <<
"Warning: Empty Hist Vector, cannot create total uncertainty" << std::endl;
2126 const TH1* HistTemplate = HistVec.at(0).first;
2131 for(
unsigned int i = 0; i < HistVec.size(); ++i ) {
2133 const TH1* nominal = HistVec.at(i).first;
2134 const TH1* error = HistVec.at(i).second;
2146 std::vector<double> TotalBinContent( numBins, 0.0);
2147 std::vector<double> HistErrorsSqr( numBins, 0.0);
2149 Int_t binNumber = 0;
2152 for(
Int_t i_bins = 0; i_bins < numBins; ++i_bins) {
2159 for(
unsigned int i_hist = 0; i_hist < numHists; ++i_hist ) {
2161 const TH1* nominal = HistVec.at(i_hist).first;
2162 const TH1* error = HistVec.at(i_hist).second;
2177 if( histError != histError ) {
2179 <<
" bin error for bin " << binNumber
2180 <<
" is NAN. Not using error!!";
2184 TotalBinContent.at(i_bins) += histValue;
2185 HistErrorsSqr.at(i_bins) += histError*histError;
2197 for(
Int_t i = 0; i < numBins; ++i) {
2205 double ErrorsSqr = HistErrorsSqr.at(i);
2206 double TotalVal = TotalBinContent.at(i);
2208 if( TotalVal <= 0 ) {
2210 <<
" is <= 0. Setting error to 0"
2217 double RelativeError =
sqrt(ErrorsSqr) / TotalVal;
2221 if( RelativeError != RelativeError ) {
2223 <<
" HistErrorsSqr: " << ErrorsSqr
2224 <<
" TotalVal: " << TotalVal;
2238 <<
" Error = " <<
sqrt(ErrorsSqr)
2239 <<
" CentralVal = " << TotalVal
2240 <<
" RelativeError = " << RelativeError <<
"\n";
2244 return std::unique_ptr<TH1>(ErrorHist);
2280 if( numBins != numParams ) {
2281 std::cout <<
"Error: In createStatConstraintTerms, encountered bad number of bins" << std::endl;
2282 std::cout <<
"Given histogram with " << numBins <<
" bins,"
2283 <<
" but require exactly " << numParams << std::endl;
2287 Int_t TH1BinNumber = 0;
2299 <<
". Type of constraint: " <<
type << std::endl;
2303 const double sigmaRel = uncertHist->
GetBinContent(TH1BinNumber);
2307 if( sigmaRel <= 0 ){
2310 <<
" because sigma = " << sigmaRel
2312 <<
" (TH1 bin number = " << TH1BinNumber <<
")"
2314 gamma.setConstant(
true);
2319 gamma.setMax( 1 + 5*sigmaRel );
2323 std::string constrName = string(
gamma.GetName()) +
"_constraint";
2324 std::string nomName = string(
"nom_") +
gamma.GetName();
2325 std::string sigmaName = string(
gamma.GetName()) +
"_sigma";
2326 std::string poisMeanName = string(
gamma.GetName()) +
"_poisMean";
2334 RooConstVar constrSigma( sigmaName.c_str(), sigmaName.c_str(), sigmaRel );
2337 RooRealVar constrNom(nomName.c_str(), nomName.c_str(), 1.0,0,10);
2342 constrNom,
gamma, constrSigma );
2348 gamma.setError(sigmaRel);
2351 double tau = 1/sigmaRel/sigmaRel;
2354 RooRealVar constrNom(nomName.c_str(), nomName.c_str(), tau);
2359 std::string scalingName = string(
gamma.GetName()) +
"_tau";
2360 RooConstVar poissonScaling( scalingName.c_str(), scalingName.c_str(), tau);
2368 RooPoisson pois(constrName.c_str(), constrName.c_str(), constrNom, constrMean);
2372 if (std::string(
gamma.GetName()).find(
"gamma_stat") != std::string::npos) {
2375 gamma.setError(sigmaRel);
2380 std::cout <<
"Error: Did not recognize Stat Error constraint term type: "
2381 <<
type <<
" for : " << paramHist.
GetName() << std::endl;
2388 if( sigmaRel < minSigma ) {
2390 <<
" and is < " << minSigma
2391 <<
". Setting: " <<
gamma.GetName() <<
" to constant"
2393 gamma.setConstant(
true);
2396 constraintTermNames.push_back( constrName );
2397 ConstraintTerms.
add( *
proto->pdf(constrName) );
2404 if( ! globalSet->
contains(*nomVarInWorkspace) ) {
2405 globalSet->
add( *nomVarInWorkspace );
2410 return ConstraintTerms;
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char cname
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
A class which maps the current values of a RooRealVar (or a set of RooRealVars) to one of a number of...
const RooArgSet * get(Int_t masterIdx) const
void setConstant(bool constant)
const RooArgList & paramList() const
static RooArgList createParamSet(RooWorkspace &w, const std::string &, const RooArgList &Vars)
Create the list of RooRealVar parameters which represent the height of the histogram bins.
void setShape(TH1 *shape)
The PiecewiseInterpolation is a class that can morph distributions into each other,...
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
void Print(Option_t *options=0) const override
Print the object to the defaultPrintStream().
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
bool contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
void Print(Option_t *options=0) const override
This method must be overridden when a class wants to print itself.
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Storage_t::size_type size() const
RooAbsData is the common abstract base class for binned and unbinned datasets.
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
void setConstant(bool value=true)
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
virtual void forceNumInt(bool flag=true)
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
RooArgList is a container object that can hold multiple RooAbsArg objects.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
RooBinWidthFunction is a class that returns the bin width (or volume) given a RooHistFunc.
Class RooBinning is an implements RooAbsBinning in terms of an array of boundary values,...
RooCategory is an object to represent discrete states.
bool setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
RooConstVar represent a constant real-valued object.
The RooDataHist is a container class to hold N-dimensional binned data.
RooDataSet is a container class to hold unbinned data.
void append(RooDataSet &data)
Add all data points of given data set to this data set.
void add(const RooArgSet &row, double weight=1.0, double weightError=0) override
Add one ore more rows of data.
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Switches the message service to a different level while the instance is alive.
RooHistFunc implements a real-valued function sampled from a multidimensional histogram.
static RooMsgService & instance()
Return reference to singleton instance.
bool isActive(const RooAbsArg *self, RooFit::MsgTopic facility, RooFit::MsgLevel level)
Check if logging is active for given object/topic/RooFit::MsgLevel combination.
void setNoRounding(bool flag=true)
Switch off/on rounding of x to the nearest integer.
A RooProduct represents the product of a given set of RooAbsReal objects.
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
RooRealVar represents a variable that can be changed from the outside.
void setVal(double value) override
Set value of variable to 'value'.
void setBins(Int_t nBins, const char *name=0)
Create a uniform binning under name 'name' for this variable.
void setMin(const char *name, double value)
Set minimum of name range to given value.
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
static void SetPrintLevel(int level)
set print level (static function)
static RooAbsData * GenerateAsimovData(const RooAbsPdf &pdf, const RooArgSet &observables)
generate the asimov data for the observables (not the global ones) need to deal with the case of a si...
TODO Here, we are missing some documentation.
void ConfigureWorkspace(RooWorkspace *)
This class encapsulates all information for the statistical interpretation of one experiment.
std::vector< RooStats::HistFactory::Data > & GetAdditionalData()
retrieve vector of additional data objects
void Print(std::ostream &=std::cout)
HistFactory::StatErrorConfig & GetStatErrorConfig()
get information about threshold for statistical uncertainties and constraint term
RooStats::HistFactory::Data & GetData()
get data object
bool CheckHistograms() const
std::vector< RooStats::HistFactory::Sample > & GetSamples()
get vector of samples for this channel
std::string GetName() const
get name of channel
void setAllInterpCodes(int code)
This class provides helper functions for creating likelihood models from histograms.
std::unique_ptr< RooProduct > CreateNormFactor(RooWorkspace *proto, std::string &channel, std::string &sigmaEpsilon, Sample &sample, bool doRatio)
static void EditSyst(RooWorkspace *proto, const char *pdfNameChar, std::map< std::string, double > gammaSyst, std::map< std::string, double > uniformSyst, std::map< std::string, double > logNormSyst, std::map< std::string, double > noSyst)
RooWorkspace * MakeSingleChannelModel(Measurement &measurement, Channel &channel)
RooWorkspace * MakeSingleChannelWorkspace(Measurement &measurement, Channel &channel)
void GuessObsNameVec(const TH1 *hist)
RooArgList createStatConstraintTerms(RooWorkspace *proto, std::vector< std::string > &constraintTerms, ParamHistFunc ¶mHist, const TH1 *uncertHist, Constraint::Type type, double minSigma)
std::vector< std::string > fPreprocessFunctions
void SetFunctionsToPreprocess(std::vector< std::string > lines)
std::vector< std::string > fObsNameVec
RooHistFunc * MakeExpectedHistFunc(const TH1 *hist, RooWorkspace *proto, std::string prefix, const RooArgList &observables) const
Create the nominal hist function from hist, and register it in the workspace.
TH1 * MakeAbsolUncertaintyHist(const std::string &Name, const TH1 *Hist)
std::map< std::string, double > fParamValues
static void ConfigureWorkspaceForMeasurement(const std::string &ModelName, RooWorkspace *ws_single, Measurement &measurement)
void MakeTotalExpected(RooWorkspace *proto, const std::string &totName, const std::vector< RooProduct * > &sampleScaleFactors, std::vector< std::vector< RooAbsArg * > > &sampleHistFuncs) const
void AddPoissonTerms(RooWorkspace *proto, std::string prefix, std::string obsPrefix, std::string expPrefix, int lowBin, int highBin, std::vector< std::string > &likelihoodTermNames)
void ConfigureHistFactoryDataset(RooDataSet &obsData, TH1 const &nominal, RooWorkspace &proto, std::vector< std::string > const &obsNameVec)
static void PrintCovarianceMatrix(RooFitResult *result, RooArgSet *params, std::string filename)
std::vector< std::string > fSystToFix
RooDataSet * MergeDataSets(RooWorkspace *combined, std::vector< std::unique_ptr< RooWorkspace > > &wspace_vec, std::vector< std::string > const &channel_names, std::string const &dataSetName, RooArgList obsList, RooCategory *channelCat)
std::unique_ptr< TH1 > MakeScaledUncertaintyHist(const std::string &Name, std::vector< std::pair< const TH1 *, const TH1 * > > HistVec) const
void AddConstraintTerms(RooWorkspace *proto, Measurement &measurement, std::string prefix, std::string interpName, std::vector< OverallSys > &systList, std::vector< std::string > &likelihoodTermNames, std::vector< std::string > &totSystTermNames)
RooArgList createObservables(const TH1 *hist, RooWorkspace *proto) const
Create observables of type RooRealVar. Creates 1 to 3 observables, depending on the type of the histo...
HistoToWorkspaceFactoryFast()
RooWorkspace * MakeCombinedModel(std::vector< std::string >, std::vector< std::unique_ptr< RooWorkspace > > &)
const std::string & GetName() const
The RooStats::HistFactory::Measurement class can be used to construct a model by combining multiple R...
std::map< std::string, double > & GetGammaSyst()
std::map< std::string, double > & GetLogNormSyst()
std::map< std::string, double > & GetNoSyst()
std::vector< std::string > & GetPOIList()
get vector of PoI names
std::map< std::string, double > & GetUniformSyst()
std::vector< std::string > & GetConstantParams()
get vector of all constant parameters
std::vector< RooStats::HistFactory::Channel > & GetChannels()
std::vector< RooStats::HistFactory::Asimov > & GetAsimovDatasets()
get vector of defined Asimov Datasets
std::vector< std::string > GetPreprocessFunctions() const
Returns a list of defined preprocess function expressions.
double GetLumi()
retrieve integrated luminosity
Configuration for an un- constrained overall systematic to scale sample normalisations.
std::string GetName() const
Configuration for a constrained overall systematic to scale sample normalisations.
const std::string & GetName() const
std::string GetName() const
get name of sample
std::vector< RooStats::HistFactory::NormFactor > & GetNormFactorList()
*Un*constrained bin-by-bin variation of affected histogram.
const TH1 * GetInitialShape() const
Constrained bin-by-bin variation of affected histogram.
Constraint::Type GetConstraintType() const
const TH1 * GetErrorHist() const
double GetRelErrorThreshold() const
Constraint::Type GetConstraintType() const
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
virtual void SetObservables(const RooArgSet &set)
Specify the observables.
virtual void SetWorkspace(RooWorkspace &ws)
void GuessObsAndNuisance(const RooAbsData &data, bool printModelConfig=true)
Makes sensible guesses of observables, parameters of interest and nuisance parameters if one or multi...
virtual void SetParametersOfInterest(const RooArgSet &set)
Specify parameters of interest.
virtual void SetGlobalObservables(const RooArgSet &set)
Specify the global observables.
const RooArgSet * GetObservables() const
get RooArgSet for observables (return NULL if not existing)
virtual void SetPdf(const RooAbsPdf &pdf)
Set the Pdf, add to the the workspace if not already there.
The RooWorkspace is a persistable container for RooFit projects.
TObject * obj(RooStringView name) const
Return any type of object (RooAbsArg, RooAbsData or generic object) with given name)
RooAbsPdf * pdf(RooStringView name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
void Print(Option_t *opts=0) const override
Print contents of the workspace.
RooArgSet allVars() const
Return set with all variable objects.
bool saveSnapshot(const char *name, const char *paramNames)
Save snapshot of values and attributes (including "Constant") of given parameters.
bool import(const RooAbsArg &arg, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg())
Import a RooAbsArg object, e.g.
bool importClassCode(const char *pat="*", bool doReplace=false)
Inport code of all classes in the workspace that have a class name that matches pattern 'pat' and whi...
RooFactoryWSTool & factory()
Return instance to factory tool.
RooRealVar * var(RooStringView name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found.
RooAbsData * data(RooStringView name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found.
const RooArgSet * set(const char *name)
Return pointer to previously defined named set with given nmame If no such set is found a null pointe...
bool loadSnapshot(const char *name)
Load the values and attributes of the parameters in the snapshot saved with the given name.
bool defineSet(const char *name, const RooArgSet &aset, bool importMissing=false)
Define a named RooArgSet with given constituents.
const Double_t * GetArray() const
Class to manage histogram axis.
Bool_t IsVariableBinSize() const
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
const TArrayD * GetXbins() const
TH1 is the base class of all histogram classes in ROOT.
virtual Bool_t Multiply(TF1 *f1, Double_t c1=1)
Performs the operation:
virtual Int_t GetNbinsY() const
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
virtual Int_t GetNbinsZ() const
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
virtual Int_t GetNbinsX() const
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error Note that this resets the bin eror option to be of Normal Type and for the non-empt...
Bool_t IsBinUnderflow(Int_t bin, Int_t axis=0) const
Return true if the bin is underflow.
Bool_t IsBinOverflow(Int_t bin, Int_t axis=0) const
Return true if the bin is overflow.
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
void SetName(const char *name) override
Change the name of this histogram.
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
TObject * Clone(const char *newname="") const override
Make a complete copy of the underlying object.
const char * GetName() const override
Returns name of object.
Mother of all ROOT objects.
virtual const char * ClassName() const
Returns name of class to which the object belongs.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
TString & ReplaceAll(const TString &s1, const TString &s2)
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
RooCmdArg RecycleConflictNodes(bool flag=true)
RooCmdArg Rename(const char *suffix)
RooCmdArg Conditional(const RooArgSet &pdfSet, const RooArgSet &depSet, bool depsAreCond=false)
RooCmdArg Index(RooCategory &icat)
RooCmdArg Import(const char *state, TH1 &histo)
RooCmdArg WeightVar(const char *name, bool reinterpretAsWeight=false)
double beta(double x, double y)
Calculates the beta function.
RVec< PromoteTypes< T0, T1 > > pow(const T0 &x, const RVec< T1 > &v)
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Namespace for the RooStats classes.
static constexpr double second
static constexpr double gauss
bool binnedFitOptimization