1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
|
#include "FastSimulation/SimplifiedGeometryPropagator/interface/ParticleFilter.h"
#include "FastSimulation/SimplifiedGeometryPropagator/interface/Particle.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "vdt/vdtMath.h"
fastsim::ParticleFilter::ParticleFilter(const edm::ParameterSet& cfg) {
// Charged particles must have pt greater than chargedPtMin [GeV]
double chargedPtMin = cfg.getParameter<double>("chargedPtMin");
chargedPtMin2_ = chargedPtMin * chargedPtMin;
// (All) particles must have energy greater than EMin [GeV]
EMin_ = cfg.getParameter<double>("EMin");
// Allow *ALL* protons with energy > protonEMin
protonEMin_ = cfg.getParameter<double>("protonEMin");
// List of invisible particles if extension needed
// Predefined: Neutrinos, Neutralino_1
skipParticles_ = cfg.getParameter<std::vector<int>>("invisibleParticles");
// Particles must have abs(eta) < etaMax (if close enough to 0,0,0)
double etaMax = cfg.getParameter<double>("etaMax");
cos2ThetaMax_ = (vdt::fast_exp(2. * etaMax) - 1.) / (vdt::fast_exp(2. * etaMax) + 1.);
cos2ThetaMax_ *= cos2ThetaMax_;
// Particles must have vertex inside the tracker
vertexRMax2_ = 129.0 * 129.0;
vertexZMax_ = 303.353;
}
bool fastsim::ParticleFilter::accepts(const fastsim::Particle& particle) const {
int absPdgId = abs(particle.pdgId());
// skip invisible particles
if (absPdgId == 12 || absPdgId == 14 || absPdgId == 16 || absPdgId == 1000022) {
return false;
}
// keep all high-energy protons
else if (absPdgId == 2212 && particle.momentum().E() >= protonEMin_) {
return true;
}
// cut on eta if the origin vertex is close to the beam
else if (particle.position().Perp2() < 25. &&
particle.momentum().Pz() * particle.momentum().Pz() / particle.momentum().P2() > cos2ThetaMax_) {
return false;
}
// possible to extend list of invisible particles
for (unsigned InvIdx = 0; InvIdx < skipParticles_.size(); InvIdx++) {
if (absPdgId == abs(skipParticles_.at(InvIdx))) {
return false;
}
}
// particles must have vertex in volume of tracker
return acceptsVtx(particle.position()) && acceptsEn(particle);
//return (acceptsVtx(particle.position()) || particle.momentum().Pz()*particle.momentum().Pz()/particle.momentum().P2() > (vdt::fast_exp(2.*3.0)-1.) / (vdt::fast_exp(2.*3.0)+1.)*(vdt::fast_exp(2.*3.0)-1.) / (vdt::fast_exp(2.*3.0)+1.)) && acceptsEn(particle);
}
bool fastsim::ParticleFilter::acceptsEn(const fastsim::Particle& particle) const {
int absPdgId = abs(particle.pdgId());
// keep all high-energy protons
if (absPdgId == 2212 && particle.momentum().E() >= protonEMin_) {
return true;
}
// cut on the energy
else if (particle.momentum().E() < EMin_) {
return false;
}
// cut on pt of charged particles
else if (particle.charge() != 0 && particle.momentum().Perp2() < chargedPtMin2_) {
return false;
}
return true;
}
bool fastsim::ParticleFilter::acceptsVtx(const math::XYZTLorentzVector& originVertex) const {
// origin vertex is within the tracker volume
return (originVertex.Perp2() < vertexRMax2_ && fabs(originVertex.Z()) < vertexZMax_);
}
|