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

Diff Legend

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