70 |
|
return vec; |
71 |
|
} |
72 |
|
|
73 |
+ |
float deltaPhi( float phi1, float phi2) { |
74 |
+ |
float result = phi1 - phi2; |
75 |
+ |
while (result > M_PI) result -= 2*M_PI; |
76 |
+ |
while (result <= -M_PI) result += 2*M_PI; |
77 |
+ |
return result; |
78 |
+ |
} |
79 |
+ |
|
80 |
|
// useful functions |
81 |
|
float deltaR( const TLorentzVector& v1, const TLorentzVector& v2 ) { |
82 |
|
// deltaR = sqrt ( deltaEta^2 + deltaPhi^2 ) |
83 |
< |
return sqrt(pow(v1.Eta() - v2.Eta(), 2) + pow(v1.Phi() - v2.Phi(), 2) ); |
83 |
> |
return sqrt(pow(v1.Eta() - v2.Eta(), 2) + pow(deltaPhi(v1.Phi(),v2.Phi()), 2) ); |
84 |
|
} |
85 |
|
|
86 |
|
float effectiveAreaElectron( float eta ) { |
260 |
|
tree->Branch("ht", &ht, "ht/F"); |
261 |
|
tree->Branch("nVertex", &nVertex, "nVertex/I"); |
262 |
|
tree->Branch("pu_weight", &pu_weight, "pu_weight/F"); |
263 |
< |
|
263 |
> |
tree->Branch("genElectron", &genElectron); |
264 |
> |
tree->Branch("genPhoton", &genPhoton); |
265 |
|
|
266 |
|
for (unsigned long jentry=0; jentry < processNEvents; ++jentry) { |
267 |
|
if ( loggingVerbosity>1 || jentry%reportEvery==0 ) |
273 |
|
jet.clear(); |
274 |
|
electron.clear(); |
275 |
|
muon.clear(); |
276 |
+ |
genElectron.clear(); |
277 |
+ |
genPhoton.clear(); |
278 |
|
ht = 0; |
279 |
|
|
280 |
|
// weights |
298 |
|
|
299 |
|
for(std::vector<susy::Photon>::iterator it = photonVector.begin(); |
300 |
|
it != photonVector.end(); ++it ) { |
301 |
< |
if( !(it->isEB() || it->isEE()) && skim ) |
301 |
> |
if( !(it->isEE() || it->isEB()) && it->momentum.Pt()<20 && it->isEBEtaGap() && it->isEBPhiGap() && it->isEERingGap() && it->isEEDeeGap() && it->isEBEEGap() && skim ) |
302 |
|
continue; |
303 |
|
tree::Photon thisphoton; |
294 |
– |
thisphoton.pt = getPtFromMatchedJet( *it, jetVector, loggingVerbosity ); |
304 |
|
|
305 |
|
thisphoton.chargedIso = chargedHadronIso_corrected(*it, event->rho25); |
306 |
|
thisphoton.neutralIso = neutralHadronIso_corrected(*it, event->rho25); |
307 |
|
thisphoton.photonIso = photonIso_corrected(*it, event->rho25); |
308 |
|
|
309 |
< |
bool loose_photon_barrel = thisphoton.pt>20 |
301 |
< |
&& it->isEB() |
302 |
< |
&& it->passelectronveto |
309 |
> |
bool loose_photon_barrel = it->isEB() |
310 |
|
&& it->hadTowOverEm<0.05 |
311 |
|
&& it->sigmaIetaIeta<0.012 |
312 |
|
&& thisphoton.chargedIso<2.6 |
313 |
|
&& thisphoton.neutralIso<3.5+0.04*thisphoton.pt |
314 |
|
&& thisphoton.photonIso<1.3+0.005*thisphoton.pt; |
315 |
< |
bool loose_photon_endcap = thisphoton.pt > 20 |
316 |
< |
&& it->isEE() |
310 |
< |
&& it->passelectronveto |
315 |
> |
|
316 |
> |
bool loose_photon_endcap = it->isEE() |
317 |
|
&& it->hadTowOverEm<0.05 |
318 |
|
&& it->sigmaIetaIeta<0.034 |
319 |
|
&& thisphoton.chargedIso<2.3 |
320 |
|
&& thisphoton.neutralIso<2.9+0.04*thisphoton.pt; |
321 |
< |
if(!(loose_photon_endcap || loose_photon_barrel || thisphoton.pt > 75 ) && skim ) |
321 |
> |
|
322 |
> |
thisphoton.ptJet = getPtFromMatchedJet( *it, jetVector, loggingVerbosity ); |
323 |
> |
if(!(loose_photon_endcap || loose_photon_barrel || thisphoton.ptJet > 75 ) && skim ) |
324 |
|
continue; |
325 |
+ |
thisphoton.pt = it->momentum.Pt(); |
326 |
|
thisphoton.eta = it->momentum.Eta(); |
327 |
|
thisphoton.phi = it->momentum.Phi(); |
328 |
|
thisphoton.r9 = it->r9; |
356 |
|
scale = it->jecScaleFactors.find("L2L3")->second; |
357 |
|
TLorentzVector corrP4 = scale * it->momentum; |
358 |
|
|
359 |
< |
if(std::abs(corrP4.Eta()) > 3.0 && skim ) continue; |
359 |
> |
// Calculate HT. |
360 |
> |
// The definiton differs from the saved jet, since trigger is described better |
361 |
> |
if( std::abs( corrP4.Eta() ) < 3 && corrP4.Pt() > 40 ) |
362 |
> |
ht += thisjet.pt; |
363 |
> |
|
364 |
> |
if(std::abs(corrP4.Eta()) > 2.6 && skim ) continue; |
365 |
|
if(corrP4.Pt() < 30 && skim ) continue; |
366 |
|
thisjet.pt = corrP4.Pt(); |
367 |
|
thisjet.eta = corrP4.Eta(); |
383 |
|
std::cout << " p_T, jet = " << thisjet.pt << std::endl; |
384 |
|
|
385 |
|
jet.push_back( thisjet ); |
372 |
– |
ht += thisjet.pt; |
386 |
|
}// for jet |
387 |
|
if( jet.size() < 2 && skim ) |
388 |
|
continue; |
390 |
|
if( loggingVerbosity > 1 ) |
391 |
|
std::cout << "Found " << jet.size() << " jets" << std::endl; |
392 |
|
|
393 |
+ |
if( ht < 450 && skim) |
394 |
+ |
continue; |
395 |
+ |
|
396 |
+ |
|
397 |
|
|
398 |
|
// met |
399 |
|
std::map<TString, susy::MET>::iterator met_it = event->metMap.find("pfMet"); |
420 |
|
// use veto electrons |
421 |
|
if( it->momentum.Pt() < 20 || it->momentum.Pt() > 1e6 ) |
422 |
|
continue; // spike rejection |
423 |
< |
float iso = ( it->chargedHadronIso + max(it->neutralHadronIso+it->photonIso |
407 |
< |
- effectiveAreaElectron(it->momentum.Eta())*event->rho25, (Float_t)0. ) |
423 |
> |
float iso = ( it->chargedHadronIso + max(it->neutralHadronIso+it->photonIso - effectiveAreaElectron(it->momentum.Eta())*event->rho25, (float)0. ) |
424 |
|
) / it->momentum.Pt(); |
425 |
|
float d0 = d0correction( *it, *event ); |
426 |
|
float dZ = std::abs( dZcorrection( *it, *event ) ); |
475 |
|
// vertices |
476 |
|
nVertex = event->vertices.size(); |
477 |
|
|
478 |
< |
if( ht < 450 && skim) |
479 |
< |
continue; |
480 |
< |
|
478 |
> |
tree::Particle thisGenParticle; |
479 |
> |
for( std::vector<susy::Particle>::iterator it = event->genParticles.begin(); it != event->genParticles.end(); ++it ) { |
480 |
> |
if( it->momentum.Pt() < 20 ) continue; |
481 |
> |
thisGenParticle.pt = it->momentum.Pt(); |
482 |
> |
thisGenParticle.eta = it->momentum.Eta(); |
483 |
> |
thisGenParticle.phi = it->momentum.Phi(); |
484 |
> |
switch( std::abs(it->pdgId) ) { |
485 |
> |
case 22: // photon |
486 |
> |
genPhoton.push_back( thisGenParticle ); |
487 |
> |
break; |
488 |
> |
case 11: // electron |
489 |
> |
// Demand a W boson as mother particle of electron |
490 |
> |
if( abs(event->genParticles[it->motherIndex].pdgId) == 24 ) |
491 |
> |
genElectron.push_back( thisGenParticle ); |
492 |
> |
break; |
493 |
> |
} |
494 |
> |
} |
495 |
|
|
496 |
|
tree->Fill(); |
497 |
|
} // for jentry |