| 183 |
|
|
| 184 |
|
} |
| 185 |
|
|
| 186 |
< |
vector<float> compute_one_upper_limit(float mceff,float mcefferr, int ibin, string mcjzb, string plotfilename, bool doexpected, int flipped) { |
| 186 |
> |
vector<float> compute_one_upper_limit(float mceff,float mcefferr, int ibin, string mcjzb, string plotfilename, bool doexpected, int flipped, bool doasymptotic=false) { |
| 187 |
|
float sigma95=-9.9,sigma95A=-9.9; |
| 188 |
|
/* |
| 189 |
|
USAGE OF ROOSTATS_CL95 |
| 229 |
|
if(flipped==0) dout << "Calling limit capsule instead of calling : CL95(" << luminosity << "," << lumiuncert*luminosity << "," << mceff << "," << mcefferr << "," << Npred[ibin] << "," << Nprederr[ibin] << "," << Nobs[ibin] << "," << false << "," << nuisancemodel<< ") " << endl; |
| 230 |
|
if(flipped>0) dout << "Calling limit capsule instead of calling : CL95(" << luminosity << "," << lumiuncert*luminosity << "," << mceff << "," << mcefferr << "," << flippedNpred[ibin] << "," << flippedNprederr[ibin] << "," << flippedNobs[ibin] << "," << false << "," << nuisancemodel<< ") " << endl; |
| 231 |
|
|
| 232 |
< |
stringstream command; |
| 233 |
< |
if(flipped==0) command << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/Limits/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << Npred[ibin] << " " << Nprederr[ibin] << " " << Nobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected; |
| 234 |
< |
if(flipped>0) command << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/Limits/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << flippedNpred[ibin] << " " << flippedNprederr[ibin] << " " << flippedNobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected; |
| 235 |
< |
dout << command.str() << endl; |
| 236 |
< |
|
| 237 |
< |
int retval = 256; |
| 238 |
< |
int attempts=0; |
| 239 |
< |
while(!(retval==0||attempts>=3)) {//try up to 3 times |
| 232 |
> |
|
| 233 |
> |
// write_warning(__FUNCTION__,"Creating special case here"); |
| 234 |
> |
// time_t msec = time(NULL)*1000; |
| 235 |
> |
// stringstream TimeCapsuleDirectory; |
| 236 |
> |
// TimeCapsuleDirectory << PlottingSetup::cbafbasedir << "/exchange/TimeCapsuleDir" << msec; |
| 237 |
> |
// gSystem->Exec(((string)"mkdir -p "+TimeCapsuleDirectory.str()).c_str()); |
| 238 |
> |
// gSystem->Exec(((string)"cp "+PlottingSetup::cbafbasedir+"/DistributedModelCalculations/Limits/TimedLimitCapsule.exec "+TimeCapsuleDirectory.str()).c_str()); |
| 239 |
> |
|
| 240 |
> |
stringstream command; |
| 241 |
> |
if(flipped==0) command << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/Limits/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << Npred[ibin] << " " << Nprederr[ibin] << " " << Nobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected << " " << doasymptotic; |
| 242 |
> |
if(flipped>0) command << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/Limits/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << flippedNpred[ibin] << " " << flippedNprederr[ibin] << " " << flippedNobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected<< " " << doasymptotic; |
| 243 |
> |
// if(flipped==0) command << TimeCapsuleDirectory.str() <<"/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << Npred[ibin] << " " << Nprederr[ibin] << " " << Nobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected << " " << doasymptotic; |
| 244 |
> |
// if(flipped>0) command << TimeCapsuleDirectory.str() <<"/TimedLimitCapsule.exec " << repname.str() << " " << luminosity << " " << luminosity*lumiuncert << " " << mceff << " " << mcefferr << " " << flippedNpred[ibin] << " " << flippedNprederr[ibin] << " " << flippedNobs[ibin] << " " << -1 << " " << PlottingSetup::basedirectory << "/" << plotfilename << " " << doexpected<< " " << doasymptotic; |
| 245 |
> |
|
| 246 |
> |
dout << command.str() << endl; |
| 247 |
> |
|
| 248 |
> |
if(doasymptotic) write_warning(__FUNCTION__, "DOING ASYMPTOTIC LIMIT!"); |
| 249 |
> |
|
| 250 |
> |
int retval = 256; |
| 251 |
> |
int attempts=0; |
| 252 |
> |
while(!(retval==0||attempts>=3)) {//try up to 3 times |
| 253 |
|
attempts++; |
| 254 |
|
dout << "Starting limit calculation (TimedLimitCapsule) now : Attempt " << attempts << endl; |
| 255 |
|
retval=gSystem->Exec(command.str().c_str()); |
| 268 |
|
remove(repname.str().c_str()); |
| 269 |
|
sigma95=limres.observed; |
| 270 |
|
|
| 271 |
+ |
// write_warning(__FUNCTION__,"Deleting special case here"); |
| 272 |
+ |
// gSystem->Exec(((string)"ls -ltrh "+TimeCapsuleDirectory.str()).c_str()); |
| 273 |
+ |
// gSystem->Exec(((string)"rm "+TimeCapsuleDirectory.str()+"/TimedLimitCapsule.exec").c_str()); |
| 274 |
+ |
// gSystem->Exec(((string)"ls -ltrh "+TimeCapsuleDirectory.str()).c_str()); |
| 275 |
|
|
| 276 |
|
///------------------------------------------ < /NEW > ---------------------------------------------------------- |
| 277 |
|
vector<float> sigmas; |
| 290 |
|
}//end of mc efficiency is ok |
| 291 |
|
} |
| 292 |
|
|
| 293 |
< |
void compute_upper_limits_from_counting_experiment(vector<vector<float> > uncertainties,vector<float> jzbcuts, string mcjzb, bool doexpected, int flipped) { |
| 293 |
> |
vector<float> compute_upper_limits_from_counting_experiment(vector<vector<float> > uncertainties,vector<float> jzbcuts, string mcjzb, bool doexpected, int flipped) { |
| 294 |
|
dout << "Doing counting experiment ... " << endl; |
| 295 |
|
vector<vector<string> > limits; |
| 296 |
|
vector<vector<float> > vlimits; |
| 402 |
|
} |
| 403 |
|
allresults.Print(); |
| 404 |
|
|
| 405 |
+ |
//--------------------------------------------- |
| 406 |
+ |
|
| 407 |
+ |
vector<float> lowestULs; |
| 408 |
+ |
for(int isample=0;isample<signalsamples.collection.size();isample++) { |
| 409 |
+ |
float lowestUL=-1; |
| 410 |
+ |
for(int icut=0;icut<jzbcuts.size();icut++) { |
| 411 |
+ |
float currUL=Round((vlimits[isample][2*icut]),3); |
| 412 |
+ |
if(currUL>0) { |
| 413 |
+ |
if(lowestUL<0) lowestUL=currUL; |
| 414 |
+ |
if(currUL<lowestUL) lowestUL=currUL; |
| 415 |
+ |
} |
| 416 |
+ |
lowestULs.push_back(lowestUL); |
| 417 |
+ |
} |
| 418 |
+ |
} |
| 419 |
|
|
| 420 |
+ |
//--------------------------------------------- |
| 421 |
+ |
return lowestULs; |
| 422 |
|
} |
| 423 |
|
|
| 424 |
|
|
