ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/Plotting/Modules/SUSYScan.C
(Generate patch)

Comparing UserCode/cbrown/AnalysisFramework/Plotting/Modules/SUSYScan.C (file contents):
Revision 1.36 by buchmann, Thu Sep 22 14:41:00 2011 UTC vs.
Revision 1.37 by buchmann, Mon Sep 26 15:43:37 2011 UTC

# Line 24 | Line 24 | using namespace std;
24   using namespace PlottingSetup;
25  
26   void susy_scan_axis_labeling(TH2F *histo) {
27 <  histo->GetXaxis()->SetTitle("#Chi_{2}^{0}-LSP");
27 >  histo->GetXaxis()->SetTitle("m_{#Chi_{2}^{0}}-m_{LSP}");
28    histo->GetXaxis()->CenterTitle();
29    histo->GetYaxis()->SetTitle("m_{#tilde{q}}");
30    histo->GetYaxis()->CenterTitle();
# Line 48 | Line 48 | void* compute_one_upper_limit_wrapper(vo
48    args->sigmas=compute_one_upper_limit(args->mceff,args->toterr,args->ibin,args->mcjzb,args->plotfilename,true);
49    std::cout << "Thread to compute limits has finished" << std::endl;
50    isThreadActive=false;
51 <  return NULL; // oder in C++: return 0;// Damit kann man Werte zurückgeben
51 >  return NULL;
52   }
53  
54   void do_limit_wrapper(float mceff,float toterr,int ibin,string mcjzb,vector<float> &sigmas,string plotfilename) {
# Line 110 | Line 110 | float get_xs(float mglu, float mlsp, str
110    stringstream addcut;
111    addcut << "(TMath::Abs("<<massgluname<<"-"<<mglu<<")<5)&&(TMath::Abs("<<massLSPname<<"-"<<mlsp<<")<5)";
112    cout << "About to calculate the process efficiencies " << endl;
113 <  for(int iproc=1;iproc<=10;iproc++) {
113 > /*  for(int iproc=1;iproc<=10;iproc++) {
114      float process_xs = GetXSecForPointAndChannel(mglu,mlsp,xsec,iproc);
115      stringstream addcutplus;
116      addcutplus<<addcut.str()<<"&&(pfJetGoodNum=="<<iproc<<")";
# Line 123 | Line 123 | float get_xs(float mglu, float mlsp, str
123      float weight=0;
124      if(nprocessevents>0) weight=(process_xs)/(nprocessevents);//*luminosity;
125      weightedpointxs+=weight*nselectedprocessevents;
126 – //    */
126    }
127    return weightedpointxs;
128 + */
129 + return 1.0;
130   }
131  
132   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) {
# Line 135 | Line 136 | void establish_SUSY_limits(string mcjzb,
136      dout << "Running on the GRID (this is job " << jobnumber << "/" << njobs << ") for a jzb cut at " << jzb_cut[ibin] << endl;
137    } else {
138      dout << "Running locally " << endl;
139 <    dout << "This will take a really, really long time - if  you want to see the results THIS week try running the LimitWorkerScript on the grid (DistributedModelCalculation/Limits/)" << endl;
139 >    dout << "This will take a really, really long time - if  you want to see the results within hours instead of weeks try running the worker script on the grid (DistributedModelCalculation/Limits/)" << endl;
140    }
141  
142    string massgluname="MassGlu";
# Line 181 | Line 182 | void establish_SUSY_limits(string mcjzb,
182    for(int mglu=mglustart;mglu<=mgluend;mglu+=mglustep) {
183      for (int mlsp=mLSPstart;mlsp<=mLSPend&&mlsp<=mglu;mlsp+=mLSPstep) Npoints++;
184    }
185 < //  TH2F *mceff     = (TH2F*)fsyst->Get((prefix+"efficiencymap"+any2string(jzb_cut[ibin])).c_str());
186 < write_warning(__FUNCTION__,"Currently the efficiencymap name was switched to Pablo's convention. This NEEDS to be switched BACK!");TH2F *mceff     = (TH2F*)fsyst->Get(("efficiency_jzbdiff"+any2string(jzb_cut[ibin])).c_str());
185 >  TH2F *mceff     = (TH2F*)fsyst->Get((prefix+"efficiencymap"+any2string(jzb_cut[ibin])).c_str());
186 > //write_warning(__FUNCTION__,"Currently the efficiencymap name was switched to Pablo's convention. This NEEDS to be switched BACK!");TH2F *mceff     = (TH2F*)fsyst->Get(("efficiency_jzbdiff"+any2string(jzb_cut[ibin])).c_str());
187    TH2F *fullerr   = (TH2F*)fsyst->Get((prefix+"systotmap"+any2string(jzb_cut[ibin])).c_str());
188    TH2F *NEvents   = (TH2F*)fsyst->Get((prefix+"Neventsmap"+any2string(jzb_cut[ibin])).c_str());
189  
190 <  TH2F *limitmap  = new TH2F((prefix+"limitmap"+any2string(jzb_cut[ibin])).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
190 >  TH2F *limitmap  = new TH2F((prefix+"limitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//observed limit
191 >  TH2F *explimitmap  = new TH2F((prefix+"explimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit
192 >  TH2F *exp1plimitmap  = new TH2F((prefix+"exp1plimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit + 1 sigma
193 >  TH2F *exp2plimitmap  = new TH2F((prefix+"exp2plimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit + 2 sigma
194 >  TH2F *exp1mlimitmap  = new TH2F((prefix+"exp1mlimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 1 sigma
195 >  TH2F *exp2mlimitmap  = new TH2F((prefix+"exp2mlimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 2 sigma
196 >  
197    TH2F *exclmap  = new TH2F((prefix+"exclusionmap"+any2string(jzb_cut[ibin])).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
198    TH2F *xsmap  = new TH2F((prefix+"crosssectionmap"+any2string(jzb_cut[ibin])).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
199    
# Line 196 | Line 203 | write_warning(__FUNCTION__,"Currently th
203      delete limcanvas;
204      return;
205    }
206 +  
207 +  bool doexpected=true;
208  
209    int ipoint=-1;
210    for(int mglu=mglustart;mglu<=mgluend;mglu+=mglustep) {
# Line 217 | Line 226 | write_warning(__FUNCTION__,"Currently th
226        //do_limit_wrapper(currmceff,currtoterr,ibin,mcjzb,sigmas,plotfilename); // no threading right now.
227        sigmas=compute_one_upper_limit(currmceff,currtoterr,ibin,mcjzb,plotfilename,true);
228        if(sigmas[0]>0) limitmap->Fill(mglu,mlsp,sigmas[0]); //anything else is an error code
229 +      if(sigmas.size()>1) {
230 +        explimitmap->Fill(mglu,mlsp,sigmas[1]);
231 +        exp1plimitmap->Fill(mglu,mlsp,sigmas[2]);
232 +        exp1mlimitmap->Fill(mglu,mlsp,sigmas[3]);
233 +        exp2plimitmap->Fill(mglu,mlsp,sigmas[4]);
234 +        exp2mlimitmap->Fill(mglu,mlsp,sigmas[5]);
235 +      } else {
236 +        write_warning(__FUNCTION__,"Watch out, limit results did not contain expected limits!");
237 +      }
238 +
239        dout << "An upper limit has been added for this point ( " << massgluname << " = " << mglu << " and " << massLSPname << " = " << mlsp << " ) at " << sigmas[0] << endl;
240        if(1) {
241            dout << "Computing exclusion status" << endl;
# Line 236 | Line 255 | write_warning(__FUNCTION__,"Currently th
255    TFile *outputfile=new TFile(("output/DistributedLimitsFromSystematics_job"+string(any2string(jobnumber))+"_of_"+string(any2string(njobs))+".root").c_str(),"UPDATE");//needs to be "UPDATE" as we can get to this point for different JZB cuts and don't want to erase the previous data :-)
256    outputfile->cd();
257    limitmap->Write();
258 +  if(doexpected) {
259 +    explimitmap->Write();
260 +    exp1plimitmap->Write();
261 +    exp1mlimitmap->Write();
262 +    exp2plimitmap->Write();
263 +    exp2mlimitmap->Write();
264 +  }
265    if(ismSUGRA) {
266      exclmap->Write();
267      xsmap->Write();
# Line 308 | Line 334 | cout << "FILE USED (name): " << (scansam
334    
335    gStyle->SetPalette(100, MyPalette);
336  
337 <  TH2F *limitmap  = new TH2F((prefix+"limitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
338 <  TH2F *exclmap  = new TH2F((prefix+"exclusionmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
339 <  TH2F *sysjesmap = new TH2F((prefix+"sysjes"+any2string(jzbSel)).c_str(),   "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
340 <  TH2F *sysjsumap = new TH2F((prefix+"sysjsu"+any2string(jzbSel)).c_str(),   "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
341 <  TH2F *sysresmap = new TH2F((prefix+"sysresmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
342 <  TH2F *efficiencymap = new TH2F((prefix+"efficiencymap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
343 <  TH2F *Neventsmap    = new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"",   (mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
344 <  TH2F *ipointmap    = new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"",   (mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
345 <  TH2F *syspdfmap = new TH2F((prefix+"syspdfmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
346 <  TH2F *systotmap = new TH2F((prefix+"systotmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
347 <  TH2F *sysstatmap = new TH2F((prefix+"sysstatmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
337 >  TH2F *limitmap  = new TH2F((prefix+"limitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//observed limit
338 >  TH2F *explimitmap  = new TH2F((prefix+"explimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit
339 >  TH2F *exp1plimitmap  = new TH2F((prefix+"exp1plimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit + 1 sigma
340 >  TH2F *exp2plimitmap  = new TH2F((prefix+"exp2plimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit + 2 sigma
341 >  TH2F *exp1mlimitmap  = new TH2F((prefix+"exp1mlimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 1 sigma
342 >  TH2F *exp2mlimitmap  = new TH2F((prefix+"exp2mlimitmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);//expected limit - 2 sigma
343 >  
344 >  TH2F *exclmap  =       new TH2F((prefix+"exclusionmap"+any2string(jzbSel)).c_str(), "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
345 >  TH2F *sysjesmap =      new TH2F((prefix+"sysjes"+any2string(jzbSel)).c_str(),   "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
346 >  TH2F *sysjsumap =      new TH2F((prefix+"sysjsu"+any2string(jzbSel)).c_str(),   "",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
347 >  TH2F *sysresmap =      new TH2F((prefix+"sysresmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
348 >  TH2F *efficiencymap =  new TH2F((prefix+"efficiencymap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
349 >  TH2F *noscefficiencymap = new TH2F((prefix+"nosc_efficiencymap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
350 >  TH2F *Neventsmap    =  new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"",   (mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
351 >  TH2F *ipointmap    =   new TH2F((prefix+"Neventsmap"+any2string(jzbSel)).c_str(),"",   (mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
352 >  TH2F *syspdfmap =      new TH2F((prefix+"syspdfmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
353 >  TH2F *systotmap =      new TH2F((prefix+"systotmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
354 >  TH2F *sysstatmap =     new TH2F((prefix+"sysstatmap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
355  
356 <  TH2F *timemap = new TH2F((prefix+"timemap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
356 >  TH2F *timemap =        new TH2F((prefix+"timemap"+any2string(jzbSel)).c_str(),"",(mgluend-mglustart)/mglustep+1,mglustart-0.5*mglustep,mgluend+0.5*mglustep,(mLSPend-mLSPstart)/mLSPstep+1,mLSPstart-0.5*mLSPstep,mLSPend+0.5*mLSPstep);
357    
358   write_warning(__FUNCTION__,"CURRENTLY SWITCHING AUTOMATIZED MODE OFF!");automatized=false;
359  
# Line 329 | Line 362 | write_warning(__FUNCTION__,"CURRENTLY SW
362  
363    TCanvas *limcanvas = new TCanvas("limcanvas","Limit canvas");
364  
365 +  
366 +  bool doexpected=true;
367 +  
368    int Npoints=0;
369    for(int mglu=mglustart;mglu<=mgluend;mglu+=mglustep) {
370      for (int mlsp=mLSPstart;mlsp<=mLSPend&&mlsp<=mglu;mlsp+=mLSPstep) Npoints++;
# Line 355 | Line 391 | write_warning(__FUNCTION__,"CURRENTLY SW
391      dout << "OK! This point ("<<ipoint<<") with configuration ("<<massgluname<<"="<<mglu<<" , "<<massLSPname<<"="<<mlsp << ") contains " << nevents << " and will therefore not be skipped."<< endl;
392        }
393        if(nevents!=0&&efficiencyonly) {
394 <          MCefficiency((scansample.collection)[0].events,result,resulterr,mcjzb,requireZ,nevents,addcut.str(),-1);
394 >          Value effwosigcont = MCefficiency((scansample.collection)[0].events,result,resulterr,mcjzb,requireZ,nevents,addcut.str(),-1);
395            efficiencymap->Fill(mglu,mlsp,result);
396 +          noscefficiencymap->Fill(mglu,mlsp,effwosigcont.getValue());
397 +          noscefficiencymap->SetBinError(mglu,mlsp,effwosigcont.getError());
398            finish = clock();
399            timemap->Fill(mglu,mlsp,((float(finish)-float(start))/CLOCKS_PER_SEC));
400        }
# Line 372 | Line 410 | write_warning(__FUNCTION__,"CURRENTLY SW
410        float sys_jes  = systematics[0][5]; // Jet Energy Scale
411        float sys_jsu  = systematics[0][6]; // JZB scale uncertainty
412        float sys_res  = systematics[0][7]; // resolution
413 +      float mcwoscef = systematics[0][8]; // efficiency without signal contamination
414 +      float mcwoscefr= systematics[0][9]; // error on efficiency without signal contamination
415        float sys_pdf   = 0;
416 <      if(systematics[0].size()>9) sys_pdf = systematics[0][9]; // PDF
416 >      if(systematics[0].size()>10) sys_pdf = systematics[0][10]; // PDF
417        efficiencymap->Fill(mglu,mlsp,mceff);
418 <      efficiencymap->SetBinError(efficiencymap->FindBin(mglu,mlsp),mcefferr);//blublu
418 >      efficiencymap->SetBinError(efficiencymap->FindBin(mglu,mlsp),mcefferr);
419 >      noscefficiencymap->Fill(mglu,mlsp,mcwoscef);
420 >      noscefficiencymap->SetBinError(efficiencymap->FindBin(mglu,mlsp),mcwoscefr);
421      
422        if(mceff!=mceff||toterr!=toterr||mceff<0) {
423      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;
# Line 389 | Line 431 | write_warning(__FUNCTION__,"CURRENTLY SW
431        cout << "back in " << __FUNCTION__ << endl;
432        if(sigmas[0]>-0.5) { // negative sigmas are the error signature of do_limit_wrapper, so we want to exclude them.
433          limitmap->Fill(mglu,mlsp,sigmas[0]);
434 +        if(sigmas.size()>1) {
435 +          explimitmap->Fill(mglu,mlsp,sigmas[1]);
436 +          exp1plimitmap->Fill(mglu,mlsp,sigmas[2]);
437 +          exp1mlimitmap->Fill(mglu,mlsp,sigmas[3]);
438 +          exp2plimitmap->Fill(mglu,mlsp,sigmas[4]);
439 +          exp2mlimitmap->Fill(mglu,mlsp,sigmas[5]);
440 +        } else {
441 +        write_warning(__FUNCTION__,"Watch out, limit results did not contain expected limits!");
442 +      }
443 +
444          sysjesmap->Fill(mglu,mlsp,sys_jes);
445          sysjsumap->Fill(mglu,mlsp,sys_jsu);
446          sysresmap->Fill(mglu,mlsp,sys_res);
# Line 439 | Line 491 | write_warning(__FUNCTION__,"CURRENTLY SW
491        if(systematicsonly) outputfile=new TFile(("output/DistributedSystematics_job"+string(any2string(jobnumber))+"_of_"+string(any2string(njobs))+".root").c_str(),"UPDATE");//needs to be "UPDATE" as we can get to this point for different JZB cuts and don't want to erase the previous data :-)
492      else outputfile=new TFile(("output/DistributedLimits_job"+string(any2string(jobnumber))+"_of_"+string(any2string(njobs))+".root").c_str(),"UPDATE");//needs to be "UPDATE" as we can get to this point for
493        limitmap->Write();
494 +      if(doexpected) {
495 +        explimitmap->Write();
496 +        exp1plimitmap->Write();
497 +        exp1mlimitmap->Write();
498 +        exp2plimitmap->Write();
499 +        exp2mlimitmap->Write();
500 +      }
501        sysjesmap->Write();
502        sysjsumap->Write();
503        sysresmap->Write();
504        efficiencymap->Write();
505 +      noscefficiencymap->Write();
506        syspdfmap->Write();
507        systotmap->Write();
508        sysstatmap->Write();
# Line 485 | Line 545 | write_warning(__FUNCTION__,"CURRENTLY SW
545        sysjsumap->Write();
546        sysresmap->Write();
547        efficiencymap->Write();
548 +      noscefficiencymap->Write();
549        Neventsmap->Write();
550        ipointmap->Write();
551        syspdfmap->Write();
# Line 517 | Line 578 | write_warning(__FUNCTION__,"CURRENTLY SW
578        TFile *outputfile=new TFile(("output/DistributedSystematics_job"+string(any2string(jobnumber))+"_of_"+string(any2string(njobs))+".root").c_str(),"UPDATE");
579        ipointmap->Write();
580        Neventsmap->Write();
581 +      noscefficiencymap->Write();
582        efficiencymap->Write();
583        timemap->Write();
584        outputfile->Close();
# Line 528 | Line 590 | write_warning(__FUNCTION__,"CURRENTLY SW
590    delete sysjesmap;
591    delete sysjsumap;
592    delete sysresmap;
593 +  delete noscefficiencymap;
594    delete efficiencymap;
595    delete Neventsmap;
596    delete ipointmap;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines