| 94 |
|
|
| 95 |
|
int srmcpretries=0; |
| 96 |
|
|
| 97 |
< |
void load_scan_sample(int a, int b, int &scanfileindex,int scantype,bool isretry=false) { |
| 97 |
> |
void load_scan_sample(int a, int b, int &scanfileindex,int scantype,int mglu, int mlsp, bool isretry=false) { |
| 98 |
|
|
| 99 |
|
// There is no need to define your sample here. That is done in Setup.C where you define the loading directory! |
| 100 |
|
|
| 101 |
< |
dout << "Going to load file with Xzone=" << a << " and Yzone=" << b << endl; |
| 101 |
> |
dout << "Going to load file with Xzone=" << a << " and Yzone=" << b << " for scantype (mSugra: " << (bool)(scantype==mSUGRA) << " , SMS: " << (bool)(scantype==SMS) << " , GMSB: " << (bool)(scantype==GMSB) << ")" << endl; |
| 102 |
|
stringstream filetoload; |
| 103 |
|
char hostname[1023]; |
| 104 |
|
gethostname(hostname,1023); |
| 105 |
< |
|
| 105 |
> |
string samplename; |
| 106 |
|
if((Contains(hostname,"t3ui")||Contains(hostname,"t3wn"))) { // local case |
| 107 |
|
filetoload << "/shome/buchmann/ntuples/"<<PlottingSetup::ScanSampleDirectory; |
| 108 |
< |
if(scantype==mSUGRA) filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 109 |
< |
if(scantype==SMS) filetoload << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 110 |
< |
if(scantype==GMSB) filetoload << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 111 |
< |
if(!Contains(((scansample.collection)[(scansample.collection).size()-1]).filename,"_"+any2string(a)+"_"+any2string(b))) { |
| 108 |
> |
if(scantype==mSUGRA) { |
| 109 |
> |
//filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 110 |
> |
filetoload.str(""); |
| 111 |
> |
filetoload << "/shome/lbaeni/jzb/" << PlottingSetup::ScanSampleDirectory << "/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root"; |
| 112 |
> |
samplename="mSUGRA"; |
| 113 |
> |
// filetoload << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 114 |
> |
} |
| 115 |
> |
if(scantype==SMS) { |
| 116 |
> |
filetoload << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 117 |
> |
samplename="SMS"; |
| 118 |
> |
} |
| 119 |
> |
if(scantype==GMSB) { |
| 120 |
> |
filetoload << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 121 |
> |
samplename="GMSB"; |
| 122 |
> |
} |
| 123 |
> |
if(scansample.collection.size()<1||!Contains(((scansample.collection)[(scansample.collection).size()-1]).filename,"_"+any2string(a)+"_"+any2string(b))) { |
| 124 |
|
dout << "The last sample is NOT the same one as the current one, possibly popping off last one and adding the new one." << endl; |
| 125 |
< |
if((scansample.collection).size()>1) { |
| 125 |
> |
if((scansample.collection).size()>=1) { |
| 126 |
|
scansample.RemoveLastSample(); |
| 127 |
|
} |
| 128 |
|
scanfileindex=(scansample.collection).size(); |
| 129 |
|
//New: Loading file when necessary, not before (avoiding high memory usage and startup times) |
| 130 |
< |
if(scanfileindex!=0) scansample.AddSample(filetoload.str(),"scansample",1,1,false,true,scanfileindex,kRed); |
| 130 |
> |
|
| 131 |
> |
scansample.AddSample(filetoload.str(),samplename,1,1,false,true,scanfileindex,kRed); |
| 132 |
> |
dout << "Just added the following file: " << filetoload.str() << endl; |
| 133 |
|
} else { |
| 134 |
|
dout << "Last sample is the same as the current one. Recycling it." << endl; |
| 135 |
|
scanfileindex=(scansample.collection).size()-1; |
| 137 |
|
dout << " Going to use the following file: " << filetoload.str() << endl; |
| 138 |
|
} // end of t3 case |
| 139 |
|
else { |
| 140 |
+ |
//CRAB case |
| 141 |
|
write_info(__FUNCTION__,"Hello, CRAB! Might need to load files ..."); |
| 142 |
|
PlottingSetup::limitpatience=PlottingSetup::limitpatienceCRAB; |
| 143 |
|
stringstream copyfile; |
| 144 |
+ |
gSystem->Exec("mkdir -p local_storage"); |
| 145 |
|
copyfile << "lcg-cp -b -T srmv2 srm://storage01.lcg.cscs.ch:8443/srm/managerv2?SFN=/pnfs/lcg.cscs.ch/cms/trivcat/store/user/buchmann/bussola/" << PlottingSetup::ScanSampleDirectory; |
| 146 |
< |
if(scantype==mSUGRA) copyfile << "/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root "; |
| 147 |
< |
if(scantype==SMS) copyfile << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root "; |
| 148 |
< |
if(scantype==GMSB) copyfile << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root "; |
| 146 |
> |
if(scantype==mSUGRA) { |
| 147 |
> |
copyfile<< "/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root local_storage/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root "; |
| 148 |
> |
samplename="mSUGRA"; |
| 149 |
> |
} |
| 150 |
> |
if(scantype==SMS) { |
| 151 |
> |
copyfile << "/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "local_storage/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root "; |
| 152 |
> |
samplename="SMS"; |
| 153 |
> |
} |
| 154 |
> |
if(scantype==GMSB) { |
| 155 |
> |
copyfile << "/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root " << "local_storage/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root "; |
| 156 |
> |
samplename="GMSB"; |
| 157 |
> |
} |
| 158 |
|
stringstream newfilename; |
| 159 |
< |
if(scantype==mSUGRA) newfilename << "mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 160 |
< |
if(scantype==SMS) newfilename << "SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 161 |
< |
if(scantype==GMSB) newfilename << "GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 159 |
> |
// if(scantype==mSUGRA) newfilename << "local_storage/mSUGRA_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 160 |
> |
if(scantype==mSUGRA) newfilename << "local_storage/mSUGRA_M0_" << mglu << "__M12_" << mlsp << ".root"; |
| 161 |
> |
if(scantype==SMS) newfilename << "local_storage/SMS_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 162 |
> |
if(scantype==GMSB) newfilename << "local_storage/GMSB_clean_splitup_" << any2string(a) << "_" << any2string(b) << ".root"; |
| 163 |
|
|
| 164 |
|
if(!Contains(((scansample.collection)[(scansample.collection).size()-1]).filename,"_"+any2string(a)+"_"+any2string(b))||(a==0&&b==0)&&!initialized_t2) { |
| 165 |
|
dout << "The last sample is NOT the same one as the current one, possibly popping off last one and downloading as well as adding the new one." << endl; |
| 166 |
< |
if((scansample.collection).size()>1) { |
| 166 |
> |
while((scansample.collection).size()>0) { |
| 167 |
|
scansample.RemoveLastSample(); |
| 168 |
|
} |
| 169 |
+ |
cout << "New scanfileindex stems from scansample size: " << (scansample.collection).size() << endl; |
| 170 |
|
scanfileindex=(scansample.collection).size(); |
| 171 |
|
//New: Loading file when necessary, not before (avoiding high memory usage and startup times) |
| 172 |
< |
if(scanfileindex!=0||isretry) { |
| 173 |
< |
delete_any_cached_scans(); |
| 174 |
< |
dout << "Going to download the scan file with the following copy command: " << copyfile.str() << endl; |
| 175 |
< |
gSystem->Exec(copyfile.str().c_str()); |
| 172 |
> |
if(1) { |
| 173 |
> |
// delete_any_cached_scans(); |
| 174 |
> |
if(!doesROOTFileExist(newfilename.str())) { |
| 175 |
> |
dout << "Going to download the scan file with the following copy command: " << copyfile.str() << endl; |
| 176 |
> |
gSystem->Exec(copyfile.str().c_str()); // download it if it hasn't been downloaded before |
| 177 |
> |
} |
| 178 |
|
if(doesROOTFileExist(newfilename.str())) { |
| 179 |
|
initialized_t2=true; |
| 180 |
< |
scansample.AddSample(newfilename.str(),"scansample",1,1,false,true,scanfileindex,kRed); |
| 180 |
> |
cout << "Going to add the following sample: " << newfilename.str() << endl; |
| 181 |
> |
gSystem->Exec(("ls -ltrh "+newfilename.str()).c_str()); |
| 182 |
> |
scansample.AddSample(newfilename.str(),samplename,1,1,false,true,scanfileindex,kRed); |
| 183 |
|
} else { |
| 184 |
|
srmcpretries++; |
| 185 |
|
if(srmcpretries<5) { |
| 186 |
|
dout << "The file could not be loaded correctly - retrying!" << endl; |
| 187 |
|
sleep(5); |
| 188 |
< |
load_scan_sample(a,b,scanfileindex,scantype,true); |
| 188 |
> |
load_scan_sample(a,b,scanfileindex,scantype,mglu,mlsp,true); |
| 189 |
|
} else { |
| 190 |
|
dout << "Have tried 5 times to load this sample. Giving up now and failing the program execution" << endl; |
| 191 |
< |
assert(0); |
| 191 |
> |
assert(srmcpretries<5); |
| 192 |
|
dout << "The command " << copyfile.str() << " failed to copy the file to the local storage (saving it as " << newfilename.str() << ")" << endl; |
| 193 |
|
} |
| 194 |
|
} |
| 199 |
|
} |
| 200 |
|
|
| 201 |
|
} |
| 202 |
+ |
|
| 203 |
+ |
|
| 204 |
+ |
// WATCH OUT THIS LINE NEEDS TO BE REMOVED |
| 205 |
+ |
scanfileindex=0; |
| 206 |
|
} |
| 207 |
|
|
| 208 |
< |
float get_xs(float &altxs, float mglu, float mlsp, string massgluname, string massLSPname, map < pair<float, float>, map<string, float> > &xsec, string mcjzb, bool requireZ) { |
| 208 |
> |
/*float get_xs(float &altxs, float mglu, float mlsp, string massgluname, string massLSPname, map < pair<float, float>, map<string, float> > &xsec, string mcjzb, bool requireZ) { |
| 209 |
|
altxs=0; |
| 210 |
< |
for(int iproc=1;iproc<=10;iproc++) { |
| 176 |
< |
float process_xs = GetXSecForPointAndChannel(mglu,mlsp,xsec,iproc); |
| 177 |
< |
altxs+=process_xs; |
| 178 |
< |
} |
| 210 |
> |
for(int iproc=1;iproc<=10;iproc++) altxs+=GetXSecForPointAndChannel(mglu,mlsp,xsec,iproc); |
| 211 |
|
return altxs; |
| 212 |
< |
} |
| 212 |
> |
}*/ |
| 213 |
|
|
| 214 |
|
void establish_SUSY_limits(string mcjzb,string datajzb,vector<float> jzb_cut,bool requireZ, float peakerror, TFile *fsyst, int ibin,float njobs=-1, float jobnumber=-1) { |
| 215 |
|
bool runninglocally=true; |
| 236 |
|
while ((key = (TKey*)nextkey())) |
| 237 |
|
{ |
| 238 |
|
TObject *obj = key->ReadObj(); |
| 239 |
< |
if(Contains((string)(obj->GetName()),"mSUGRA")) scantype=mSUGRA; |
| 239 |
> |
if(Contains((string)(obj->GetName()),"SUGRA")) scantype=mSUGRA; |
| 240 |
|
if(Contains((string)(obj->GetName()),"GMSB")) scantype=GMSB; |
| 241 |
|
} |
| 242 |
|
|
| 281 |
|
TH2F *fullerr = (TH2F*)fsyst->Get((prefix+"systotmap"+any2string(jzb_cut[ibin])).c_str()); |
| 282 |
|
TH2F *NEvents = (TH2F*)fsyst->Get((prefix+"Neventsmap"+any2string(jzb_cut[ibin])).c_str()); |
| 283 |
|
TH2F *lflip = (TH2F*)fsyst->Get((prefix+"flipmap"+any2string(jzb_cut[ibin])).c_str()); |
| 284 |
+ |
TH2F *XSmap = (TH2F*)fsyst->Get((prefix+"XS"+any2string(jzb_cut[ibin])).c_str()); |
| 285 |
|
|
| 286 |
|
TH2F *limitmap = new TH2F((prefix+"limitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//observed limit |
| 287 |
|
TH2F *explimitmap = new TH2F((prefix+"explimitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit |
| 291 |
|
TH2F *exp2mlimitmap = new TH2F((prefix+"exp2mlimitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 2 sigma |
| 292 |
|
|
| 293 |
|
TH2F *exclmap = new TH2F((prefix+"exclusionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 261 |
– |
TH2F *xsmap = new TH2F((prefix+"crosssectionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 262 |
– |
TH2F *totxsmap = new TH2F((prefix+"absolutecrosssectionmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 294 |
|
TH2F *flipmap = new TH2F((prefix+"limitflipmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 295 |
< |
|
| 295 |
> |
TH2F *asymptoticmap = new TH2F((prefix+"asymptoticmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 296 |
|
|
| 297 |
|
|
| 298 |
|
if(!fullerr || !mceff || !NEvents) { |
| 324 |
|
continue; |
| 325 |
|
} |
| 326 |
|
vector<float> sigmas; |
| 327 |
< |
sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped); |
| 328 |
< |
|
| 327 |
> |
if(scantype!=mSUGRA) sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped); |
| 328 |
> |
else { |
| 329 |
> |
//this is a bit trickier - what we want to do is compute the limit fast, and if it is in [0.5 xs, 2xs], compute it again but "correctly"! |
| 330 |
> |
vector<float> asigmas; |
| 331 |
> |
asigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped,true); // asymptotic limit first |
| 332 |
> |
float strength=asigmas[0]/XSmap->GetBinContent(GlobalBin); |
| 333 |
> |
if(strength>0.5&&strength<2) { |
| 334 |
> |
sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true, flipped); // asymptotic limit first |
| 335 |
> |
asymptoticmap->SetBinContent(GlobalBin,0); |
| 336 |
> |
} |
| 337 |
> |
else { |
| 338 |
> |
sigmas=asigmas; |
| 339 |
> |
asymptoticmap->SetBinContent(GlobalBin,1); |
| 340 |
> |
} |
| 341 |
> |
exclmap->SetBinContent(GlobalBin,strength); |
| 342 |
> |
} |
| 343 |
|
|
| 344 |
|
if(sigmas[0]>0) limitmap->SetBinContent(GlobalBin,sigmas[0]); //anything else is an error code |
| 345 |
|
if(sigmas.size()>1) { |
| 351 |
|
} |
| 352 |
|
|
| 353 |
|
dout << "An upper limit has been added for this point ( " << massgluname << " = " << mglu << " and " << massLSPname << " = " << mlsp << " ) at " << sigmas[0] << endl; |
| 309 |
– |
if(scantype==mSUGRA) { // for SMS this is a bit easier at the moment - we have a reference XS file which we use when plotting |
| 310 |
– |
dout << "Computing exclusion status" << endl; |
| 311 |
– |
float rel_limit=0; |
| 312 |
– |
float totxs; |
| 313 |
– |
float xs=get_xs(totxs,mglu,mlsp,massgluname,massLSPname,xsec,mcjzb,requireZ); |
| 314 |
– |
if(xs>0) rel_limit=sigmas[0]/xs; |
| 315 |
– |
exclmap->SetBinContent(GlobalBin,rel_limit); |
| 316 |
– |
xsmap->SetBinContent(GlobalBin,xs); |
| 317 |
– |
totxsmap->SetBinContent(GlobalBin,totxs); |
| 318 |
– |
} |
| 354 |
|
} |
| 355 |
|
} |
| 356 |
|
|
| 359 |
|
outputfile->cd(); |
| 360 |
|
limitmap->Write(); |
| 361 |
|
flipmap->Write(); |
| 362 |
+ |
asymptoticmap->Write(); |
| 363 |
|
if(doexpected) { |
| 364 |
|
explimitmap->Write(); |
| 365 |
|
exp1plimitmap->Write(); |
| 367 |
|
exp2plimitmap->Write(); |
| 368 |
|
exp2mlimitmap->Write(); |
| 369 |
|
} |
| 370 |
< |
if(scantype==mSUGRA) { |
| 335 |
< |
exclmap->Write(); |
| 336 |
< |
xsmap->Write(); |
| 337 |
< |
totxsmap->Write(); |
| 338 |
< |
} |
| 370 |
> |
if(scantype==mSUGRA) exclmap->Write(); |
| 371 |
|
outputfile->Close(); |
| 372 |
|
delete limcanvas; |
| 373 |
|
} |
| 391 |
|
bool runninglocally=true; |
| 392 |
|
if(njobs>-1&&jobnumber>-1) { |
| 393 |
|
runninglocally=false; |
| 394 |
+ |
dout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl; |
| 395 |
|
dout << "Running on the GRID (this is job " << jobnumber << "/" << njobs << ") for a jzb cut at " << jzb_cut[ibin] << endl; |
| 396 |
|
} else { |
| 397 |
|
dout << "Running locally " << endl; |
| 409 |
|
string prefix="SMS_"; |
| 410 |
|
// up to here, everything is set up for SMS; now we need to switch stuff around if we're dealing with an mSUGRA scan! |
| 411 |
|
int scantype=SMS; |
| 412 |
< |
if(Contains((scansample.collection)[0].samplename,"mSUGRA")) scantype=mSUGRA; |
| 412 |
> |
if(Contains((scansample.collection)[0].samplename,"SUGRA")) scantype=mSUGRA; |
| 413 |
|
if(Contains((scansample.collection)[0].samplename,"GMSB")) scantype=GMSB; |
| 414 |
|
|
| 415 |
|
if(scantype==mSUGRA) { |
| 451 |
|
TH2F *exp1mlimitmap = new TH2F((prefix+"exp1mlimitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 1 sigma |
| 452 |
|
TH2F *exp2mlimitmap = new TH2F((prefix+"exp2mlimitmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 2 sigma |
| 453 |
|
|
| 454 |
< |
TH2F *exclmap = new TH2F((prefix+"exclusionmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 455 |
< |
TH2F *sysjesmap = new TH2F((prefix+"sysjes"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 456 |
< |
TH2F *sysjsumap = new TH2F((prefix+"sysjsu"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 457 |
< |
TH2F *sysresmap = new TH2F((prefix+"sysresmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 458 |
< |
TH2F *efficiencymap = new TH2F((prefix+"efficiencymap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 459 |
< |
TH2F *efficiencymapmet = new TH2F((prefix+"efficiencymapmet"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 460 |
< |
TH2F *efficiencymapAcc = new TH2F((prefix+"efficiencymapAcc"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 461 |
< |
TH2F *efficiencymapID = new TH2F((prefix+"efficiencymapID"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 454 |
> |
TH2F *exclmap = new TH2F((prefix+"exclusionmap"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 455 |
> |
TH2F *sysjesmap = new TH2F((prefix+"sysjes"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 456 |
> |
TH2F *sysjsumap = new TH2F((prefix+"sysjsu"+any2string(jzbSel)).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 457 |
> |
TH2F *sysresmap = new TH2F((prefix+"sysresmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 458 |
> |
TH2F *efficiencymap = new TH2F((prefix+"efficiencymap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 459 |
> |
TH2F *efficiencymapmet = new TH2F((prefix+"efficiencymapmet"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 460 |
> |
TH2F *efficiencymapAcc = new TH2F((prefix+"efficiencymapAcc"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 461 |
> |
TH2F *efficiencymapID = new TH2F((prefix+"efficiencymapID"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 462 |
|
TH2F *efficiencymapJets = new TH2F((prefix+"efficiencymapJets"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 463 |
|
TH2F *efficiencymapSignal = new TH2F((prefix+"efficiencymapSignal"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 464 |
+ |
TH2F *efficiencymapContam = new TH2F((prefix+"efficiencymapContam"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend- |
| 465 |
+ |
mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 466 |
|
TH2F *noscefficiencymap = new TH2F((prefix+"noscefficiencymap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 467 |
< |
TH2F *Neventsmap = new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"", (int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 468 |
< |
TH2F *ipointmap = new TH2F((prefix+"ipointmap"+any2string(jzbSel)).c_str(),"", (int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 469 |
< |
TH2F *syspdfmap = new TH2F((prefix+"syspdfmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 470 |
< |
TH2F *systotmap = new TH2F((prefix+"systotmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 471 |
< |
TH2F *sysstatmap = new TH2F((prefix+"sysstatmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 467 |
> |
TH2F *Neventsmap = new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"", (int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 468 |
> |
TH2F *ipointmap = new TH2F((prefix+"ipointmap"+any2string(jzbSel)).c_str(),"", (int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 469 |
> |
TH2F *syspdfmap = new TH2F((prefix+"syspdfmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 470 |
> |
TH2F *systotmap = new TH2F((prefix+"systotmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 471 |
> |
TH2F *sysstatmap = new TH2F((prefix+"sysstatmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 472 |
> |
TH2F *XSmap = new TH2F((prefix+"XS"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 473 |
|
|
| 474 |
< |
TH2F *imposedxmap = new TH2F((prefix+"imposedxmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 475 |
< |
TH2F *realxmap = new TH2F((prefix+"realxmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 474 |
> |
TH2F *imposedxmap; |
| 475 |
> |
TH2F *realxmap; |
| 476 |
|
|
| 477 |
< |
TH2F *timemap = new TH2F((prefix+"timemap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 477 |
> |
TH2F *centermap; |
| 478 |
> |
TH2F *widthmap; |
| 479 |
> |
|
| 480 |
> |
TH2F *timemap = new TH2F((prefix+"timemap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 481 |
> |
|
| 482 |
> |
TH2F *flipmap = new TH2F((prefix+"flipmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 483 |
> |
|
| 484 |
> |
TH2F *nZmap = new TH2F((prefix+"nZmap"+any2string(jzb_cut[ibin])).c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 485 |
|
|
| 486 |
< |
TH2F *flipmap = new TH2F((prefix+"flipmap"+any2string(jzbSel)).c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 486 |
> |
if(ibin==0) { |
| 487 |
> |
centermap = new TH2F((prefix+"centermap").c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 488 |
> |
widthmap = new TH2F((prefix+"widthmap").c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 489 |
> |
|
| 490 |
> |
imposedxmap = new TH2F((prefix+"imposedxmap").c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 491 |
> |
realxmap = new TH2F((prefix+"realxmap").c_str(),"",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,int((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 492 |
> |
} |
| 493 |
> |
|
| 494 |
> |
vector<int> susy_Zdecay_origin; |
| 495 |
> |
vector<TH2F*> h_susy_Zdecay_origin; |
| 496 |
> |
// we only generate these histograms IF we're dealing with the lowest cut! (we don't need them n times, they're the same for all of them) |
| 497 |
> |
for(int i=1000022;i<1000039;i++) { |
| 498 |
> |
if(ibin>0) continue; |
| 499 |
> |
susy_Zdecay_origin.push_back(i); |
| 500 |
> |
TH2F *susy_dec = new TH2F((prefix+"SUSY_Decay_From_"+any2string(i)+"_to_Z").c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 501 |
> |
h_susy_Zdecay_origin.push_back(susy_dec); |
| 502 |
> |
} |
| 503 |
> |
for(int i=2000006;i<2000016;i++) { |
| 504 |
> |
if(ibin>0) continue; |
| 505 |
> |
susy_Zdecay_origin.push_back(i); |
| 506 |
> |
TH2F *susy_dec = new TH2F((prefix+"SUSY_Decay_From_"+any2string(i)+"_to_Z").c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 507 |
> |
h_susy_Zdecay_origin.push_back(susy_dec); |
| 508 |
> |
} |
| 509 |
|
|
| 510 |
+ |
TH2F *PrimaryZSource = new TH2F((prefix+"PrimaryZSource").c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 511 |
+ |
TH2F *SecondaryZSource = new TH2F((prefix+"SecondaryZSource").c_str(), "",(int)((mgluend-mglustart)/mglustep+1),mglustart-0.5*mglustep,mgluend+0.5*mglustep,(int)((mLSPend-mLSPstart)/mLSPstep+1),mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep); |
| 512 |
+ |
|
| 513 |
|
write_warning(__FUNCTION__,"CURRENTLY SWITCHING AUTOMATIZED MODE OFF!");automatized=false; |
| 514 |
|
|
| 515 |
|
float rightmargin=gStyle->GetPadRightMargin(); |
| 524 |
|
for(int mglu=(int)mglustart;mglu<=mgluend;mglu+=(int)mglustep) { |
| 525 |
|
for (int mlsp=(int)mLSPstart;mlsp<=mLSPend&&(scantype==mSUGRA||mlsp<=mglu);mlsp+=(int)mLSPstep) Npoints++; |
| 526 |
|
} |
| 527 |
< |
|
| 527 |
> |
|
| 528 |
> |
map < pair<float, float>, map<string, float> > xsec; |
| 529 |
> |
if(scantype==mSUGRA) { |
| 530 |
> |
xsec=getXsec(PlottingSetup::mSUGRAxsFile); |
| 531 |
> |
} |
| 532 |
> |
|
| 533 |
|
int ipoint=-1; |
| 534 |
|
for(int mglu=(int)mglustart;mglu<=mgluend;mglu+=(int)mglustep) { |
| 535 |
|
for (int mlsp=(int)mLSPstart;mlsp<=mLSPend&&(scantype==mSUGRA||mlsp<=mglu);mlsp+=(int)mLSPstep) |
| 544 |
|
int a = int((PlottingSetup::ScanXzones*mlsp)/(m12trees+1)); |
| 545 |
|
int b = int((PlottingSetup::ScanYzones*mglu)/(m0trees+1)); |
| 546 |
|
int scanfileindex=PlottingSetup::ScanYzones*a+b; |
| 547 |
< |
load_scan_sample(a,b,scanfileindex,scantype); |
| 547 |
> |
load_scan_sample(a,b,scanfileindex,scantype,mglu,mlsp); |
| 548 |
|
|
| 549 |
|
clock_t start,finish; |
| 550 |
|
start = clock(); // starting the clock to measure how long the computation takes! |
| 554 |
|
float nevents = (scansample.collection)[scanfileindex].events->GetSelectedRows(); |
| 555 |
|
vector<vector<float> > systematics; |
| 556 |
|
if(nevents<10) { |
| 557 |
< |
dout << "This point ("<<ipoint<<") with configuration ("<<massgluname<<"="<<mglu<<" , "<<massLSPname<<"="<<mlsp << ") does not contain enough events and will be skipped."<< endl; |
| 557 |
> |
dout << "This point ("<<ipoint<<" / " << Npoints << ") with configuration ("<<massgluname<<"="<<mglu<<" , "<<massLSPname<<"="<<mlsp << ") does not contain enough events and will be \033[1;31m skipped \033[0m."<< endl; |
| 558 |
|
continue; |
| 559 |
|
} else { |
| 560 |
< |
dout << "OK! This point ("<<ipoint<<") with configuration ("<<massgluname<<"="<<mglu<<" , "<<massLSPname<<"="<<mlsp << ") contains " << nevents << " and will therefore not be skipped. " << endl; |
| 561 |
< |
write_analysis_type(PlottingSetup::RestrictToMassPeak); |
| 560 |
> |
dout << "OK! This point ("<<ipoint<<") with configuration ("<<massgluname<<"="<<mglu<<" , "<<massLSPname<<"="<<mlsp << ") contains " << nevents << " and will therefore \033[1;32m not be skipped \033[0m. " << endl; |
| 561 |
> |
if(PlottingSetup::RestrictToMassPeak==true) cout << "Analysis type: on-peak" << endl; |
| 562 |
> |
if(PlottingSetup::RestrictToMassPeak==false) cout << "Analysis type: offpeak" << endl; |
| 563 |
> |
} |
| 564 |
> |
float hitrecord=0; |
| 565 |
> |
int hitrecordholder=0; |
| 566 |
> |
float hitsecondplace=0; |
| 567 |
> |
int hitsecondplaceholder=0; |
| 568 |
> |
for(int imo=0;imo<susy_Zdecay_origin.size()&&ibin==0;imo++) { |
| 569 |
> |
if(ibin>0) continue; |
| 570 |
> |
int MotherParticlePDGid = susy_Zdecay_origin[imo]; |
| 571 |
> |
(scansample.collection)[scanfileindex].events->Draw("eventNum",("("+addcut.str()+")&&(abs(SourceOfZ-"+any2string(MotherParticlePDGid)+")<0.5)").c_str(),"goff"); |
| 572 |
> |
float hits = (scansample.collection)[scanfileindex].events->GetSelectedRows(); |
| 573 |
> |
hits/=nevents; |
| 574 |
> |
h_susy_Zdecay_origin[imo]->SetBinContent(GlobalBin,hits); |
| 575 |
> |
if(hits>hitrecord) { |
| 576 |
> |
hitsecondplace=hitrecord; |
| 577 |
> |
hitsecondplaceholder=hitrecordholder; |
| 578 |
> |
hitrecord=hits; |
| 579 |
> |
hitrecordholder=imo; |
| 580 |
> |
} |
| 581 |
> |
if(hits<hitrecord&&hits>hitsecondplace) { |
| 582 |
> |
hitsecondplace=hits; |
| 583 |
> |
hitsecondplaceholder=imo; |
| 584 |
> |
} |
| 585 |
> |
} |
| 586 |
> |
if(ibin==0) { |
| 587 |
> |
PrimaryZSource->SetBinContent(GlobalBin,hitrecordholder); |
| 588 |
> |
SecondaryZSource->SetBinContent(GlobalBin,hitsecondplaceholder); |
| 589 |
> |
} |
| 590 |
> |
|
| 591 |
> |
if(ibin==0) { |
| 592 |
> |
TH1F *realxhisto = new TH1F("realxhisto","realxhisto",1000,0,1); |
| 593 |
> |
(scansample.collection)[scanfileindex].events->Draw("realx>>realxhisto",addcut.str().c_str(),"goff"); |
| 594 |
> |
realxmap->SetBinContent(GlobalBin,realxhisto->GetMean()); |
| 595 |
> |
realxmap->SetBinError(GlobalBin,realxhisto->GetRMS()); |
| 596 |
> |
(scansample.collection)[scanfileindex].events->Draw("imposedx>>realxhisto",addcut.str().c_str(),"goff"); |
| 597 |
> |
imposedxmap->SetBinContent(GlobalBin,realxhisto->GetMean()); |
| 598 |
> |
imposedxmap->SetBinError(GlobalBin,realxhisto->GetRMS()); |
| 599 |
> |
stringstream specialexpression; |
| 600 |
> |
specialexpression << mcjzb << ">>jzbhisto"; |
| 601 |
> |
TH1F *jzbhisto = new TH1F("jzbhisto","jzbhisto",100,-500,500); |
| 602 |
> |
(scansample.collection)[scanfileindex].events->Draw(specialexpression.str().c_str(),addcut.str().c_str(),"goff"); |
| 603 |
> |
centermap->SetBinContent(GlobalBin,jzbhisto->GetMean()); |
| 604 |
> |
widthmap->SetBinContent(GlobalBin,jzbhisto->GetRMS()); |
| 605 |
> |
delete realxhisto; |
| 606 |
> |
delete jzbhisto; |
| 607 |
|
} |
| 490 |
– |
TH1F *realxhisto = new TH1F("realxhisto","realxhisto",1000,0,1); |
| 491 |
– |
(scansample.collection)[scanfileindex].events->Draw("xSMS>>realxhisto",addcut.str().c_str(),"goff"); |
| 492 |
– |
realxmap->SetBinContent(GlobalBin,realxhisto->GetMean()); |
| 493 |
– |
realxmap->SetBinError(GlobalBin,realxhisto->GetRMS()); |
| 494 |
– |
float impx=0.0; |
| 495 |
– |
if(Contains(((scansample.collection)[0].filename),"T5zz")) impx=0.5; |
| 496 |
– |
if(Contains(((scansample.collection)[0].filename),"T5zzl")) impx=0.75; |
| 497 |
– |
if(Contains(((scansample.collection)[0].filename),"T5zzh")) impx=0.25; |
| 498 |
– |
imposedxmap->SetBinContent(GlobalBin,impx); |
| 499 |
– |
delete realxhisto; |
| 608 |
|
|
| 609 |
|
|
| 610 |
+ |
if(scantype==mSUGRA) { |
| 611 |
+ |
float absxs=0; |
| 612 |
+ |
for(int i=0;i<12;i++) absxs+=GetXSecForPointAndChannel(mglu,mlsp,xsec,i); |
| 613 |
+ |
XSmap->SetBinContent(GlobalBin,absxs); |
| 614 |
+ |
} |
| 615 |
+ |
|
| 616 |
|
if(nevents!=0&&(efficiencyonly||systematicsonly)) { |
| 617 |
< |
Value effwosigcont = MCefficiency((scansample.collection)[scanfileindex].events,result,resulterr,flipped,mcjzb,requireZ,(int)nevents,addcut.str(),-1); |
| 618 |
< |
if(result<0) { |
| 617 |
> |
Value effwosigcont = MCefficiency((scansample.collection)[scanfileindex].events,result,resulterr,flipped,mcjzb,requireZ,(int)nevents,scantype,xsec,addcut.str(),-1); |
| 618 |
> |
if(result<0&&allowflipping) { |
| 619 |
|
write_info(__FUNCTION__,"Found negative efficiency (signal distribution populating JZB<0 more than JZB>0) - flipping analysis!"); |
| 620 |
|
flipped=1; |
| 621 |
< |
MCefficiency((scansample.collection)[scanfileindex].events,result,resulterr,flipped,mcjzb,requireZ,(int)nevents,addcut.str(),-1); |
| 621 |
> |
MCefficiency((scansample.collection)[scanfileindex].events,result,resulterr,flipped,mcjzb,requireZ,(int)nevents,scantype,xsec,addcut.str(),-1); |
| 622 |
|
} |
| 623 |
|
efficiencymap->SetBinContent(GlobalBin,result); |
| 624 |
|
efficiencymap->SetBinError(GlobalBin,resulterr); |
| 625 |
|
flipmap->SetBinContent(GlobalBin,flipped); |
| 626 |
|
noscefficiencymap->SetBinContent(GlobalBin,effwosigcont.getValue()); |
| 627 |
|
noscefficiencymap->SetBinError(GlobalBin,effwosigcont.getError()); |
| 628 |
+ |
efficiencymapContam->SetBinContent(GlobalBin,effwosigcont.getValue()-result); |
| 629 |
+ |
efficiencymapContam->SetBinError(GlobalBin,TMath::Sqrt(effwosigcont.getError()+resulterr)); |
| 630 |
|
//Partial efficiencies |
| 631 |
|
float resultAcc, resultID, resultJets, resultSignal, resultmet; |
| 632 |
|
float resulterrAcc, resulterrID, resulterrJets, resulterrSignal, resulterrmet; |
| 651 |
|
ipointmap->SetBinContent(GlobalBin,ipoint); |
| 652 |
|
if(efficiencyonly) continue; |
| 653 |
|
|
| 654 |
< |
do_systematics_for_one_file((scansample.collection)[scanfileindex].events,(int)nevents,"SUSY SCAN", systematics,flipped, mcjzb,datajzb,peakerror,requireZ, addcut.str(),true);//mSUGRA is now always true here because we always want to compute PDF systematics |
| 654 |
> |
do_systematics_for_one_file((scansample.collection)[scanfileindex].events,(int)nevents,"SUSY SCAN", systematics,flipped, xsec, mcjzb,datajzb,peakerror,requireZ, addcut.str(),true,scantype); |
| 655 |
|
float JZBcutat = systematics[0][0]; |
| 656 |
|
float mceff = systematics[0][1]; |
| 657 |
|
float mcefferr = systematics[0][2];//MC stat error |
| 663 |
|
float mcwoscefr= systematics[0][9]; // error on efficiency without signal contamination |
| 664 |
|
float sys_pdf = 0; |
| 665 |
|
if(systematics[0].size()>10) sys_pdf = systematics[0][10]; // PDF |
| 666 |
< |
|
| 666 |
> |
cout << "Going to store: total error is: " << toterr << " and the relative tot error, i.e. total error / mc efficiency is " << toterr/mceff << endl; |
| 667 |
|
if(mceff!=mceff||toterr!=toterr||mceff<0 && (!systematicsonly&&!efficiencyonly)) { |
| 668 |
|
dout << "Limits can't be calculated for this configuration (mglu="<<mglu<<" , mlsp="<<mlsp << ") as either the efficiency or its error are not positive numbers! (mceff="<<mceff<<" and toterr="<<toterr<<")"<< endl; |
| 669 |
|
continue; |
| 745 |
|
sysjsumap->Write(); |
| 746 |
|
sysresmap->Write(); |
| 747 |
|
efficiencymap->Write(); |
| 748 |
+ |
if(ibin==0) centermap->Write(); |
| 749 |
+ |
if(ibin==0) widthmap->Write(); |
| 750 |
|
flipmap->Write(); |
| 751 |
|
efficiencymapmet->Write(); |
| 752 |
|
efficiencymapAcc->Write(); |
| 754 |
|
efficiencymapJets->Write(); |
| 755 |
|
efficiencymapSignal->Write(); |
| 756 |
|
noscefficiencymap->Write(); |
| 757 |
+ |
efficiencymapContam->Write(); |
| 758 |
|
syspdfmap->Write(); |
| 759 |
|
systotmap->Write(); |
| 760 |
|
sysstatmap->Write(); |
| 761 |
< |
imposedxmap->Write(); |
| 762 |
< |
realxmap->Write(); |
| 761 |
> |
XSmap->Write(); |
| 762 |
> |
if(ibin==0) imposedxmap->Write(); |
| 763 |
> |
if(ibin==0) realxmap->Write(); |
| 764 |
|
Neventsmap->Write(); |
| 765 |
|
ipointmap->Write(); |
| 766 |
|
timemap->Write(); |
| 767 |
+ |
for(int i=0;i<susy_Zdecay_origin.size();i++) { |
| 768 |
+ |
h_susy_Zdecay_origin[i]->Write(); |
| 769 |
+ |
} |
| 770 |
|
outputfile->Close(); |
| 771 |
|
} |
| 772 |
|
} |
| 806 |
|
sysjsumap->Write(); |
| 807 |
|
sysresmap->Write(); |
| 808 |
|
flipmap->Write(); |
| 809 |
+ |
if(ibin==0) centermap->Write(); |
| 810 |
+ |
if(ibin==0) widthmap->Write(); |
| 811 |
|
efficiencymap->Write(); |
| 812 |
|
efficiencymapmet->Write(); |
| 813 |
|
efficiencymapAcc->Write(); |
| 815 |
|
efficiencymapJets->Write(); |
| 816 |
|
efficiencymapSignal->Write(); |
| 817 |
|
noscefficiencymap->Write(); |
| 818 |
+ |
efficiencymapContam->Write(); |
| 819 |
|
Neventsmap->Write(); |
| 820 |
|
ipointmap->Write(); |
| 821 |
|
syspdfmap->Write(); |
| 822 |
|
systotmap->Write(); |
| 823 |
|
sysstatmap->Write(); |
| 824 |
< |
imposedxmap->Write(); |
| 825 |
< |
realxmap->Write(); |
| 824 |
> |
XSmap->Write(); |
| 825 |
> |
if(ibin==0) imposedxmap->Write(); |
| 826 |
> |
if(ibin==0) realxmap->Write(); |
| 827 |
|
timemap->Write(); |
| 828 |
+ |
for(int i=0;i<susy_Zdecay_origin.size();i++) { |
| 829 |
+ |
h_susy_Zdecay_origin[i]->Write(); |
| 830 |
+ |
} |
| 831 |
|
outputfile->Close(); |
| 832 |
|
} |
| 833 |
|
}//end of systematics only |
| 893 |
|
ipointmap->Write(); |
| 894 |
|
Neventsmap->Write(); |
| 895 |
|
noscefficiencymap->Write(); |
| 896 |
+ |
efficiencymapContam->Write(); |
| 897 |
|
efficiencymap->Write(); |
| 898 |
+ |
if(ibin==0) centermap->Write(); |
| 899 |
+ |
if(ibin==0) widthmap->Write(); |
| 900 |
|
flipmap->Write(); |
| 901 |
|
efficiencymapmet->Write(); |
| 902 |
|
efficiencymapAcc->Write(); |
| 908 |
|
} |
| 909 |
|
}//end of efficiencies only |
| 910 |
|
|
| 911 |
+ |
for(int i=0;i<susy_Zdecay_origin.size();i++) { |
| 912 |
+ |
delete h_susy_Zdecay_origin[i]; |
| 913 |
+ |
} |
| 914 |
+ |
|
| 915 |
+ |
delete PrimaryZSource; |
| 916 |
+ |
delete SecondaryZSource; |
| 917 |
+ |
|
| 918 |
+ |
|
| 919 |
|
delete limitmap; |
| 920 |
|
delete exclmap; |
| 921 |
|
delete sysjesmap; |
| 922 |
|
delete sysjsumap; |
| 923 |
|
delete sysresmap; |
| 924 |
|
delete noscefficiencymap; |
| 925 |
+ |
delete efficiencymapContam; |
| 926 |
|
delete efficiencymap; |
| 927 |
|
delete efficiencymapmet; |
| 928 |
|
delete efficiencymapAcc; |
| 934 |
|
delete syspdfmap; |
| 935 |
|
delete systotmap; |
| 936 |
|
delete sysstatmap; |
| 937 |
+ |
delete XSmap; |
| 938 |
|
delete limcanvas; |
| 939 |
|
} |
| 940 |
|
|
