ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/NonMCBackground/src/plotSignalFakes.cc
(Generate patch)

Comparing UserCode/MitHzz4l/NonMCBackground/src/plotSignalFakes.cc (file contents):
Revision 1.2 by dkralph, Thu Jul 19 11:15:32 2012 UTC vs.
Revision 1.4 by dkralph, Tue Jul 31 16:46:22 2012 UTC

# Line 3 | Line 3
3   #include <iomanip>
4   // #include <pair.h>
5  
6 #include "RooHistPdf.h"
7 #include "RooGlobalFunc.h"
8 #include "RooPlot.h"
9 #include "RooRealVar.h"
10 #include "RooDataHist.h"
11 #include "RooGenericPdf.h"
12 #include "RooDataSet.h"
13 #include "RooAddPdf.h"
14 #include "RooFormulaVar.h"
15 #include "RooFitResult.h"
16 #include "RooCurve.h"
17 #include "RooBinning.h"
18 #include "RooExponential.h"
19 #include "RooExtendPdf.h"
20
6   #include "TCanvas.h"
7   #include "TChain.h"
8   #include "TString.h"
9   #include "TStyle.h"
10   #include "TH1D.h"
11 + #include "TNtuple.h"
12  
13   #include "CPlot.h"
14   #include "FOArgs.h"
# Line 31 | Line 17
17   #include "Various.h"
18   #include "CommonDefs.h"
19   #include "PlotHeaders.h"
20 + #include "TriggerUtils.h"
21 + #include "JetDefs.h"
22 + #include "JetInfoStruct.h"
23 +
24 + #include "TMVA/Reader.h"
25 + #include "TMVA/Tools.h"
26 + #include "TMVA/Config.h"
27 + #include "TMVA/MethodBDT.h"
28  
29   #ifndef CMSSW_BASE
30   #define CMSSW_BASE "../../"
# Line 40 | Line 34 | using namespace std;
34   using namespace RooFit;
35   using namespace mithep;
36  
37 + TH1D* hpu_2011;
38 + TH1D* hpu_2012;
39 + //----------------------------------------------------------------------------------------
40 + class EvtsOfInterest
41 + {
42 + public:
43 +  EvtsOfInterest(TString config="");
44 +  bool has(unsigned run, unsigned evt);
45 +  vector<unsigned> runs,evts;
46 + };
47 + EvtsOfInterest::EvtsOfInterest(TString config)
48 + {
49 +  if(config!="") {
50 +    ifstream ifs(config);
51 +    assert(ifs.is_open());
52 +    string line;
53 +    while(getline(ifs,line)) {
54 +      stringstream ss(line);
55 +      unsigned run,evt;
56 +      ss >> run >> evt;
57 +      runs.push_back(run);
58 +      evts.push_back(evt);
59 +    }
60 +    ifs.close();
61 +  }
62 + }
63 + bool EvtsOfInterest::has(unsigned run, unsigned evt)
64 + {
65 +  for(unsigned i=0; i<runs.size(); i++) {
66 +    if(runs[i]==run && evts[i]==evt)
67 +      return true;
68 +  }
69 +  return false;
70 + }
71 + //----------------------------------------------------------------------------------------
72   TCanvas *can;
73 + //----------------------------------------------------------------------------------------
74 + class lepFrs
75 + {
76 + public:
77 +  double fwgt_3,fwgt_4;
78 +  double fwgt_lo_3,fwgt_lo_4;
79 +  double fwgt_hi_3,fwgt_hi_4;
80 + };
81 + //----------------------------------------------------------------------------------------
82 + class dr_struct
83 + {
84 + public:
85 +  dr_struct();
86 +  void check();
87 +  float l1j1,l2j1,l3j1,l4j1;
88 +  float l1j2,l2j2,l3j2,l4j2;
89 +  float l1j3,l2j3,l3j3,l4j3;
90 +  float l1j4,l2j4,l3j4,l4j4;
91 + };
92 + dr_struct::dr_struct()
93 + {  
94 +  l1j1 = l2j1 = l3j1 = l4j1 = l1j2 = l2j2 = l3j2 = l4j2 = l1j3 = l2j3 = l3j3 = l4j3 = l1j4 = l2j4 = l3j4 = l4j4 = -1;
95 + }
96 + void dr_struct::check()
97 + {
98 +  if(fabs(l1j1) < 0.2) { cout << "l1j1 too small: " << l1j1 << endl; }
99 +  if(fabs(l2j1) < 0.2) { cout << "l2j1 too small: " << l2j1 << endl; }
100 +  if(fabs(l3j1) < 0.2) { cout << "l3j1 too small: " << l3j1 << endl; }
101 +  if(fabs(l4j1) < 0.2) { cout << "l4j1 too small: " << l4j1 << endl; }
102 +  if(fabs(l1j2) < 0.2) { cout << "l1j2 too small: " << l1j2 << endl; }
103 +  if(fabs(l2j2) < 0.2) { cout << "l2j2 too small: " << l2j2 << endl; }
104 +  if(fabs(l3j2) < 0.2) { cout << "l3j2 too small: " << l3j2 << endl; }
105 +  if(fabs(l4j2) < 0.2) { cout << "l4j2 too small: " << l4j2 << endl; }
106 +  if(fabs(l1j3) < 0.2) { cout << "l1j3 too small: " << l1j3 << endl; }
107 +  if(fabs(l2j3) < 0.2) { cout << "l2j3 too small: " << l2j3 << endl; }
108 +  if(fabs(l3j3) < 0.2) { cout << "l3j3 too small: " << l3j3 << endl; }
109 +  if(fabs(l4j3) < 0.2) { cout << "l4j3 too small: " << l4j3 << endl; }
110 +  if(fabs(l1j4) < 0.2) { cout << "l1j4 too small: " << l1j4 << endl; }
111 +  if(fabs(l2j4) < 0.2) { cout << "l2j4 too small: " << l2j4 << endl; }
112 +  if(fabs(l3j4) < 0.2) { cout << "l3j4 too small: " << l3j4 << endl; }
113 +  if(fabs(l4j4) < 0.2) { cout << "l4j4 too small: " << l4j4 << endl; }
114 + }  
115  
116 + bool findGoodJets(vector<SimpleLepton> &goodJets, filestuff *fs,
117 +                  SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4);
118   void fillHist(CSample *cs, TString channel, TString type, TString var, double val, double wgt, double wgt_lo=0, double wgt_hi=0);
119 + void fillAllHists(FOFlags ctrl, CSample *cs, TString channel, TString type, filestuff *fs, KinematicsStruct kine,
120 +                  SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4,
121 +                  double wgt, double wgt_lo=0, double wgt_hi=0);
122 + void fillAllJetHists(FOFlags ctrl,  CSample *cs, TString channel, TString type, filestuff *fs, KinematicsStruct kine,
123 +                      SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4,
124 +                      FusionMva &fusion, vector<SimpleLepton> &goodJets, JetInfoStruct ji,
125 +                      double wgt, double wgt_lo=0, double wgt_hi=0);
126   void makeHTML(FOFlags &ctrl, TString type, TString plotLabel);
127   map<TString,map<TString,TH1D*>* > init_hists(FOFlags &ctrl, TString str="");
128 <
128 > bool passHlt(FOFlags &ctrl, TrigInfo ti, InfoStruct *info, unsigned lep1matchBits,
129 >             unsigned lep2matchBits, unsigned lep3matchBits,
130 >             unsigned lep4matchBits);
131 > dr_struct fill_dr_struct(vector<SimpleLepton> jets, SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4);
132 > void init_cuts(vector<TString> &cutstrs, map<TString,int> &cutvec);
133 > void fillLeptonFakeWeights(lepFrs &fwgts, FR_struct *fr, SimpleLepton lep3, SimpleLepton lep4);
134 > void resetCutVect(map<TString,int> &cutvec) {
135 >  for(map<TString,int>::iterator it=cutvec.begin(); it!=cutvec.end(); it++)
136 >    cutvec[(*it).first] = 0;
137 > }
138 > map<TString,TH1D*> setHistSet(CSample *cs, TString type, TString var);
139 > void shiftOverflows(map<TString,TH1D*> &hset);
140 > void scaleHists(map<TString,TH1D*> &hset);
141 > void fillFakeTuple(CSample *cs, filestuff *fs, EventData evtdata, unsigned ientry, double fillweight);
142   //----------------------------------------------------------------------------------------
143   int main(int argc, char** argv)
144   {
145 <  double lumi = 5000;//10093;
145 >  initPUWeights();
146 >  vector<TString> cutstrs;
147 >  map<TString,int> cutvec;
148 >  init_cuts(cutstrs, cutvec);
149 >
150 >  EvtsOfInterest eoiAdish("/tmp/adish-uniq");
151 >  EvtsOfInterest eoiMe("/tmp/me-uniq");
152  
153 +  // arguments...
154    FOFlags ctrl;
155    parse_foargs( argc, argv, ctrl );
156    ctrl.dump();
157 +  bool makeFakeTuples = (ctrl.extraArgs.Contains("makeFakeTuples"));
158 +  FusionMva fusion(ctrl.uncert);//"./weights/againstZZ-fusion_BDTG.weights.xml");
159 +  bool makeJetTuple=false; // jet tuple for vbf mva training
160 +  assert(ctrl.faketype=="SR"   || // plot events in signal region, if requested including ntuples of the extrapolation from 2P2F and 3P1F regions
161 +         ctrl.faketype=="2P2F" || // plot 2P2F region, and if requested write out ntuple which gives extrapolation to SR (ntuple would then be read by "SR", above)
162 +         ctrl.faketype=="3P1F");  // plot 3P1F region (including the extrapolation to there from the 2P2F region), and if requested write out ntuple with SR prediction
163 +  if(ctrl.heavyFlavor) assert(ctrl.faketype=="2P2F");
164 +  if(ctrl.ssof) assert(ctrl.faketype=="2P2F");
165  
166    can = new TCanvas("can","can");
167  
# Line 64 | Line 172 | int main(int argc, char** argv)
172    string xspath = (string(cmsswpath)+"/MitPhysics/data/xs.dat");
173    SimpleTable xstab(xspath.c_str());
174  
175 +  TrigInfo ti;
176 +  initAnalysisTriggers(ti);
177 +
178    vector<CSample*> samplev;
179    TString ntupledir(""),label(""),plotLabel(""),jsonFile("");
180    int puTarget;
181    readConfigFile(ctrl.config, ntupledir, label, plotLabel, jsonFile, puTarget, samplev, &ctrl, init_hists);
182  
183 <  vector<pair<unsigned,unsigned> > runEvtv; // vector to veto duplicate events
183 >  KinematicsStruct kine; // NOTE: this is not the same as fs->kine -- the one in fs is filled with whichever choice of leptons was made in applyZPlusX
184 >  JetInfoStruct ji;
185 >
186    UInt_t minRun=999999999,maxRun=0;
187    for(unsigned ics=0; ics<samplev.size(); ics++) {
188      CSample *cs = samplev[ics];
189      cout << cs->name << endl;
190 +
191 +    // output for vbf mva training
192 +    unsigned evtVar;
193 +    AngleTuple *jettuple;
194 +    if(makeJetTuple) {
195 +      jettuple = new AngleTuple( (const char*)(TString(cs->name)+"_pass.root").Data(), (const char*)"zznt");
196 +      jettuple->makeKinematicsBranch(kine);
197 +      jettuple->makeJetInfoBranch(ji);
198 +      const char *str1 = "evt";
199 +      const char *str2 = "evt/i";
200 +      jettuple->makeBranch(str1,(void *)(&evtVar),str2);
201 +    }
202 +
203 +    vector<pair<unsigned,unsigned> > runEvtv; // vector to veto duplicate events
204 +    resetCutVect(cutvec);
205      for(unsigned ifs=0; ifs<(cs->fsv).size(); ifs++) {
206        filestuff *fs = (cs->fsv)[ifs];
207 <      cout << "\t" << fs->fname_ << endl;
207 >      if(makeFakeTuples && fs->dataset_!="fakes") { // write a copy of the input trees for events that will be used for the SR fake prediction
208 >        TString fakefilename(fs->fname_);
209 >        fakefilename.ReplaceAll(".root","-fakes.root");
210 >        fs->makeOutputFile(fakefilename);
211 >      }
212 >      cout << "\t" << fs->fname_; cout.flush();
213        FR_struct *fr = (fs->era_==2011) ? &fr2011 : &fr2012;
214        unsigned nDuplSkipped=0;
215        for(unsigned ientry=0; ientry<fs->getentries("FOtree"); ientry++) {
216 <        fs->getentry(ientry,"","FOtree");
216 >        fs->getentry(ientry,"FOs","FOtree");
217 >        fs->getentry(ientry,"Zleptons","FOtree");
218          fs->getentry(ientry,"info","zznt");
219  
220 <        if(fs->isdata_) {
220 >        // if(eoiAdish.has(fs->info->run, fs->info->evt) &&
221 >        //    !eoiMe.has(fs->info->run, fs->info->evt)) {
222 >        //   cout << "found an event: " << fs->info->run << " " << fs->info->evt << endl;
223 >        //   ctrl.debug = true;
224 >        // } else
225 >        //   ctrl.debug = false;
226 >
227 >        cutvec["start"] = fs->total_entries_;
228 >        if(fs->isdata_) {// && !(fs->dataset_=="fakes")) {
229 >          // WARNING: when I write the fakeTuple from here, I fill an event *twice* if it has two fake z2 candidates, so event numbers are duplicated.
230 >          // in plotH4l I *don't* want to remove these, but here I *do*, because in this file I again loop over all z2 candidates
231            setMinMaxRun(fs->info->run, minRun, maxRun);
232            bool dupl = takeCareOfDuplicateEvents(fs->info->run, fs->info->evt, runEvtv, nDuplSkipped);
233 <          if(dupl) continue;
234 <          pair<unsigned,unsigned> rl(fs->info->run, fs->info->lumi);
91 <          if(!fs->rlrm_.HasRunLumi(rl)) continue;
233 >          if(dupl) { if(ctrl.debug) cout << " failing dupl" << endl; continue; }
234 >          if(!fs->rlrm_.HasRunLumi(pair<unsigned,unsigned> (fs->info->run, fs->info->lumi))) { if(ctrl.debug) cout << " failing rlrm" << endl; continue; }
235          }
236 +        cutvec["rlrmAndDupl"] += 1;
237  
238          double wgt=1;
239          if(!fs->isdata_) {
240 <          double xsWgt = lumi*xstab.Get(fs->dataset_)/fs->total_entries_;
241 <          // double puWgt = weightTrue2012(fs->info->npu);
242 <          wgt = xsWgt;//*puWgt;
240 >          double xsWgt = fs->lumi_*xstab.Get(fs->dataset_)/fs->total_entries_;
241 >          double puWgt = getPUWeight(fs->era_,fs->info->npu);
242 >          if(fs->era_==2012)
243 >            puWgt = weightTrue2012(fs->info->npu);
244 >          wgt = xsWgt*puWgt;
245 >        }
246 >        if(fs->dataset_=="fakes") { // if we're reading in the fake ntuple that we wrote out earlier, use the event weight that was stored (this should be the same as what we'd get if we calculated it below)
247 >          fs->getentry(ientry,"weights","zznt");
248 >          wgt = fs->weights->w;
249          }
250  
251          unsigned npass = fs->passingL->size();
252          unsigned nfail = fs->failingL->size();
253  
254 +        // look for a z1
255          pair<int,int> best_z_indices;
256 <        findZ1(fs->passingL,best_z_indices,40);
256 >        double best_mz1 = findZ1(fs->passingL,best_z_indices,40);
257 >        if(best_mz1<=0) {
258 >          cout << "WARNING: best z1 mass: " << best_mz1 << endl;
259 >        }
260 >        if(ctrl.heavyFlavor && best_mz1<60) { if(ctrl.debug) cout << " failing mz1" << endl; continue; }
261          EventData evtdata;
262          evtdata.Z1leptons.push_back((*fs->passingL)[best_z_indices.first]);
263          evtdata.Z1leptons.push_back((*fs->passingL)[best_z_indices.second]);
# Line 110 | Line 265 | int main(int argc, char** argv)
265          // Look for a z2
266          vector<pair<SimpleLepton,SimpleLepton> > z2cands;
267          vector<TString> types;
268 <        findZ2Candidates(fs->failingL, z2cands, types, true);
269 <        
268 >        TString signFlavor;
269 >        if(ctrl.heavyFlavor) signFlavor = "";
270 >        else if(ctrl.ssof)   signFlavor = "SS-SF";
271 >        else                 signFlavor = "OS-SF";
272 >        double minMz2 = ctrl.heavyFlavor ? 0 : 12;
273 >        int hiPtPassingZ2 = findZ2CandidatesPassFail(fs->passingL, z2cands, types, signFlavor, minMz2, "pass",&best_z_indices); // signal region
274 >        int hiPtFailingZ2 = findZ2CandidatesPassFail(fs->failingL, z2cands, types, signFlavor, minMz2, "fail"); // 2P+2F
275 >        int hiPtPfZ2      = findZ2CandidatesPassFail(fs->passingL, z2cands, types, signFlavor, minMz2, "PF", &best_z_indices, fs->failingL); // 3P+1F
276 >
277 >        if(ctrl.debug) {
278 >          cout << "z2 cands: " << endl;
279 >          for(unsigned iz2=0; iz2<z2cands.size(); iz2++)
280 >            cout << " ---" << types[iz2] << "'ing z2" << endl;
281 >        }
282 >
283          for(unsigned iz2=0; iz2<z2cands.size(); iz2++) {
284            SimpleLepton lep1 = evtdata.Z1leptons[0];
285            SimpleLepton lep2 = evtdata.Z1leptons[1];
286            SimpleLepton lep3 = z2cands[iz2].first;
287            SimpleLepton lep4 = z2cands[iz2].second;
120            
121          TString channel = getChannel(lep1,lep2,lep3,lep4);
288            evtdata.Z2leptons.clear();
289            evtdata.Z2leptons.push_back(lep3);
290            evtdata.Z2leptons.push_back(lep4);
125          KinematicsStruct kine;
291            fillKinematics(evtdata,kine);
292 +          lepFrs fwgts;
293 +          fillLeptonFakeWeights(fwgts, fr, lep3, lep4);
294 +          TString channel = getChannel(lep1,lep2,lep3,lep4);
295 +          // bool hltPass = fs->isdata_ ? passHlt(ctrl,ti,fs->info,lep1.scID,lep2.scID,lep3.scID,lep4.scID) : true;
296  
297 <          if(fabs(lep3.ip3dSig) >= 4 || fabs(lep4.ip3dSig) >= 4)        continue; // |sip| < 4 is applied in ZPlusX to lep1, lep2
298 <          if(!finalLeptonPtReqs(lep1,lep2,lep3,lep4))                   continue;
299 <          if(!leptonDrReqs(lep1,lep2,lep3,lep4))                        continue;
300 <          if(!resonanceKilling(lep1,lep2,lep3,lep4))                    continue;
301 <          if(kine.m4l < 100)                                            continue;
302 <
303 <          if(types[iz2]=="pass") {
304 <            fillHist( cs, channel, "obs", "run"      , fs->info->run  , wgt);
305 <            fillHist( cs, channel, "obs", "mZ1"      , kine.mZ1  ,      wgt);                
306 <            fillHist( cs, channel, "obs", "mZ2"      , kine.mZ2  ,      wgt);                
307 <            fillHist( cs, channel, "obs", "mZ2_lo"   , kine.mZ2  ,      wgt);                
308 <            fillHist( cs, channel, "obs", "m4l"      , kine.m4l  ,      wgt);                
309 <            fillHist( cs, channel, "obs", "Z1pt"     , kine.Z1pt ,      wgt);                
310 <            fillHist( cs, channel, "obs", "ZZpt"     , kine.ZZpt ,      wgt);                
311 <            fillHist( cs, channel, "obs", "ip3ds_l1" , lep1.ip3dSig ,   wgt);
312 <            fillHist( cs, channel, "obs", "ip3ds_l2" , lep2.ip3dSig ,   wgt);
313 <            fillHist( cs, channel, "obs", "ip3ds_l3" , lep3.ip3dSig ,   wgt);
314 <            fillHist( cs, channel, "obs", "ip3ds_l4" , lep4.ip3dSig ,   wgt);
315 <            fillHist( cs, channel, "obs", "pt_l1"    , lep1.vec.Pt() ,  wgt);
316 <            fillHist( cs, channel, "obs", "pt_l2"    , lep2.vec.Pt() ,  wgt);
317 <            fillHist( cs, channel, "obs", "pt_l3"    , lep3.vec.Pt() ,  wgt);
318 <            fillHist( cs, channel, "obs", "pt_l4"    , lep4.vec.Pt() ,  wgt);
319 <            fillHist( cs, channel, "obs", "eta_l1"   , lep1.vec.Eta() , wgt);
320 <            fillHist( cs, channel, "obs", "eta_l2"   , lep2.vec.Eta() , wgt);
321 <            fillHist( cs, channel, "obs", "eta_l3"   , lep3.vec.Eta() , wgt);
322 <            fillHist( cs, channel, "obs", "eta_l4"   , lep4.vec.Eta() , wgt);
323 <            fillHist( cs, channel, "obs", "dR_l3l4_lo",dr(lep3,lep4),   wgt);
324 <            fillHist( cs, channel, "obs", "dR_l3l4",   dr(lep3,lep4),   wgt);
325 <          } else if(types[iz2]=="fail") {
326 <            double fwgt_1    = get_fake_weight("",lep3,*fr);
327 <            double fwgt_lo_1 = get_fake_weight("lo",lep3,*fr);
328 <            double fwgt_hi_1 = get_fake_weight("hi",lep3,*fr);
329 <            double fwgt_2    = get_fake_weight("",lep4,*fr);
330 <            double fwgt_lo_2 = get_fake_weight("lo",lep4,*fr);
331 <            double fwgt_hi_2 = get_fake_weight("hi",lep4,*fr);
332 <
333 <            double centr = wgt*fwgt_1*fwgt_2;
334 <            double lo    = wgt*fwgt_lo_1*fwgt_lo_2;
335 <            double hi    = wgt*fwgt_hi_1*fwgt_hi_2;
336 <            fillHist( cs, channel, "pred", "run",       fs->info->run,  centr, lo, hi);
337 <            fillHist( cs, channel, "pred", "mZ1",       kine.mZ1,       centr, lo, hi);
338 <            fillHist( cs, channel, "pred", "mZ2",       kine.mZ2,       centr, lo, hi);
339 <            fillHist( cs, channel, "pred", "mZ2_lo",    kine.mZ2,       centr, lo, hi);
340 <            fillHist( cs, channel, "pred", "m4l",       kine.m4l,       centr, lo, hi);
341 <            fillHist( cs, channel, "pred", "Z1pt",      kine.Z1pt,      centr, lo, hi);
342 <            fillHist( cs, channel, "pred", "ZZpt",      kine.ZZpt,      centr, lo, hi);
343 <            fillHist( cs, channel, "pred", "ip3ds_l1",  lep1.ip3dSig,   centr, lo, hi);
344 <            fillHist( cs, channel, "pred", "ip3ds_l2",  lep2.ip3dSig,   centr, lo, hi);
345 <            fillHist( cs, channel, "pred", "ip3ds_l3",  lep3.ip3dSig,   centr, lo, hi);
346 <            fillHist( cs, channel, "pred", "ip3ds_l4",  lep4.ip3dSig,   centr, lo, hi);
347 <            fillHist( cs, channel, "pred", "pt_l1",     lep1.vec.Pt(),  centr, lo, hi);
348 <            fillHist( cs, channel, "pred", "pt_l2",     lep2.vec.Pt(),  centr, lo, hi);
349 <            fillHist( cs, channel, "pred", "pt_l3",     lep3.vec.Pt(),  centr, lo, hi);
350 <            fillHist( cs, channel, "pred", "pt_l4",     lep4.vec.Pt(),  centr, lo, hi);
351 <            fillHist( cs, channel, "pred", "eta_l1",    lep1.vec.Eta(), centr, lo, hi);
352 <            fillHist( cs, channel, "pred", "eta_l2",    lep2.vec.Eta(), centr, lo, hi);
353 <            fillHist( cs, channel, "pred", "eta_l3",    lep3.vec.Eta(), centr, lo, hi);
354 <            fillHist( cs, channel, "pred", "eta_l4",    lep4.vec.Eta(), centr, lo, hi);
355 <            fillHist( cs, channel, "pred", "dR_l3l4_lo",dr(lep3,lep4),  centr, lo, hi);
356 <            fillHist( cs, channel, "pred", "dR_l3l4",   dr(lep3,lep4),  centr, lo, hi);
297 >          if(ctrl.heavyFlavor) {
298 >            if(fabs(lep1.dz) > 0.01 || fabs(lep2.dz) > 0.01)            { cout << " failing XXX" << endl; continue; }
299 >            if( !(fabs(lep3.ip3dSig) > 2 && fabs(lep4.ip3dSig) > 8) &&
300 >                !(fabs(lep4.ip3dSig) > 2 && fabs(lep3.ip3dSig) > 8)  )  { cout << " failing XXX" << endl; continue; }
301 >            if( !(fabs(lep3.d0) > .002 && fabs(lep4.d0) > .01) &&
302 >                !(fabs(lep4.d0) > .002 && fabs(lep3.d0) > .01)  )       { cout << " failing XXX" << endl; continue; }
303 >            if(fabs(lep3.d0) > 0.1 || fabs(lep4.d0) > 0.1)              { cout << " failing XXX" << endl; continue; }
304 >            if(fabs(lep3.dz) > 0.2 || fabs(lep4.dz) > 0.2)              { cout << " failing XXX" << endl; continue; }
305 >          } else {
306 >            if(fabs(lep3.ip3dSig) >= 4 || fabs(lep4.ip3dSig) >= 4)      { if(ctrl.debug) cout << " failing sip" << endl; continue; } // |sip| < 4 is applied to lep1, lep2 in ZPlusX
307 >            if(!leptonDrReqs(lep1,lep2,lep3,lep4))                      { if(ctrl.debug) cout << " failing dr" << endl; continue; }
308 >            if(kine.m4l < 100)                                          { if(ctrl.debug) cout << " failing m4l 100" << endl; continue; }
309 >            if(!resonanceKilling(lep1,lep2,lep3,lep4))                  { if(ctrl.debug) cout << " failing reskill" << endl; continue; }
310 >          }
311 >          if(!finalLeptonPtReqs(lep1,lep2,lep3,lep4))                   { if(ctrl.debug) cout << " failing lepton pt " << endl; continue; }
312 >
313 >          cutvec["4lsele"] += 1;
314 >
315 >          //
316 >          // jets!
317 >          //
318 >          bool doJets=false;
319 >          vector<SimpleLepton> goodJets;
320 >          if(doJets) {
321 >            fs->getentry(ientry,"jets","FOtree");
322 >            if(fs->jets->size() > 1) cutvec["twoJets"] += 1;
323 >            bool hazJetz = findGoodJets(goodJets, fs, lep1, lep2, lep3, lep4);
324 >            if(goodJets.size() > 1) cutvec["twoJetsAfter"] += 1;
325 >            if(!hazJetz) { cout << " failing XXX" << endl; continue; }
326 >            ControlFlags ctrlTmp;
327 >            fillJetInfo(goodJets,ji,ctrlTmp);
328 >            if(makeJetTuple) {
329 >              evtVar = fs->info->evt;
330 >              jettuple->Fill();
331 >            }
332 >            fusion.setValues(ji,kine);
333 >          }
334 >
335 >          cutvec["filling"] += 1;
336 >
337 >          if(types[iz2]=="pass" || ((fs->dataset_=="fakes") && (types[iz2]=="fail"))) { // for the fakes ntuple, we plot the "failing" z2s as passing ones, i.e. put them in hObs, but with the proper weight
338 >            // llll events
339 >            if((fs->dataset_!="fakes") && (iz2!=hiPtPassingZ2)) { cout << " failing XXX" << endl; continue; } // for non-fakes, we only want the  highest-pt "passing" z2 (for the fakes these are by design failing z2s, and we want all of 'em)
340 >            cutvec["fillPass"] += 1;
341 >            fillAllHists( ctrl, cs, channel, "obs", fs, kine, lep1, lep2, lep3, lep4, wgt);
342 >            if(doJets) fillAllJetHists( ctrl, cs, channel, "obs", fs, kine, lep1, lep2, lep3, lep4, fusion, goodJets, ji, wgt);
343 >
344 >          } else if(types[iz2]=="fail") { // lljj events
345 >            if(iz2!=hiPtFailingZ2 && ctrl.plotWholeSample) { cout << " failing XXX" << endl; continue; }
346 >            double centr=0,lo=0,hi=0;
347 >            if(ctrl.faketype=="2P2F") {
348 >              if(ctrl.plotWholeSample) { // make plots without weighting with the FR
349 >                centr = wgt;
350 >              } else {
351 >                centr = wgt*fwgts.fwgt_3*fwgts.fwgt_4;
352 >                lo    = wgt*fwgts.fwgt_lo_3*fwgts.fwgt_lo_4;
353 >                hi    = wgt*fwgts.fwgt_hi_3*fwgts.fwgt_hi_4;
354 >              }
355 >            } else if(ctrl.faketype=="3P1F") { // apply the weights for extrapolation from 2P2F region to 3P1F region
356 >              centr = wgt*(fwgts.fwgt_3    + fwgts.fwgt_4);
357 >              lo    = wgt*(fwgts.fwgt_lo_3 + fwgts.fwgt_lo_4);
358 >              hi    = wgt*(fwgts.fwgt_hi_3 + fwgts.fwgt_hi_4);
359 >            }
360 >
361 >            // runfile << fixed << setprecision(0) << fs->info->run << " " << fixed << setprecision(0) << fs->info->evt << endl;
362 >
363 >            if(makeFakeTuples && fs->dataset_!="fakes") {
364 >              double fillwgt;
365 >              bool dofill=true;
366 >              if(ctrl.faketype=="2P2F") {
367 >                assert(fs->isdata_);
368 >                fillwgt = centr;
369 >              } else if(ctrl.faketype=="3P1F") {
370 >                fillwgt = - wgt*(fwgts.fwgt_3*fwgts.fwgt_4 + fwgts.fwgt_4*fwgts.fwgt_3);
371 >                dofill = (cs->name=="data" || cs->name=="2012" || cs->name=="2011");
372 >              } else assert(0);
373 >              fillFakeTuple(cs,fs,evtdata,ientry,fillwgt);
374 >            }
375 >            cutvec["fillFail"] += 1;
376 >            fillAllHists( ctrl, cs, channel, "pred", fs, kine, lep1, lep2, lep3, lep4, centr, lo, hi);
377 >            if(doJets) fillAllJetHists( ctrl, cs, channel, "pred", fs, kine, lep1, lep2, lep3, lep4, fusion, goodJets, ji, centr, lo, hi);
378 >
379 >          } else if(types[iz2]=="PF") { // lllj events
380 >            if(fs->dataset_!="fakes" && iz2!=hiPtPfZ2) continue; // shouldn't have this when I extrapolate to signal region
381 >
382 >            if(makeFakeTuples && fs->dataset_!="fakes" && ctrl.faketype=="3P1F") {
383 >              assert(lep3.isLoose || lep4.isLoose);
384 >              assert(!(lep3.isLoose && lep4.isLoose));
385 >              double fwgt = lep3.isLoose ? fwgts.fwgt_4 : fwgts.fwgt_3;
386 >              double fillwgt;// only fill for data or ZZ
387 >              if(fs->isdata_)           fillwgt = wgt*fwgt;     // set the weight to be the fake weight (maybe check that this is the same weight as I read in from the fakes dataset?)
388 >              else if(cs->name=="zz")   fillwgt = -wgt*fwgt;    // subtract off the zz contribution in the 3P1F region
389 >              else assert(0);
390 >              fillFakeTuple(cs,fs,evtdata,ientry,fillwgt);
391 >            }
392 >            
393 >            cutvec["fillPF"] += 1;
394 >            fillAllHists( ctrl, cs, channel, "PF", fs, kine, lep1, lep2, lep3, lep4, wgt);
395 >            if(doJets) fillAllJetHists( ctrl, cs, channel, "PF", fs, kine, lep1, lep2, lep3, lep4, fusion, goodJets, ji, wgt);
396            } else assert(0);
397          }
398        }
399 <      cout << "\t\tWARNING: skipped " << nDuplSkipped << " duplicate events" << endl;
399 >      cout << "\t\t" << setw(7) << nDuplSkipped << " duplicate events skipped" << endl;
400 >      if(makeFakeTuples) {
401 >        fs->outFotuple->getFile()->cd();
402 >        fs->outFotuple->getTree()->Write();
403 >        fs->outtuple->WriteClose();
404 >      }
405      }
406 +    for(unsigned icut=0; icut<cutstrs.size(); icut++) {
407 +      cout << setw(22) << cutstrs[icut] << setw(22) << cutvec[cutstrs[icut]] << endl;
408 +    }
409 +    if(makeJetTuple)
410 +      jettuple->WriteClose();
411    }
412    cout << "run range: " << setw(12) << minRun << setw(12) << maxRun << endl;
413 +
414 +  // plot!
415    assert(samplev.size()>0);
416    vector<TString> typev;
417    typev.push_back("all");
418 <  typev.push_back("4e");
419 <  typev.push_back("4m");
420 <  typev.push_back("2e2m");
418 >
419 >  bool plotAllCats=false;
420 >  if(plotAllCats) {
421 >    typev.push_back("4e");
422 >    typev.push_back("4m");
423 >    typev.push_back("2e2m");
424 >  }
425 >
426 >  // loop over channels
427    for(unsigned itype=0; itype<typev.size(); itype++) {
428      TString type(typev[itype]);
429      gSystem->mkdir(ctrl.outdir+"/"+plotLabel+"/"+type,true);
430      TFile runHistFile(ctrl.outdir+"/"+plotLabel+"/"+type+"/runs.root","recreate");
431 +    // loop over variables
432      map<TString,TH1D*>::iterator it_v;
433      for(it_v=(*(samplev[0]->hists)[type+"_obs"]).begin(); it_v!=(*(samplev[0]->hists)[type+"_obs"]).end(); it_v++) {
434        TString var((*it_v).first);
435        CPlot cplot(var,"",(*(samplev[0]->hists)[type+"_obs"])[var]->GetXaxis()->GetTitle(),"events",ctrl.outdir+"/"+plotLabel+"/"+type+"/plots");
436 +      double fakeCounter=0;
437        double ymax=0;
438 +      // loop over csamples
439        for(unsigned isam=0; isam<samplev.size(); isam++) {
440          CSample *cs = samplev[isam];
441 <        TH1D *hObs = (*(cs->hists)[type+"_obs"])[var];
442 <        TH1D *hPred = (*(cs->hists)[type+"_pred"])[var];
443 <        TH1D *hPred_lo = (*(cs->hists)[type+"_pred_lo"])[var];
444 <        TH1D *hPred_hi = (*(cs->hists)[type+"_pred_hi"])[var];
445 <
446 <        if(cs->isdata) {
447 <          double ratio = integrateHist(hObs) / integrateHist(hPred);
448 <          cplot.AddHist1D(hObs,cs->legend+": "+integral_str(hObs,0)+", ratio: "+f_to_a(ratio),"E");
449 <          cplot.AddHist1D(hPred,"FR predic: "+integral_str(hPred,2),"hist",kRed);
450 <          cplot.AddHist1D(hPred_lo,"stat lo: "+integral_str(hPred_lo,2),"hist",kRed,kDashed);
451 <          cplot.AddHist1D(hPred_hi,"stat hi: "+integral_str(hPred_hi,2),"hist",kRed,kDashed);
452 <          ymax = max(ymax,histMax(hObs,hPred,hPred_lo,hPred_hi));
453 <        } else {
454 <          cplot.AddToStack(hPred,cs->legend+" pred.: "+integral_str(hPred,2),cs->color);
455 <          ymax = max(ymax,hPred->GetMaximum());
456 <          // cplot.AddHist1D(hObs,cs->legend+": "+integral_str(hObs),"Ehist",kBlue);//cs->color);
441 >        map<TString,TH1D*> hs = setHistSet(cs,type,var);
442 >        if(var=="met")  //var=="mjj" || var=="etaProd" || var=="ncj")
443 >          shiftOverflows(hs);
444 >        bool plotShapeOnly = false;
445 >        if(plotShapeOnly)
446 >          scaleHists(hs);
447 >        
448 >        if(cs->isdata && !(cs->name=="fakes")) {
449 >          if(ctrl.faketype=="SR" || (ctrl.heavyFlavor && !ctrl.plotWholeSample)) {
450 >            cplot.AddHist1D(hs["obs"],cs->legend+": "+integral_str(hs["obs"],0),"E");
451 >            ymax = max(ymax,hs["obs"]->GetMaximum());
452 >          }
453 >          if(ctrl.faketype=="2P2F") {
454 >            if(ctrl.plotWholeSample) {
455 >              cplot.AddHist1D(hs["pred"],cs->legend+": "+integral_str(hs["pred"],0),"E");
456 >              ymax = max(ymax,hs["pred"]->GetMaximum());
457 >            } else {
458 >              cplot.AddHist1D(hs["pred"],   "FR predic: "+integral_str(hs["pred"],3),   "hist",kRed);
459 >              cplot.AddHist1D(hs["pred_lo"],  "stat lo: "+integral_str(hs["pred_lo"],2),"hist",kRed,kDashed);
460 >              cplot.AddHist1D(hs["pred_hi"],  "stat hi: "+integral_str(hs["pred_hi"],2),"hist",kRed,kDashed);
461 >              ymax = max(ymax,histMax(hs["pred"],hs["pred_lo"],hs["pred_hi"]));
462 >            }
463 >          }
464 >          if(ctrl.faketype=="3P1F") {
465 >            double ratio = 0;//integrateHist(hs["PF"]) / integrateHist(hs["pred"]);
466 >            cplot.AddHist1D(hs["PF"],cs->legend+": "+integral_str(hs["PF"],0)+", ratio: "+f_to_a(ratio),"E",cs->color);
467 >            if(!ctrl.plotWholeSample) {
468 >              cplot.AddToStack(hs["pred"],cs->legend+"(2P2F extrap.): "+integral_str(hs["pred"],1), kRed);
469 >              ymax = max(ymax,histMax(hs["PF"],cplot.GetStack(),hs["pred"]));
470 >            } else {
471 >              ymax = max(ymax,hs["PF"]->GetMaximum());
472 >            }
473 >            fakeCounter = integrateHist(hs["PF"]) - integrateHist(hs["pred"]);
474 >          }
475 >        } else { // mc
476 >          TH1D *mchist = 0;
477 >          if(ctrl.faketype=="SR")   mchist = hs["obs"];
478 >          if(ctrl.faketype=="2P2F") {
479 >            if(ctrl.plotWholeSample)  mchist = hs["pred"];
480 >            else                      mchist = hs["obs"];
481 >          }
482 >          if(ctrl.faketype=="3P1F")   mchist = hs["PF"];
483 >          assert(mchist);
484 >          bool stackMc = true;
485 >          if(stackMc) cplot.AddToStack(mchist,cs->legend+": "+integral_str(mchist,2),cs->color);
486 >          else        cplot.AddHist1D(mchist,cs->legend +": "+integral_str(mchist,3),"Ehist",cs->color);//cs->color);
487 >          if(cs->name=="zz") fakeCounter -= integrateHist(hs["PF"]);
488 >          ymax = max(ymax,histMax(mchist,cplot.GetStack()));
489          }
490          if(cs->isdata && var=="run") {
491 <          assert(hObs->GetBinContent(0)==0 && hObs->GetBinContent(hObs->GetXaxis()->GetNbins()+1)==0);
492 <          hObs->Write("runs");
491 >          assert(hs["obs"]->GetBinContent(0)==0 && hs["obs"]->GetBinContent(hs["obs"]->GetXaxis()->GetNbins()+1)==0);
492 >          hs["obs"]->Write("runs");
493          }
494 <        if(var.Contains("ip3d")) {
234 <          double maxVal = max(hObs->GetMaximum(), hPred->GetMaximum());
235 <          double isoMin = 0.0005*maxVal;
236 <          double isoMax = 1.8*maxVal;
237 <          cplot.SetYRange(isoMin,isoMax);
238 <          // cplot.SetLogy();
239 <        } else {
494 >        if(!var.Contains("fusionMVA"))
495            cplot.SetYRange(0,1.2*ymax);
496 <        }
496 >        if(ctrl.faketype=="3P1F") cplot.AddTextBox("diff: "+f_to_a(fakeCounter),.7,.5,.8,.6);
497        }
498        cplot.Draw(can,true,"png");
499 +      if(var.Contains("fusionMVA")) {
500 +        cplot.SetLogy();
501 +        cplot.Draw(can,true,"png",0,var+"-log");
502 +        cplot.SetLogy(false);
503 +      }
504      }
505      runHistFile.Close();
506      makeHTML(ctrl,type,plotLabel);
# Line 254 | Line 514 | map<TString,map<TString,TH1D*>* > init_h
514    map<TString,TH1D*> *h_4e_pred         = new map<TString,TH1D*>;      hists["4e_pred"]         = h_4e_pred;
515    map<TString,TH1D*> *h_4e_pred_lo      = new map<TString,TH1D*>;      hists["4e_pred_lo"]      = h_4e_pred_lo;
516    map<TString,TH1D*> *h_4e_pred_hi      = new map<TString,TH1D*>;      hists["4e_pred_hi"]      = h_4e_pred_hi;
517 +  map<TString,TH1D*> *h_4e_PF           = new map<TString,TH1D*>;      hists["4e_PF"]           = h_4e_PF;
518 +  map<TString,TH1D*> *h_4e_PF_lo        = new map<TString,TH1D*>;      hists["4e_PF_lo"]        = h_4e_PF_lo;
519 +  map<TString,TH1D*> *h_4e_PF_hi        = new map<TString,TH1D*>;      hists["4e_PF_hi"]        = h_4e_PF_hi;
520    map<TString,TH1D*> *h_4m_obs          = new map<TString,TH1D*>;      hists["4m_obs"]          = h_4m_obs;
521    map<TString,TH1D*> *h_4m_pred         = new map<TString,TH1D*>;      hists["4m_pred"]         = h_4m_pred;
522    map<TString,TH1D*> *h_4m_pred_lo      = new map<TString,TH1D*>;      hists["4m_pred_lo"]      = h_4m_pred_lo;
523    map<TString,TH1D*> *h_4m_pred_hi      = new map<TString,TH1D*>;      hists["4m_pred_hi"]      = h_4m_pred_hi;
524 +  map<TString,TH1D*> *h_4m_PF           = new map<TString,TH1D*>;      hists["4m_PF"]           = h_4m_PF;
525 +  map<TString,TH1D*> *h_4m_PF_lo        = new map<TString,TH1D*>;      hists["4m_PF_lo"]        = h_4m_PF_lo;
526 +  map<TString,TH1D*> *h_4m_PF_hi        = new map<TString,TH1D*>;      hists["4m_PF_hi"]        = h_4m_PF_hi;
527    map<TString,TH1D*> *h_2e2m_obs        = new map<TString,TH1D*>;      hists["2e2m_obs"]        = h_2e2m_obs;
528    map<TString,TH1D*> *h_2e2m_pred       = new map<TString,TH1D*>;      hists["2e2m_pred"]       = h_2e2m_pred;
529    map<TString,TH1D*> *h_2e2m_pred_lo    = new map<TString,TH1D*>;      hists["2e2m_pred_lo"]    = h_2e2m_pred_lo;
530    map<TString,TH1D*> *h_2e2m_pred_hi    = new map<TString,TH1D*>;      hists["2e2m_pred_hi"]    = h_2e2m_pred_hi;
531 +  map<TString,TH1D*> *h_2e2m_PF         = new map<TString,TH1D*>;      hists["2e2m_PF"]         = h_2e2m_PF;
532 +  map<TString,TH1D*> *h_2e2m_PF_lo      = new map<TString,TH1D*>;      hists["2e2m_PF_lo"]      = h_2e2m_PF_lo;
533 +  map<TString,TH1D*> *h_2e2m_PF_hi      = new map<TString,TH1D*>;      hists["2e2m_PF_hi"]      = h_2e2m_PF_hi;
534    map<TString,TH1D*> *h_all_obs         = new map<TString,TH1D*>;      hists["all_obs"]         = h_all_obs;
535    map<TString,TH1D*> *h_all_pred        = new map<TString,TH1D*>;      hists["all_pred"]        = h_all_pred;
536    map<TString,TH1D*> *h_all_pred_lo     = new map<TString,TH1D*>;      hists["all_pred_lo"]     = h_all_pred_lo;
537    map<TString,TH1D*> *h_all_pred_hi     = new map<TString,TH1D*>;      hists["all_pred_hi"]     = h_all_pred_hi;
538 +  map<TString,TH1D*> *h_all_PF          = new map<TString,TH1D*>;      hists["all_PF"]          = h_all_PF;
539 +  map<TString,TH1D*> *h_all_PF_lo       = new map<TString,TH1D*>;      hists["all_PF_lo"]       = h_all_PF_lo;
540 +  map<TString,TH1D*> *h_all_PF_hi       = new map<TString,TH1D*>;      hists["all_PF_hi"]       = h_all_PF_hi;
541    map<TString,map<TString,TH1D*>* >::iterator it_h;
542  
543    for(it_h=hists.begin(); it_h!=hists.end(); it_h++) {
544 +    (*((*it_h).second))["npv"]          = new TH1D(TString("npv") +"_"+(*it_h).first+str,";#bf{N PV};",         50,-0.5,49.5);      (*((*it_h).second))["npv"]->Sumw2();
545      (*((*it_h).second))["run"]          = new TH1D(TString("run") +"_"+(*it_h).first+str,";#bf{run};",          50,160800,196535);   (*((*it_h).second))["run"]->Sumw2();
546 <    (*((*it_h).second))["mZ1"]          = new TH1D(TString("mZ1") +"_"+(*it_h).first+str,";#bf{mZ1 [GeV]};",    20,60,120);          (*((*it_h).second))["mZ1"]->Sumw2();
546 >    (*((*it_h).second))["mZ1"]          = new TH1D(TString("mZ1") +"_"+(*it_h).first+str,";#bf{mZ1 [GeV]};",    20,40,120);          (*((*it_h).second))["mZ1"]->Sumw2();
547      (*((*it_h).second))["mZ2"]          = new TH1D(TString("mZ2") +"_"+(*it_h).first+str,";#bf{mZ2 [GeV]};",    30,0,115);           (*((*it_h).second))["mZ2"]->Sumw2();
548      (*((*it_h).second))["mZ2_lo"]       = new TH1D(TString("mZ2_lo") +"_"+(*it_h).first+str,";#bf{mZ2 [GeV]};", 30,0,10);           (*((*it_h).second))["mZ2_lo"]->Sumw2();
549 <    (*((*it_h).second))["m4l"]          = new TH1D(TString("m4l") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};",    20,50,500);          (*((*it_h).second))["m4l"]->Sumw2();
549 >    (*((*it_h).second))["m4l_lo"]       = new TH1D(TString("m4l_lo") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};", 35,70,120);          (*((*it_h).second))["m4l_lo"]->Sumw2();
550 >    (*((*it_h).second))["m4l_med"]      = new TH1D(TString("m4l_med") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};",35,100,180);          (*((*it_h).second))["m4l_med"]->Sumw2();
551 >    (*((*it_h).second))["m4l"]          = new TH1D(TString("m4l") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};",    35,100,600);          (*((*it_h).second))["m4l"]->Sumw2();
552      (*((*it_h).second))["Z1pt"]         = new TH1D(TString("Z1pt") +"_"+(*it_h).first+str,";#bf{Z1pt [GeV]};",  20,0,200);           (*((*it_h).second))["Z1pt"]->Sumw2();
553      (*((*it_h).second))["ZZpt"]         = new TH1D(TString("ZZpt") +"_"+(*it_h).first+str,";#bf{ZZpt [GeV]};",  30,0,200);           (*((*it_h).second))["ZZpt"]->Sumw2();
554 +    (*((*it_h).second))["met"]          = new TH1D(TString("met") +"_"+(*it_h).first+str,";#bf{met [GeV]};",    30,0,100);           (*((*it_h).second))["met"]->Sumw2();
555  
556      (*((*it_h).second))["pt_l1"]        = new TH1D(TString("pt_l1") +"_"+(*it_h).first+str,";#bf{l1 p_{T} [GeV]};",     37,0,120);  (*((*it_h).second))["pt_l1"]->Sumw2();
557      (*((*it_h).second))["pt_l2"]        = new TH1D(TString("pt_l2") +"_"+(*it_h).first+str,";#bf{l2 p_{T} [GeV]};",     37,0,65);  (*((*it_h).second))["pt_l2"]->Sumw2();
# Line 287 | Line 563 | map<TString,map<TString,TH1D*>* > init_h
563      (*((*it_h).second))["eta_l3"]       = new TH1D(TString("eta_l3")+"_"+(*it_h).first+str,";#bf{l3 #eta};",            25,-2.5,2.5);     (*((*it_h).second))["eta_l3"]->Sumw2();
564      (*((*it_h).second))["eta_l4"]       = new TH1D(TString("eta_l4")+"_"+(*it_h).first+str,";#bf{l4 #eta};",            25,-2.5,2.5);     (*((*it_h).second))["eta_l4"]->Sumw2();
565  
566 <    (*((*it_h).second))["ip3ds_l1"]     = new TH1D(TString("ip3ds_l1") +"_"+(*it_h).first+str,";#bf{ip3ds l1};",30,-5,5);   (*((*it_h).second))["ip3ds_l1"]->Sumw2();
567 <    (*((*it_h).second))["ip3ds_l2"]     = new TH1D(TString("ip3ds_l2") +"_"+(*it_h).first+str,";#bf{ip3ds l2};",30,-5,5);   (*((*it_h).second))["ip3ds_l2"]->Sumw2();
568 <    (*((*it_h).second))["ip3ds_l3"]     = new TH1D(TString("ip3ds_l3") +"_"+(*it_h).first+str,";#bf{ip3ds l3};",30,-10,10);   (*((*it_h).second))["ip3ds_l3"]->Sumw2();
569 <    (*((*it_h).second))["ip3ds_l4"]     = new TH1D(TString("ip3ds_l4") +"_"+(*it_h).first+str,";#bf{ip3ds l4};",30,-10,10);   (*((*it_h).second))["ip3ds_l4"]->Sumw2();
566 >    (*((*it_h).second))["ip3ds_l1"]     = new TH1D(TString("ip3ds_l1") +"_"+(*it_h).first+str,";#bf{ip3ds l1};",37,-3,3);       (*((*it_h).second))["ip3ds_l1"]->Sumw2();
567 >    (*((*it_h).second))["ip3ds_l2"]     = new TH1D(TString("ip3ds_l2") +"_"+(*it_h).first+str,";#bf{ip3ds l2};",37,-3,3);       (*((*it_h).second))["ip3ds_l2"]->Sumw2();
568 >    (*((*it_h).second))["ip3ds_l3_lo"]  = new TH1D(TString("ip3ds_l3_lo") +"_"+(*it_h).first+str,";#bf{ip3ds l3_lo};",25,-6,6); (*((*it_h).second))["ip3ds_l3_lo"]->Sumw2();
569 >    (*((*it_h).second))["ip3ds_l4_lo"]  = new TH1D(TString("ip3ds_l4_lo") +"_"+(*it_h).first+str,";#bf{ip3ds l4_lo};",25,-6,6); (*((*it_h).second))["ip3ds_l4_lo"]->Sumw2();
570 >    (*((*it_h).second))["ip3ds_l3"]     = new TH1D(TString("ip3ds_l3") +"_"+(*it_h).first+str,";#bf{ip3ds l3};",55,-40,40);     (*((*it_h).second))["ip3ds_l3"]->Sumw2();
571 >    (*((*it_h).second))["ip3ds_l4"]     = new TH1D(TString("ip3ds_l4") +"_"+(*it_h).first+str,";#bf{ip3ds l4};",55,-40,40);     (*((*it_h).second))["ip3ds_l4"]->Sumw2();
572 >
573 >    (*((*it_h).second))["d0_l1"]        = new TH1D(TString("d0_l1") +"_"+(*it_h).first+str,";#bf{d0 l1};",50,-.008,.008);      (*((*it_h).second))["d0_l1"]->Sumw2();
574 >    (*((*it_h).second))["d0_l2"]        = new TH1D(TString("d0_l2") +"_"+(*it_h).first+str,";#bf{d0 l2};",50,-.008,.008);      (*((*it_h).second))["d0_l2"]->Sumw2();
575 >    (*((*it_h).second))["d0_l3_lo"]     = new TH1D(TString("d0_l3_lo") +"_"+(*it_h).first+str,";#bf{d0 l3_lo};",25,-.03,.03);  (*((*it_h).second))["d0_l3_lo"]->Sumw2();
576 >    (*((*it_h).second))["d0_l4_lo"]     = new TH1D(TString("d0_l4_lo") +"_"+(*it_h).first+str,";#bf{d0 l4_lo};",25,-.03,.03);  (*((*it_h).second))["d0_l4_lo"]->Sumw2();
577 >    (*((*it_h).second))["d0_l3"]        = new TH1D(TString("d0_l3") +"_"+(*it_h).first+str,";#bf{d0 l3};",50,-.1,.1);          (*((*it_h).second))["d0_l3"]->Sumw2();
578 >    (*((*it_h).second))["d0_l4"]        = new TH1D(TString("d0_l4") +"_"+(*it_h).first+str,";#bf{d0 l4};",50,-.1,.1);          (*((*it_h).second))["d0_l4"]->Sumw2();
579 >
580 >    (*((*it_h).second))["dz_l1"]        = new TH1D(TString("dz_l1") +"_"+(*it_h).first+str,";#bf{dz l1};",50,-.01,.01);   (*((*it_h).second))["dz_l1"]->Sumw2();
581 >    (*((*it_h).second))["dz_l2"]        = new TH1D(TString("dz_l2") +"_"+(*it_h).first+str,";#bf{dz l2};",50,-.01,.01);   (*((*it_h).second))["dz_l2"]->Sumw2();
582 >    (*((*it_h).second))["dz_l3"]        = new TH1D(TString("dz_l3") +"_"+(*it_h).first+str,";#bf{dz l3};",50,-.2,.2);   (*((*it_h).second))["dz_l3"]->Sumw2();
583 >    (*((*it_h).second))["dz_l4"]        = new TH1D(TString("dz_l4") +"_"+(*it_h).first+str,";#bf{dz l4};",50,-.2,.2);   (*((*it_h).second))["dz_l4"]->Sumw2();
584  
585      (*((*it_h).second))["dR_l3l4_lo"]   = new TH1D(TString("dR_l3l4_lo") +"_"+(*it_h).first+str,";#bf{#Delta R l3,l4};",75,0,.2);    (*((*it_h).second))["dR_l3l4_lo"]->Sumw2();
586      (*((*it_h).second))["dR_l3l4"]      = new TH1D(TString("dR_l3l4") +"_"+(*it_h).first+str,";#bf{#Delta R l3,l4};",75,0,6);        (*((*it_h).second))["dR_l3l4"]->Sumw2();
587 +
588 +    // jets
589 +    (*((*it_h).second))["pt_allJet"]    = new TH1D(TString("pt_allJet") +"_"+(*it_h).first+str,";#bf{jet 1 p_{T} [GeV]};",37,0,400);  (*((*it_h).second))["pt_allJet"]->Sumw2();
590 +    (*((*it_h).second))["eta_allJet"]   = new TH1D(TString("eta_allJet") +"_"+(*it_h).first+str,";#bf{jet 1 p_{T} [GeV]};",37,-4.8,4.8);  (*((*it_h).second))["eta_allJet"]->Sumw2();
591 +    (*((*it_h).second))["pt_j1"]        = new TH1D(TString("pt_j1") +"_"+(*it_h).first+str,";#bf{jet 1 p_{T} [GeV]};",  37,0,400);  (*((*it_h).second))["pt_j1"]->Sumw2();
592 +    (*((*it_h).second))["pt_j2"]        = new TH1D(TString("pt_j2") +"_"+(*it_h).first+str,";#bf{jet 2 p_{T} [GeV]};",  37,0,125);  (*((*it_h).second))["pt_j2"]->Sumw2();
593 +    (*((*it_h).second))["pt_j3"]        = new TH1D(TString("pt_j3") +"_"+(*it_h).first+str,";#bf{jet 3 p_{T} [GeV]};",  37,0,125);  (*((*it_h).second))["pt_j3"]->Sumw2();
594 +    (*((*it_h).second))["pt_j4"]        = new TH1D(TString("pt_j4") +"_"+(*it_h).first+str,";#bf{jet 4 p_{T} [GeV]};",  37,0,125);  (*((*it_h).second))["pt_j4"]->Sumw2();
595 +    (*((*it_h).second))["eta_j1"]       = new TH1D(TString("eta_j1")+"_"+(*it_h).first+str,";#bf{jet 1 #eta};",         55,-4.8,4.8);     (*((*it_h).second))["eta_j1"]->Sumw2();
596 +    (*((*it_h).second))["eta_j2"]       = new TH1D(TString("eta_j2")+"_"+(*it_h).first+str,";#bf{jet 2 #eta};",         55,-4.8,4.8);     (*((*it_h).second))["eta_j2"]->Sumw2();
597 +    (*((*it_h).second))["eta_j3"]       = new TH1D(TString("eta_j3")+"_"+(*it_h).first+str,";#bf{jet 3 #eta};",         55,-4.8,4.8);     (*((*it_h).second))["eta_j3"]->Sumw2();
598 +    (*((*it_h).second))["eta_j4"]       = new TH1D(TString("eta_j4")+"_"+(*it_h).first+str,";#bf{jet 4 #eta};",         55,-4.8,4.8);     (*((*it_h).second))["eta_j4"]->Sumw2();
599 +    (*((*it_h).second))["phi_j1"]       = new TH1D(TString("phi_j1")+"_"+(*it_h).first+str,";#bf{jet 1 #phi};",         25,-3.2,3.2);     (*((*it_h).second))["phi_j1"]->Sumw2();
600 +    (*((*it_h).second))["phi_j2"]       = new TH1D(TString("phi_j2")+"_"+(*it_h).first+str,";#bf{jet 2 #phi};",         25,-3.2,3.2);     (*((*it_h).second))["phi_j2"]->Sumw2();
601 +    (*((*it_h).second))["phi_j3"]       = new TH1D(TString("phi_j3")+"_"+(*it_h).first+str,";#bf{jet 3 #phi};",         25,-3.2,3.2);     (*((*it_h).second))["phi_j3"]->Sumw2();
602 +    (*((*it_h).second))["phi_j4"]       = new TH1D(TString("phi_j4")+"_"+(*it_h).first+str,";#bf{jet 4 #phi};",         25,-3.2,3.2);     (*((*it_h).second))["phi_j4"]->Sumw2();
603 +    (*((*it_h).second))["idmva_j1"]     = new TH1D(TString("idmva_j1")+"_"+(*it_h).first+str,";#bf{jet 1 MVA};",        25,-1.1,1.1);     (*((*it_h).second))["idmva_j1"]->Sumw2();
604 +    (*((*it_h).second))["idmva_j2"]     = new TH1D(TString("idmva_j2")+"_"+(*it_h).first+str,";#bf{jet 2 MVA};",        25,-1.1,1.1);     (*((*it_h).second))["idmva_j2"]->Sumw2();
605 +    (*((*it_h).second))["idmva_j3"]     = new TH1D(TString("idmva_j3")+"_"+(*it_h).first+str,";#bf{jet 3 MVA};",        25,-1.1,1.1);     (*((*it_h).second))["idmva_j3"]->Sumw2();
606 +    (*((*it_h).second))["idmva_j4"]     = new TH1D(TString("idmva_j4")+"_"+(*it_h).first+str,";#bf{jet 4 MVA};",        25,-1.1,1.1);     (*((*it_h).second))["idmva_j4"]->Sumw2();
607 +    (*((*it_h).second))["mjj"]          = new TH1D(TString("mjj")+"_"+(*it_h).first+str,";#bf{m_{jj}};",                25,10,1550);     (*((*it_h).second))["mjj"]->Sumw2();
608 +    (*((*it_h).second))["dEta"]         = new TH1D(TString("dEta")+"_"+(*it_h).first+str,";#bf{#Delta #eta};",          25,-6,6);       (*((*it_h).second))["dEta"]->Sumw2();
609 +    (*((*it_h).second))["etaProd"]      = new TH1D(TString("etaProd")+"_"+(*it_h).first+str,";#bf{#eta_{1} * #eta_{2}};",25,-12,12);      (*((*it_h).second))["etaProd"]->Sumw2();
610 +    (*((*it_h).second))["njets"]        = new TH1D(TString("njets")+"_"+(*it_h).first+str,";#bf{N Jets};",      7,-0.5,6.5);      (*((*it_h).second))["njets"]->Sumw2();
611 +
612 +    (*((*it_h).second))["dR_l1j1"]      = new TH1D(TString("dR_l1j1") +"_"+(*it_h).first+str,";#bf{#Delta R l1,jet 1};",75,0,6);        (*((*it_h).second))["dR_l1j1"]->Sumw2();
613 +    (*((*it_h).second))["dR_l2j1"]      = new TH1D(TString("dR_l2j1") +"_"+(*it_h).first+str,";#bf{#Delta R l2,jet 1};",75,0,6);        (*((*it_h).second))["dR_l2j1"]->Sumw2();
614 +    (*((*it_h).second))["dR_l3j1"]      = new TH1D(TString("dR_l3j1") +"_"+(*it_h).first+str,";#bf{#Delta R l3,jet 1};",75,0,6);        (*((*it_h).second))["dR_l3j1"]->Sumw2();
615 +    (*((*it_h).second))["dR_l4j1"]      = new TH1D(TString("dR_l4j1") +"_"+(*it_h).first+str,";#bf{#Delta R l4,jet 1};",75,0,6);        (*((*it_h).second))["dR_l4j1"]->Sumw2();
616 +    (*((*it_h).second))["dR_l1j2"]      = new TH1D(TString("dR_l1j2") +"_"+(*it_h).first+str,";#bf{#Delta R l1,jet 2};",75,0,6);        (*((*it_h).second))["dR_l1j2"]->Sumw2();
617 +    (*((*it_h).second))["dR_l2j2"]      = new TH1D(TString("dR_l2j2") +"_"+(*it_h).first+str,";#bf{#Delta R l2,jet 2};",75,0,6);        (*((*it_h).second))["dR_l2j2"]->Sumw2();
618 +    (*((*it_h).second))["dR_l3j2"]      = new TH1D(TString("dR_l3j2") +"_"+(*it_h).first+str,";#bf{#Delta R l3,jet 2};",75,0,6);        (*((*it_h).second))["dR_l3j2"]->Sumw2();
619 +    (*((*it_h).second))["dR_l4j2"]      = new TH1D(TString("dR_l4j2") +"_"+(*it_h).first+str,";#bf{#Delta R l4,jet 2};",75,0,6);        (*((*it_h).second))["dR_l4j2"]->Sumw2();
620 +    (*((*it_h).second))["dR_l1j3"]      = new TH1D(TString("dR_l1j3") +"_"+(*it_h).first+str,";#bf{#Delta R l1,jet 3};",75,0,6);        (*((*it_h).second))["dR_l1j3"]->Sumw2();
621 +    (*((*it_h).second))["dR_l2j3"]      = new TH1D(TString("dR_l2j3") +"_"+(*it_h).first+str,";#bf{#Delta R l2,jet 3};",75,0,6);        (*((*it_h).second))["dR_l2j3"]->Sumw2();
622 +    (*((*it_h).second))["dR_l3j3"]      = new TH1D(TString("dR_l3j3") +"_"+(*it_h).first+str,";#bf{#Delta R l3,jet 3};",75,0,6);        (*((*it_h).second))["dR_l3j3"]->Sumw2();
623 +    (*((*it_h).second))["dR_l4j3"]      = new TH1D(TString("dR_l4j3") +"_"+(*it_h).first+str,";#bf{#Delta R l4,jet 3};",75,0,6);        (*((*it_h).second))["dR_l4j3"]->Sumw2();
624 +    (*((*it_h).second))["dR_l1j4"]      = new TH1D(TString("dR_l1j4") +"_"+(*it_h).first+str,";#bf{#Delta R l1,jet 4};",75,0,6);        (*((*it_h).second))["dR_l1j4"]->Sumw2();
625 +    (*((*it_h).second))["dR_l2j4"]      = new TH1D(TString("dR_l2j4") +"_"+(*it_h).first+str,";#bf{#Delta R l2,jet 4};",75,0,6);        (*((*it_h).second))["dR_l2j4"]->Sumw2();
626 +    (*((*it_h).second))["dR_l3j4"]      = new TH1D(TString("dR_l3j4") +"_"+(*it_h).first+str,";#bf{#Delta R l3,jet 4};",75,0,6);        (*((*it_h).second))["dR_l3j4"]->Sumw2();
627 +    (*((*it_h).second))["dR_l4j4"]      = new TH1D(TString("dR_l4j4") +"_"+(*it_h).first+str,";#bf{#Delta R l4,jet 4};",75,0,6);        (*((*it_h).second))["dR_l4j4"]->Sumw2();
628 +
629 +    (*((*it_h).second))["fusionMVA"]            = new TH1D(TString("fusionMVA")+"_"+(*it_h).first+str,";#bf{MVA output};",        100,-1,1);      (*((*it_h).second))["fusionMVA"]->Sumw2();
630 +    (*((*it_h).second))["fusionMVA_lo"]         = new TH1D(TString("fusionMVA_lo")+"_"+(*it_h).first+str,";#bf{MVA output};",        100,-1,-.83);      (*((*it_h).second))["fusionMVA_lo"]->Sumw2();
631 +    (*((*it_h).second))["fusionMVA_med"]                = new TH1D(TString("fusionMVA_med")+"_"+(*it_h).first+str,";#bf{MVA output};",        100,-.5,1);      (*((*it_h).second))["fusionMVA_med"]->Sumw2();
632 +    (*((*it_h).second))["fusionMVA_hi"]         = new TH1D(TString("fusionMVA_hi")+"_"+(*it_h).first+str,";#bf{MVA output};",        100,.75,1);      (*((*it_h).second))["fusionMVA_hi"]->Sumw2();
633 +    (*((*it_h).second))["m4l_jet_lo"]           = new TH1D(TString("m4l_jet_lo") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};",     45,100,180);          (*((*it_h).second))["m4l_jet_lo"]->Sumw2();
634 +    (*((*it_h).second))["m4l_jet"]              = new TH1D(TString("m4l_jet") +"_"+(*it_h).first+str,";#bf{m4l [GeV]};",        45,100,600);          (*((*it_h).second))["m4l_jet"]->Sumw2();
635    }
298  
636    return hists;
637   }
638   //--------------------------------------------------------------------------------------------------
639   void makeHTML(FOFlags &ctrl, TString type, TString plotLabel)
640   {
641 <  TString title("Signal region: "+type);
642 <  ofstream htmlfile;
643 <  char htmlfname[100];
644 <  sprintf(htmlfname,"%s/plots.html",TString(ctrl.outdir+"/"+plotLabel+"/"+type).Data());
645 <  htmlfile.open(htmlfname);
641 >  TString tmpString(plotLabel);
642 >  tmpString.ReplaceAll("/"," ");
643 >  TString title(ctrl.faketype+", "+tmpString+", "+type);
644 >  TString htmlfname(ctrl.outdir+"/"+plotLabel+"/"+type+"/plots.html");
645 >  ofstream htmlfile(htmlfname);
646  
647    htmlfile << "<!DOCTYPE html" << endl;
648    htmlfile << "    PUBLIC \"-//W3C//DTD HTML 3.2//EN\">" << endl;
649    htmlfile << "<html>" << endl;
650  
651    htmlfile << "<head><title>"+title+"</title></head>" << endl;
652 <  htmlfile << "<body bgcolor=\"EEEEEE\">" << endl;
652 >  htmlfile << "<body bgcolor=\"000000\">" << endl;
653    htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">"+title+"</h3>" << endl;
654  
655 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">boson kinematics</h3>" << endl;
656    htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
657  
658    htmlfile << "<tr>" << endl;
659 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l_lo.png\"><img src=\"plots/m4l_lo.png\" alt=\"plots/m4l_lo.png\" width=\"100%\"></a></td>" << endl;
660 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l_med.png\"><img src=\"plots/m4l_med.png\" alt=\"plots/m4l_med.png\" width=\"100%\"></a></td>" << endl;
661 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l.png\"><img src=\"plots/m4l.png\" alt=\"plots/m4l.png\" width=\"100%\"></a></td>" << endl;
662 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
663 +  htmlfile << "</tr>" << endl;
664 +
665 +  htmlfile << "<tr>" << endl;
666    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/mZ1.png\"><img src=\"plots/mZ1.png\" alt=\"plots/mZ1.png\" width=\"100%\"></a></td>" << endl;
667    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/mZ2_lo.png\"><img src=\"plots/mZ2_lo.png\" alt=\"plots/mZ2_lo.png\" width=\"100%\"></a></td>" << endl;
668    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/mZ2.png\"><img src=\"plots/mZ2.png\" alt=\"plots/mZ2.png\" width=\"100%\"></a></td>" << endl;
324  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l.png\"><img src=\"plots/m4l.png\" alt=\"plots/m4l.png\" width=\"100%\"></a></td>" << endl;
669    htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
670    htmlfile << "</tr>" << endl;
671  
672    htmlfile << "<tr>" << endl;
673    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/Z1pt.png\"><img src=\"plots/Z1pt.png\" alt=\"plots/Z1pt.png\" width=\"100%\"></a></td>" << endl;
674    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ZZpt.png\"><img src=\"plots/ZZpt.png\" alt=\"plots/ZZpt.png\" width=\"100%\"></a></td>" << endl;
675 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/met.png\"><img src=\"plots/met.png\" alt=\"plots/met.png\" width=\"100%\"></a></td>" << endl;
676    htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
677    htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
678    htmlfile << "</tr>" << endl;
679  
680 +  htmlfile << "</table>" << endl;
681 +  htmlfile << "<hr />" << endl;
682 +  htmlfile << "control: " << endl;
683 +  htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
684 +
685    htmlfile << "<tr>" << endl;
686 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l1.png\"><img src=\"plots/ip3ds_l1.png\" alt=\"plots/ip3ds_l1.png\" width=\"100%\"></a></td>" << endl;
687 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l2.png\"><img src=\"plots/ip3ds_l2.png\" alt=\"plots/ip3ds_l2.png\" width=\"100%\"></a></td>" << endl;
688 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l3.png\"><img src=\"plots/ip3ds_l3.png\" alt=\"plots/ip3ds_l3.png\" width=\"100%\"></a></td>" << endl;
689 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l4.png\"><img src=\"plots/ip3ds_l4.png\" alt=\"plots/ip3ds_l4.png\" width=\"100%\"></a></td>" << endl;
686 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/run.png\"><img src=\"plots/run.png\" alt=\"plots/run.png\" width=\"100%\"></a></td>" << endl;
687 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/npv.png\"><img src=\"plots/npv.png\" alt=\"plots/npv.png\" width=\"100%\"></a></td>" << endl;
688 >  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
689 >  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
690    htmlfile << "</tr>" << endl;
691  
692 +  htmlfile << "</table>" << endl;
693 +
694 +  htmlfile << "</body>" << endl;
695 +  htmlfile << "</html>" << endl;
696 +  htmlfile.close();
697 +
698 +  TString leptonfname(htmlfname);
699 +  leptonfname.ReplaceAll("plots.html","leptonplots.html");
700 +  htmlfile.open(leptonfname);
701 +  cout << "opening: " << leptonfname << endl;
702 +
703 +  htmlfile << "<!DOCTYPE html" << endl;
704 +  htmlfile << "    PUBLIC \"-//W3C//DTD HTML 3.2//EN\">" << endl;
705 +  htmlfile << "<html>" << endl;
706 +
707 +  htmlfile << "<head><title>"+title+"</title></head>" << endl;
708 +  htmlfile << "<body bgcolor=\"000000\">" << endl;
709 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">"+title+"</h3>" << endl;
710 +
711 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">lepton plots</h3>" << endl;
712 +  htmlfile << "<hr />" << endl;
713 +  htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
714 +  htmlfile << "</table>" << endl;
715 +
716 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">Z1 leptons:</h3>" << endl;
717 +  htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
718 +
719    htmlfile << "<tr>" << endl;
720    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l1.png\"><img src=\"plots/pt_l1.png\" alt=\"plots/pt_l1.png\" width=\"100%\"></a></td>" << endl;
721    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l2.png\"><img src=\"plots/pt_l2.png\" alt=\"plots/pt_l2.png\" width=\"100%\"></a></td>" << endl;
722 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l3.png\"><img src=\"plots/pt_l3.png\" alt=\"plots/pt_l3.png\" width=\"100%\"></a></td>" << endl;
723 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l4.png\"><img src=\"plots/pt_l4.png\" alt=\"plots/pt_l4.png\" width=\"100%\"></a></td>" << endl;
722 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l1.png\"><img src=\"plots/eta_l1.png\" alt=\"plots/eta_l1.png\" width=\"100%\"></a></td>" << endl;
723 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l2.png\"><img src=\"plots/eta_l2.png\" alt=\"plots/eta_l2.png\" width=\"100%\"></a></td>" << endl;
724    htmlfile << "</tr>" << endl;
725  
726    htmlfile << "<tr>" << endl;
727 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3l4_lo.png\"><img src=\"plots/dR_l3l4_lo.png\" alt=\"plots/dR_l3l4_lo.png\" width=\"100%\"></a></td>" << endl;
728 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3l4.png\"><img src=\"plots/dR_l3l4.png\" alt=\"plots/dR_l3l4.png\" width=\"100%\"></a></td>" << endl;
729 <  htmlfile << "<td width=\"25%\"><a><</a></td>" << endl;
730 <  htmlfile << "<td width=\"25%\"><a><</a></td>" << endl;
727 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dz_l1.png\"><img src=\"plots/dz_l1.png\" alt=\"plots/dz_l1.png\" width=\"100%\"></a></td>" << endl;
728 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dz_l2.png\"><img src=\"plots/dz_l2.png\" alt=\"plots/dz_l2.png\" width=\"100%\"></a></td>" << endl;
729 >  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
730 >  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
731    htmlfile << "</tr>" << endl;
732  
733    htmlfile << "<tr>" << endl;
734 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l1.png\"><img src=\"plots/eta_l1.png\" alt=\"plots/eta_l1.png\" width=\"100%\"></a></td>" << endl;
735 <  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l2.png\"><img src=\"plots/eta_l2.png\" alt=\"plots/eta_l2.png\" width=\"100%\"></a></td>" << endl;
734 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l1.png\"><img src=\"plots/ip3ds_l1.png\" alt=\"plots/ip3ds_l1.png\" width=\"100%\"></a></td>" << endl;
735 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l2.png\"><img src=\"plots/ip3ds_l2.png\" alt=\"plots/ip3ds_l2.png\" width=\"100%\"></a></td>" << endl;
736 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l1.png\"><img src=\"plots/d0_l1.png\" alt=\"plots/d0_l1.png\" width=\"100%\"></a></td>" << endl;
737 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l2.png\"><img src=\"plots/d0_l2.png\" alt=\"plots/d0_l2.png\" width=\"100%\"></a></td>" << endl;
738 >  htmlfile << "</tr>" << endl;
739 >
740 >  // htmlfile << "<tr>" << endl;
741 >  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3l4_lo.png\"><img src=\"plots/dR_l3l4_lo.png\" alt=\"plots/dR_l3l4_lo.png\" width=\"100%\"></a></td>" << endl;
742 >  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3l4.png\"><img src=\"plots/dR_l3l4.png\" alt=\"plots/dR_l3l4.png\" width=\"100%\"></a></td>" << endl;
743 >  // htmlfile << "<td width=\"25%\"><a><</a></td>" << endl;
744 >  // htmlfile << "<td width=\"25%\"><a><</a></td>" << endl;
745 >  // htmlfile << "</tr>" << endl;
746 >
747 >  htmlfile << "</table>" << endl;
748 >  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">extra leptons (l3 and l4): </h3>" << endl;
749 >  htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
750 >
751 >  htmlfile << "<tr>" << endl;
752 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l3.png\"><img src=\"plots/pt_l3.png\" alt=\"plots/pt_l3.png\" width=\"100%\"></a></td>" << endl;
753 >  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_l4.png\"><img src=\"plots/pt_l4.png\" alt=\"plots/pt_l4.png\" width=\"100%\"></a></td>" << endl;
754    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l3.png\"><img src=\"plots/eta_l3.png\" alt=\"plots/eta_l3.png\" width=\"100%\"></a></td>" << endl;
755    htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_l4.png\"><img src=\"plots/eta_l4.png\" alt=\"plots/eta_l4.png\" width=\"100%\"></a></td>" << endl;
756    htmlfile << "</tr>" << endl;
757 +
758 +  htmlfile << "<tr>" << endl;
759 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l3_lo.png\"><img src=\"plots/ip3ds_l3_lo.png\" alt=\"plots/ip3ds_l3_lo.png\" width=\"100%\"></a></td>" << endl;
760 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l4_lo.png\"><img src=\"plots/ip3ds_l4_lo.png\" alt=\"plots/ip3ds_l4_lo.png\" width=\"100%\"></a></td>" << endl;
761 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l3.png\"><img src=\"plots/ip3ds_l3.png\" alt=\"plots/ip3ds_l3.png\" width=\"100%\"></a></td>" << endl;
762 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/ip3ds_l4.png\"><img src=\"plots/ip3ds_l4.png\" alt=\"plots/ip3ds_l4.png\" width=\"100%\"></a></td>" << endl;
763 +  htmlfile << "</tr>" << endl;
764 +
765 +  htmlfile << "<tr>" << endl;
766 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l3_lo.png\"><img src=\"plots/d0_l3_lo.png\" alt=\"plots/d0_l3_lo.png\" width=\"100%\"></a></td>" << endl;
767 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l4_lo.png\"><img src=\"plots/d0_l4_lo.png\" alt=\"plots/d0_l4_lo.png\" width=\"100%\"></a></td>" << endl;
768 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l3.png\"><img src=\"plots/d0_l3.png\" alt=\"plots/d0_l3.png\" width=\"100%\"></a></td>" << endl;
769 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/d0_l4.png\"><img src=\"plots/d0_l4.png\" alt=\"plots/d0_l4.png\" width=\"100%\"></a></td>" << endl;
770 +  htmlfile << "</tr>" << endl;
771 +
772 +  htmlfile << "<tr>" << endl;
773 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dz_l3.png\"><img src=\"plots/dz_l3.png\" alt=\"plots/dz_l3.png\" width=\"100%\"></a></td>" << endl;
774 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dz_l4.png\"><img src=\"plots/dz_l4.png\" alt=\"plots/dz_l4.png\" width=\"100%\"></a></td>" << endl;
775 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
776 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
777 +  htmlfile << "</tr>" << endl;
778 +
779 +  htmlfile << "</table>" << endl;
780 +
781 +  htmlfile << "</body>" << endl;
782 +  htmlfile << "</html>" << endl;
783 +  htmlfile.close();
784 +
785 +  TString jetfname(htmlfname);
786 +  jetfname.ReplaceAll("plots.html","jetplots.html");
787 +  htmlfile.open(jetfname);
788 +
789 +  htmlfile << "<!DOCTYPE html" << endl;
790 +  htmlfile << "    PUBLIC \"-//W3C//DTD HTML 3.2//EN\">" << endl;
791 +  htmlfile << "<html>" << endl;
792 +
793 +  htmlfile << "<head><title>"+title+"</title></head>" << endl;
794 +  htmlfile << "<body bgcolor=\"000000\">" << endl;
795 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">"+title+"</h3>" << endl;
796 +
797 +  htmlfile << "<h3 style=\"text-align:left; color:DD6600;\">jet plots</h3>" << endl;
798 +  htmlfile << "<hr />" << endl;
799 +  htmlfile << "<table border=\"0\" cellspacing=\"5\" width=\"100%\">" << endl;  
800 +
801 +  htmlfile << "<tr>" << endl;
802 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_j1.png\"><img src=\"plots/pt_j1.png\" alt=\"plots/pt_j1.png\" width=\"100%\"></a></td>" << endl;
803 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_j2.png\"><img src=\"plots/pt_j2.png\" alt=\"plots/pt_j2.png\" width=\"100%\"></a></td>" << endl;
804 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_j2.png\"><img src=\"plots/pt_j2.png\" alt=\"plots/pt_j2.png\" width=\"100%\"></a></td>" << endl;
805 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/pt_j4.png\"><img src=\"plots/pt_j4.png\" alt=\"plots/pt_j4.png\" width=\"100%\"></a></td>" << endl;
806 +  htmlfile << "</tr>" << endl;
807 +
808 +  htmlfile << "<tr>" << endl;
809 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_j1.png\"><img src=\"plots/eta_j1.png\" alt=\"plots/eta_j1.png\" width=\"100%\"></a></td>" << endl;
810 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_j2.png\"><img src=\"plots/eta_j2.png\" alt=\"plots/eta_j2.png\" width=\"100%\"></a></td>" << endl;
811 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_j3.png\"><img src=\"plots/eta_j3.png\" alt=\"plots/eta_j3.png\" width=\"100%\"></a></td>" << endl;
812 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/eta_j4.png\"><img src=\"plots/eta_j4.png\" alt=\"plots/eta_j4.png\" width=\"100%\"></a></td>" << endl;
813 +  htmlfile << "</tr>" << endl;
814 +
815 +  htmlfile << "<tr>" << endl;
816 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/idmva_j1.png\"><img src=\"plots/idmva_j1.png\" alt=\"plots/idmva_j1.png\" width=\"100%\"></a></td>" << endl;
817 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/idmva_j2.png\"><img src=\"plots/idmva_j2.png\" alt=\"plots/idmva_j2.png\" width=\"100%\"></a></td>" << endl;
818 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/idmva_j3.png\"><img src=\"plots/idmva_j3.png\" alt=\"plots/idmva_j3.png\" width=\"100%\"></a></td>" << endl;
819 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/idmva_j4.png\"><img src=\"plots/idmva_j4.png\" alt=\"plots/idmva_j4.png\" width=\"100%\"></a></td>" << endl;
820 +  htmlfile << "</tr>" << endl;
821 +
822 +  htmlfile << "<tr>" << endl;
823 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/phi_j1.png\"><img src=\"plots/phi_j1.png\" alt=\"plots/phi_j1.png\" width=\"100%\"></a></td>" << endl;
824 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/phi_j2.png\"><img src=\"plots/phi_j2.png\" alt=\"plots/phi_j2.png\" width=\"100%\"></a></td>" << endl;
825 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/phi_j3.png\"><img src=\"plots/phi_j3.png\" alt=\"plots/phi_j3.png\" width=\"100%\"></a></td>" << endl;
826 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/phi_j4.png\"><img src=\"plots/phi_j4.png\" alt=\"plots/phi_j4.png\" width=\"100%\"></a></td>" << endl;
827 +  htmlfile << "</tr>" << endl;
828 +
829 +  htmlfile << "<tr>" << endl;
830 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/mjj.png\"><img src=\"plots/mjj.png\" alt=\"plots/mjj.png\" width=\"100%\"></a></td>" << endl;
831 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dEta.png\"><img src=\"plots/dEta.png\" alt=\"plots/dEta.png\" width=\"100%\"></a></td>" << endl;
832 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/etaProd.png\"><img src=\"plots/etaProd.png\" alt=\"plots/etaProd.png\" width=\"100%\"></a></td>" << endl;
833 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/njets.png\"><img src=\"plots/njets.png\" alt=\"plots/njets.png\" width=\"100%\"></a></td>" << endl;
834 +  htmlfile << "</tr>" << endl;
835 +
836 +  // htmlfile << "<tr>" << endl;
837 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l1j1.png\"><img src=\"plots/dR_l1j1.png\" alt=\"plots/dR_l1j1.png\" width=\"100%\"></a></td>" << endl;
838 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l2j1.png\"><img src=\"plots/dR_l2j1.png\" alt=\"plots/dR_l2j1.png\" width=\"100%\"></a></td>" << endl;
839 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3j1.png\"><img src=\"plots/dR_l3j1.png\" alt=\"plots/dR_l3j1.png\" width=\"100%\"></a></td>" << endl;
840 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l4j1.png\"><img src=\"plots/dR_l4j1.png\" alt=\"plots/dR_l4j1.png\" width=\"100%\"></a></td>" << endl;
841 +  // htmlfile << "</tr>" << endl;
842 +
843 +  // htmlfile << "<tr>" << endl;
844 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l1j2.png\"><img src=\"plots/dR_l1j2.png\" alt=\"plots/dR_l1j2.png\" width=\"100%\"></a></td>" << endl;
845 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l2j2.png\"><img src=\"plots/dR_l2j2.png\" alt=\"plots/dR_l2j2.png\" width=\"100%\"></a></td>" << endl;
846 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3j2.png\"><img src=\"plots/dR_l3j2.png\" alt=\"plots/dR_l3j2.png\" width=\"100%\"></a></td>" << endl;
847 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l4j2.png\"><img src=\"plots/dR_l4j2.png\" alt=\"plots/dR_l4j2.png\" width=\"100%\"></a></td>" << endl;
848 +  // htmlfile << "</tr>" << endl;
849 +
850 +  // htmlfile << "<tr>" << endl;
851 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l1j3.png\"><img src=\"plots/dR_l1j3.png\" alt=\"plots/dR_l1j3.png\" width=\"100%\"></a></td>" << endl;
852 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l2j3.png\"><img src=\"plots/dR_l2j3.png\" alt=\"plots/dR_l2j3.png\" width=\"100%\"></a></td>" << endl;
853 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3j3.png\"><img src=\"plots/dR_l3j3.png\" alt=\"plots/dR_l3j3.png\" width=\"100%\"></a></td>" << endl;
854 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l4j3.png\"><img src=\"plots/dR_l4j3.png\" alt=\"plots/dR_l4j3.png\" width=\"100%\"></a></td>" << endl;
855 +  // htmlfile << "</tr>" << endl;
856 +
857 +  // htmlfile << "<tr>" << endl;
858 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l1j4.png\"><img src=\"plots/dR_l1j4.png\" alt=\"plots/dR_l1j4.png\" width=\"100%\"></a></td>" << endl;
859 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l2j4.png\"><img src=\"plots/dR_l2j4.png\" alt=\"plots/dR_l2j4.png\" width=\"100%\"></a></td>" << endl;
860 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l3j4.png\"><img src=\"plots/dR_l3j4.png\" alt=\"plots/dR_l3j4.png\" width=\"100%\"></a></td>" << endl;
861 +  // htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/dR_l4j4.png\"><img src=\"plots/dR_l4j4.png\" alt=\"plots/dR_l4j4.png\" width=\"100%\"></a></td>" << endl;
862 +  // htmlfile << "</tr>" << endl;
863 +
864 +  htmlfile << "<tr>" << endl;
865 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA.png\"><img src=\"plots/fusionMVA.png\" alt=\"plots/fusionMVA.png\" width=\"100%\"></a></td>" << endl;
866 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA_lo.png\"><img src=\"plots/fusionMVA_lo.png\" alt=\"plots/fusionMVA_lo.png\" width=\"100%\"></a></td>" << endl;
867 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA_hi.png\"><img src=\"plots/fusionMVA_hi.png\" alt=\"plots/fusionMVA_hi.png\" width=\"100%\"></a></td>" << endl;
868 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
869 +  htmlfile << "</tr>" << endl;
870 +
871 +  htmlfile << "<tr>" << endl;
872 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA-log.png\"><img src=\"plots/fusionMVA-log.png\" alt=\"plots/fusionMVA-log.png\" width=\"100%\"></a></td>" << endl;
873 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA_lo-log.png\"><img src=\"plots/fusionMVA_lo-log.png\" alt=\"plots/fusionMVA_lo-log.png\" width=\"100%\"></a></td>" << endl;
874 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/fusionMVA_hi-log.png\"><img src=\"plots/fusionMVA_hi-log.png\" alt=\"plots/fusionMVA_hi-log.png\" width=\"100%\"></a></td>" << endl;
875 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
876 +  htmlfile << "</tr>" << endl;
877 +
878 +  htmlfile << "<tr>" << endl;
879 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l_jet_lo.png\"><img src=\"plots/m4l_jet_lo.png\" alt=\"plots/m4l_jet_lo.png\" width=\"100%\"></a></td>" << endl;
880 +  htmlfile << "<td width=\"25%\"><a target=\"_blank\" href=\"plots/m4l_jet.png\"><img src=\"plots/m4l_jet.png\" alt=\"plots/m4l_jet.png\" width=\"100%\"></a></td>" << endl;
881 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
882 +  htmlfile << "<td width=\"25%\"><a></a></td>" << endl;
883 +  htmlfile << "</tr>" << endl;
884 +
885 +
886    htmlfile << "</table>" << endl;
887  
888    htmlfile << "<hr />" << endl;
# Line 372 | Line 896 | void makeHTML(FOFlags &ctrl, TString typ
896   //----------------------------------------------------------------------------------------
897   void fillHist(CSample *cs, TString channel, TString type, TString var, double val, double wgt, double wgt_lo, double wgt_hi)
898   {
375  // // name = channel+"_"+type;
376  // // name = "all_"+type+"_lo";
377  // // name = "all_"+type+"_hi";
378  // // name = channel+"_"+type+"_lo";
379  // // name = channel+"_"+type+"_hi";
380
381  // cout << endl << "all_"+type;
382  // (*(cs->hists)["all_"+type])[var]->Fill(        val,   wgt);
383  // if(channel!="") {
384  //   cout << endl << channel+"_"+type;
385  //   (*(cs->hists)[channel+"_"+type])[var]->Fill(   val,   wgt);
386  // }
387  // if(type=="pred") {
388  //   cout << endl << "all_"+type+"_lo";
389  //   (*(cs->hists)["all_"+type+"_lo"])[var]->Fill(   val,   wgt_lo);
390  //   cout << endl << "all_"+type+"_hi";
391  //   (*(cs->hists)["all_"+type+"_hi"])[var]->Fill(   val,   wgt_hi);
392  //   if(channel!="") {
393  //     cout << endl << channel+"_"+type+"_lo";
394  //     (*(cs->hists)[channel+"_"+type+"_lo"])[var]->Fill(   val,   wgt_lo);
395  //     cout << endl << channel+"_"+type+"_hi" << endl;
396  //     (*(cs->hists)[channel+"_"+type+"_hi"])[var]->Fill(   val,   wgt_hi);
397  //   }
398  // }
899    (*(cs->hists)["all_"+type])[var]->Fill(           val,   wgt);
900    if(channel!="")
901       (*(cs->hists)[channel+"_"+type])[var]->Fill(   val,   wgt);
902 <  if(type=="pred") {
902 >  if(type=="pred" || type=="PF") {
903      (*(cs->hists)["all_"+type+"_lo"])[var]->Fill(   val,   wgt_lo);
904      (*(cs->hists)["all_"+type+"_hi"])[var]->Fill(   val,   wgt_hi);
905      if(channel!="") {
# Line 408 | Line 908 | void fillHist(CSample *cs, TString chann
908      }
909    }
910   }
911 + //----------------------------------------------------------------------------------------
912 + void fillAllHists(FOFlags ctrl, CSample *cs, TString channel, TString type, filestuff *fs, KinematicsStruct kine,
913 +                  SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4,
914 +                  double wgt, double wgt_lo, double wgt_hi)
915 + {
916 +  fillHist( cs, channel, type, "npv"         , fs->info->npv  , wgt, wgt_lo, wgt_hi);
917 +  fillHist( cs, channel, type, "run"         , fs->info->run  , wgt, wgt_lo, wgt_hi);
918 +  fillHist( cs, channel, type, "mZ1"         , kine.mZ1       , wgt, wgt_lo, wgt_hi);                
919 +  fillHist( cs, channel, type, "mZ2"         , kine.mZ2       , wgt, wgt_lo, wgt_hi);                
920 +  fillHist( cs, channel, type, "mZ2_lo"      , kine.mZ2       , wgt, wgt_lo, wgt_hi);                
921 +  fillHist( cs, channel, type, "m4l_lo"      , kine.m4l       , wgt, wgt_lo, wgt_hi);                
922 +  fillHist( cs, channel, type, "m4l_med"     , kine.m4l       , wgt, wgt_lo, wgt_hi);                
923 +  fillHist( cs, channel, type, "m4l"         , kine.m4l       , wgt, wgt_lo, wgt_hi);                
924 +  fillHist( cs, channel, type, "Z1pt"        , kine.Z1pt      , wgt, wgt_lo, wgt_hi);                
925 +  fillHist( cs, channel, type, "ZZpt"        , kine.ZZpt      , wgt, wgt_lo, wgt_hi);                
926 +  fillHist( cs, channel, type, "met"         , fs->info->met  , wgt, wgt_lo, wgt_hi);                
927 +  fillHist( cs, channel, type, "ip3ds_l1"    , lep1.ip3dSig   , wgt, wgt_lo, wgt_hi);
928 +  fillHist( cs, channel, type, "ip3ds_l2"    , lep2.ip3dSig   , wgt, wgt_lo, wgt_hi);
929 +  fillHist( cs, channel, type, "ip3ds_l3"    , lep3.ip3dSig   , wgt, wgt_lo, wgt_hi);
930 +  fillHist( cs, channel, type, "ip3ds_l4"    , lep4.ip3dSig   , wgt, wgt_lo, wgt_hi);
931 +  fillHist( cs, channel, type, "ip3ds_l3_lo" , lep3.ip3dSig   , wgt, wgt_lo, wgt_hi);
932 +  fillHist( cs, channel, type, "ip3ds_l4_lo" , lep4.ip3dSig   , wgt, wgt_lo, wgt_hi);
933 +  fillHist( cs, channel, type, "d0_l1"       , lep1.d0        , wgt, wgt_lo, wgt_hi);
934 +  fillHist( cs, channel, type, "d0_l2"       , lep2.d0        , wgt, wgt_lo, wgt_hi);
935 +  fillHist( cs, channel, type, "d0_l3"       , lep3.d0        , wgt, wgt_lo, wgt_hi);
936 +  fillHist( cs, channel, type, "d0_l4"       , lep4.d0        , wgt, wgt_lo, wgt_hi);
937 +  fillHist( cs, channel, type, "d0_l3_lo"    , lep3.d0        , wgt, wgt_lo, wgt_hi);
938 +  fillHist( cs, channel, type, "d0_l4_lo"    , lep4.d0        , wgt, wgt_lo, wgt_hi);
939 +  fillHist( cs, channel, type, "dz_l1"       , lep1.dz        , wgt, wgt_lo, wgt_hi);
940 +  fillHist( cs, channel, type, "dz_l2"       , lep2.dz        , wgt, wgt_lo, wgt_hi);
941 +  fillHist( cs, channel, type, "dz_l3"       , lep3.dz        , wgt, wgt_lo, wgt_hi);
942 +  fillHist( cs, channel, type, "dz_l4"       , lep4.dz        , wgt, wgt_lo, wgt_hi);
943 +  fillHist( cs, channel, type, "pt_l1"       , lep1.vec.Pt()  , wgt, wgt_lo, wgt_hi);
944 +  fillHist( cs, channel, type, "pt_l2"       , lep2.vec.Pt()  , wgt, wgt_lo, wgt_hi);
945 +  fillHist( cs, channel, type, "pt_l3"       , lep3.vec.Pt()  , wgt, wgt_lo, wgt_hi);
946 +  fillHist( cs, channel, type, "pt_l4"       , lep4.vec.Pt()  , wgt, wgt_lo, wgt_hi);
947 +  fillHist( cs, channel, type, "eta_l1"      , lep1.vec.Eta() , wgt, wgt_lo, wgt_hi);
948 +  fillHist( cs, channel, type, "eta_l2"      , lep2.vec.Eta() , wgt, wgt_lo, wgt_hi);
949 +  fillHist( cs, channel, type, "eta_l3"      , lep3.vec.Eta() , wgt, wgt_lo, wgt_hi);
950 +  fillHist( cs, channel, type, "eta_l4"      , lep4.vec.Eta() , wgt, wgt_lo, wgt_hi);
951 +  fillHist( cs, channel, type, "dR_l3l4_lo"  , dr(lep3,lep4)  , wgt, wgt_lo, wgt_hi);
952 +  fillHist( cs, channel, type, "dR_l3l4"     , dr(lep3,lep4)  , wgt, wgt_lo, wgt_hi);
953 + }
954 + //----------------------------------------------------------------------------------------
955 + void fillAllJetHists( FOFlags ctrl, CSample *cs, TString channel, TString type, filestuff *fs, KinematicsStruct kine,
956 +                      SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4,
957 +                      FusionMva &fusion, vector<SimpleLepton> &goodJets, JetInfoStruct ji,
958 +                      double wgt, double wgt_lo, double wgt_hi)
959 + {
960 +  double fMVAval = (ctrl.uncert=="") ? 0 : fusion.reader->EvaluateMVA("BDTG");
961 +  dr_struct drs = fill_dr_struct(goodJets,lep1,lep2,lep3,lep4);
962 +  drs.check();
963 +  fillHist( cs, channel, type, "pt_j1"  ,     ji.ptJet1,        wgt, wgt_lo, wgt_hi);
964 +  fillHist( cs, channel, type, "pt_j2"  ,     ji.ptJet2,        wgt, wgt_lo, wgt_hi);
965 +  fillHist( cs, channel, type, "pt_j3"  ,     ji.ptJet3,        wgt, wgt_lo, wgt_hi);
966 +  fillHist( cs, channel, type, "pt_j4"  ,     ji.ptJet4,        wgt, wgt_lo, wgt_hi);
967 +  fillHist( cs, channel, type, "eta_j1" ,     ji.etaJet1,       wgt, wgt_lo, wgt_hi);
968 +  fillHist( cs, channel, type, "eta_j2" ,     ji.etaJet2,       wgt, wgt_lo, wgt_hi);
969 +  fillHist( cs, channel, type, "eta_j3" ,     ji.etaJet3,       wgt, wgt_lo, wgt_hi);
970 +  fillHist( cs, channel, type, "eta_j4" ,     ji.etaJet4,       wgt, wgt_lo, wgt_hi);
971 +  fillHist( cs, channel, type, "phi_j1" ,     ji.phiJet1,       wgt, wgt_lo, wgt_hi);
972 +  fillHist( cs, channel, type, "phi_j2" ,     ji.phiJet2,       wgt, wgt_lo, wgt_hi);
973 +  fillHist( cs, channel, type, "phi_j3" ,     ji.phiJet3,       wgt, wgt_lo, wgt_hi);
974 +  fillHist( cs, channel, type, "phi_j4" ,     ji.phiJet4,       wgt, wgt_lo, wgt_hi);
975 +  fillHist( cs, channel, type, "idmva_j1",    ji.mvaJet1,       wgt, wgt_lo, wgt_hi);
976 +  fillHist( cs, channel, type, "idmva_j2",    ji.mvaJet2,       wgt, wgt_lo, wgt_hi);
977 +  fillHist( cs, channel, type, "idmva_j3",    ji.mvaJet3,       wgt, wgt_lo, wgt_hi);
978 +  fillHist( cs, channel, type, "idmva_j4",    ji.mvaJet4,       wgt, wgt_lo, wgt_hi);
979 +  fillHist( cs, channel, type, "mjj"    ,     ji.mjj,           wgt, wgt_lo, wgt_hi);
980 +  fillHist( cs, channel, type, "dEta"   ,     ji.dEta,          wgt, wgt_lo, wgt_hi);
981 +  fillHist( cs, channel, type, "etaProd",     ji.etaProd,       wgt, wgt_lo, wgt_hi);
982 +  fillHist( cs, channel, type, "njets"    ,   goodJets.size(),  wgt, wgt_lo, wgt_hi);
983 +  fillHist( cs, channel, type, "dR_l1j1",     drs.l1j1,         wgt, wgt_lo, wgt_hi);
984 +  fillHist( cs, channel, type, "dR_l2j1",     drs.l2j1,         wgt, wgt_lo, wgt_hi);
985 +  fillHist( cs, channel, type, "dR_l3j1",     drs.l3j1,         wgt, wgt_lo, wgt_hi);
986 +  fillHist( cs, channel, type, "dR_l4j1",     drs.l4j1,         wgt, wgt_lo, wgt_hi);
987 +  fillHist( cs, channel, type, "dR_l1j2",     drs.l1j2,         wgt, wgt_lo, wgt_hi);
988 +  fillHist( cs, channel, type, "dR_l2j2",     drs.l2j2,         wgt, wgt_lo, wgt_hi);
989 +  fillHist( cs, channel, type, "dR_l3j2",     drs.l3j2,         wgt, wgt_lo, wgt_hi);
990 +  fillHist( cs, channel, type, "dR_l4j2",     drs.l4j2,         wgt, wgt_lo, wgt_hi);
991 +  fillHist( cs, channel, type, "dR_l1j3",     drs.l1j3,         wgt, wgt_lo, wgt_hi);
992 +  fillHist( cs, channel, type, "dR_l2j3",     drs.l2j3,         wgt, wgt_lo, wgt_hi);
993 +  fillHist( cs, channel, type, "dR_l3j3",     drs.l3j3,         wgt, wgt_lo, wgt_hi);
994 +  fillHist( cs, channel, type, "dR_l4j3",     drs.l4j3,         wgt, wgt_lo, wgt_hi);
995 +  fillHist( cs, channel, type, "dR_l1j4",     drs.l1j4,         wgt, wgt_lo, wgt_hi);
996 +  fillHist( cs, channel, type, "dR_l2j4",     drs.l2j4,         wgt, wgt_lo, wgt_hi);
997 +  fillHist( cs, channel, type, "dR_l3j4",     drs.l3j4,         wgt, wgt_lo, wgt_hi);
998 +  fillHist( cs, channel, type, "dR_l4j4",     drs.l4j4,         wgt, wgt_lo, wgt_hi);
999 +
1000 +  fillHist( cs, channel, type, "fusionMVA_lo", fMVAval,         wgt, wgt_lo, wgt_hi);
1001 +  fillHist( cs, channel, type, "fusionMVA_med",fMVAval,         wgt, wgt_lo, wgt_hi);
1002 +  fillHist( cs, channel, type, "fusionMVA_hi", fMVAval,         wgt, wgt_lo, wgt_hi);
1003 +  fillHist( cs, channel, type, "fusionMVA",    fMVAval,         wgt, wgt_lo, wgt_hi);
1004 +  if(fMVAval > -.3) {
1005 +    fillHist( cs, channel, type, "m4l_jet_lo",  kine.m4l  ,     wgt, wgt_lo, wgt_hi);                
1006 +    fillHist( cs, channel, type, "m4l_jet",     kine.m4l  ,     wgt, wgt_lo, wgt_hi);
1007 +  }
1008 + }
1009 + //----------------------------------------------------------------------------------------
1010 + bool passHlt(FOFlags &ctrl, TrigInfo ti, InfoStruct *info, unsigned lep1matchBits,
1011 +             unsigned lep2matchBits, unsigned lep3matchBits,
1012 +             unsigned lep4matchBits)
1013 + {
1014 +  bool pass=false;
1015 +  bool onePassedButWasOutside=false;
1016 +  bitset<30> firedBits(info->status);
1017 +  if(ctrl.debug) cout << setw(12) << info->run << setw(12) << info->evt << setw(65) << firedBits <<endl;
1018 +  for(unsigned ibit=0; ibit<ti.trigBits.size(); ibit++) {
1019 +    bool outsideRunRange=false;
1020 +    if(info->run < ti.minRuns[ibit] || info->run > ti.maxRuns[ibit]) {
1021 +      outsideRunRange = true;
1022 +    }
1023 +    int baconBit = ti.trigBits[ibit];
1024 +    bool passBit = firedBits.test(ibit);
1025 +    if(ctrl.debug) cout << "   ibit: " << setw(4) << ibit << " (bacon bit: " << baconBit << ") fired: " << passBit << endl;
1026 +    if(ctrl.debug && passBit && outsideRunRange) cout << "        would have passed but it's outside run range (" << ti.minRuns[ibit] << "," << ti.maxRuns[ibit] << ")" << endl;
1027 +    if(passBit && outsideRunRange) {
1028 +       onePassedButWasOutside = true;
1029 +    }
1030 +    if(passBit && !outsideRunRange) pass = true;
1031 +  }
1032 +  if(!pass && onePassedButWasOutside) cout << "WARNING: would have passed, but one was outside its run range" << endl;
1033 +  return pass;
1034 + }
1035 + //----------------------------------------------------------------------------------------
1036 + void fillLeptonFakeWeights(lepFrs &fwgts, FR_struct *fr, SimpleLepton lep3, SimpleLepton lep4)
1037 + {
1038 +  fwgts.fwgt_3    = get_fake_weight("",lep3,*fr);
1039 +  fwgts.fwgt_lo_3 = get_fake_weight("lo",lep3,*fr);
1040 +  fwgts.fwgt_hi_3 = get_fake_weight("hi",lep3,*fr);
1041 +  fwgts.fwgt_4    = get_fake_weight("",lep4,*fr);
1042 +  fwgts.fwgt_lo_4 = get_fake_weight("lo",lep4,*fr);
1043 +  fwgts.fwgt_hi_4 = get_fake_weight("hi",lep4,*fr);
1044 + }
1045 + //----------------------------------------------------------------------------------------
1046 + dr_struct fill_dr_struct(vector<SimpleLepton> jets, SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4)
1047 + {
1048 +  dr_struct drs;
1049 +  if(jets.size() > 0) {
1050 +    drs.l1j1 = dr(lep1,jets[0]);
1051 +    drs.l2j1 = dr(lep2,jets[0]);
1052 +    drs.l3j1 = dr(lep3,jets[0]);
1053 +    drs.l4j1 = dr(lep4,jets[0]);
1054 +  }
1055 +  if(jets.size() > 1) {
1056 +    drs.l1j2 = dr(lep1,jets[1]);
1057 +    drs.l2j2 = dr(lep2,jets[1]);
1058 +    drs.l3j2 = dr(lep3,jets[1]);
1059 +    drs.l4j2 = dr(lep4,jets[1]);
1060 +  }
1061 +  if(jets.size() > 2) {
1062 +    drs.l1j3 = dr(lep1,jets[2]);
1063 +    drs.l2j3 = dr(lep2,jets[2]);
1064 +    drs.l3j3 = dr(lep3,jets[2]);
1065 +    drs.l4j3 = dr(lep4,jets[2]);
1066 +  }
1067 +  if(jets.size() > 3) {
1068 +    drs.l1j4 = dr(lep1,jets[3]);
1069 +    drs.l2j4 = dr(lep2,jets[3]);
1070 +    drs.l3j4 = dr(lep3,jets[3]);
1071 +    drs.l4j4 = dr(lep4,jets[3]);
1072 +  }
1073 +  return drs;
1074 + }
1075 + //----------------------------------------------------------------------------------------
1076 + void init_cuts( vector<TString> &cutstrs, map<TString,int> &cutvec)
1077 + {
1078 +  cutstrs.push_back("start"  );    cutvec["start"]  = 0;
1079 +  cutstrs.push_back("rlrmAndDupl"  );    cutvec["rlrmAndDupl"]  = 0;
1080 +  // cutstrs.push_back("skim"  );     cutvec["skim"]  = 0;
1081 +  cutstrs.push_back("4lsele"  );     cutvec["4lsele"]  = 0;
1082 +  cutstrs.push_back("twoJets"  );  cutvec["twoJets"] = 0;
1083 +  cutstrs.push_back("twoJetsAfter"  );  cutvec["twoJetsAfter"] = 0;
1084 +  cutstrs.push_back("filling"  );  cutvec["filling"] = 0;
1085 +  cutstrs.push_back("fillPass"  );  cutvec["fillPass"] = 0;
1086 +  cutstrs.push_back("fillFail"  );  cutvec["fillFail"] = 0;
1087 +  cutstrs.push_back("fillPF"  );  cutvec["fillPF"] = 0;
1088 + }
1089 + //----------------------------------------------------------------------------------------
1090 + bool findGoodJets(vector<SimpleLepton> &goodJets, filestuff *fs,
1091 +                  SimpleLepton lep1, SimpleLepton lep2, SimpleLepton lep3, SimpleLepton lep4)
1092 + {
1093 +  int nHighPtJets=0,nLowPtJets=0;
1094 +  for(unsigned ijet=0; ijet<fs->jets->size(); ijet++) {
1095 +    SimpleLepton *jet = &((*fs->jets)[ijet]);
1096 +    if(dr(*jet,lep1) < 0.2)             continue;
1097 +    if(dr(*jet,lep2) < 0.2)             continue;
1098 +    if(dr(*jet,lep3) < 0.2)             continue;
1099 +    if(dr(*jet,lep4) < 0.2)             continue;
1100 +    if(!(jet->status.passPre()))        continue;
1101 +    if(!jet->isLoose)                   continue;
1102 +    if(fabs(jet->vec.Eta()) > 4.5)      continue;
1103 +    goodJets.push_back(*jet);
1104 +    if(jet->vec.Pt() > 40)              nHighPtJets++;
1105 +    if(jet->vec.Pt() > 20)              nLowPtJets++;
1106 +  }
1107 +  sort( goodJets.begin(), goodJets.end(), SimpleLepton::lep_pt_sort ); // note: they're not sorted, 'cause jet corrections were applied in applyZPlusX
1108 +  if((goodJets.size()>0 && nHighPtJets<1) || // require pt > 40,20 if the jets are there
1109 +     (goodJets.size()>1 && nLowPtJets<2))
1110 +    return false;
1111 +  else
1112 +    return true;
1113 + }
1114 + //----------------------------------------------------------------------------------------
1115 + map<TString,TH1D*> setHistSet(CSample *cs, TString type, TString var)
1116 + {
1117 +  map<TString,TH1D*> hset;
1118 +  hset["obs"]            = (*(cs->hists)[type+"_obs"])[var];
1119 +  hset["pred"]           = (*(cs->hists)[type+"_pred"])[var];
1120 +  hset["pred_lo"]        = (*(cs->hists)[type+"_pred_lo"])[var];
1121 +  hset["pred_hi"]        = (*(cs->hists)[type+"_pred_hi"])[var];
1122 +  hset["PF"]             = (*(cs->hists)[type+"_PF"])[var];
1123 +  hset["PF_lo"]          = (*(cs->hists)[type+"_PF_lo"])[var];
1124 +  hset["PF_hi"]          = (*(cs->hists)[type+"_PF_hi"])[var];
1125 +  return hset;
1126 + }
1127 + //----------------------------------------------------------------------------------------
1128 + void shiftOverflows(map<TString,TH1D*> &hset)
1129 + {
1130 +  for(map<TString,TH1D*>::iterator it=hset.begin(); it!=hset.end(); it++)
1131 +    shiftOverflows((*it).second);
1132 + }
1133 + //----------------------------------------------------------------------------------------
1134 + void scaleHists(map<TString,TH1D*> &hset)
1135 + {
1136 +  for(map<TString,TH1D*>::iterator it=hset.begin(); it!=hset.end(); it++)
1137 +    (*it).second->Scale(1./integrateHist((*it).second));
1138 + }
1139 + //----------------------------------------------------------------------------------------
1140 + void fillFakeTuple(CSample *cs, filestuff *fs, EventData evtdata, unsigned ientry, double fillweight)
1141 + {
1142 +  fs->getentry(ientry,"","zznt"); // make sure we've got this entry for all the branches
1143 +  fillKinematics(evtdata,*fs->kine); // reminder: kine has the fake information, while fs->kine *had*, until this line, whatever ZPlusX put in it
1144 +  fs->weights->w = fillweight;
1145 +  fs->outtuple->Fill();
1146 +  fs->outFotuple->Fill();
1147 + }
1148 +
1149 +          // double mjj=0,dEta=0,etaProd=0,nCentralJets=0;
1150 +          // if(jet1 && jet2) {
1151 +          //   cutvec["twoJets"] += 1;
1152 +          //   if(ctrl.debug) cout << "looping through " << goodJets.size() << " good jets" << endl;
1153 +          //   if(ctrl.debug) {cout << "jet 1 (" << jet1 << "): "; jet1->print();}
1154 +          //   if(ctrl.debug) {cout << "jet 2 (" << jet2 << "): "; jet2->print();}
1155 +          //   for(unsigned ijet=0; ijet<goodJets.size(); ijet++) {
1156 +          //     SimpleLepton *jet = goodJets[ijet];
1157 +          //     if(jet == jet1) {
1158 +          //    if(ctrl.debug) cout << "  " << jet << " this is jet 1" << endl;
1159 +          //    continue;
1160 +          //     }
1161 +          //     if(jet == jet2) {
1162 +          //    if(ctrl.debug) cout << "  " << jet << " this is jet 2" << endl;
1163 +          //    continue;
1164 +          //     }
1165 +          //     double eta = jet->vec.Eta();
1166 +          //     double eta1 = jet1->vec.Eta();
1167 +          //     double eta2 = jet2->vec.Eta();
1168 +          //     if(eta1 > eta2) {
1169 +          //    if(eta > eta2 && eta < eta1)
1170 +          //      nCentralJets++;
1171 +          //     } else if(eta2 > eta1) {
1172 +          //    if(eta > eta1 && eta < eta2)
1173 +          //      nCentralJets++;
1174 +          //     }
1175 +          //     if(ctrl.debug) cout << "  nCentralJets: " << nCentralJets << endl;
1176 +          //   }
1177 +          //   double mjjMin = 400;
1178 +          //   double dEtaMin = 4;
1179 +          //   TLorentzVector dijet(jet1->vec + jet2->vec);
1180 +          //   mjj = dijet.M();
1181 +          //   if(ctrl.debug) cout << "setting mjj: " << mjj << "(using " << jet1->vec.Pt() << " " << jet2->vec.Pt() << endl;
1182 +          //   dEta = jet1->vec.Eta() - jet2->vec.Eta();
1183 +          //   etaProd = jet1->vec.Eta()*jet2->vec.Eta();
1184 +          //   bool isVbf = (mjj > mjjMin) && (fabs(dEta) > dEtaMin) && (etaProd < 0) && (nCentralJets==0);
1185 +          // }
1186 +

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines