ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/MitHzz4l/Selection/src/Selection.cc

Comparing UserCode/MitHzz4l/Selection/src/Selection.cc (file contents):
Revision 1.1.1.1 by khahn, Thu Sep 8 13:33:18 2011 UTC vs.
Revision 1.25 by dkralph, Tue Oct 23 11:39:52 2012 UTC

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

Diff Legend

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