Skip to content

Commit

Permalink
Merge pull request #24754 from andreh7/2018-09-em-enrichment-tracked-…
Browse files Browse the repository at this point in the history 
…parameters

made double em enrichment filter parameters tracked
  • Loading branch information
@cmsbuild
cmsbuild committed Oct 11, 2018
2 parents 497653c + 0a5a575 commit f559cdf
Show file tree
Hide file tree
Showing 3 changed files with 165 additions and 65 deletions.
46 changes: 23 additions & 23 deletions Configuration/Generator/python/doubleEMEnrichingHepMCfilter_cfi.py
View file
Expand Up @@ -2,39 +2,39 @@

doubleEMenrichingHepMCfilterParams = cms.PSet(
# seed thresholds
PtSeedThr = cms.untracked.double(5.0),
EtaSeedThr = cms.untracked.double(2.8),
PtSeedThr = cms.double(5.0),
EtaSeedThr = cms.double(2.8),

# photon thresholds
PtGammaThr = cms.untracked.double(0.0),
EtaGammaThr = cms.untracked.double(2.8),
PtGammaThr = cms.double(0.0),
EtaGammaThr = cms.double(2.8),

# electron threshold
PtElThr = cms.untracked.double(2.0),
EtaElThr = cms.untracked.double(2.8),
PtElThr = cms.double(2.0),
EtaElThr = cms.double(2.8),

dRSeedMax = cms.untracked.double(0.0),
dPhiSeedMax = cms.untracked.double(0.2),
dEtaSeedMax = cms.untracked.double(0.12),
dRNarrowCone = cms.untracked.double(0.02),
dRSeedMax = cms.double(0.0),
dPhiSeedMax = cms.double(0.2),
dEtaSeedMax = cms.double(0.12),
dRNarrowCone = cms.double(0.02),

PtTkThr = cms.untracked.double(1.6),
EtaTkThr = cms.untracked.double(2.2),
dRTkMax = cms.untracked.double(0.2),
PtTkThr = cms.double(1.6),
EtaTkThr = cms.double(2.2),
dRTkMax = cms.double(0.2),

PtMinCandidate1 = cms.untracked.double(15.0),
PtMinCandidate2 = cms.untracked.double(15.0),
EtaMaxCandidate = cms.untracked.double(3.0),
PtMinCandidate1 = cms.double(15.0),
PtMinCandidate2 = cms.double(15.0),
EtaMaxCandidate = cms.double(3.0),

NTkConeMax = cms.untracked.int32(2),
NTkConeSum = cms.untracked.int32(4),
NTkConeMax = cms.int32(2),
NTkConeSum = cms.int32(4),

# mass 80..Inf GeV
InvMassMin = cms.untracked.double(80.0),
InvMassMax = cms.untracked.double(14000.0),
InvMassMin = cms.double(80.0),
InvMassMax = cms.double(14000.0),

EnergyCut = cms.untracked.double(1.0),
EnergyCut = cms.double(1.0),

AcceptPrompts = cms.untracked.bool(True),
PromptPtThreshold = cms.untracked.double(15.0),
AcceptPrompts = cms.bool(True),
PromptPtThreshold = cms.double(15.0),
)
65 changes: 61 additions & 4 deletions GeneratorInterface/Core/interface/PythiaHepMCFilterGammaGamma.h
View file
Expand Up @@ -32,14 +32,71 @@ class PythiaHepMCFilterGammaGamma : public BaseHepMCFilter {
private:

const edm::EDGetTokenT<edm::HepMCProduct> token_;
const int maxEvents;

const double ptSeedThr, etaSeedThr, ptGammaThr, etaGammaThr, ptTkThr, etaTkThr;
const double ptElThr, etaElThr, dRTkMax, dRSeedMax, dPhiSeedMax, dEtaSeedMax, dRNarrowCone, pTMinCandidate1, pTMinCandidate2, etaMaxCandidate;
//----------
// filter parameters
//----------

/** minimum pt and maximum absolute eta for electron and photon seeds */
const double ptSeedThr, etaSeedThr;

/** minimum pt and maximum absolute eta for photons to be added
to seeds to form candidates (see also ptElThr, etaElThr) */
const double ptGammaThr, etaGammaThr;

/** minimum pt and maximum absolute eta for charged (stable) particles
to be counted in the isolation cone */
const double ptTkThr, etaTkThr;

/** minimum pt and maximum absolute eta for electrons to be added
to seeds to form candidates (see also ptGammaThr, etaGammaThr) */
const double ptElThr, etaElThr;

/** delta R of cone around candidates in which charged tracks are counted */
const double dRTkMax;

/** delta R of cone around seeds in which other electrons/photons are
added to seeds to form candidates */
const double dRSeedMax;

/** maximum difference in phi and eta for which other electrons/photons
are added to seeds to form candidates.
Note that electrons/photons are accepted if they are within the cone
specified by dRSeedMax or if they are within the rectangular region
specified by (dPhiSeedMax, dEtaSeedMax). */
const double dPhiSeedMax, dEtaSeedMax;

/** this parameter is effectively unused */
const double dRNarrowCone;

/** minimum pt for leading and subleading candidate */
const double pTMinCandidate1, pTMinCandidate2;

/** maximum absolute eta for candidates */
const double etaMaxCandidate;

/** invariant mass range for mass of a pair of candidates */
const double invMassMin, invMassMax;

/** minimum energy for both candidates */
const double energyCut;
const int nTkConeMax, nTkConeSum;

/** maximum number of charged particles in the isolation cone
around each candidate */
const int nTkConeMax;

/** maximum number of charged particles summed over both
cones of a pair of candidates */
const int nTkConeSum;

/** if true, prompt seeds (electrons/photons with no mother
or only ancestors of the same type) will be considered
as having zero charged tracks in the isolation cone */
const bool acceptPrompts;

/** minimum pt for prompt seed particles to be considered (only
effective if acceptPrompts is true) */
const double promptPtThreshold;

};
Expand Down
119 changes: 81 additions & 38 deletions GeneratorInterface/Core/src/PythiaHepMCFilterGammaGamma.cc
View file
Expand Up @@ -14,33 +14,30 @@ using namespace HepMC;


PythiaHepMCFilterGammaGamma::PythiaHepMCFilterGammaGamma(const edm::ParameterSet& iConfig) :
maxEvents(iConfig.getUntrackedParameter<int>("maxEvents", 0)),
ptSeedThr(iConfig.getUntrackedParameter<double>("PtSeedThr")),
etaSeedThr(iConfig.getUntrackedParameter<double>("EtaSeedThr")),
ptGammaThr(iConfig.getUntrackedParameter<double>("PtGammaThr")),
etaGammaThr(iConfig.getUntrackedParameter<double>("EtaGammaThr")),
ptTkThr(iConfig.getUntrackedParameter<double>("PtTkThr")),
etaTkThr(iConfig.getUntrackedParameter<double>("EtaTkThr")),
ptElThr(iConfig.getUntrackedParameter<double>("PtElThr")),
etaElThr(iConfig.getUntrackedParameter<double>("EtaElThr")),
dRTkMax(iConfig.getUntrackedParameter<double>("dRTkMax")),
dRSeedMax(iConfig.getUntrackedParameter<double>("dRSeedMax")),
dPhiSeedMax(iConfig.getUntrackedParameter<double>("dPhiSeedMax")),
dEtaSeedMax(iConfig.getUntrackedParameter<double>("dEtaSeedMax")),
dRNarrowCone(iConfig.getUntrackedParameter<double>("dRNarrowCone")),
pTMinCandidate1(iConfig.getUntrackedParameter<double>("PtMinCandidate1")),
pTMinCandidate2(iConfig.getUntrackedParameter<double>("PtMinCandidate2")),
etaMaxCandidate(iConfig.getUntrackedParameter<double>("EtaMaxCandidate")),
invMassMin(iConfig.getUntrackedParameter<double>("InvMassMin")),
invMassMax(iConfig.getUntrackedParameter<double>("InvMassMax")),
energyCut(iConfig.getUntrackedParameter<double>("EnergyCut")),
nTkConeMax(iConfig.getUntrackedParameter<int>("NTkConeMax")),
nTkConeSum(iConfig.getUntrackedParameter<int>("NTkConeSum")),
acceptPrompts(iConfig.getUntrackedParameter<bool>("AcceptPrompts")),
promptPtThreshold(iConfig.getUntrackedParameter<double>("PromptPtThreshold")) {

if (maxEvents != 0) edm::LogInfo("PythiaFilterGammaGamma::PythiaFilterGammaGamma") << "WARNING, ignoring unsuported option, maxEvents = " << maxEvents << endl;

ptSeedThr(iConfig.getParameter<double>("PtSeedThr")),
etaSeedThr(iConfig.getParameter<double>("EtaSeedThr")),
ptGammaThr(iConfig.getParameter<double>("PtGammaThr")),
etaGammaThr(iConfig.getParameter<double>("EtaGammaThr")),
ptTkThr(iConfig.getParameter<double>("PtTkThr")),
etaTkThr(iConfig.getParameter<double>("EtaTkThr")),
ptElThr(iConfig.getParameter<double>("PtElThr")),
etaElThr(iConfig.getParameter<double>("EtaElThr")),
dRTkMax(iConfig.getParameter<double>("dRTkMax")),
dRSeedMax(iConfig.getParameter<double>("dRSeedMax")),
dPhiSeedMax(iConfig.getParameter<double>("dPhiSeedMax")),
dEtaSeedMax(iConfig.getParameter<double>("dEtaSeedMax")),
dRNarrowCone(iConfig.getParameter<double>("dRNarrowCone")),
pTMinCandidate1(iConfig.getParameter<double>("PtMinCandidate1")),
pTMinCandidate2(iConfig.getParameter<double>("PtMinCandidate2")),
etaMaxCandidate(iConfig.getParameter<double>("EtaMaxCandidate")),
invMassMin(iConfig.getParameter<double>("InvMassMin")),
invMassMax(iConfig.getParameter<double>("InvMassMax")),
energyCut(iConfig.getParameter<double>("EnergyCut")),
nTkConeMax(iConfig.getParameter<int>("NTkConeMax")),
nTkConeSum(iConfig.getParameter<int>("NTkConeSum")),
acceptPrompts(iConfig.getParameter<bool>("AcceptPrompts")),
promptPtThreshold(iConfig.getParameter<double>("PromptPtThreshold")) {

}

PythiaHepMCFilterGammaGamma::~PythiaHepMCFilterGammaGamma()
Expand All @@ -51,10 +48,22 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {

bool accepted = false;

std::vector<const GenParticle*> seeds, egamma, stable;
// electron/photon seeds
std::vector<const GenParticle*> seeds;

// other electrons/photons to be added to seeds
// to form candidates
std::vector<const GenParticle*> egamma;

// charged tracks to be taken into account in the isolation cones
// around candidates
std::vector<const GenParticle*> stable;

std::vector<const GenParticle*>::const_iterator itPart, itStable, itEn;

// Loop on egamma
//----------
// 1. find electron/photon seeds
//----------
for(HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin(); p != myGenEvent->particles_end(); ++p) {

if (
Expand All @@ -74,11 +83,11 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {
if ((*p)->status() == 1) {

// save charged stable tracks
if (abs((*p)->pdg_id()) == 211 ||
abs((*p)->pdg_id()) == 321 ||
abs((*p)->pdg_id()) == 11 ||
abs((*p)->pdg_id()) == 13 ||
abs((*p)->pdg_id()) == 15) {
if (abs((*p)->pdg_id()) == 211 || // charged pion
abs((*p)->pdg_id()) == 321 || // charged kaon
abs((*p)->pdg_id()) == 11 || // electron
abs((*p)->pdg_id()) == 13 || // muon
abs((*p)->pdg_id()) == 15) { // tau
// check if it passes the cut
if ((*p)->momentum().perp() > ptTkThr &&
fabs((*p)->momentum().eta()) < etaTkThr) {
Expand All @@ -102,8 +111,24 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {

if (seeds.size() < 2) return accepted;

//----------
// 2. loop over seeds to build candidates
//
// (adding nearby electrons/photons
// to the seed electrons/photons to obtain the total
// electromagnetic energy)
//----------

// number of tracks around each of the candidates
std::vector<int> nTracks;
std::vector<TLorentzVector> candidate, candidateNarrow, candidateSeed;

// the candidates (four momenta) formed from the
// seed electrons/photons and nearby electrons/photons
std::vector<TLorentzVector> candidate;

// these are filled but then not used afterwards (could be removed)
std::vector<TLorentzVector> candidateNarrow, candidateSeed;

std::vector<const GenParticle*>::iterator itSeed;

const GenParticle* mom;
Expand All @@ -123,7 +148,8 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {
double DEta = (*itEn)->momentum().eta()-(*itSeed)->momentum().eta();
if(DPhi<0) DPhi=-DPhi;
if(DEta<0) DEta=-DEta;


// accept if within cone or within rectangular region around seed
if (DR < dRSeedMax || (DPhi<dPhiSeedMax&&DEta<dEtaSeedMax)) {
energy += temp1;
}
Expand All @@ -132,13 +158,16 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {
}
}

// number of stable charged particles found within dRTkMax
// around candidate
int counter = 0;

if ( energy.Et() != 0. ) {
if (fabs(energy.Eta()) < etaMaxCandidate) {

temp2.SetXYZM(energy.Px(), energy.Py(), energy.Pz(), 0);


// count number of stable particles within cone around candidate
for(itStable = stable.begin(); itStable != stable.end(); ++itStable) {
temp1.SetXYZM((*itStable)->momentum().px(), (*itStable)->momentum().py(), (*itStable)->momentum().pz(), 0);
double DR = temp1.DeltaR(temp2);
Expand All @@ -149,6 +178,8 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {

if ((*itSeed)->momentum().perp()>promptPtThreshold)
{
// check if *itSeed is a prompt particle

bool isPrompt=true;
this_id = (*itSeed)->pdg_id();
mom = (*itSeed);
Expand All @@ -167,7 +198,7 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {

if (mom->status() == 2 && std::abs(first_different_id)>100) isPrompt=false;


// ignore charged particles around prompt particles
if(isPrompt) counter=0;
}
}
Expand All @@ -184,6 +215,12 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {

TLorentzVector minvMin, minvMax;

//----------
// 3. perform further checks on candidates
//
// (energy, charged isolation requirements etc.)
//----------

int i1, i2;
for(unsigned int i=0; i<candidate.size()-1; ++i) {

Expand All @@ -192,12 +229,16 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {
if (fabs(candidate[i].Eta()) > etaMaxCandidate) continue;

for(unsigned int j=i+1; j<candidate.size(); ++j) {

// check features of second candidate alone
if (candidate[j].Energy() < energyCut) continue;
if(nTracks[j]>nTkConeMax) continue;
if (fabs(candidate[j].Eta()) > etaMaxCandidate) continue;

// check requirement on sum of tracks in both isolation cones
if (nTracks[i] + nTracks[j] > nTkConeSum) continue;

// swap candidates to have pt[i1] >= pt[i2]
if (candidate[i].Pt() > candidate[j].Pt()) {
i1 = i;
i2 = j;
Expand All @@ -207,8 +248,10 @@ bool PythiaHepMCFilterGammaGamma::filter(const HepMC::GenEvent* myGenEvent) {
i2 = i;
}

// require minimum pt on leading and subleading candidate
if (candidate[i1].Pt() < pTMinCandidate1 || candidate[i2].Pt() < pTMinCandidate2) continue;

// apply requirements on candidate pair mass
minvMin = candidate[i] + candidate[j];
if (minvMin.M() < invMassMin) continue;

Expand Down

0 comments on commit f559cdf

Please sign in to comment.