| 1 |
#include <iostream>
|
| 2 |
#include <vector>
|
| 3 |
#include <sys/stat.h>
|
| 4 |
#include <fstream>
|
| 5 |
|
| 6 |
#include "../../Plotting/Modules/LimitDroplet.C"
|
| 7 |
#include "../../Plotting/Modules/GeneralToolBox.C"
|
| 8 |
#include "../../Plotting/Modules/SampleClass.C"
|
| 9 |
#include "../../Plotting/Modules/Setup.C"
|
| 10 |
|
| 11 |
#include <TSystem.h>
|
| 12 |
#include <TROOT.h>
|
| 13 |
|
| 14 |
/*
|
| 15 |
|
| 16 |
This capsule computes the limits; the advantage is that in this way the capsule can fail (segfault/...) without problems, while the main program is not affected.
|
| 17 |
|
| 18 |
*/
|
| 19 |
|
| 20 |
|
| 21 |
bool isThreadActive=false;
|
| 22 |
vector<string> all_args;
|
| 23 |
void* do_limit_computation( void *arg );
|
| 24 |
|
| 25 |
void* hello_thread( void *arg ) {
|
| 26 |
printf( "hello " );
|
| 27 |
return( 0 );
|
| 28 |
}
|
| 29 |
|
| 30 |
string threadoutputfile;
|
| 31 |
|
| 32 |
|
| 33 |
bool log_is_freaking_out() {
|
| 34 |
//this function basically checks if the last ten lines are mostly "WARNING:Generation -- skip: -nan" (i.e. problematic)
|
| 35 |
// it returns true if there are too many such errors, false if everything is ok.
|
| 36 |
return false;
|
| 37 |
/* ifstream ifs(threadoutputfile.c_str());
|
| 38 |
string temp;
|
| 39 |
int nlines=0;
|
| 40 |
while( getline( ifs, temp ) ) nlines++;
|
| 41 |
|
| 42 |
cout << "Going to read from file : " << threadoutputfile << " containing " << nlines << " lines " << endl;
|
| 43 |
int nwarnings=0;
|
| 44 |
int iline=0;
|
| 45 |
ifstream ifs2(threadoutputfile.c_str());
|
| 46 |
while( getline( ifs2, temp ) ) {
|
| 47 |
iline++;
|
| 48 |
if(nlines-iline>20) continue;
|
| 49 |
if(Contains(temp,"WARNING:Generation")&&Contains(temp,"nan")) nwarnings++;
|
| 50 |
}
|
| 51 |
cout << "Have found " << nwarnings << " warnings in the log ... " << endl;
|
| 52 |
|
| 53 |
|
| 54 |
if(nwarnings>10) return true;
|
| 55 |
else return false;
|
| 56 |
*/
|
| 57 |
}
|
| 58 |
|
| 59 |
bool CapsuleWrapper() {
|
| 60 |
pthread_t limitthread;
|
| 61 |
long t=0;
|
| 62 |
//pthread_create( &limitthread, NULL, (void*) do_limit_computation, NULL);
|
| 63 |
pthread_create( &limitthread, NULL, do_limit_computation, (void*) NULL );
|
| 64 |
int counter=0;
|
| 65 |
int counterinterval=10;
|
| 66 |
sleep(1); //waiting a second for the process to become active
|
| 67 |
while(counter<PlottingSetup::limitpatience*60 && isThreadActive) {
|
| 68 |
std::cout << "Limits are being calculated; Checking round " << counter/counterinterval << " ( corresponds to " << seconds_to_time(counter) << " ) , patience will end in " << seconds_to_time(60*PlottingSetup::limitpatience-counter) << std::endl;
|
| 69 |
counter+=counterinterval;
|
| 70 |
if(log_is_freaking_out()) counter=PlottingSetup::limitpatience*60+1000;
|
| 71 |
sleep(counterinterval);
|
| 72 |
}
|
| 73 |
gROOT->ProcessLine(".>");
|
| 74 |
std::cout << "\e[1;34m Limit calculation is finishing."<< "\e[0m" << std::endl;
|
| 75 |
remove(threadoutputfile.c_str());
|
| 76 |
|
| 77 |
if(!isThreadActive) {
|
| 78 |
write_info(__FUNCTION__,"Thread finished sucessfully");
|
| 79 |
pthread_join( limitthread, NULL );
|
| 80 |
return true;
|
| 81 |
} else {
|
| 82 |
pthread_cancel(limitthread);
|
| 83 |
write_error(__FUNCTION__,"DID NOT TERMINATE IN TIME - ABORTED!");
|
| 84 |
return false;
|
| 85 |
}
|
| 86 |
}
|
| 87 |
|
| 88 |
|
| 89 |
///***********************************************************************************************************************************************
|
| 90 |
|
| 91 |
|
| 92 |
void* do_limit_computation( void *arg ) {
|
| 93 |
|
| 94 |
isThreadActive=true;
|
| 95 |
string reportname=all_args[1];
|
| 96 |
float luminosity=atof(all_args[2].c_str());
|
| 97 |
float lumiuncert=atof(all_args[3].c_str());
|
| 98 |
float mceff=atof(all_args[4].c_str());
|
| 99 |
float mcefferr=atof(all_args[5].c_str());
|
| 100 |
float Npred=atof(all_args[6].c_str());
|
| 101 |
float Nprederr=atof(all_args[7].c_str());
|
| 102 |
float Nobs=atof(all_args[8].c_str());
|
| 103 |
int JZBcutused=atoi(all_args[9].c_str());
|
| 104 |
string plotname=all_args[10];
|
| 105 |
int doexpected=atoi(all_args[11].c_str());
|
| 106 |
int doasymptotic=false;
|
| 107 |
if(all_args.size()>12) doasymptotic=atoi(all_args[12].c_str());
|
| 108 |
|
| 109 |
if(doasymptotic) write_warning(__FUNCTION__,"Doing asymptotic limit");
|
| 110 |
threadoutputfile=(plotname+"__threadoutput.txt");
|
| 111 |
|
| 112 |
if(doasymptotic) dout << "\e[0;34m " << endl;
|
| 113 |
else dout << "\e[0;31m " << endl;
|
| 114 |
dout << " " << endl;
|
| 115 |
dout << " .+8OOOO? " << endl;
|
| 116 |
dout << " ..OOOOOOOO$ZZZZZI " << endl;
|
| 117 |
dout << " ...=:: :ZZZZZZZZZZZZO? " << endl;
|
| 118 |
dout << " ..~:~ ,ZZZZZZZZZZOO$. " << endl;
|
| 119 |
dout << " +:, :ZZZZZZOOOOO$. " << endl;
|
| 120 |
dout << " =: ,ZZZOOOZZ$Z$I. " << endl;
|
| 121 |
if(doasymptotic) dout << " :~ blue +??ZZZZZZZZ$. " << endl;
|
| 122 |
else dout << " :~ red +??ZZZZZZZZ$. " << endl;
|
| 123 |
dout << " :~ pill :ZZZZZZZ$7, " << endl;
|
| 124 |
dout << " :: =ZZZZZ8$?: " << endl;
|
| 125 |
dout << " ,: ~ZZ8ZI " << endl;
|
| 126 |
dout << " :::,, ,,::::~=8Z7+~ " << endl;
|
| 127 |
dout << " :::::: " << endl;
|
| 128 |
dout << " " << endl;
|
| 129 |
dout << " \e[00m" << endl;
|
| 130 |
|
| 131 |
dout << "Limit capsule : " << endl;
|
| 132 |
dout << "Initialized with the following parameters: " << endl;
|
| 133 |
dout << "Compute expected limits? " << doexpected << " (will be calculated anyway) " << endl;
|
| 134 |
dout << "Luminosity: " << luminosity << " +/- " << lumiuncert << endl;
|
| 135 |
dout << "Efficiency: " << mceff << " +/- " << mcefferr << endl;
|
| 136 |
dout << "Observed: " << Nobs << endl;
|
| 137 |
dout << "Predicted: " << Npred << " +/- " << Nprederr << endl;
|
| 138 |
dout << "Asymptotic limit? " << doasymptotic << endl;
|
| 139 |
|
| 140 |
|
| 141 |
/*
|
| 142 |
// TODO:
|
| 143 |
- *everything* needs to have absolute paths
|
| 144 |
- need to create a directory for every process (otherwise ws.root's will be overwritten!) --> Good possibility: use the reportname and add some time specific string.
|
| 145 |
*/
|
| 146 |
char currentpath[1024];
|
| 147 |
char *path = getcwd(currentpath,1024);
|
| 148 |
extract_cbaf_dir(currentpath);
|
| 149 |
stringstream loadroostats;
|
| 150 |
if(!doasymptotic) loadroostats << ".L " << PlottingSetup::cbafbasedir << "/Plotting/Modules/external/roostats_cl95.C+";
|
| 151 |
else loadroostats << ".L " << PlottingSetup::cbafbasedir << "/Plotting/Modules/external/roostats_cl95_ASYMPTOTIC.C+";
|
| 152 |
stringstream loadLimitDroplet;
|
| 153 |
loadLimitDroplet << ".L " << PlottingSetup::cbafbasedir << "/Plotting/Modules/LimitDroplet.C";
|
| 154 |
|
| 155 |
cout << "The output of this thread is written to : " << threadoutputfile << "; this Limit Capsule will monitor the output for signs of crahes and if necessary abort. After termination, the log will be deleted to avoid unnecessary clutter." << endl;
|
| 156 |
|
| 157 |
gROOT->ProcessLine("#include <exception>");
|
| 158 |
// gROOT->ProcessLine(((string)".> "+threadoutputfile).c_str());
|
| 159 |
gROOT->ProcessLine(loadroostats.str().c_str());
|
| 160 |
gROOT->ProcessLine(loadLimitDroplet.str().c_str());
|
| 161 |
|
| 162 |
stringstream allargs95;
|
| 163 |
allargs95 << "("<<luminosity<<","<<lumiuncert<<","<<mceff<<","<<mcefferr<<","<<Npred<<","<<Nprederr<<","<<Nobs<<","<<false<<","<<PlottingSetup::nuisancemodel<<",\""<<PlottingSetup::limitmethod<<"\",\""<<plotname<<"\"";
|
| 164 |
stringstream allargslm;
|
| 165 |
allargslm << "("<<luminosity<<","<<lumiuncert<<","<<mceff<<","<<mcefferr<<","<<Npred<<","<<Nprederr<<","<<Nobs<<","<<false<<","<<PlottingSetup::nuisancemodel<<",\""<<PlottingSetup::limitmethod<<"\",\"" << plotname << "\",0";
|
| 166 |
|
| 167 |
stringstream obscommand;
|
| 168 |
// obscommand <<"float observed=roostats_cl95"<<allargs95.str()<<");";
|
| 169 |
// gROOT->ProcessLine(obscommand.str().c_str());
|
| 170 |
stringstream expcommand;
|
| 171 |
expcommand <<"LimitResult limit = roostats_limit"<<allargslm.str()<<");";
|
| 172 |
dout << expcommand.str() << endl;
|
| 173 |
gROOT->ProcessLine(expcommand.str().c_str());
|
| 174 |
gROOT->ProcessLine("LimitDroplet limres;");
|
| 175 |
stringstream setjzb;
|
| 176 |
setjzb<<"limres.JZB="<<JZBcutused<<";";
|
| 177 |
gROOT->ProcessLine(setjzb.str().c_str());
|
| 178 |
gROOT->ProcessLine("limres.observed=limit.GetObservedLimit();");
|
| 179 |
gROOT->ProcessLine("limres.upper68=limit.GetOneSigmaHighRange();");
|
| 180 |
gROOT->ProcessLine("limres.lower68=limit.GetOneSigmaLowRange();");
|
| 181 |
gROOT->ProcessLine("limres.expected=limit.GetExpectedLimit();");
|
| 182 |
gROOT->ProcessLine("limres.upper95= limit.GetTwoSigmaHighRange();");
|
| 183 |
gROOT->ProcessLine("limres.lower95=limit.GetTwoSigmaLowRange();");
|
| 184 |
|
| 185 |
gROOT->ProcessLine(".>");
|
| 186 |
stringstream savecommand;
|
| 187 |
savecommand << "limres.saveDroplet(\""<<reportname<<"\");";
|
| 188 |
gROOT->ProcessLine(savecommand.str().c_str());
|
| 189 |
gROOT->ProcessLine(".q");
|
| 190 |
|
| 191 |
isThreadActive=false;
|
| 192 |
|
| 193 |
}
|
| 194 |
|
| 195 |
int main(int nargs, char* args[]) {
|
| 196 |
/* we obtain lumi etc. directly from the template; the only things passed are the following:
|
| 197 |
- report filename [1]
|
| 198 |
- luminosity [2]
|
| 199 |
- lumi uncert [3]
|
| 200 |
- MC efficiency [4]
|
| 201 |
- MC efficiency error [5]
|
| 202 |
- Npred [6]
|
| 203 |
- Nprederr [7]
|
| 204 |
- Nobs [8]
|
| 205 |
- JZB cut [9]
|
| 206 |
- plot name [10]
|
| 207 |
*/
|
| 208 |
|
| 209 |
if(nargs<12) {
|
| 210 |
dout << "NOT ENOUGH ARGUMENTS!" << endl;
|
| 211 |
dout << "You're supposed to provide: " << endl;
|
| 212 |
dout << "LUMI LUMIUNCERT MCEFF MCEFFERR NPRED NPREDERR NOBS JZBCUT PLOTNAME DOEXPECTED" << endl;
|
| 213 |
return -1;
|
| 214 |
}
|
| 215 |
|
| 216 |
for(int i=0;i<nargs;i++) all_args.push_back((string)args[i]);
|
| 217 |
|
| 218 |
if(CapsuleWrapper()) {
|
| 219 |
cout << "Timed limit capsule finished successfully." << endl;
|
| 220 |
return 0;
|
| 221 |
} else {
|
| 222 |
cout << "Timed limit capsule had to be aborted." << endl;
|
| 223 |
return 1;
|
| 224 |
}
|
| 225 |
return 1;
|
| 226 |
}
|
