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

Comparing UserCode/MitHzz4l/Selection/src/Selection.cc (file contents):
Revision 1.6 by khahn, Thu Oct 13 14:20:55 2011 UTC vs.
Revision 1.19 by khahn, Mon Feb 13 09:47:12 2012 UTC

# Line 1 | Line 1
1 < #include "Selection.h"
2 < #include "PassHLT.h"
3 < #include "HZZCiCElectronSelection.h"
4 < #include "HZZLikelihoodElectronSelection.h"
5 < #include "HZZBDTElectronSelection.h"
6 < #include "RunLumiRangeMap.h"
1 > #include "SelectionStatus.h"
2 > #include "EventData.h"
3 > #include "SimpleLepton.h"
4 > #include "EfficiencyWeightsInterface.h"
5  
6 < RunLumiRangeMap rlrm;
6 > #include "HZZBDTElectronSelection.h"
7 > #include "IsolationSelection.h"
8 > #include "PassHLT.h"
9 > #include "Selection.h"
10  
11 < void initRunLumiRangeMap() {
12 <  rlrm.AddJSONFile(std::string("./data/Cert_136033-149442_7TeV_Apr21ReReco_Collisions10_JSON.txt"));
12 <  rlrm.AddJSONFile(std::string("./data/Cert_160404-173244_7TeV_PromptReco_Collisions11_JSON_v2.txt"));
13 <  rlrm.AddJSONFile(std::string("./data/Cert_160404-163869_7TeV_May10ReReco_Collisions11_JSON_v3.txt"));  
14 <  rlrm.AddJSONFile(std::string("./data/Cert_170249-172619_7TeV_ReReco5Aug_Collisions11_JSON.txt"));  
15 < };
16 <
17 < unsigned fails_HZZ4L_selection(ControlFlags &ctrl,           // input control
18 <                               mithep::TEventInfo *info,     // input event inof
19 <                               TClonesArray *electronArr,    // input electrons
20 <                               TClonesArray *muonArr,        // input muons
21 <                               double eventweight,           // weight
22 <                               TTree * passtuple ) {
23 <
24 <  fails_HZZ4L_selection( ctrl, info, electronArr, muonArr, eventweight, passtuple, NULL );
25 <
26 < };
27 <
28 < unsigned fails_HZZ4L_selection(ControlFlags &ctrl,           // input control
29 <                               mithep::TEventInfo *info,     // input event inof
30 <                               TClonesArray *electronArr,    // input electrons
31 <                               TClonesArray *muonArr,        // input muons
32 <                               double eventweight,           // weight
33 <                               LabVectors  *l ) {
34 <
35 <  fails_HZZ4L_selection( ctrl, info, electronArr, muonArr, eventweight, NULL, l );
36 <
37 < };
38 <
39 <
40 < unsigned fails_HZZ4L_selection(ControlFlags &ctrl,           // input control
41 <                               mithep::TEventInfo *info,     // input event inof
42 <                               TClonesArray *electronArr,    // input electrons
43 <                               TClonesArray *muonArr,        // input muons
44 <                               double eventweight,           // weight
45 <                               TTree * passtuple,
46 <                               LabVectors * l) {       // output ntuple
11 > #include "ExternData.h"
12 > #include "SelectionDefs.h"
13  
48  unsigned evtfail = 0x0;
14  
15 + // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
16 + // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
17 + // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
18 + EventData apply_HZZ4L_selection(ControlFlags &ctrl,           // input control
19 +                                mithep::TEventInfo *info,     // input event info
20 +                                TClonesArray *electronArr,    // input electrons
21 +                                SelectionStatus (*ElectronPreSelector)( ControlFlags &, const mithep::TElectron*),
22 +                                SelectionStatus (*ElectronIDSelector)( ControlFlags &, const mithep::TElectron*),
23 +                                SelectionStatus (*ElectronIsoSelector)( ControlFlags &, const mithep::TElectron*),
24 +                                TClonesArray *muonArr,         // input muons
25 +                                SelectionStatus (*MuonPreSelector)( ControlFlags &, const mithep::TMuon*),
26 +                                SelectionStatus (*MuonIDSelector)( ControlFlags &, const mithep::TMuon*),
27 +                                SelectionStatus (*MuonIsoSelector)( ControlFlags &, const mithep::TMuon*) )
28 +  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
29 + {      
30  
31 +  EventData ret;
32 +  unsigned evtfail = 0x0;
33 +  TRandom3 r;
34 +  
35    if( ctrl.debug ) {
36      cout << "Run: " << info->runNum
37           << "\tEvt: " << info->evtNum
# Line 60 | Line 44 | unsigned fails_HZZ4L_selection(ControlFl
44      RunLumiRangeMap::RunLumiPairType rl(info->runNum, info->lumiSec);      
45      if( !(rlrm.HasRunLumi(rl)) )  {
46        if( ctrl.debug ) cout << "\tfails JSON" << endl;
47 <      evtfail |= (1<<EVTFAIL_JSON);
48 <      return evtfail;
47 >      ret.status.setStatus(0);
48 >      return ret;
49      }
50    }
51    
52    
69
70  
71  
53    //********************************************************
54    // Trigger
55    //********************************************************
56    if( !ctrl.mc ) {
57 <    //  if( !(passHLT(info->triggerBits, info->runNum, channel) )  ) {
57 >    //if( !(passHLT(info->triggerBits, info->runNum, channel) )  ) {
58      if( !(passHLT(info->triggerBits, info->runNum, 999) )  ) {
59 +      if( ctrl.debug ) cout << "\tfails trigger" << endl;
60        evtfail |= (1<<EVTFAIL_TRIGGER);
61 <      return evtfail;
61 >      ret.status.setStatus(0);
62 >      return ret;
63      }  
64 <  } else {
82 <    if( !(passHLTMC(info->triggerBits)) ) {
83 <      evtfail |= (1<<EVTFAIL_TRIGGER);
84 <      return evtfail;
85 <    }      
86 <    //    cout << "MC trigger bits: " << hex << info->triggerBits << dec << endl;
87 <  }
88 <
64 >  }
65    if( ctrl.debug ) {
66      cout << "presel nlep: " << muonArr->GetEntries() + electronArr->GetEntries()
67           << "\tnmuon: "    << muonArr->GetEntries()
# Line 101 | Line 77 | unsigned fails_HZZ4L_selection(ControlFl
77    //    
78    if( ctrl.debug ) cout << "\tnMuons: " << muonArr->GetEntries() << endl;
79    //----------------------------------------------------
80 <  for(Int_t i=0; i<muonArr->GetEntries(); i++) {
81 <    const mithep::TMuon *mu = (mithep::TMuon*)((*muonArr)[i]);      
82 <    unsigned muonfail = passMuonSelectionZZ(mu);
83 <    if( ctrl.debug ) {
84 <      cout << "muon:: pt: " << mu->pt
85 <           << "\teta: " << mu->eta
86 <           << "\tmask: 0x" << hex << muonfail << dec
87 <           << endl;
88 <    }
89 <    if ( !muonfail ) {
90 <      SimpleLepton tmplep;
91 <      tmplep.vec.SetPtEtaPhiM(mu->pt,
92 <                              mu->eta,
93 <                              mu->phi,
94 <                              105.658369e-3);
95 <      tmplep.type    = 13;
96 <      tmplep.index   = i;
97 <      tmplep.charge  = mu->q;
98 <      tmplep.isoTrk  = mu->trkIso03;
99 <      tmplep.isoEcal = mu->emIso03;
100 <      tmplep.isoHcal = mu->hadIso03;
101 <      tmplep.ip3dSig = mu->ip3dSig;
102 <      tmplep.is4l    = false;
103 <      tmplep.isEB    = (fabs(mu->eta) < 1.479 ? 1 : 0 );
104 <      lepvec.push_back(tmplep);
105 <      if( ctrl.debug ) { cout << "muon passes ... " << endl;}
106 <    }
107 <  }
80 >  for(Int_t i=0; i<muonArr->GetEntries(); i++)
81 >    {
82 >      const mithep::TMuon *mu = (mithep::TMuon*)((*muonArr)[i]);      
83 >      
84 >      SelectionStatus musel;
85 >      if(ctrl.debug) cout << "musel.status  before anything: " << musel.getStatus() << endl;
86 >      musel |= (*MuonPreSelector)(ctrl,mu);
87 >      if(ctrl.debug) cout << "musel.status  after presel: " << musel.getStatus() << endl;
88 >      musel |= (*MuonIDSelector)(ctrl,mu);
89 >      if(ctrl.debug) cout << "musel.status  after ID: " << musel.getStatus() << endl;
90 >      musel |= (*MuonIsoSelector)(ctrl,mu);
91 >      if(ctrl.debug) cout << "musel.status  after iso: " << musel.getStatus() << endl;
92 >
93 >      if( ctrl.debug ) {
94 >        cout << "muon:: pt: " << mu->pt
95 >             << "\teta: " << mu->eta
96 >             << "\tisorel: " <<  mu->pfIso03/mu->pt
97 >             << "\tstatus: " << hex << musel.getStatus() << dec
98 >             << endl;
99 >      }
100 >      
101 >      if ( musel.pass() ) {
102 >
103 >        SimpleLepton tmplep;
104 >        float pt = mu->pt;
105 >        tmplep.vecorig->SetPtEtaPhiM(pt,
106 >                                    mu->eta,
107 >                                    mu->phi,
108 >                                    MUON_MASS);
109 >        
110 >        if( ctrl.do_escale_up ) {
111 >          pt=scale_smear_muon_Up(pt, 1,  r);
112 >        }
113 >        if( ctrl.do_escale_down ) {
114 >          pt=scale_smear_muon_Down(pt, 1,  r);
115 >        }
116 >        
117 >        tmplep.vec->SetPtEtaPhiM(pt,
118 >                                mu->eta,
119 >                                mu->phi,
120 >                                MUON_MASS);
121 >        
122 >        tmplep.type    = 13;
123 >        tmplep.index   = i;
124 >        tmplep.charge  = mu->q;
125 >        tmplep.isoTrk  = mu->trkIso03;
126 >        tmplep.isoEcal = mu->emIso03;
127 >        tmplep.isoHcal = mu->hadIso03;
128 >        tmplep.isoPF03 = mu->pfIso03;
129 >        tmplep.isoPF04 = mu->pfIso04;
130 >        tmplep.ip3dSig = mu->ip3dSig;
131 >        tmplep.is4l    = false;
132 >        tmplep.isEB    = (fabs(mu->eta) < 1.479 ? 1 : 0 );
133 >        tmplep.isTight = musel.tight();
134 >        tmplep.isLoose = musel.loose();
135 >        lepvec.push_back(tmplep);
136 >        if( ctrl.debug ) { cout << "muon passes ... " << endl;}
137 >      }
138 >      //  }
139 >    }    
140 >  
141    
142 <  if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
142 >    
143 >    //
144 >    if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
145 >    // --------------------------------------------------------------------------------
146 >    for(Int_t i=0; i<electronArr->GetEntries(); i++)
147 >      {
148 >        const mithep::TElectron *ele = (mithep::TElectron*)((*electronArr)[i]);
149 >        
150 >        Bool_t isMuonOverlap = kFALSE;
151 >        for (int k=0; k<lepvec.size(); ++k) {
152 >          TVector3 tmplep;
153 >          tmplep.SetPtEtaPhi(ele->pt, ele->eta, ele->phi);
154 >          if ( lepvec[k].isLoose && lepvec[k].type == 13 && lepvec[k].vec->Vect().DrEtaPhi(tmplep) < 0.1 ) {
155 >            if( ctrl.debug ) cout << "-----> isMuonOverlap! " << endl;
156 >            isMuonOverlap = kTRUE;
157 >            break;
158 >          }
159 >        }
160 >        
161 >        SelectionStatus elesel;
162 >        if( ctrl.debug ) cout << "--> status before anything: " << hex << elesel.getStatus() << dec << endl;
163 >        elesel |= (*ElectronPreSelector)(ctrl,ele);
164 >        if( ctrl.debug ) cout << "--> status after presel: " << hex << elesel.getStatus() << dec << endl;
165 >        elesel |= (*ElectronIDSelector)(ctrl,ele);
166 >        if( ctrl.debug ) cout << "--> status after ID: " << hex << elesel.getStatus() << dec << endl;
167 >        elesel |= (*ElectronIsoSelector)(ctrl,ele);
168 >        if( ctrl.debug ) cout << "--> status after iso: " << hex << elesel.getStatus() << dec << endl;
169 >        
170 >        if( ctrl.debug ){
171 >          cout << "\tscEt: " << ele->scEt
172 >               << "\tscEta: " << ele->scEta
173 >               << "\tncluster: " << ele->ncluster
174 >               << "\tisorel: " <<  ele->pfIso04/ele->pt
175 >               << "\tstatus: " << hex << elesel.getStatus() << dec
176 >               << endl;
177 >        }
178  
179 <  //----------------------------------------------------
180 <  for(Int_t i=0; i<electronArr->GetEntries(); i++) {
181 <    const mithep::TElectron *ele = (mithep::TElectron*)((*electronArr)[i]);
179 >        if ( elesel.pass() && !isMuonOverlap )
180 >          {
181 >            SimpleLepton tmplep;
182 >            
183 >            float pt = ele->pt;
184 >            tmplep.vecorig->SetPtEtaPhiM( pt,
185 >                                          ele->eta,
186 >                                          ele->phi,
187 >                                          ELECTRON_MASS );
188 >            
189 >            if( ctrl.do_escale ) {
190 >              pt=scale_smear_electron(pt, ele->isEB, r);
191 >            }
192 >            if( ctrl.do_escale_up ) {
193 >              pt=scale_smear_electron_Up(pt, ele->isEB,  r);
194 >            }
195 >            if( ctrl.do_escale_down ) {
196 >              pt=scale_smear_electron_Down(pt, ele->isEB,  r);
197 >            }
198  
199 <    Bool_t isMuonOverlap = kFALSE;
200 <    for (int k=0; k<lepvec.size(); ++k) {
201 <      TVector3 tmplep;
202 <      tmplep.SetPtEtaPhi(ele->pt, ele->eta, ele->phi);
203 <      if ( lepvec[k].type == 13 && lepvec[k].vec.Vect().DrEtaPhi(tmplep) < 0.1 ) {
204 <        if( ctrl.debug ) cout << "-----> isMuonOverlap! " << endl;
205 <        isMuonOverlap = kTRUE;
206 <        break;
199 >            
200 >            tmplep.vec->SetPtEtaPhiM( pt,
201 >                                      ele->eta,
202 >                                      ele->phi,
203 >                                      ELECTRON_MASS );
204 >            
205 >            tmplep.type    = 11;
206 >            tmplep.index   = i;
207 >            tmplep.charge  = ele->q;
208 >            tmplep.isoTrk  = ele->trkIso03;
209 >            tmplep.isoEcal = ele->emIso03;
210 >            tmplep.isoHcal = ele->hadIso03;
211 >            tmplep.isoPF03 = ele->pfIso03;
212 >            tmplep.isoPF04 = ele->pfIso04;
213 >            tmplep.ip3dSig = ele->ip3dSig;
214 >            tmplep.is4l    = false;
215 >            tmplep.isEB    = ele->isEB;
216 >            tmplep.scID    = ele->scID;
217 >            tmplep.isTight = elesel.tight();
218 >            tmplep.isLoose = elesel.loose();
219 >            lepvec.push_back(tmplep);
220 >            if( ctrl.debug ) { cout << "\telectron passes ... " << endl; }
221 >          }
222        }
223 +    
224 +    
225 +    //********************************************************
226 +    // Dump Stuff
227 +    //********************************************************
228 +    sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
229 +    int nmu=0, nele=0;
230 +    for( int i=0; i<lepvec.size(); i++ ) {
231 +      if(ctrl.debug) cout << "lepvec :: index: " << i
232 +                          << "\tpt: " << lepvec[i].vec->Pt()
233 +                          << "\ttype: " << lepvec[i].type
234 +                          << endl;
235 +      if( abs(lepvec[i].type) == 11 ) nele++;
236 +      else nmu++;
237      }
238 <
239 <    unsigned FAIL=0;
240 <    CICStruct ciccuts = getCiCCuts(ctrl.eleSeleScheme);
241 <    unsigned  failsCIC=0;
242 <    if(ctrl.eleSele=="cic") {
154 <      failsCIC = failsCicSelection(ctrl, ele, ciccuts, ctrl.kinematics);
155 <      FAIL = failsCIC;
156 <    }
157 <    LikStruct likcuts;
158 <    unsigned failsLike=0;
159 <    if(ctrl.eleSele=="lik") {
160 <      likcuts = getLikCuts(ctrl.eleSeleScheme);
161 <      failsLike = failsLikelihoodSelection(ele, likcuts, ctrl.kinematics);
162 <      FAIL = failsLike;
238 >    if( ctrl.debug ) {
239 >      cout << "postsel nlep: " << lepvec.size()
240 >           << "\tnmuon: " << nmu
241 >           << "\tnelectron: " << nele
242 >           << endl;
243      }
244 <    unsigned  failsBDT=0;
245 <    if(ctrl.eleSele=="bdt") {
246 <      failsBDT = failsBDTSelection(ele,ctrl.kinematics,ctrl.eleSeleScheme);
247 <      FAIL = failsBDT;
244 >    
245 >    
246 >    //******************************************************************************
247 >    // Z1 Selection
248 >    //******************************************************************************
249 >    int Z1LeptonPlusIndex = -1;
250 >    int Z1LeptonMinusIndex = -1;
251 >    double BestZ1Mass = -999;
252 >    if( ctrl.debug ) { cout << "looking for a Z1 ..." << endl; }
253 >    for(int i = 0; i < lepvec.size(); ++i) {
254 >      if( !(lepvec[i].isLoose) ) continue;
255 >      for(int j = i+1; j < lepvec.size(); ++j) {
256 >        if( !(lepvec[j].isLoose) ) continue;
257 >        if( ctrl.debug ) { cout << "\tconsidering leptons " << i << " & " << j << endl; }
258 >        if (!(lepvec[i].vec->Pt() > 20.0 || lepvec[j].vec->Pt() > 20.0)) continue;
259 >        if( ctrl.debug ) { cout << "\tat least one is > 20 GeV" << endl; }
260 >        if (!(lepvec[i].vec->Pt() > 10.0 && lepvec[j].vec->Pt() > 10.0)) continue;
261 >        if( ctrl.debug ) { cout << "\tthe other  is > 10 GeV" << endl; }
262 >        if (lepvec[i].charge == lepvec[j].charge) continue;          
263 >        if( ctrl.debug ) { cout << "\tthey're opposite charge" << endl; }
264 >        if (fabs(lepvec[i].type) != fabs(lepvec[j].type)) continue;  
265 >        if( ctrl.debug ) { cout << "\tthey're same flavor" << endl; }
266 >        
267 >        //Make Z1 hypothesis
268 >        TLorentzVector *leptonPlus, *leptonMinus;
269 >        if ( lepvec[i].charge > 0 ) {
270 >          leptonPlus  = lepvec[i].vec;
271 >          leptonMinus = lepvec[j].vec;
272 >        } else {
273 >          leptonPlus  = lepvec[j].vec;
274 >          leptonMinus = lepvec[i].vec;
275 >        }
276 >        
277 >        float tmpZ1Mass = (*leptonPlus + *leptonMinus).M();
278 >        if( ctrl.debug ) cout << "Z1 selection, tmpZ1Mass: " << tmpZ1Mass << endl;
279 >        if( tmpZ1Mass > 60 ) {
280 >          if (fabs(tmpZ1Mass - Z_MASS) < fabs(BestZ1Mass - Z_MASS)) {
281 >            BestZ1Mass = tmpZ1Mass;
282 >            if( ctrl.debug ) cout << "Z1 selection, new BestZ1Mass: " << BestZ1Mass
283 >                                  << "\tdM: " << fabs(BestZ1Mass - Z_MASS)
284 >                                  << endl;
285 >            if (lepvec[i].charge > 0) {
286 >              Z1LeptonPlusIndex = i;
287 >              Z1LeptonMinusIndex = j;
288 >            } else {
289 >              Z1LeptonPlusIndex = j;
290 >              Z1LeptonMinusIndex = i;
291 >            }
292 >          }
293 >        }
294 >      }
295      }
296 <
297 <    if( ctrl.debug ){
298 <      cout << "CIC category: " << cicCategory(ele)
299 <           << "\tlikelihood: " << ele->likelihood
300 <           << "\tFAIL:  0x"     << hex << FAIL      << dec
301 <           << "\tfailsCIC:  0x" << hex << failsCIC  << dec
175 <           << "\tfailsLike: 0x" << hex << failsLike << dec
176 <           << "\tfailsBDT:  0x" << hex << failsBDT << dec
177 <           << "\tscEt: " << ele->scEt
178 <           << "\tscEta: " << ele->scEta
179 <           << "\tncluster: " << ele->ncluster
180 <           << endl;
296 >    // stop if no Z1 candidate is found
297 >    if( BestZ1Mass < 0 ) {
298 >      evtfail |= (1<<EVTFAIL_Z1);
299 >      //ret.status = evtfail;
300 >      ret.status.setStatus(0);
301 >      return ret;
302      }
303 <    if ( !FAIL && !isMuonOverlap ) {
304 <      SimpleLepton tmplep;
305 <      tmplep.vec.SetPtEtaPhiM( ele->pt,
306 <                               ele->eta,
307 <                               ele->phi,
308 <                               0.51099892e-3 );
309 <      tmplep.type    = 11;
310 <      tmplep.index   = i;
311 <      tmplep.charge  = ele->q;
312 <      tmplep.isoTrk  = ele->trkIso03;
313 <      tmplep.isoEcal = ele->emIso03;
314 <      tmplep.isoHcal = ele->hadIso03;
315 <      tmplep.ip3dSig = ele->ip3dSig;
316 <      tmplep.is4l    = false;
317 <      tmplep.isEB    = ele->isEB;
197 <      lepvec.push_back(tmplep);
198 <      if( ctrl.debug ) { cout << "\telectron passes ... " << endl; }
303 >    if( ctrl.debug ) cout << "\tgot a Z1 ... run: " << info->runNum << "\tevt: " << info->evtNum  << endl;
304 >    if( ctrl.debug ) cout << "\tZ1 plusindex: " << Z1LeptonPlusIndex << "\tminusindex: " << Z1LeptonMinusIndex << endl;
305 >    TLorentzVector Z1LeptonPlus  = *(lepvec[Z1LeptonPlusIndex].vec);
306 >    TLorentzVector Z1LeptonMinus = *(lepvec[Z1LeptonMinusIndex].vec);
307 >    TLorentzVector Z1Candidate   =  Z1LeptonPlus + Z1LeptonMinus;
308 >    
309 >    
310 >    //******************************************************************************
311 >    // Z1 + l
312 >    //******************************************************************************
313 >    if( lepvec.size() < 3 ) {
314 >      evtfail |= (1<<EVTFAIL_Z1_PLUSL);
315 >      //ret.status = evtfail;
316 >      ret.status.setStatus(0);
317 >      return ret;
318      }
319 <  }
320 <  
321 <  sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
322 <
323 <  int nmu=0, nele=0;
324 <  for( int i=0; i<lepvec.size(); i++ ) {
325 <    if( abs(lepvec[i].type) == 11 ) nele++;
326 <    else nmu++;
327 <  }
328 <  if( ctrl.debug ) {
329 <    cout << "postsel nlep: " << lepvec.size()
330 <         << "\tnmuon: " << nmu
331 <         << "\tnelectron: " << nele
332 <         << endl;
214 <  }
215 <
216 <  //******************************************************************************
217 <  //Z1 Selection
218 <  //******************************************************************************
219 <  int Z1LeptonPlusIndex = -1;
220 <  int Z1LeptonMinusIndex = -1;
221 <  double BestZ1Mass = -999;
222 <  if( ctrl.debug ) { cout << "looking for a Z1 ..." << endl; }
223 <  for(int i = 0; i < lepvec.size(); ++i) {
224 <    for(int j = i+1; j < lepvec.size(); ++j) {
225 <      if( ctrl.debug ) { cout << "\tconsidering leptons " << i << " & " << j << endl; }
226 <      if (!(lepvec[i].vec.Pt() > 20.0 || lepvec[j].vec.Pt() > 20.0)) continue;
227 <      if( ctrl.debug ) { cout << "\tat least one is > 20 GeV" << endl; }
228 <      if (!(lepvec[i].vec.Pt() > 10.0 && lepvec[j].vec.Pt() > 10.0)) continue;
229 <      if( ctrl.debug ) { cout << "\tthe other  is > 10 GeV" << endl; }
230 <      if (lepvec[i].charge == lepvec[j].charge) continue;          
231 <      if( ctrl.debug ) { cout << "\tthey're opposite charge" << endl; }
232 <      if (fabs(lepvec[i].type) != fabs(lepvec[j].type)) continue;  
233 <      if( ctrl.debug ) { cout << "\tthey're same flavor" << endl; }
319 >    
320 >    //******************************************************************************
321 >    // 4l/Z2 Selection
322 >    //******************************************************************************
323 >    Int_t Z2LeptonPlusIndex = -1;
324 >    Int_t Z2LeptonMinusIndex = -1;
325 >    Double_t BestZ2Mass = -1;
326 >    if( ctrl.debug ) cout << "looking for a Z2 ... out of " << lepvec.size() << " leptons" <<endl;
327 >    for(int i = 0; i < lepvec.size(); ++i) {
328 >      
329 >      if( ctrl.debug)  cout << "i: " << i
330 >                            << "\tpt: " << lepvec[i].vec->Pt()
331 >                            << "\ttype: " << lepvec[i].type
332 >                            << endl;
333        
334 <      //Make Z1 hypothesis
335 <      TLorentzVector leptonPlus, leptonMinus;
336 <      if ( lepvec[i].charge > 0 ) {
337 <        leptonPlus  = lepvec[i].vec;
239 <        leptonMinus = lepvec[j].vec;
240 <      } else {
241 <        leptonPlus  = lepvec[j].vec;
242 <        leptonMinus = lepvec[i].vec;
334 >      if( ctrl.eleSeleScheme == "mediumloose" &&
335 >          !(lepvec[i].isTight) ) {
336 >        if( ctrl.debug)  cout << "it's not tight, skipping ... " << endl;
337 >        continue;
338        }
339        
340 <      float tmpZ1Mass = (leptonPlus+leptonMinus).M();
341 <      if( ctrl.debug ) cout << "Z1 selection, tmpZ1Mass: " << tmpZ1Mass << endl;
342 <      if( tmpZ1Mass > 60 ) {
343 <        if (fabs(tmpZ1Mass - 91.1876) < fabs(BestZ1Mass - 91.1876)) {
344 <          BestZ1Mass = tmpZ1Mass;
345 <          if( ctrl.debug ) cout << "Z1 selection, new BestZ1Mass: " << BestZ1Mass
346 <               << "\tdM: " << fabs(BestZ1Mass - 91.1876)
347 <               << endl;
340 >      for(int j = i+1; j < lepvec.size(); ++j) {
341 >        if( ctrl.debug)  cout << "\t\tj: " << j
342 >                              << "\tpt: " << lepvec[j].vec->Pt()
343 >                              << "\ttype: " << lepvec[j].type
344 >                              << endl;
345 >        
346 >        if( ctrl.eleSeleScheme == "mediumloose" &&
347 >            !(lepvec[j].isTight) ) {
348 >          if( ctrl.debug)  cout << "it's not tight, skipping ... " << endl;
349 >          continue;
350 >        }
351 >        
352 >        
353 >        if (i == Z1LeptonPlusIndex || i == Z1LeptonMinusIndex) {
354 >          if( ctrl.debug)  cout << "\ti matches a Z1 index, skipping ..." << endl;
355 >          continue; //skip Z1 leptons
356 >        }
357 >        if (j == Z1LeptonPlusIndex || j == Z1LeptonMinusIndex) {
358 >          if( ctrl.debug)  cout << "\tj matches a Z1 index, skipping ..." << endl;
359 >          continue; //skip Z1 leptons
360 >        }
361 >        if (lepvec[i].charge == lepvec[j].charge) {
362 >          if( ctrl.debug)  cout << "\ti and j are same sign, skipping ..." << endl;
363 >          continue;         //require opp sign
364 >        }
365 >        if (fabs(lepvec[i].type) != fabs(lepvec[j].type)) {
366 >          if( ctrl.debug) cout << "\ti and j are not same flavor, skipping ..." << endl;
367 >          continue; //require same flavor
368 >        }
369 >        
370 >        
371 >        //Make Z2 hypothesis
372 >        TLorentzVector *leptonPlus, *leptonMinus;
373 >        
374 >        if (lepvec[i].charge > 0 ) {
375 >          leptonPlus  = lepvec[i].vec;
376 >          leptonMinus = lepvec[j].vec;
377 >        } else {
378 >          leptonPlus  = lepvec[j].vec;
379 >          leptonMinus = lepvec[i].vec;
380 >        }
381 >        
382 >        TLorentzVector dilepton = *leptonPlus + *leptonMinus;
383 >        TLorentzVector fourLepton = Z1Candidate + dilepton;
384 >        
385 >        if( ctrl.debug ) cout << "dilepton.M() : " << dilepton.M() << endl;
386 >        if( ctrl.debug ) cout << "fourLepton.M() : " << fourLepton.M() << endl;
387 >        
388 >        if (!(dilepton.M() > 12.0)) continue;
389 >        if (!(fourLepton.M() > 100.0)) continue;
390 >        
391 >        //for 4e and 4mu, require at least 1 of the other opp sign lepton pairs have mass > 12
392 >        if (fabs(lepvec[i].type) == fabs(lepvec[Z1LeptonPlusIndex].type)) {
393 >          TLorentzVector pair1 = Z1LeptonPlus + *leptonMinus;
394 >          TLorentzVector pair2 = Z1LeptonMinus + *leptonPlus;
395 >          if( ctrl.debug ) cout << "pair1: " << pair1.M() << "\tpair2: "<< pair2.M() << endl;
396 >          if (!(pair1.M() > 12 || pair2.M() > 12)) continue;
397 >        }
398 >      
399 >      
400 >        //Disambiguiation is done by choosing the pair with the largest ptMax and largest ptMin
401 >        if (Z2LeptonPlusIndex < 0) {
402            if (lepvec[i].charge > 0) {
403 <            Z1LeptonPlusIndex = i;
404 <            Z1LeptonMinusIndex = j;
403 >            Z2LeptonPlusIndex = i;
404 >            Z2LeptonMinusIndex = j;
405            } else {
406 <            Z1LeptonPlusIndex = j;
407 <            Z1LeptonMinusIndex = i;
406 >            Z2LeptonPlusIndex = j;
407 >            Z2LeptonMinusIndex = i;
408 >          }
409 >        } else {
410 >          Double_t BestPairPtMax = lepvec[Z2LeptonPlusIndex].vec->Pt();              
411 >          Double_t BestPairPtMin = lepvec[Z2LeptonMinusIndex].vec->Pt();
412 >          if (lepvec[Z2LeptonMinusIndex].vec->Pt() > BestPairPtMax) {
413 >            BestPairPtMax = lepvec[Z2LeptonMinusIndex].vec->Pt();
414 >            BestPairPtMin = lepvec[Z2LeptonPlusIndex].vec->Pt();
415            }
260        }
261      }
262    }
263  }
264  // stop if no Z1 candidate is found
265  if( BestZ1Mass < 0 ) {
266    evtfail |= (1<<EVTFAIL_Z1);
267    return evtfail;
268  }
269  if( ctrl.debug ) cout << "\tgot a Z1 ... run: " << info->runNum << "\tevt: " << info->evtNum  << endl;
270  if( ctrl.debug ) cout << "\tZ1 plusindex: " << Z1LeptonPlusIndex << "\tminusindex: " << Z1LeptonMinusIndex << endl;
271  TLorentzVector Z1LeptonPlus  = lepvec[Z1LeptonPlusIndex].vec;
272  TLorentzVector Z1LeptonMinus = lepvec[Z1LeptonMinusIndex].vec;
273  TLorentzVector Z1Candidate   =  Z1LeptonPlus + Z1LeptonMinus;
274  if( l != NULL ) {
275    l->vecz1 = Z1Candidate;
276    l->vecl1p = Z1LeptonPlus;
277    l->vecl1m = Z1LeptonMinus;
278  }
279  
280  //******************************************************************************
281  // Z1 + l
282  //******************************************************************************
283  if( lepvec.size() < 3 ) {
284    evtfail |= (1<<EVTFAIL_Z1_PLUSL);
285    return evtfail;
286  }
416          
417 <  //******************************************************************************
418 <  // 4l/Z2 Selection
419 <  //******************************************************************************
420 <  Int_t Z2LeptonPlusIndex = -1;
421 <  Int_t Z2LeptonMinusIndex = -1;
422 <  Double_t BestZ2Mass = -1;
423 <  if( ctrl.debug ) cout << "looking for a Z2 ... out of " << lepvec.size() << " leptons" <<endl;
424 <  for(int i = 0; i < lepvec.size(); ++i) {
425 <    for(int j = i+1; j < lepvec.size(); ++j) {
426 <      //            cout << "i: " << i << "\tj: " << j << endl;
427 <      if (i == Z1LeptonPlusIndex || i == Z1LeptonMinusIndex) {
299 <        //            cout << "\ti matches a Z1 index, skipping ..." << endl;
300 <        continue; //skip Z1 leptons
301 <      }
302 <      if (j == Z1LeptonPlusIndex || j == Z1LeptonMinusIndex) {
303 <        //            cout << "\tj matches a Z1 index, skipping ..." << endl;
304 <              continue; //skip Z1 leptons
305 <            }
306 <            if (lepvec[i].charge == lepvec[j].charge) {
307 <              //              cout << "\ti and j are same sign, skipping ..." << endl;
308 <              continue;         //require opp sign
309 <            }
310 <            if (fabs(lepvec[i].type) != fabs(lepvec[j].type)) {
311 <              //              cout << "\ti and j are not same flavor, skipping ..." << endl;
312 <              continue; //require same flavor
313 <            }
314 <            
315 <            
316 <            //Make Z2 hypothesis
317 <            TLorentzVector leptonPlus, leptonMinus;
318 <            
319 <            if (lepvec[i].charge > 0 ) {
320 <              leptonPlus  = lepvec[i].vec;
321 <              leptonMinus = lepvec[j].vec;
322 <            } else {
323 <              leptonPlus  = lepvec[j].vec;
324 <              leptonMinus = lepvec[i].vec;
325 <            }
326 <
327 <            TLorentzVector dilepton = leptonPlus+leptonMinus;
328 <            TLorentzVector fourLepton = Z1Candidate + dilepton;
329 <
330 <            if( ctrl.debug ) cout << "dilepton.M() : " << dilepton.M() << endl;
331 <            if( ctrl.debug ) cout << "fourLepton.M() : " << fourLepton.M() << endl;
332 <
333 <            if (!(dilepton.M() > 12.0)) continue;
334 <            if (!(fourLepton.M() > 100.0)) continue;
335 <
336 <            
337 <            //for 4e and 4mu, require at least 1 of the other opp sign lepton pairs have mass > 12
338 <            if (fabs(lepvec[i].type) == fabs(lepvec[Z1LeptonPlusIndex].type)) {
339 <              TLorentzVector pair1 = Z1LeptonPlus+leptonMinus;
340 <              TLorentzVector pair2 = Z1LeptonMinus+leptonPlus;
341 <              if( ctrl.debug ) cout << "pair1: " << pair1.M() << "\tpair2: "<< pair2.M() << endl;
342 <              if (!(pair1.M() > 12 || pair2.M() > 12)) continue;
343 <            }
344 <            
345 <            //Disambiguiation is done by choosing the pair with the largest ptMax and largest ptMin
346 <            if (Z2LeptonPlusIndex < 0) {
347 <              if (lepvec[i].charge > 0) {
348 <                Z2LeptonPlusIndex = i;
349 <                Z2LeptonMinusIndex = j;
350 <              } else {
351 <                Z2LeptonPlusIndex = j;
352 <                Z2LeptonMinusIndex = i;
353 <              }
354 <            } else {
355 <              Double_t BestPairPtMax = lepvec[Z2LeptonPlusIndex].vec.Pt();              
356 <              Double_t BestPairPtMin = lepvec[Z2LeptonMinusIndex].vec.Pt();
357 <              if (lepvec[Z2LeptonMinusIndex].vec.Pt() > BestPairPtMax) {
358 <                BestPairPtMax = lepvec[Z2LeptonMinusIndex].vec.Pt();
359 <                BestPairPtMin = lepvec[Z2LeptonPlusIndex].vec.Pt();
360 <              }
361 <
362 <              Double_t CurrentPairPtMax = lepvec[i].vec.Pt();              
363 <              Double_t CurrentPairPtMin = lepvec[j].vec.Pt();
364 <              if (lepvec[j].vec.Pt() > CurrentPairPtMax) {
365 <                CurrentPairPtMax = lepvec[j].vec.Pt();
366 <                CurrentPairPtMin = lepvec[i].vec.Pt();
367 <              }
368 <
369 <              if (CurrentPairPtMax > BestPairPtMax) {
370 <                if (lepvec[i].charge > 0) {
371 <                  Z2LeptonPlusIndex = i;
372 <                  Z2LeptonMinusIndex = j;
373 <                } else {
374 <                  Z2LeptonPlusIndex = j;
375 <                  Z2LeptonMinusIndex = i;
376 <                }
377 <              } else if (CurrentPairPtMax  == BestPairPtMax) {
378 <                if (CurrentPairPtMin > BestPairPtMin) {
379 <                  if (lepvec[i].charge > 0) {
380 <                    Z2LeptonPlusIndex = i;
381 <                    Z2LeptonMinusIndex = j;
382 <                  } else {
383 <                    Z2LeptonPlusIndex = j;
384 <                    Z2LeptonMinusIndex = i;
385 <                  }                  
386 <                }
387 <              }
388 <            }            
389 <          }
390 <        }
391 <
392 <        // stop if no Z2 candidate is found
393 <        if (Z2LeptonPlusIndex == -1) {
394 <          evtfail |= ( 1<<EVTFAIL_4L );
395 <          return evtfail;
396 <          //      h_evtfail->Fill( evtfail );
397 <          //      cout << "evtfail: " << hex << evtfail << dec << endl;
398 <          //      continue;
399 <        }
400 <        if( ctrl.debug ) cout << "\tgot a Z2 ..." << endl;
401 <        if( ctrl.debug ) cout << "\tZ2 plusindex: " << Z2LeptonPlusIndex
402 <                              << "\tminusindex: " << Z2LeptonMinusIndex << endl;
403 <        TLorentzVector  Z2LeptonPlus  = lepvec[Z2LeptonPlusIndex].vec;
404 <        TLorentzVector  Z2LeptonMinus = lepvec[Z2LeptonMinusIndex].vec;
405 <        TLorentzVector  Z2Candidate   = Z2LeptonPlus+Z2LeptonMinus;
406 <        TLorentzVector  ZZSystem      = Z1Candidate + Z2Candidate;
407 <        if( l != NULL ) {
408 <          l->vecz2 = Z2Candidate;
409 <          l->vecl2p = Z2LeptonPlus;
410 <          l->vecl2m = Z2LeptonMinus;
411 <          l->vec4l = ZZSystem;
412 <        }
413 <        lepvec[Z1LeptonPlusIndex].is4l = true;
414 <        lepvec[Z1LeptonMinusIndex].is4l = true;
415 <        lepvec[Z2LeptonPlusIndex].is4l = true;
416 <        lepvec[Z2LeptonMinusIndex].is4l = true;
417 <
418 <        //***************************************************************
419 <        // Isolation
420 <        //***************************************************************
421 <        bool failiso=false;
422 <
423 <        /*
424 <        int i,j;
425 <        i=Z1LeptonPlusIndex;
426 <        j=Z1LeptonMinusIndex;
427 <        float RIso1 = (lepvec[i].isoTrk+lepvec[i].isoEcal+lepvec[i].isoHcal)/lepvec[i].vec.Pt();
428 <        float RIso2 = (lepvec[j].isoTrk+lepvec[j].isoEcal+lepvec[j].isoHcal)/lepvec[j].vec.Pt();
429 <        float comboIso12 = RIso1 + RIso2;
430 <        i=Z2LeptonPlusIndex;
431 <        j=Z2LeptonMinusIndex;
432 <        float RIso3 = (lepvec[i].isoTrk+lepvec[i].isoEcal+lepvec[i].isoHcal)/lepvec[i].vec.Pt();
433 <        float RIso4 = (lepvec[j].isoTrk+lepvec[j].isoEcal+lepvec[j].isoHcal)/lepvec[j].vec.Pt();
434 <        float comboIso34 = RIso3 + RIso4;
435 <        if( comboIso12 > 0.35 || comboIso34 > 0.35 ) {
436 <              failiso = true;
437 <        }
438 <        */
439 <
440 <        float rho = info->rho;
441 <        for( int i=0; i<lepvec.size(); i++ ) {
442 <          if( !(lepvec[i].is4l) ) continue;
443 <            float effArea_ecal_i, effArea_hcal_i;
444 <            if( lepvec[i].isEB ) {
445 <              if( lepvec[i].type == 11 ) {
446 <                effArea_ecal_i = 0.101;
447 <                effArea_hcal_i = 0.021;
448 <              } else {
449 <                effArea_ecal_i = 0.074;
450 <                effArea_hcal_i = 0.022;
451 <              }
417 >          Double_t CurrentPairPtMax = lepvec[i].vec->Pt();              
418 >          Double_t CurrentPairPtMin = lepvec[j].vec->Pt();
419 >          if (lepvec[j].vec->Pt() > CurrentPairPtMax) {
420 >            CurrentPairPtMax = lepvec[j].vec->Pt();
421 >            CurrentPairPtMin = lepvec[i].vec->Pt();
422 >          }
423 >        
424 >          if (CurrentPairPtMax > BestPairPtMax) {
425 >            if (lepvec[i].charge > 0) {
426 >              Z2LeptonPlusIndex = i;
427 >              Z2LeptonMinusIndex = j;
428              } else {
429 <              if( lepvec[i].type == 11 ) {
430 <                effArea_ecal_i = 0.046;
455 <                effArea_hcal_i = 0.040;
456 <              } else {
457 <                effArea_ecal_i = 0.045;
458 <                effArea_hcal_i = 0.030;
459 <              }
429 >              Z2LeptonPlusIndex = j;
430 >              Z2LeptonMinusIndex = i;
431              }
432 <            float isoEcal_corr_i = lepvec[i].isoEcal - (effArea_ecal_i*rho);
433 <            float isoHcal_corr_i = lepvec[i].isoHcal - (effArea_hcal_i*rho);
434 <            for( int j=i+1; j<lepvec.size(); j++ ) {
435 <              if( !(lepvec[j].is4l) ) continue;
436 <              float effArea_ecal_j, effArea_hcal_j;
466 <              if( lepvec[j].isEB ) {
467 <                if( lepvec[j].type == 11 ) {
468 <                  effArea_ecal_j = 0.101;
469 <                  effArea_hcal_j = 0.021;
470 <                } else {
471 <                  effArea_ecal_j = 0.074;
472 <                  effArea_hcal_j = 0.022;
473 <                }
432 >          } else if (CurrentPairPtMax  == BestPairPtMax) {
433 >            if (CurrentPairPtMin > BestPairPtMin) {
434 >              if (lepvec[i].charge > 0) {
435 >                Z2LeptonPlusIndex = i;
436 >                Z2LeptonMinusIndex = j;
437                } else {
438 <                if( lepvec[j].type == 11 ) {
439 <                  effArea_ecal_j = 0.046;
440 <                  effArea_hcal_j = 0.040;
478 <                } else {
479 <                  effArea_ecal_j = 0.045;
480 <                  effArea_hcal_j = 0.030;
481 <                }
482 <              }
483 <              float isoEcal_corr_j = lepvec[j].isoEcal - (effArea_ecal_j*rho);
484 <              float isoHcal_corr_j = lepvec[j].isoHcal - (effArea_hcal_j*rho);
485 <              float RIso_i = (lepvec[i].isoTrk+isoEcal_corr_i+isoHcal_corr_i)/lepvec[i].vec.Pt();
486 <              float RIso_j = (lepvec[j].isoTrk+isoEcal_corr_j+isoHcal_corr_j)/lepvec[j].vec.Pt();          
487 <              float comboIso = RIso_i + RIso_j;
488 <              if( info->evtNum == 1038911933 ) {
489 <                float tmpdR = lepvec[i].vec.DrEtaPhi(lepvec[j].vec);
490 <                cout << "i: " << i
491 <                     << "\tdR: " << tmpdR
492 <                     << "\trho: " << rho
493 <                     << "\tRIso_i: " << RIso_i  
494 <                     << "\ttkrel: " << lepvec[i].isoTrk/lepvec[i].vec.Pt()
495 <                     << "\tecalrel: " << lepvec[i].isoEcal/lepvec[i].vec.Pt()
496 <                     << "\tecalrelcor: " << isoEcal_corr_i/lepvec[i].vec.Pt()
497 <                     << "\thcalrel: " << lepvec[i].isoHcal/lepvec[i].vec.Pt()
498 <                     << "\thcalrelcor: " << isoHcal_corr_i/lepvec[i].vec.Pt()
499 <                     << "\tpt_i: " << lepvec[i].vec.Pt()
500 <                     << "\tj: " << j
501 <                     << "\tRIso_j: " << RIso_j  
502 <                     << "\ttkrel: "  << lepvec[j].isoTrk/lepvec[j].vec.Pt()
503 <                     << "\tecalrel: " << lepvec[j].isoEcal/lepvec[j].vec.Pt()
504 <                     << "\tecalrelcor: " << isoEcal_corr_j/lepvec[j].vec.Pt()
505 <                     << "\thcalrel: " << lepvec[j].isoHcal/lepvec[j].vec.Pt()
506 <                     << "\thcalrelcor: " << isoHcal_corr_j/lepvec[j].vec.Pt()
507 <                     << "\tpt_j: " << lepvec[j].vec.Pt()
508 <                     << "\tcombo: " << comboIso
509 <                     << endl;
510 <                cout.flush();
511 <              }
512 <              if( comboIso > 0.35 ) {
513 <                if( ctrl.debug ) cout << "combo failing for indices: " << i << "," << j << endl;
514 <                failiso = true;
515 <                //            break;
516 <              }
438 >                Z2LeptonPlusIndex = j;
439 >                Z2LeptonMinusIndex = i;
440 >              }                  
441              }
518        }
519        if( failiso ) {
520          evtfail |= ( 1<<EVTFAIL_ISOLATION );
521          return evtfail;
522          //h_evtfail->Fill( evtfail, eventweight );
523          //      h_evtfail->Fill( evtfail );
524          //      cout << "evtfail: " << hex << evtfail << dec << endl;
525          //      continue;
526        }
527
528        //***************************************************************
529        // IP significance
530        //***************************************************************
531        bool failip = false;
532        for( int i=0; i<lepvec.size(); i++ ) {
533          if( !(lepvec[i].is4l) ) continue;
534          if( lepvec[i].ip3dSig > 4 ) {
535            failip=true;
536            break;
442            }
443 <        }
444 <        if( failip ) {
445 <          evtfail |= (1<<EVTFAIL_IP );
446 <          return evtfail;
447 <          //h_evtfail->Fill( evtfail, eventweight );
448 <          //      h_evtfail->Fill( evtfail );
449 <          //      cout << "evtfail: " << hex << evtfail << dec << endl;
450 <          //      continue;
451 <        }
452 <
453 <        //***************************************************************
454 <        // remaining kinematic cuts
455 <        //***************************************************************
456 <        double Z2massCut=0;
457 <        if      ( ctrl.kinematics == "loose" ) Z2massCut = 12;
458 <        else if ( ctrl.kinematics == "tight" ) Z2massCut = 20;
459 <        else { cout <<  "error! kinematic tightness not defined!" << endl; assert(0); }
460 <
461 <        if ( Z1Candidate.M() > 120        ||
462 <             Z2Candidate.M() < Z2massCut  ||
463 <             Z2Candidate.M() > 120        ||
464 <             !(lepvec[Z1LeptonPlusIndex].vec.Pt() > 20.0 || lepvec[Z1LeptonMinusIndex].vec.Pt() > 20.0) ||
465 <             !(lepvec[Z1LeptonPlusIndex].vec.Pt() > 10.0 && lepvec[Z1LeptonMinusIndex].vec.Pt() > 10.0)
561 <             ) {
562 <          evtfail |= (1<<EVTFAIL_KINEMATICS );
563 <          return evtfail;
564 <          //h_evtfail->Fill( evtfail, eventweight );
565 <          //      h_evtfail->Fill( evtfail );
566 <          //      cout << "evtfail: " << hex << evtfail << dec << endl;
567 <          //      continue;
568 <        }
569 <        
570 <        int channel;
571 <        if( lepvec[Z1LeptonMinusIndex].type == 11 && lepvec[Z2LeptonMinusIndex].type == 11 ) channel=0;
572 <        if( lepvec[Z1LeptonMinusIndex].type == 13 && lepvec[Z2LeptonMinusIndex].type == 13 ) channel=1;
573 <        if( (lepvec[Z1LeptonMinusIndex].type == 11 && lepvec[Z2LeptonMinusIndex].type == 13) ||
574 <            (lepvec[Z1LeptonMinusIndex].type == 13 && lepvec[Z2LeptonMinusIndex].type == 11)) channel=2;
575 <        
576 <
577 <
578 <        if( passtuple != NULL ) {
579 <          unsigned run   = info->runNum;
580 <          unsigned evt   = info->evtNum;
581 <          unsigned lumi  = info->lumiSec;
582 <          unsigned chan  = channel;
583 <          float mZ1      = Z1Candidate.M() ;
584 <          float mZ2      = Z2Candidate.M() ;
585 <          float m4l      = ZZSystem.M() ;
586 <          float pt4l     = ZZSystem.Pt() ;
587 <          passtuple->SetBranchAddress("run",  &run);
588 <          passtuple->SetBranchAddress("evt",  &evt);
589 <          passtuple->SetBranchAddress("lumi", &lumi);
590 <          passtuple->SetBranchAddress("mZ1",  &mZ1);
591 <          passtuple->SetBranchAddress("mZ2",  &mZ2);
592 <          passtuple->SetBranchAddress("m4l",  &m4l);
593 <          passtuple->SetBranchAddress("pt4l", &pt4l);
594 <          passtuple->SetBranchAddress("w",    &eventweight);
595 <          passtuple->Fill( );
596 <        }
597 <
598 <        if( ctrl.debug ) cout  << "run: " << info->runNum  
599 <              << "\tevt: " << info->evtNum
600 <              << "\tZ1channel: " << lepvec[Z1LeptonMinusIndex].type
601 <              << "\tZ2channel: " << lepvec[Z2LeptonMinusIndex].type
602 <              << "\tmZ1: " << Z1Candidate.M()
603 <              << "\tmZ2: " << Z2Candidate.M()
604 <              << "\tm4l: " << ZZSystem.M()
605 <              << "\tevtfail: " << hex << evtfail << dec
606 <              << "\ttrigbits: " << hex << info->triggerBits << dec
607 <          //          << "\ttree: " << inputFiles[q][f]
608 <              << endl;
443 >        }            
444 >      }
445 >    }
446 >  
447 >    // stop if no Z2 candidate is found
448 >    if (Z2LeptonPlusIndex == -1) {
449 >      evtfail |= ( 1<<EVTFAIL_4L );
450 >      //      ret.status = evtfail;
451 >      ret.status.setStatus(0);
452 >      return ret;
453 >    }
454 >    if( ctrl.debug ) cout << "\tgot a Z2 ..." << endl;
455 >    if( ctrl.debug ) cout << "\tZ2 plusindex: " << Z2LeptonPlusIndex
456 >                          << "\tminusindex: " << Z2LeptonMinusIndex << endl;
457 >    TLorentzVector  Z2LeptonPlus  = *(lepvec[Z2LeptonPlusIndex].vec);
458 >    TLorentzVector  Z2LeptonMinus = *(lepvec[Z2LeptonMinusIndex].vec);
459 >    TLorentzVector  Z2Candidate   = Z2LeptonPlus+Z2LeptonMinus;
460 >    TLorentzVector  ZZSystem      = Z1Candidate + Z2Candidate;
461 >    lepvec[Z1LeptonPlusIndex].is4l = true;
462 >    lepvec[Z1LeptonMinusIndex].is4l = true;
463 >    lepvec[Z2LeptonPlusIndex].is4l = true;
464 >    lepvec[Z2LeptonMinusIndex].is4l = true;
465 >  
466  
467 <        return evtfail;
468 <          
469 < }
467 >
468 >  
469 >    //***************************************************************
470 >    // remaining kinematic cuts
471 >    //***************************************************************
472 >    double Z2massCut=0;
473 >    if      ( ctrl.kinematics == "loose" ) Z2massCut = 12;
474 >    else if ( ctrl.kinematics == "tight" ) Z2massCut = 20;
475 >    else { cout <<  "error! kinematic tightness not defined!" << endl; assert(0); }
476 >  
477 >    if ( Z1Candidate.M() > 120        ||
478 >         Z2Candidate.M() < Z2massCut  ||
479 >         Z2Candidate.M() > 120        ||
480 >         !(lepvec[Z1LeptonPlusIndex].vec->Pt() > 20.0 || lepvec[Z1LeptonMinusIndex].vec->Pt() > 20.0) ||
481 >         !(lepvec[Z1LeptonPlusIndex].vec->Pt() > 10.0 && lepvec[Z1LeptonMinusIndex].vec->Pt() > 10.0)
482 >         ) {
483 >      evtfail |= (1<<EVTFAIL_KINEMATICS );
484 >      //      ret.status = evtfail;
485 >      ret.status.setStatus(0);
486 >      return ret;
487 >    }
488 >  
489 >    unsigned channel;
490 >    if( lepvec[Z1LeptonMinusIndex].type == 11 && lepvec[Z2LeptonMinusIndex].type == 11 ) channel=0;
491 >    if( lepvec[Z1LeptonMinusIndex].type == 13 && lepvec[Z2LeptonMinusIndex].type == 13 ) channel=1;
492 >    if( (lepvec[Z1LeptonMinusIndex].type == 11 && lepvec[Z2LeptonMinusIndex].type == 13) ||
493 >        (lepvec[Z1LeptonMinusIndex].type == 13 && lepvec[Z2LeptonMinusIndex].type == 11)) channel=2;
494 >  
495 >  
496 >  
497 >  
498 >    if( ctrl.debug ) cout  << "run: " << info->runNum  
499 >                           << "\tevt: " << info->evtNum
500 >                           << "\tZ1channel: " << lepvec[Z1LeptonMinusIndex].type
501 >                           << "\tZ2channel: " << lepvec[Z2LeptonMinusIndex].type
502 >                           << "\tmZ1: " << Z1Candidate.M()
503 >                           << "\tmZ2: " << Z2Candidate.M()
504 >                           << "\tm4l: " << ZZSystem.M()
505 >                           << "\tevtfail: " << hex << evtfail << dec
506 >                           << "\ttrigbits: " << hex << info->triggerBits << dec
507 >      //              << "\ttree: " << inputFiles[q][f]
508 >                           << endl;
509 >  
510  
511  
512 +    //***************************************************************
513 +    // finish
514 +    //***************************************************************
515 +
516 +    if( !evtfail ) {
517 +      ret.status.setStatus(SelectionStatus::EVTPASS);
518 +      ret.Z1leptons.push_back(lepvec[Z1LeptonMinusIndex]);
519 +      ret.Z1leptons.push_back(lepvec[Z1LeptonPlusIndex]);
520 +      ret.Z2leptons.push_back(lepvec[Z2LeptonMinusIndex]);
521 +      ret.Z2leptons.push_back(lepvec[Z2LeptonPlusIndex]);
522 +    }
523  
524 +    return ret;
525 + }
526  
527  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines