85 #define alpha_Low "-5" 86 #define alpha_High "5" 87 #define NoHistConst_Low "0" 88 #define NoHistConst_High "2000" 98 namespace HistFactory{
100 HistoToWorkspaceFactoryFast::HistoToWorkspaceFactoryFast() :
101 fNomLumi(1.0), fLumiError(0),
102 fLowBin(0), fHighBin(0)
109 fSystToFix( measurement.GetConstantParams() ),
112 fLumiError( measurement.GetLumi()*measurement.GetLumiRelErr() ),
113 fLowBin( measurement.GetBinLow() ),
114 fHighBin( measurement.GetBinHigh() ) {
129 if( proto_config == NULL ) {
130 std::cout <<
"Error: Did not find 'ModelConfig' object in file: " << ws_single->
GetName()
135 std::vector<std::string> poi_list = measurement.
GetPOIList();
136 if( poi_list.size()==0 ) {
137 std::cout <<
"Warining: No Parametetrs of interest are set" << std::endl;
140 cout <<
"Setting Parameter(s) of Interest as: ";
141 for(
unsigned int i = 0; i < poi_list.size(); ++i) {
142 cout << poi_list.at(i) <<
" ";
147 for(
unsigned int i = 0; i < poi_list.size(); ++i ) {
148 std::string poi_name = poi_list.at(i);
154 std::cout <<
"WARNING: Can't find parameter of interest: " << poi_name
155 <<
" in Workspace. Not setting in ModelConfig." << std::endl;
162 std::string NewModelName =
"newSimPdf";
176 proto_config->
SetPdf( *ws_single->
pdf(
"newSimPdf" ) );
183 std::cout <<
"Error: Failed to find dataset: " << expData
184 <<
" in workspace" << std::endl;
187 if(poi_list.size()!=0){
199 if( !pdf ) pdf = ws_single->
pdf( ModelName.c_str() );
200 const RooArgSet* observables = ws_single->
set(
"observables");
203 std::string SnapShotName =
"NominalParamValues";
210 std::string AsimovName = asimov.
GetName();
212 std::cout <<
"Generating additional Asimov Dataset: " << AsimovName << std::endl;
217 std::cout <<
"Importing Asimov dataset" << std::endl;
218 bool failure = ws_single->
import(*asimov_dataset,
Rename(AsimovName.c_str()));
220 std::cout <<
"Error: Failed to import Asimov dataset: " << AsimovName
222 delete asimov_dataset;
231 delete asimov_dataset;
253 string ch_name = channel.
GetName();
257 if( ws_single == NULL ) {
258 std::cout <<
"Error: Failed to make Single-Channel workspace for channel: " << ch_name
259 <<
" and measurement: " << measurement.
GetName() << std::endl;
287 vector<RooWorkspace*> channel_workspaces;
288 vector<string> channel_names;
290 for(
unsigned int chanItr = 0; chanItr < measurement.
GetChannels().size(); ++chanItr ) {
295 std::cout <<
"MakeModelAndMeasurementsFast: Channel: " << channel.
GetName()
296 <<
" has uninitialized histogram pointers" << std::endl;
300 string ch_name = channel.
GetName();
301 channel_names.push_back(ch_name);
306 channel_workspaces.push_back(ws_single);
320 for (vector<RooWorkspace*>::iterator iter = channel_workspaces.begin() ; iter != channel_workspaces.end() ; ++iter) {
330 string prefix,
string productPrefix,
333 cout <<
"processing hist " << hist->
GetName() << endl;
335 cout <<
"hist is empty" << endl;
345 unsigned int histndim(1);
346 std::string classname = hist->
ClassName();
347 if (classname.find(
"TH1")==0) { histndim=1; }
348 else if (classname.find(
"TH2")==0) { histndim=2; }
349 else if (classname.find(
"TH3")==0) { histndim=3; }
354 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
355 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
356 if ( !proto->
var(itr->c_str()) ) {
358 if (idx==0) { axis = hist->
GetXaxis(); }
359 if (idx==1) { axis = hist->
GetYaxis(); }
360 if (idx==2) { axis = hist->
GetZaxis(); }
366 proto->
var(itr->c_str())->setBins(nbins);
368 observables.add( *proto->
var(itr->c_str()) );
377 proto->
factory((
"prod:"+productPrefix+
"("+prefix+
"_nominal,"+systTerm+
")").c_str() );
390 for(
Int_t i=lowBin; i<highBin; ++i){
391 std::stringstream str;
398 for(
int i=lowBin; i<highBin; ++i){
399 for(
int j=0; j<highBin-lowBin; ++j){
400 if(i==j) { Cov(i,j) =
sqrt(mean(i)); }
401 else { Cov(i,j) = 0; }
408 floating, mean, Cov);
410 proto->
import(constraint);
412 constraintTermNames.push_back(constraint.GetName());
416 std::vector<HistoSys> histoSysList,
417 string prefix,
string productPrefix,
419 vector<string>& constraintTermNames){
429 unsigned int histndim(1);
430 std::string classname = nominal->
ClassName();
431 if (classname.find(
"TH1")==0) { histndim=1; }
432 else if (classname.find(
"TH2")==0) { histndim=2; }
433 else if (classname.find(
"TH3")==0) { histndim=3; }
439 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
440 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
441 if ( !proto->
var(itr->c_str()) ) {
443 if (idx==0) { axis = nominal->
GetXaxis(); }
444 else if (idx==1) { axis = nominal->
GetYaxis(); }
445 else if (idx==2) { axis = nominal->
GetZaxis(); }
447 std::cout <<
"Error: Too many observables. " 448 <<
"HistFactory only accepts up to 3 observables (3d) " 457 proto->
var(itr->c_str())->setBins(nbins);
459 observables.add( *proto->
var(itr->c_str()) );
471 for(
unsigned int j=0; j<histoSysList.size(); ++j){
472 std::stringstream str;
475 HistoSys& histoSys = histoSysList.at(j);
476 string histoSysName = histoSys.
GetName();
484 string command=(
"Gaussian::alpha_"+histoSysName+
"Constraint(alpha_"+histoSysName+
",nom_alpha_"+histoSysName+
"[0.,-10,10],1.)");
485 cout << command << endl;
486 constraintTermNames.push_back( proto->
factory( command.c_str() )->
GetName() );
488 const_cast<RooArgSet*
>(proto->
set(
"globalObservables"))->add(*proto->
var((
"nom_alpha_"+histoSysName).c_str()));
495 vector<double> low, high;
498 for(
unsigned int j=0; j<histoSysList.size(); ++j){
499 std::stringstream str;
502 HistoSys& histoSys = histoSysList.at(j);
507 lowSet.
add(*lowFunc);
508 highSet.
add(*highFunc);
514 interp.setAllInterpCodes(4);
517 interp.setBinIntegrator(observableSet);
518 interp.forceNumInt();
523 proto->
factory((
"prod:"+productPrefix+
"("+prefix+
","+systTerm+
")").c_str() );
529 string overallNorm_times_sigmaEpsilon ;
533 vector<string> normFactorNames, rangeNames;
535 if(normList.size() > 0){
537 for(vector<NormFactor>::iterator itr = normList.begin(); itr != normList.end(); ++itr){
542 if(!prodNames.empty()) prodNames +=
",";
544 varname = norm.
GetName() +
"_" + channel;
553 std::stringstream range;
556 if( proto->
obj(varname.c_str()) == NULL) {
557 cout <<
"making normFactor: " << norm.
GetName() << endl;
559 proto->
factory((varname + range.str()).c_str());
565 cout <<
"WARNING: Const attribute to <NormFactor> tag is deprecated, will ignore." <<
566 " Instead, add \n\t<ParamSetting Const=\"True\">" << varname <<
"</ParamSetting>\n" <<
567 " to your top-level XML's <Measurment> entry" << endl;
570 rangeNames.push_back(range.str());
571 normFactorNames.push_back(varname);
574 overallNorm_times_sigmaEpsilon = sample.
GetName() +
"_" + channel +
"_overallNorm_x_sigma_epsilon";
575 proto->
factory((
"prod::" + overallNorm_times_sigmaEpsilon +
"(" + prodNames +
"," + sigmaEpsilon +
")").c_str());
578 unsigned int rangeIndex=0;
579 for( vector<string>::iterator nit = normFactorNames.begin(); nit!=normFactorNames.end(); ++nit){
580 if( count (normFactorNames.begin(), normFactorNames.end(), *nit) > 1 ){
581 cout <<
"WARNING: <NormFactor Name =\""<<*nit<<
"\"> is duplicated for <Sample Name=\"" 582 << sample.
GetName() <<
"\">, but only one factor will be included. \n Instead, define something like" 583 <<
"\n\t<Function Name=\""<<*nit<<
"Squared\" Expresion=\""<<*nit<<
"*"<<*nit<<
"\" Var=\""<<*nit<<rangeNames.at(rangeIndex)
584 <<
"\"> \nin your top-level XML's <Measurment> entry and use <NormFactor Name=\""<<*nit<<
"Squared\" in your channel XML file."<< endl;
589 if(!overallNorm_times_sigmaEpsilon.empty())
590 return overallNorm_times_sigmaEpsilon;
597 std::vector<OverallSys>& systList,
598 vector<string>& constraintTermNames,
599 vector<string>& totSystTermNames) {
605 totSystTermNames.push_back(prefix);
608 vector<double> lowVec, highVec;
610 std::map<std::string, double>::iterator itconstr;
611 for(
unsigned int i = 0; i < systList.size(); ++i) {
614 std::string strname = sys.
GetName();
615 const char *
name = strname.c_str();
619 std::cout <<
"HistoToWorkspaceFast::AddConstraintTerm - skip systematic " << sys.
GetName() << std::endl;
626 std::cout <<
"HistoToWorkspaceFast::AddConstraintTerm - zero uncertainty assigned - skip systematic " << sys.
GetName() << std::endl;
629 double tauVal = 1./(relerr*relerr);
630 double sqtau = 1./relerr;
643 alphaOfBeta->
Print(
"t");
645 constraintTermNames.push_back(gamma->
GetName());
649 const_cast<RooArgSet*
>(proto->
set(
"globalObservables"))->add(*yvar);
652 params.
add(*alphaOfBeta);
653 std::cout <<
"Added a gamma constraint for " << name << std::endl;
662 double gaussSigma = 1;
665 std::cout <<
"Added a uniform constraint for " << name <<
" as a gaussian constraint with a very large sigma " << std::endl;
676 constraintTermNames.push_back( gausConstraint->
GetName() );
677 proto->
var((
"nom_" + prefix + sys.
GetName()).c_str())->setConstant();
678 const_cast<RooArgSet*
>(proto->
set(
"globalObservables"))->add(*nomAlpha);
694 double tauVal = 1./relerr;
695 std::string tauName =
"tau_" + sys.
GetName();
697 double kappaVal = 1. + relerr;
698 std::string kappaName =
"kappa_" + sys.
GetName();
700 const char * alphaName = alpha->
GetName();
702 std::string alphaOfBetaName =
"alphaOfBeta_" + sys.
GetName();
704 tauName.c_str(),kappaName.c_str(),alphaName,
705 tauName.c_str(),kappaName.c_str(),alphaName ) );
707 std::cout <<
"Added a log-normal constraint for " << name << std::endl;
708 alphaOfBeta->
Print(
"t");
709 params.add(*alphaOfBeta);
714 double low = sys.
GetLow();
716 lowVec.push_back(low);
717 highVec.push_back(high);
721 if(systList.size() > 0){
725 assert( params.getSize() > 0);
726 assert(
int(lowVec.size()) == params.getSize() );
728 FlexibleInterpVar interp( (interpName).c_str(),
"", params, 1., lowVec, highVec);
729 interp.setAllInterpCodes(4);
746 vector<string>& syst_x_expectedPrefixNames,
747 vector<string>& normByNames){
757 double binWidth(1.0);
758 std::string obsNameVecStr;
759 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
761 std::string obsName = *itr;
762 binWidth *= proto->
var(obsName.c_str())->numBins()/(proto->
var(obsName.c_str())->getMax() - proto->
var(obsName.c_str())->getMin()) ;
763 if (obsNameVecStr.size()>0) { obsNameVecStr +=
"_"; }
764 obsNameVecStr += obsName;
768 for(
unsigned int j=0; j<syst_x_expectedPrefixNames.size();++j){
769 std::stringstream str;
773 command=string(
Form(
"binWidth_%s_%d[%e]",obsNameVecStr.c_str(),j,binWidth));
774 proto->
factory(command.c_str());
776 coeffList+=prepend+
"binWidth_"+obsNameVecStr+str.str();
778 command=
"prod::L_x_"+syst_x_expectedPrefixNames.at(j)+
"("+normByNames.at(j)+
","+syst_x_expectedPrefixNames.at(j)+
")";
780 proto->
factory(command.c_str());
781 shapeList+=prepend+
"L_x_"+syst_x_expectedPrefixNames.at(j);
790 proto->
defineSet(
"coefList",coeffList.c_str());
791 proto->
defineSet(
"shapeList",shapeList.c_str());
795 tot.specialIntegratorConfig(
kTRUE)->method2D().setLabel(
"RooBinIntegrator") ;
796 tot.specialIntegratorConfig(
kTRUE)->methodND().setLabel(
"RooBinIntegrator") ;
800 tot.setAttribute(
"GenerateBinned");
825 vector<string>& likelihoodTermNames){
830 for(
Int_t i=lowBin; i<highBin; ++i){
831 std::stringstream str;
834 string command(
"Poisson::"+prefix+str.str()+
"("+obsPrefix+str.str()+
","+expPrefix+str.str()+
",1)");
838 cout <<
"Poisson Term " << command << endl;
842 likelihoodTermNames.push_back( temp->GetName() );
851 TTree*
tree =
new TTree();
855 for(
Int_t i=lowBin; i<highBin; ++i){
856 std::stringstream str;
859 cout <<
"expected number of events called: " << expPrefix << endl;
865 cout <<
"setting obs"+str.str()+
" to expected = " << exp->
getVal() <<
" check: " << obs->
getVal() << endl;
868 obsForTree[i] = exp->
getVal();
869 tree->Branch((obsPrefix+str.str()).c_str(), obsForTree+i ,(obsPrefix+str.str()+
"/D").c_str());
872 cout <<
"problem retrieving obs or exp " << obsPrefix+str.str() << obs <<
" " << expPrefix+str.str() << exp << endl;
879 delete [] obsForTree;
890 map<string,double> gammaSyst,
891 map<string,double> uniformSyst,
892 map<string,double> logNormSyst,
893 map<string,double> noSyst) {
894 string pdfName(pdfNameChar);
897 if( combined_config==NULL ) {
898 std::cout <<
"Error: Failed to find object 'ModelConfig' in workspace: " 899 << proto->
GetName() << std::endl;
904 string edit=
"EDIT::newSimPdf("+pdfName+
",";
906 string lastPdf=pdfName;
908 unsigned int numReplacements = 0;
909 unsigned int nskipped = 0;
910 map<string,double>::iterator it;
914 for(it=gammaSyst.begin(); it!=gammaSyst.end(); ++it) {
916 if(! proto->
var((
"alpha_"+it->first).c_str())){
923 double relativeUncertainty = it->second;
924 double scale = 1/
sqrt((1+1/
pow(relativeUncertainty,2)));
927 proto->
factory(
Form(
"beta_%s[1,0,10]",it->first.c_str()));
928 proto->
factory(
Form(
"y_%s[%f]",it->first.c_str(),1./
pow(relativeUncertainty,2))) ;
929 proto->
factory(
Form(
"theta_%s[%f]",it->first.c_str(),
pow(relativeUncertainty,2))) ;
930 proto->
factory(
Form(
"Gamma::beta_%sConstraint(beta_%s,sum::k_%s(y_%s,one[1]),theta_%s,zero[0])",
935 it->first.c_str())) ;
951 proto->
factory(
Form(
"PolyVar::alphaOfBeta_%s(beta_%s,{%f,%f})",it->first.c_str(),it->first.c_str(),-1./scale,1./scale));
964 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
966 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
976 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
977 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
981 cout <<
"Going to issue this edit command\n" << edit<< endl;
982 proto->
factory( edit.c_str() );
985 cout <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
991 for(it=uniformSyst.begin(); it!=uniformSyst.end(); ++it) {
992 cout <<
"edit for " << it->first <<
"with rel uncert = " << it->second << endl;
993 if(! proto->
var((
"alpha_"+it->first).c_str())){
994 cout <<
"systematic not there" << endl;
1001 proto->
factory(
Form(
"beta_%s[1,0,10]",it->first.c_str()));
1002 proto->
factory(
Form(
"Uniform::beta_%sConstraint(beta_%s)",it->first.c_str(),it->first.c_str()));
1003 proto->
factory(
Form(
"PolyVar::alphaOfBeta_%s(beta_%s,{-1,1})",it->first.c_str(),it->first.c_str()));
1014 cout <<
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint" << endl;
1015 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
1017 cout <<
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first << endl;
1018 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
1020 if( proto->
pdf((
"alpha_"+it->first+
"Constraint").c_str()) && proto->
var((
"alpha_"+it->first).c_str()) )
1021 cout <<
" checked they are there" << proto->
pdf((
"alpha_"+it->first+
"Constraint").c_str()) <<
" " << proto->
var((
"alpha_"+it->first).c_str()) << endl;
1023 cout <<
"NOT THERE" << endl;
1026 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
1027 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
1031 cout << edit<< endl;
1032 proto->
factory( edit.c_str() );
1035 cout <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
1045 for(it=logNormSyst.begin(); it!=logNormSyst.end(); ++it) {
1046 cout <<
"edit for " << it->first <<
"with rel uncert = " << it->second << endl;
1047 if(! proto->
var((
"alpha_"+it->first).c_str())){
1048 cout <<
"systematic not there" << endl;
1054 double relativeUncertainty = it->second;
1055 double kappa = 1+relativeUncertainty;
1059 double scale = relativeUncertainty;
1062 const char * cname = it->first.c_str();
1069 cname, cname, cname, cname)) ;
1070 proto->
factory(
TString::Format(
"PolyVar::alphaOfBeta_%s(beta_%s,{%f,%f})",cname,cname,-1./scale,1./scale));
1082 cout <<
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint" << endl;
1083 editList+=precede +
"alpha_"+it->first+
"Constraint=beta_" + it->first+
"Constraint";
1085 cout <<
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first << endl;
1086 editList+=precede +
"alpha_"+it->first+
"=alphaOfBeta_"+ it->first;
1088 if( proto->
pdf((
"alpha_"+it->first+
"Constraint").c_str()) && proto->
var((
"alpha_"+it->first).c_str()) )
1089 cout <<
" checked they are there" << proto->
pdf((
"alpha_"+it->first+
"Constraint").c_str()) <<
" " << proto->
var((
"alpha_"+it->first).c_str()) << endl;
1091 cout <<
"NOT THERE" << endl;
1094 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
1095 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
1099 cout << edit<< endl;
1100 proto->
factory( edit.c_str() );
1103 cout <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
1107 const RooArgSet * gobs = proto->
set(
"globalObservables");
1111 proto->
defineSet(
"globalObservables",gobsNew);
1119 for(it=noSyst.begin(); it!=noSyst.end(); ++it) {
1121 cout <<
"remove constraint for parameter" << it->first << endl;
1122 if(! proto->
var((
"alpha_"+it->first).c_str()) || ! proto->
pdf((
"alpha_"+it->first+
"Constraint").c_str()) ) {
1123 cout <<
"systematic not there" << endl;
1130 if ( !proto->
var(
"one") ) { proto->
factory(
"one[1.0]"); }
1134 cout <<
"alpha_"+it->first+
"Constraint=one" << endl;
1135 editList+=precede +
"alpha_"+it->first+
"Constraint=one";
1139 if(numReplacements%10 == 0 && numReplacements+nskipped!=gammaSyst.size()){
1140 edit=
"EDIT::"+lastPdf+
"_("+lastPdf+
","+editList+
")";
1144 cout << edit << endl;
1145 proto->
factory( edit.c_str() );
1147 if(!newOne) { cout <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl; }
1154 edit=
"EDIT::newSimPdf("+lastPdf+
","+editList+
")";
1155 cout << edit<< endl;
1156 proto->
factory( edit.c_str() );
1161 combined_config->
SetPdf(*newOne);
1164 cout <<
"\n\n ---------------------\n WARNING: failed to make EDIT\n\n" << endl;
1171 FILE* covFile = fopen ((filename).c_str(),
"w");
1176 fprintf(covFile,
" ") ;
1179 fprintf(covFile,
" & %s", myargi->
GetName());
1181 fprintf(covFile,
"\\\\ \\hline \n" );
1185 fprintf(covFile,
"%s", myargi->
GetName());
1190 fprintf(covFile,
" & %.2f", result->
correlation(*myargi, *myargj));
1193 fprintf(covFile,
" \\\\\n");
1206 Error(
"MakeSingleChannelWorkspace",
"The input Channel does not contain any sample - return a nullptr");
1218 string channel_name = channel.
GetName();
1226 TH1* channel_hist_template = channel.
GetSamples().at(0).GetHisto();
1229 for (
unsigned int idx=0; idx<
fObsNameVec.size(); ++idx ) {
1240 cout <<
"\n\n-------------------\nStarting to process " << channel_name <<
" channel with " <<
fObsNameVec.size() <<
" observables" << endl;
1246 auto proto_config = make_unique<ModelConfig>(
"ModelConfig",
proto);
1247 proto_config->SetWorkspace(*proto);
1252 cout <<
"will preprocess this line: " << *funcIter <<endl;
1253 proto->
factory(funcIter->c_str());
1257 RooArgSet likelihoodTerms(
"likelihoodTerms"), constraintTerms(
"constraintTerms");
1258 vector<string> likelihoodTermNames, constraintTermNames, totSystTermNames, syst_x_expectedPrefixNames, normalizationNames;
1260 vector< pair<string,string> > statNamePairs;
1261 vector< pair<TH1*,TH1*> > statHistPairs;
1262 std::string statFuncName;
1263 std::string statNodeName;
1267 string prefix, range;
1272 std::stringstream lumiStr;
1275 proto->
factory((
"Lumi"+lumiStr.str()).c_str());
1276 cout <<
"lumi str = " << lumiStr.str() << endl;
1278 std::stringstream lumiErrorStr;
1280 proto->
factory((
"Gaussian::lumiConstraint(Lumi,"+lumiErrorStr.str()+
")").c_str());
1282 proto->
defineSet(
"globalObservables",
"nominalLumi");
1284 constraintTermNames.push_back(
"lumiConstraint");
1285 cout <<
"lumi Error str = " << lumiErrorStr.str() << endl;
1292 vector<Sample>::iterator it = channel.
GetSamples().begin();
1293 for(; it!=channel.
GetSamples().end(); ++it) {
1297 string overallSystName = sample.
GetName() +
"_" + channel_name +
"_epsilon";
1299 string systSourcePrefix =
"alpha_";
1307 overallSystName =
AddNormFactor(proto, channel_name, overallSystName, sample, doRatio);
1312 string syst_x_expectedPrefix =
"";
1332 cout << sample.
GetName() +
"_" + channel_name +
" has no variation histograms " << endl;
1333 string expPrefix = sample.
GetName() +
"_" + channel_name;
1334 syst_x_expectedPrefix = sample.
GetName() +
"_" + channel_name +
"_overallSyst_x_Exp";
1342 string constraintPrefix = sample.
GetName() +
"_" + channel_name +
"_Hist_alpha";
1343 syst_x_expectedPrefix = sample.
GetName() +
"_" + channel_name +
"_overallSyst_x_HistSyst";
1348 constraintPrefix, syst_x_expectedPrefix, overallSystName,
1349 constraintTermNames);
1359 std::cout <<
"Cannot include Stat Error for histograms of more than 3 dimensions." 1368 std::cout <<
"Sample: " << sample.
GetName() <<
" to be included in Stat Error " 1369 <<
"for channel " << channel_name
1399 string UncertName = syst_x_expectedPrefix +
"_StatAbsolUncert";
1400 TH1* statErrorHist = NULL;
1404 std::cout <<
"Making Statistical Uncertainty Hist for " 1405 <<
" Channel: " << channel_name
1406 <<
" Sample: " << sample.
GetName()
1417 std::cout <<
"Using external histogram for Stat Errors for " 1418 <<
" Channel: " << channel_name
1419 <<
" Sample: " << sample.
GetName()
1421 std::cout <<
"Error Histogram: " << statErrorHist->
GetName() << std::endl;
1424 statErrorHist->
Multiply( nominal );
1425 statErrorHist->
SetName( UncertName.c_str() );
1430 statHistPairs.push_back( pair<TH1*,TH1*>(nominal, statErrorHist) );
1447 statFuncName =
"mc_stat_" + channel_name;
1449 if( paramHist == NULL ) {
1454 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1455 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
1456 observables.
add( *proto->
var(itr->c_str()) );
1461 std::string ParamSetPrefix =
"gamma_stat_" + channel_name;
1465 ParamSetPrefix.c_str(),
1467 gammaMin, gammaMax);
1469 ParamHistFunc statUncertFunc(statFuncName.c_str(), statFuncName.c_str(),
1470 observables, statFactorParams );
1481 statNodeName = sample.
GetName() +
"_" + channel_name +
"_overallSyst_x_StatUncert";
1484 RooProduct nodeWithMcStat(statNodeName.c_str(), statNodeName.c_str(),
1492 syst_x_expectedPrefix = nodeWithMcStat.GetName();
1505 std::cout <<
"Cannot include Stat Error for histograms of more than 3 dimensions." 1510 std::cout <<
"Sample: " << sample.
GetName() <<
" in channel: " << channel_name
1511 <<
" to be include a ShapeFactor." 1514 std::vector<ParamHistFunc*> paramHistFuncList;
1515 std::vector<std::string> shapeFactorNameList;
1521 std::string funcName = channel_name +
"_" + shapeFactor.
GetName() +
"_shapeFactor";
1523 if( paramHist == NULL ) {
1526 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1527 for (
int idx=0; itr!=
fObsNameVec.end(); ++itr, ++idx ) {
1528 observables.
add( *proto->
var(itr->c_str()) );
1532 std::string funcParams =
"gamma_" + shapeFactor.
GetName();
1538 observables, 0, 1000);
1541 ParamHistFunc shapeFactorFunc( funcName.c_str(), funcName.c_str(),
1542 observables, shapeFactorParams );
1547 std::cout <<
"Setting Shape Factor: " << shapeFactor.
GetName()
1548 <<
" to have initial shape from hist: " 1551 shapeFactorFunc.setShape( initialShape );
1556 std::cout <<
"Setting Shape Factor: " << shapeFactor.
GetName()
1557 <<
" to be constant" << std::endl;
1558 shapeFactorFunc.setConstant(
true);
1566 paramHistFuncList.push_back(paramHist);
1567 shapeFactorNameList.push_back(funcName);
1576 std::string shapeFactorNodeName = syst_x_expectedPrefix;
1577 for(
unsigned int i=0; i < shapeFactorNameList.size(); ++i) {
1578 shapeFactorNodeName +=
"_x_" + shapeFactorNameList.at(i);
1583 for(
unsigned int i=0; i < paramHistFuncList.size(); ++i) {
1584 nodesForProduct.
add( *paramHistFuncList.at(i) );
1589 RooProduct nodeWithShapeFactor(shapeFactorNodeName.c_str(),
1590 shapeFactorNodeName.c_str(),
1598 syst_x_expectedPrefix = nodeWithShapeFactor.GetName();
1611 std::cout <<
"Cannot include Stat Error for histograms of more than 3 dimensions." 1617 std::vector<string> ShapeSysNames;
1632 std::cout <<
"Sample: " << sample.
GetName() <<
" in channel: " << channel_name
1633 <<
" to include a ShapeSys." << std::endl;
1635 std::string funcName = channel_name +
"_" + shapeSys.
GetName() +
"_ShapeSys";
1636 ShapeSysNames.push_back( funcName );
1638 if( paramHist == NULL ) {
1644 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1646 observables.
add( *proto->
var(itr->c_str()) );
1650 std::string funcParams =
"gamma_" + shapeSys.
GetName();
1653 observables, 0, 10);
1656 ParamHistFunc shapeFactorFunc( funcName.c_str(), funcName.c_str(),
1657 observables, shapeFactorParams );
1680 Double_t minShapeUncertainty = 0.0;
1682 *paramHist, shapeErrorHist,
1684 minShapeUncertainty);
1692 std::string NodeName = syst_x_expectedPrefix;
1695 ShapeSysForNode.
add( *expFunc );
1696 for(
unsigned int i = 0; i < ShapeSysNames.size(); ++i ) {
1697 std::string ShapeSysName = ShapeSysNames.at(i);
1698 ShapeSysForNode.
add( *proto->
function(ShapeSysName.c_str()) );
1699 NodeName = NodeName +
"_x_" + ShapeSysName;
1703 RooProduct nodeWithShapeFactor(NodeName.c_str(), NodeName.c_str(), ShapeSysForNode );
1709 syst_x_expectedPrefix = nodeWithShapeFactor.GetName();
1718 syst_x_expectedPrefixNames.push_back(syst_x_expectedPrefix);
1724 normalizationNames.push_back(
"Lumi" );
1727 TString lumiParamString;
1728 lumiParamString += measurement.
GetLumi();
1729 lumiParamString.ReplaceAll(
' ', TString());
1730 normalizationNames.push_back(lumiParamString.Data());
1738 if( statHistPairs.size() > 0 ) {
1743 if( fracStatError == NULL ) {
1744 std::cout <<
"Error: Failed to make ScaledUncertaintyHist for: " 1745 << statNodeName << std::endl;
1752 std::cout <<
"About to create Constraint Terms from: " 1753 << chanStatUncertFunc->
GetName()
1754 <<
" params: " << chanStatUncertFunc->
paramList()
1765 std::cout <<
"Using Gaussian StatErrors in channel: " << channel.
GetName() << std::endl;
1768 std::cout <<
"Using Poisson StatErrors in channel: " << channel.
GetName() << std::endl;
1773 *chanStatUncertFunc, fracStatError,
1775 statRelErrorThreshold);
1779 for (
unsigned int i = 0; i < statHistPairs.size() ; ++i )
1780 delete statHistPairs[i].
second;
1782 statHistPairs.clear();
1784 delete fracStatError;
1794 syst_x_expectedPrefixNames, normalizationNames);
1795 likelihoodTermNames.push_back(channel_name+
"_model");
1799 for(
unsigned int i=0; i<systToFix.size(); ++i){
1807 if(systToFix.at(i)==
"Lumi"){
1808 auxMeas = proto->
var(
"nominalLumi");
1814 const_cast<RooArgSet*
>(proto->
set(
"globalObservables"))->
remove(*auxMeas);
1816 cout <<
"could not corresponding auxiliary measurement " 1820 cout <<
"could not find variable " << systToFix.at(i)
1821 <<
" could not set it to constant" << endl;
1827 for(
unsigned int i=0; i<constraintTermNames.size(); ++i){
1828 RooAbsArg* proto_arg = (proto->
arg(constraintTermNames[i].c_str()));
1829 if( proto_arg==NULL ) {
1830 std::cout <<
"Error: Cannot find arg set: " << constraintTermNames.at(i)
1831 <<
" in workspace: " << proto->
GetName() << std::endl;
1834 constraintTerms.
add( *proto_arg );
1837 for(
unsigned int i=0; i<likelihoodTermNames.size(); ++i){
1838 RooAbsArg* proto_arg = (proto->
arg(likelihoodTermNames[i].c_str()));
1839 if( proto_arg==NULL ) {
1840 std::cout <<
"Error: Cannot find arg set: " << likelihoodTermNames.at(i)
1841 <<
" in workspace: " << proto->
GetName() << std::endl;
1844 likelihoodTerms.add( *proto_arg );
1846 proto->
defineSet(
"constraintTerms",constraintTerms);
1847 proto->
defineSet(
"likelihoodTerms",likelihoodTerms);
1852 std::string observablesStr;
1854 std::vector<std::string>::iterator itr =
fObsNameVec.begin();
1856 observables.
add( *proto->
var(itr->c_str()) );
1857 if (!observablesStr.empty()) { observablesStr +=
","; }
1858 observablesStr += *itr;
1869 cout <<
"-----------------------------------------"<<endl;
1870 cout <<
"import model into workspace" << endl;
1872 auto model = make_unique<RooProdPdf>(
1873 (
"model_"+channel_name).c_str(),
1874 "product of Poissons accross bins for a single channel",
1875 constraintTerms,
Conditional(likelihoodTerms,observables));
1878 proto_config->SetPdf(*
model);
1879 proto_config->SetObservables(observables);
1880 proto_config->SetGlobalObservables(*proto->
set(
"globalObservables"));
1890 const char* weightName=
"weightVar";
1900 RooDataSet* asimovDataUnbinned = new RooDataSet("asimovData","",*proto->set("obsAndWeight"),weightName);
1901 for(int i=0; i<asimov_data->numEntries(); ++i){
1902 asimov_data->get(i)->Print("v");
1903 //cout << "GREPME : " << i << " " << data->weight() <<endl;
1904 asimovDataUnbinned->add( *asimov_data->get(i), asimov_data->weight() );
1906 proto->import(*asimovDataUnbinned);
1912 proto->
import(dynamic_cast<RooDataSet&>(*asimov_dataset),
Rename(
"asimovData"));
1920 std::cout <<
"Error: Data histogram for channel: " << channel.
GetName()
1921 <<
" is NULL" << std::endl;
1926 auto obsDataUnbinned = make_unique<RooDataSet>(
"obsData",
"",*proto->
set(
"obsAndWeight"),weightName);
1965 proto->import(*obsDataUnbinned);
1972 std::string dataName = data.
GetName();
1975 std::cout <<
"Error: Additional Data histogram for channel: " << channel.
GetName()
1976 <<
" with name: " << dataName <<
" is NULL" << std::endl;
1981 auto obsDataUnbinned = make_unique<RooDataSet>(dataName.c_str(), dataName.c_str(),
1982 *proto->
set(
"obsAndWeight"), weightName);
2020 proto->import(*obsDataUnbinned);
2032 std::vector<std::string> ObsNameVec) {
2038 if (ObsNameVec.empty() ) {
2039 Error(
"ConfigureHistFactoryDataset",
"Invalid input - return");
2049 for (
int i=1; i<=ax->
GetNbins(); ++i) {
2052 proto->
var( ObsNameVec[0].c_str() )->
setVal( xval );
2054 if(ObsNameVec.size()==1) {
2056 obsDataUnbinned->
add( *proto->
set(
"obsAndWeight"), fval );
2059 for(
int j=1; j<=ay->
GetNbins(); ++j) {
2061 proto->
var( ObsNameVec[1].c_str() )->
setVal( yval );
2063 if(ObsNameVec.size()==2) {
2065 obsDataUnbinned->
add( *proto->
set(
"obsAndWeight"), fval );
2068 for(
int k=1; k<=az->
GetNbins(); ++k) {
2070 proto->
var( ObsNameVec[2].c_str() )->
setVal( zval );
2072 obsDataUnbinned->
add( *proto->
set(
"obsAndWeight"), fval );
2085 unsigned int histndim(1);
2086 std::string classname = hist->
ClassName();
2087 if (classname.find(
"TH1")==0) { histndim=1; }
2088 else if (classname.find(
"TH2")==0) { histndim=2; }
2089 else if (classname.find(
"TH3")==0) { histndim=3; }
2091 for (
unsigned int idx=0; idx<histndim; ++idx ) {
2104 if (ch_names.empty() || chs.empty() ) {
2105 Error(
"MakeCombinedModel",
"Input vectors are empty - return a nullptr");
2108 if (chs.size() < ch_names.size() ) {
2109 Error(
"MakeCombinedModel",
"Input vector of workspace has an invalid size - return a nullptr");
2117 map<string, RooAbsPdf*> pdfMap;
2118 vector<RooAbsPdf*> models;
2122 for(
unsigned int i = 0; i< ch_names.size(); ++i){
2126 cout <<
"full list of observables:"<<endl;
2130 for(
unsigned int i = 0; i< ch_names.size(); ++i){
2131 string channel_name=ch_names[i];
2133 if (ss.str().empty()) ss << channel_name ;
2134 else ss <<
',' << channel_name ;
2138 if(!model) cout <<
"failed to find model for channel"<<endl;
2140 models.push_back(model);
2141 globalObs.
add(*ch->
set(
"globalObservables"));
2144 pdfMap[channel_name]=
model;
2148 cout <<
"\n\n------------------\n Entering combination" << endl;
2159 combined->
import(globalObs);
2160 combined->
defineSet(
"globalObservables",globalObs);
2166 cout <<
"-----------------------------------------"<<endl;
2167 cout <<
"create toy data for " << ss.str() << endl;
2173 combined->
factory(
"weightVar[0,-1e10,1e10]");
2174 obsList.
add(*combined->
var(
"weightVar"));
2195 if( asimov_combined ) {
2196 combined->
import( *asimov_combined,
Rename(
"asimovData"));
2199 std::cout <<
"Error: Failed to create combined asimov dataset" << std::endl;
2202 delete asimov_combined;
2205 if(chs[0]->
data(
"obsData") != NULL) {
2206 MergeDataSets(combined, chs, ch_names,
"obsData", obsList, channelCat);
2254 obsList.
add(*channelCat);
2255 combined->
defineSet(
"observables",obsList);
2260 cout <<
"\n\n----------------\n Importing combined model" << endl;
2266 std::map< std::string, double>::iterator param_itr =
fParamValues.begin();
2269 std::string paramName = param_itr->first;
2270 double paramVal = param_itr->second;
2274 temp->
setVal( paramVal );
2275 cout <<
"setting " << paramName <<
" to the value: " << paramVal << endl;
2277 cout <<
"could not find variable " << paramName <<
" could not set its value" << endl;
2281 for(
unsigned int i=0; i<
fSystToFix.size(); ++i){
2286 cout <<
"setting " <<
fSystToFix.at(i) <<
" constant" << endl;
2288 cout <<
"could not find variable " <<
fSystToFix.at(i) <<
" could not set it to constant" << endl;
2296 combined_config->
SetPdf(*simPdf);
2299 combined->
import(*combined_config,combined_config->
GetName());
2304 delete combined_config;
2312 std::vector<RooWorkspace*> wspace_vec,
2313 std::vector<std::string> channel_names,
2314 std::string dataSetName,
2323 for(
unsigned int i = 0; i< channel_names.size(); ++i){
2326 std::cout <<
"Merging data for channel " << channel_names[i].c_str() << std::endl;
2328 if( !obsDataInChannel ) {
2329 std::cout <<
"Error: Can't find DataSet: " << dataSetName
2330 <<
" in channel: " << channel_names.at(i)
2337 obsList,
Index(*channelCat),
2339 Import(channel_names[i].c_str(),*obsDataInChannel));
2341 simData->
append(*tempData);
2352 combined->
import(*simData,
Rename(dataSetName.c_str()));
2357 std::cout <<
"Error: Unable to merge observable datasets" << std::endl;
2372 TH1* ErrorHist = (
TH1*) Nominal->
Clone( Name.c_str() );
2376 Int_t binNumber = 0;
2379 for(
Int_t i_bin = 0; i_bin < numBins; ++i_bin) {
2390 if( histError != histError ) {
2391 std::cout <<
"Warning: In histogram " << Nominal->
GetName()
2392 <<
" bin error for bin " << i_bin
2393 <<
" is NAN. Not using Error!!!" 2401 if( histError < 0 ) {
2402 std::cout <<
"Warning: In histogram " << Nominal->
GetName()
2403 <<
" bin error for bin " << binNumber
2404 <<
" is < 0. Setting Error to 0" 2430 unsigned int numHists = HistVec.size();
2432 if( numHists == 0 ) {
2433 std::cout <<
"Warning: Empty Hist Vector, cannot create total uncertainty" << std::endl;
2437 TH1* HistTemplate = HistVec.at(0).first;
2442 for(
unsigned int i = 0; i < HistVec.size(); ++i ) {
2444 TH1* nominal = HistVec.at(i).first;
2445 TH1* error = HistVec.at(i).second;
2448 std::cout <<
"Error: Provided hists have unequal bins" << std::endl;
2452 std::cout <<
"Error: Provided hists have unequal bins" << std::endl;
2457 std::vector<double> TotalBinContent( numBins, 0.0);
2458 std::vector<double> HistErrorsSqr( numBins, 0.0);
2460 Int_t binNumber = 0;
2463 for(
Int_t i_bins = 0; i_bins < numBins; ++i_bins) {
2470 for(
unsigned int i_hist = 0; i_hist < numHists; ++i_hist ) {
2472 TH1* nominal = HistVec.at(i_hist).first;
2473 TH1* error = HistVec.at(i_hist).second;
2488 if( histError != histError ) {
2489 std::cout <<
"Warning: In histogram " << error->
GetName()
2490 <<
" bin error for bin " << binNumber
2491 <<
" is NAN. Not using error!!" 2497 TotalBinContent.at(i_bins) += histValue;
2498 HistErrorsSqr.at(i_bins) += histError*histError;
2506 TH1* ErrorHist = (
TH1*) HistTemplate->
Clone( Name.c_str() );
2510 for(
Int_t i = 0; i < numBins; ++i) {
2518 Double_t ErrorsSqr = HistErrorsSqr.at(i);
2519 Double_t TotalVal = TotalBinContent.at(i);
2521 if( TotalVal <= 0 ) {
2522 std::cout <<
"Warning: Sum of histograms for bin: " << binNumber
2523 <<
" is <= 0. Setting error to 0" 2534 if( RelativeError != RelativeError ) {
2535 std::cout <<
"Error: bin " << i <<
" error is NAN" << std::endl;
2536 std::cout <<
" HistErrorsSqr: " << ErrorsSqr
2537 <<
" TotalVal: " << TotalVal
2548 std::cout <<
"Making Total Uncertainty for bin " << binNumber
2549 <<
" Error = " <<
sqrt(ErrorsSqr)
2550 <<
" Val = " << TotalVal
2551 <<
" RelativeError = " << RelativeError
2593 if( numBins != numParams ) {
2594 std::cout <<
"Error: In createStatConstraintTerms, encountered bad number of bins" << std::endl;
2595 std::cout <<
"Given histogram with " << numBins <<
" bins," 2596 <<
" but require exactly " << numParams << std::endl;
2600 Int_t TH1BinNumber = 0;
2611 std::cout <<
"Creating constraint for: " << gamma.
GetName()
2612 <<
". Type of constraint: " << type << std::endl;
2621 std::cout <<
"Not creating constraint term for " 2623 <<
" because sigma = " << sigma
2625 <<
" (TH1 bin number = " << TH1BinNumber <<
")" 2632 gamma.
setMax( 1 + 5*sigma );
2637 std::string constrName = string(gamma.
GetName()) +
"_constraint";
2638 std::string nomName = string(
"nom_") + gamma.
GetName();
2639 std::string sigmaName = string(gamma.
GetName()) +
"_sigma";
2640 std::string poisMeanName = string(gamma.
GetName()) +
"_poisMean";
2653 RooRealVar constrNom(nomName.c_str(), nomName.c_str(), 1.0,0,10);
2658 constrNom,
gamma, constrSigma );
2668 RooRealVar constrNom(nomName.c_str(), nomName.c_str(),
tau);
2670 constrNom.setConstant(
true );
2673 std::string scalingName = string(gamma.
GetName()) +
"_tau";
2674 RooConstVar poissonScaling( scalingName.c_str(), scalingName.c_str(),
tau);
2677 RooProduct constrMean( poisMeanName.c_str(), poisMeanName.c_str(),
RooArgSet(gamma, poissonScaling) );
2682 RooPoisson pois(constrName.c_str(), constrName.c_str(), constrNom, constrMean);
2688 std::cout <<
"Error: Did not recognize Stat Error constraint term type: " 2689 << type <<
" for : " << paramHist.
GetName() << std::endl;
2696 if( sigma < minSigma ) {
2697 std::cout <<
"Warning: Bin " << i <<
" = " << sigma
2698 <<
" and is < " << minSigma
2699 <<
". Setting: " << gamma.
GetName() <<
" to constant" 2704 constraintTermNames.push_back( constrName );
2705 ConstraintTerms.
add( *proto->
pdf(constrName.c_str()) );
2711 RooRealVar* nomVarInWorkspace = proto->
var(nomName.c_str());
2712 if( ! globalSet->
contains(*nomVarInWorkspace) ) {
2713 globalSet->
add( *nomVarInWorkspace );
2718 return ConstraintTerms;
std::vector< RooStats::HistFactory::Asimov > & GetAsimovDatasets()
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...
virtual const char * GetName() const
Returns name of object.
Bool_t saveSnapshot(const char *name, const char *paramNames)
Save snapshot of values and attributes (including "Constant") of parameters 'params' If importValues ...
RooStats::HistFactory::StatError & GetStatError()
TIterator * createIterator(Bool_t dir=kIterForward) const
RooWorkspace * MakeCombinedModel(std::vector< std::string >, std::vector< RooWorkspace *>)
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
std::map< std::string, double > & GetNoSyst()
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Bool_t IsBinUnderflow(Int_t bin, Int_t axis=0) const
Return true if the bin is underflow.
const RooArgSet * GetObservables() const
get RooArgSet for observables (return NULL if not existing)
void AddMultiVarGaussConstraint(RooWorkspace *proto, std::string prefix, int lowBin, int highBin, std::vector< std::string > &likelihoodTermNames)
const RooArgList & paramList() const
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
std::vector< RooStats::HistFactory::NormFactor > & GetNormFactorList()
void AddPoissonTerms(RooWorkspace *proto, std::string prefix, std::string obsPrefix, std::string expPrefix, int lowBin, int highBin, std::vector< std::string > &likelihoodTermNames)
std::map< std::string, double > fParamValues
Bool_t IsBinOverflow(Int_t bin, Int_t axis=0) const
Return true if the bin is overflow.
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)
virtual void SetWorkspace(RooWorkspace &ws)
Double_t getVal(const RooArgSet *set=0) const
void GuessObsNameVec(TH1 *hist)
void GuessObsAndNuisance(const RooAbsData &data)
guesses Observables and ParametersOfInterest if not already set
Constraint::Type GetConstraintType()
RooWorkspace * MakeSingleChannelWorkspace(Measurement &measurement, Channel &channel)
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
RooDataSet * MergeDataSets(RooWorkspace *combined, std::vector< RooWorkspace *> wspace_vec, std::vector< std::string > channel_names, std::string dataSetName, RooArgList obsList, RooCategory *channelCat)
void ConfigureWorkspace(RooWorkspace *)
virtual Int_t GetNbinsZ() const
RooCmdArg Conditional(const RooArgSet &pdfSet, const RooArgSet &depSet, Bool_t depsAreCond=kFALSE)
std::vector< std::string > & GetPOIList()
RooAbsReal * function(const char *name) const
Retrieve function (RooAbsReal) with given name. Note that all RooAbsPdfs are also RooAbsReals...
bool GetNormalizeByTheory()
RooAbsArg * arg(const char *name) const
Return RooAbsArg with given name. A null pointer is returned if none is found.
Class RooRealSumPdf implements a PDF constructed from a sum of functions:
RooCmdArg RecycleConflictNodes(Bool_t flag=kTRUE)
virtual void SetObservables(const RooArgSet &set)
specify the observables
std::vector< RooStats::HistFactory::HistoSys > & GetHistoSysList()
void setMax(const char *name, Double_t value)
Set maximum of name range to given value.
static void ConfigureWorkspaceForMeasurement(const std::string &ModelName, RooWorkspace *ws_single, Measurement &measurement)
std::string AddNormFactor(RooWorkspace *proto, std::string &channel, std::string &sigmaEpsilon, Sample &sample, bool doRatio)
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
double beta(double x, double y)
Calculates the beta function.
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE)
Set value by specifying the name of the desired state If printError is set, a message will be printed...
Bool_t importClassCode(const char *pat="*", Bool_t doReplace=kFALSE)
Inport code of all classes in the workspace that have a class name that matches pattern 'pat' and whi...
std::vector< std::string > GetPreprocessFunctions()
void SetFunctionsToPreprocess(std::vector< std::string > lines)
virtual const char * ClassName() const
Returns name of class to which the object belongs.
double GetRelErrorThreshold()
RooDataSet is a container class to hold N-dimensional binned data.
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t tau
RooHistFunc implements a real-valued function sampled from a multidimensional histogram.
void ProcessExpectedHisto(TH1 *hist, RooWorkspace *proto, std::string prefix, std::string productPrefix, std::string systTerm)
double pow(double, double)
std::vector< RooStats::HistFactory::Channel > & GetChannels()
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
std::vector< std::string > fObsNameVec
virtual void SetPdf(const RooAbsPdf &pdf)
Set the Pdf, add to the the workspace if not already there.
RooConstVar represent a constant real-valued object.
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)
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
static RooArgList createParamSet(RooWorkspace &w, const std::string &, const RooArgList &Vars)
std::map< std::string, double > & GetUniformSyst()
static constexpr double second
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
void LinInterpWithConstraint(RooWorkspace *proto, TH1 *nominal, std::vector< HistoSys >, std::string prefix, std::string productPrefix, std::string systTerm, std::vector< std::string > &likelihoodTermNames)
RooRealVar represents a fundamental (non-derived) real valued object.
RooAbsData * data(const char *name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found...
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Class to manage histogram axis.
void setConstant(Bool_t value=kTRUE)
void setMin(const char *name, Double_t value)
Set minimum of name range to given value.
Bool_t defineSet(const char *name, const RooArgSet &aset, Bool_t importMissing=kFALSE)
Define a named RooArgSet with given constituents.
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...
std::vector< RooStats::HistFactory::ShapeSys > & GetShapeSysList()
virtual Bool_t Multiply(TF1 *f1, Double_t c1=1)
Performs the operation:
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
char * Form(const char *fmt,...)
RooAbsData is the common abstract base class for binned and unbinned datasets.
Bool_t loadSnapshot(const char *name)
Load the values and attributes of the parameters in the snapshot saved with the given name...
void setNoRounding(bool flag=kTRUE)
RooCmdArg Rename(const char *suffix)
RooArgSet allVars() const
Return set with all variable objects.
void SetObsToExpected(RooWorkspace *proto, std::string obsPrefix, std::string expPrefix, int lowBin, int highBin)
RooDataSet is a container class to hold unbinned data.
RooCategory represents a fundamental (non-derived) discrete value object.
std::vector< RooStats::HistFactory::OverallSys > & GetOverallSysList()
RooProduct a RooAbsReal implementation that represent the product of a given set of other RooAbsReal ...
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the 'data' argset, to the data set...
TObject * obj(const char *name) const
Return any type of object (RooAbsArg, RooAbsData or generic object) with given name) ...
std::vector< std::string > & GetConstantParams()
RooCmdArg Import(const char *state, TH1 &histo)
virtual void SetName(const char *name)
Change the name of this histogram.
RooWorkspace * MakeSingleChannelModel(Measurement &measurement, Channel &channel)
Multivariate Gaussian p.d.f.
HistoToWorkspaceFactoryFast()
RooCmdArg Index(RooCategory &icat)
static constexpr double gauss
Namespace for the RooStats classes.
HistFactory::StatErrorConfig & GetStatErrorConfig()
Constraint::Type GetConstraintType()
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Bool_t removeSet(const char *name)
Remove a named set from the workspace.
void append(RooDataSet &data)
Add all data points of given data set to this data set.
virtual ~HistoToWorkspaceFactoryFast()
RooCmdArg WeightVar(const char *name, Bool_t reinterpretAsWeight=kFALSE)
std::map< std::string, double > & GetLogNormSyst()
std::vector< RooStats::HistFactory::Sample > & GetSamples()
static void PrintCovarianceMatrix(RooFitResult *result, RooArgSet *params, std::string filename)
RooRealVar * var(const char *name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found...
std::vector< RooStats::HistFactory::Data > & GetAdditionalData()
RooFactoryWSTool & factory()
Return instance to factory tool.
std::vector< std::string > fPreprocessFunctions
TH1 * MakeScaledUncertaintyHist(const std::string &Name, std::vector< std::pair< TH1 *, TH1 *> > HistVec)
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Bool_t 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.
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...
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
void setPositiveDefinite(bool flag=true)
virtual void SetParametersOfInterest(const RooArgSet &set)
std::vector< std::string > fSystToFix
Bool_t contains(const RooAbsArg &var) const
virtual Int_t GetNbinsX() const
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
virtual void SetGlobalObservables(const RooArgSet &set)
specify the global observables
virtual const char * GetName() const
Returns name of object.
TH1 * MakeAbsolUncertaintyHist(const std::string &Name, const TH1 *Hist)
void Print(Option_t *opts=0) const
Print contents of the workspace.
RooStats::HistFactory::Data & GetData()
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset...
void ConfigureHistFactoryDataset(RooDataSet *obsData, TH1 *nominal, RooWorkspace *proto, std::vector< std::string > obsNameVec)
void MakeTotalExpected(RooWorkspace *proto, std::string totName, std::vector< std::string > &syst_x_expectedPrefixNames, std::vector< std::string > &normByNames)
RooArgList createStatConstraintTerms(RooWorkspace *proto, std::vector< std::string > &constraintTerms, ParamHistFunc ¶mHist, TH1 *uncertHist, Constraint::Type type, Double_t minSigma)
Double_t correlation(const RooAbsArg &par1, const RooAbsArg &par2) const
std::vector< RooStats::HistFactory::ShapeFactor > & GetShapeFactorList()
Bool_t isConstant() const
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
The RooWorkspace is a persistable container for RooFit projects.
virtual Int_t GetNbinsY() const
std::map< std::string, double > & GetGammaSyst()
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.