struct HypoTestOptions {
bool noSystematics = false;
double nToysRatio = 4;
double poiValue = -1;
int printLevel=0;
bool generateBinned = false;
bool useProof = false;
bool enableDetailedOutput = false;
};
HypoTestOptions optHT;
void StandardHypoTestDemo(const char* infile = "",
const char* workspaceName = "combined",
const char* modelSBName = "ModelConfig",
const char* modelBName = "",
const char* dataName = "obsData",
int calcType = 0,
int testStatType = 3,
int ntoys = 5000,
bool useNC = false,
const char * nuisPriorName = 0)
{
bool noSystematics = optHT.noSystematics;
double nToysRatio = optHT.nToysRatio;
double poiValue = optHT.poiValue;
int printLevel = optHT.printLevel;
bool generateBinned = optHT.generateBinned;
bool useProof = optHT.useProof;
bool enableDetOutput = optHT.enableDetailedOutput;
SimpleLikelihoodRatioTestStat::SetAlwaysReuseNLL(true);
ProfileLikelihoodTestStat::SetAlwaysReuseNLL(true);
RatioOfProfiledLikelihoodsTestStat::SetAlwaysReuseNLL(true);
const char* filename = "";
if (!strcmp(infile,"")) {
filename = "results/example_combined_GaussExample_model.root";
if (!fileExist) {
#ifdef _WIN32
cout << "HistFactory file cannot be generated on Windows - exit" << endl;
return;
#endif
cout <<"will run standard hist2workspace example"<<endl;
gROOT->ProcessLine(
".! prepareHistFactory .");
gROOT->ProcessLine(
".! hist2workspace config/example.xml");
cout <<"\n\n---------------------"<<endl;
cout <<"Done creating example input"<<endl;
cout <<"---------------------\n\n"<<endl;
}
}
else
filename = infile;
if(!file ){
cout <<"StandardRooStatsDemoMacro: Input file " << filename << " is not found" << endl;
return;
}
if(!w){
cout <<"workspace not found" << endl;
return;
}
if(!data || !sbModel){
cout << "data or ModelConfig was not found" <<endl;
return;
}
if (noSystematics) {
if (nuisPar && nuisPar->
getSize() > 0) {
std::cout << "StandardHypoTestInvDemo" << " - Switch off all systematics by setting them constant to their initial values" << std::endl;
}
if (bModel) {
if (bnuisPar)
}
}
if (!bModel ) {
Info(
"StandardHypoTestInvDemo",
"The background model %s does not exist",modelBName);
Info(
"StandardHypoTestInvDemo",
"Copy it from ModelConfig %s and set POI to zero",modelSBName);
bModel->
SetName(TString(modelSBName)+TString(
"B_only"));
if (!var) return;
double oldval = var->
getVal();
}
Info(
"StandardHypoTestDemo",
"Model %s has no snapshot - make one using model poi",modelSBName);
if (!var) return;
double oldval = var->
getVal();
if (poiValue > 0) var->
setVal(poiValue);
if (poiValue > 0) var->
setVal(oldval);
}
if (enableDetOutput) {
ropl->EnableDetailedOutput();
}
AsymptoticCalculator::SetPrintLevel(printLevel);
else if (calcType == 1) hypoCalc=
new HybridCalculator(*data, *sbModel, *bModel);
if (calcType == 0) {
}
if (calcType == 1) {
}
if (calcType == 2 ) {
if (testStatType != 2 && testStatType != 3)
Warning(
"StandardHypoTestDemo",
"Only the PL test statistic can be used with AsymptoticCalculator - use by default a two-sided PL");
}
if (nuisPriorName) nuisPdf = w->
pdf(nuisPriorName);
if (!nuisPdf) {
Info(
"StandardHypoTestDemo",
"No nuisance pdf given for the HybridCalculator - try to deduce pdf from the model");
else
}
if (!nuisPdf ) {
Info(
"StandardHypoTestDemo",
"No nuisance pdf given - try to use %s that is defined as a prior pdf in the B model",nuisPdf->
GetName());
}
else {
Error(
"StandardHypoTestDemo",
"Cannot run Hybrid calculator because no prior on the nuisance parameter is specified or can be derived");
return;
}
}
assert(nuisPdf);
Info(
"StandardHypoTestDemo",
"Using as nuisance Pdf ... " );
Warning(
"StandardHypoTestDemo",
"Prior nuisance does not depend on nuisance parameters. They will be smeared in their full range");
}
delete np;
}
if (sampler && (calcType == 0 || calcType == 1) ) {
if (useNC)
Warning(
"StandardHypoTestDemo",
"Pdf is extended: but number counting flag is set: ignore it ");
}
else {
if (!useNC) {
Info(
"StandardHypoTestDemo",
"Pdf is not extended: number of events to generate taken from observed data set is %d",nEvents);
}
else {
Info(
"StandardHypoTestDemo",
"using a number counting pdf");
}
}
Info(
"StandardHypoTestDemo",
"Data set is weighted, nentries = %d and sum of weights = %8.1f but toy generation is unbinned - it would be faster to set generateBinned to true\n",data->
numEntries(), data->
sumEntries());
}
if (useProof) {
}
}
delete sampler;
delete slrts;
delete ropl;
delete profll;
if (calcType != 2) {
}
else {
std::cout << "Asymptotic results " << std::endl;
}
if (calcType != 2) {
htExp.Append(htr);
double p[5];
double q[5];
for (int i = 0; i < 5; ++i) {
double sig = -2 + i;
}
for (int i = 0; i < 5; ++i) {
htExp.SetTestStatisticData( q[i] );
double sig = -2 + i;
std::cout << " Expected p -value and significance at " << sig << " sigma = "
<< htExp.NullPValue() << " significance " << htExp.Significance() << " sigma " << std::endl;
}
}
else {
for (int i = 0; i < 5; ++i) {
double sig = -2 + i;
std::cout << " Expected p -value and significance at " << sig << " sigma = "
}
}
bool writeResult = (calcType != 2);
if (enableDetOutput) {
writeResult=true;
Info(
"StandardHypoTestDemo",
"Detailed output will be written in output result file");
}
if (htr !=
NULL && writeResult) {
const char * calcTypeName = (calcType == 0) ? "Freq" : (calcType == 1) ? "Hybr" : "Asym";
TString resultFileName =
TString::Format(
"%s_HypoTest_ts%d_",calcTypeName,testStatType);
TString name = infile;
name.Replace(0, name.Last('/')+1, "");
TFile * fileOut = new TFile(resultFileName,"RECREATE");
Info(
"StandardHypoTestDemo",
"HypoTestResult has been written in the file %s",resultFileName.Data());
fileOut->Close();
}
}