CMS 3D CMS Logo

FFTJetVertexAdder.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: FFTJetProducers
4 // Class: FFTJetVertexAdder
5 //
13 //
14 // Original Author: Igor Volobouev
15 // Created: Thu Jun 21 19:19:40 CDT 2012
16 //
17 //
18 
19 #include <iostream>
20 #include "CLHEP/Random/RandGauss.h"
21 
22 // framework include files
31 
34 
37 
39 
40 #define init_param(type, varname) varname(ps.getParameter<type>(#varname))
41 
42 //
43 // class declaration
44 //
46 public:
47  explicit FFTJetVertexAdder(const edm::ParameterSet&);
48  FFTJetVertexAdder() = delete;
49  FFTJetVertexAdder(const FFTJetVertexAdder&) = delete;
51  ~FFTJetVertexAdder() override;
52 
53 protected:
54  // methods
55  void beginJob() override;
56  void produce(edm::Event&, const edm::EventSetup&) override;
57  void endJob() override;
58 
59 private:
62 
65 
67 
68  const bool useBeamSpot;
69  const bool addExistingVertices;
70 
71  const double fixedX;
72  const double fixedY;
73  const double fixedZ;
74 
75  const double sigmaX;
76  const double sigmaY;
77  const double sigmaZ;
78 
79  const double nDof;
80  const double chi2;
81  const double errX;
82  const double errY;
83  const double errZ;
84 
85  const unsigned nVerticesToMake;
86 };
87 
88 //
89 // constructors and destructor
90 //
97  init_param(double, fixedX),
98  init_param(double, fixedY),
99  init_param(double, fixedZ),
100  init_param(double, sigmaX),
101  init_param(double, sigmaY),
102  init_param(double, sigmaZ),
103  init_param(double, nDof),
104  init_param(double, chi2),
105  init_param(double, errX),
106  init_param(double, errY),
107  init_param(double, errZ),
108  init_param(unsigned, nVerticesToMake) {
109  if (useBeamSpot)
110  beamSpotToken = consumes<reco::BeamSpot>(beamSpotLabel);
112  existingVerticesToken = consumes<reco::VertexCollection>(existingVerticesLabel);
113  produces<reco::VertexCollection>(outputLabel);
114 }
115 
117 
118 // ------------ method called to produce the data ------------
121  CLHEP::RandGauss rGauss(rng->getEngine(iEvent.streamID()));
122 
123  // get PFCandidates
124  auto pOutput = std::make_unique<reco::VertexCollection>();
125 
126  double xmean = fixedX;
127  double ymean = fixedY;
128  double zmean = fixedZ;
129 
130  double xwidth = sigmaX;
131  double ywidth = sigmaY;
132  double zwidth = sigmaZ;
133 
134  if (useBeamSpot) {
135  edm::Handle<reco::BeamSpot> beamSpotHandle;
136  iEvent.getByToken(beamSpotToken, beamSpotHandle);
137  if (!beamSpotHandle.isValid())
138  throw cms::Exception("FFTJetBadConfig") << "ERROR in FFTJetVertexAdder:"
139  " could not find beam spot information"
140  << std::endl;
141 
142  xmean = beamSpotHandle->x0();
143  ymean = beamSpotHandle->y0();
144  zmean = beamSpotHandle->z0();
145 
146  xwidth = beamSpotHandle->BeamWidthX();
147  ywidth = beamSpotHandle->BeamWidthY();
148  zwidth = beamSpotHandle->sigmaZ();
149  }
150 
152  for (unsigned i = 0; i < 3; ++i)
153  for (unsigned j = 0; j < 3; ++j)
154  err[i][j] = 0.0;
155  err[0][0] = errX * errX;
156  err[1][1] = errY * errY;
157  err[2][2] = errZ * errZ;
158 
159  for (unsigned iv = 0; iv < nVerticesToMake; ++iv) {
160  const double x0 = rGauss(xmean, xwidth);
161  const double y0 = rGauss(ymean, ywidth);
162  const double z0 = rGauss(zmean, zwidth);
163  const reco::Vertex::Point position(x0, y0, z0);
164  pOutput->push_back(reco::Vertex(position, err, chi2, nDof, 0));
165  }
166 
167  if (addExistingVertices) {
168  typedef reco::VertexCollection::const_iterator IV;
169 
172  if (!vertices.isValid())
173  throw cms::Exception("FFTJetBadConfig") << "ERROR in FFTJetVertexAdder:"
174  " could not find existing collection of vertices"
175  << std::endl;
176 
177  const IV vertend(vertices->end());
178  for (IV iv = vertices->begin(); iv != vertend; ++iv)
179  pOutput->push_back(*iv);
180  }
181 
182  iEvent.put(std::move(pOutput), outputLabel);
183 }
184 
185 // ------------ method called once each job just before starting event loop
188  if (!rng.isAvailable()) {
189  throw cms::Exception("FFTJetBadConfig") << "ERROR in FFTJetVertexAdder:"
190  " failed to initialize the random number generator"
191  << std::endl;
192  }
193 }
194 
195 // ------------ method called once each job just after ending the event loop
197 
198 //define this as a plug-in
fftjetvertexadder_cfi.errZ
errZ
Definition: fftjetvertexadder_cfi.py:39
edm::RandomNumberGenerator::getEngine
virtual CLHEP::HepRandomEngine & getEngine(StreamID const &)=0
Use this engine in event methods.
Handle.h
electrons_cff.bool
bool
Definition: electrons_cff.py:366
mps_fire.i
i
Definition: mps_fire.py:428
BeamSpotPI::sigmaX
Definition: BeamSpotPayloadInspectorHelper.h:34
fftjetvertexadder_cfi.nDof
nDof
Definition: fftjetvertexadder_cfi.py:35
BeamSpotPI::sigmaZ
Definition: BeamSpotPayloadInspectorHelper.h:36
EDProducer.h
FFTJetVertexAdder::sigmaX
const double sigmaX
Definition: FFTJetVertexAdder.cc:75
FFTJetVertexAdder::outputLabel
const std::string outputLabel
Definition: FFTJetVertexAdder.cc:66
fftjetvertexadder_cfi.errY
errY
Definition: fftjetvertexadder_cfi.py:38
reco::BeamSpot::z0
double z0() const
z coordinate
Definition: BeamSpot.h:65
edm::EDGetTokenT< reco::BeamSpot >
edm
HLT enums.
Definition: AlignableModifier.h:19
RandomNumberGenerator.h
reco::Vertex::Error
math::Error< dimension >::type Error
covariance error matrix (3x3)
Definition: Vertex.h:44
gpuVertexFinder::iv
int32_t *__restrict__ iv
Definition: gpuClusterTracksDBSCAN.h:42
FFTJetVertexAdder::errX
const double errX
Definition: FFTJetVertexAdder.cc:81
HLT_FULL_cff.InputTag
InputTag
Definition: HLT_FULL_cff.py:89301
FFTJetVertexAdder::sigmaZ
const double sigmaZ
Definition: FFTJetVertexAdder.cc:77
recoTrackAccumulator_cfi.outputLabel
outputLabel
Definition: recoTrackAccumulator_cfi.py:13
fftjetvertexadder_cfi.existingVerticesLabel
existingVerticesLabel
Definition: fftjetvertexadder_cfi.py:10
FFTJetVertexAdder::fixedY
const double fixedY
Definition: FFTJetVertexAdder.cc:72
hltPixelTracks_cff.chi2
chi2
Definition: hltPixelTracks_cff.py:25
FFTJetVertexAdder
Definition: FFTJetVertexAdder.cc:45
edm::Handle< reco::BeamSpot >
edm::Service::isAvailable
bool isAvailable() const
Definition: Service.h:40
fftjetvertexadder_cfi.addExistingVertices
addExistingVertices
Definition: fftjetvertexadder_cfi.py:21
AlignmentTracksFromVertexSelector_cfi.vertices
vertices
Definition: AlignmentTracksFromVertexSelector_cfi.py:5
FFTJetVertexAdder::nVerticesToMake
const unsigned nVerticesToMake
Definition: FFTJetVertexAdder.cc:85
FFTJetVertexAdder::fixedX
const double fixedX
Definition: FFTJetVertexAdder.cc:71
reco::BeamSpot::sigmaZ
double sigmaZ() const
sigma z
Definition: BeamSpot.h:76
MakerMacros.h
fftjetvertexadder_cfi.nVerticesToMake
nVerticesToMake
Definition: fftjetvertexadder_cfi.py:42
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
BeamSpot.h
fftjetvertexadder_cfi.fixedY
fixedY
Definition: fftjetvertexadder_cfi.py:26
Service.h
fftjetvertexadder_cfi.fixedZ
fixedZ
Definition: fftjetvertexadder_cfi.py:27
PrimaryVertexMonitor_cff.beamSpotLabel
beamSpotLabel
Definition: PrimaryVertexMonitor_cff.py:9
HLTMuonOfflineAnalyzer_cfi.z0
z0
Definition: HLTMuonOfflineAnalyzer_cfi.py:98
FFTJetVertexAdder::existingVerticesToken
edm::EDGetTokenT< reco::VertexCollection > existingVerticesToken
Definition: FFTJetVertexAdder.cc:64
fftjetvertexadder_cfi.fixedX
fixedX
Definition: fftjetvertexadder_cfi.py:25
FFTJetVertexAdder::nDof
const double nDof
Definition: FFTJetVertexAdder.cc:79
BeamSpotPI::sigmaY
Definition: BeamSpotPayloadInspectorHelper.h:35
FFTJetVertexAdder::beamSpotLabel
const edm::InputTag beamSpotLabel
Definition: FFTJetVertexAdder.cc:60
FFTJetVertexAdder::errZ
const double errZ
Definition: FFTJetVertexAdder.cc:83
edm::ParameterSet
Definition: ParameterSet.h:47
FFTJetVertexAdder::operator=
FFTJetVertexAdder & operator=(const FFTJetVertexAdder &)=delete
Event.h
FFTJetVertexAdder::produce
void produce(edm::Event &, const edm::EventSetup &) override
Definition: FFTJetVertexAdder.cc:119
position
static int position[264][3]
Definition: ReadPGInfo.cc:289
FFTJetVertexAdder::existingVerticesLabel
const edm::InputTag existingVerticesLabel
Definition: FFTJetVertexAdder.cc:61
edm::Service< edm::RandomNumberGenerator >
iEvent
int iEvent
Definition: GenABIO.cc:224
init_param
#define init_param(type, varname)
Definition: FFTJetVertexAdder.cc:40
FFTJetVertexAdder::beginJob
void beginJob() override
Definition: FFTJetVertexAdder.cc:186
edm::EventSetup
Definition: EventSetup.h:58
FFTJetVertexAdder::~FFTJetVertexAdder
~FFTJetVertexAdder() override
Definition: FFTJetVertexAdder.cc:116
reco::BeamSpot::BeamWidthX
double BeamWidthX() const
beam width X
Definition: BeamSpot.h:82
submitPVResolutionJobs.err
err
Definition: submitPVResolutionJobs.py:85
FFTJetVertexAdder::fixedZ
const double fixedZ
Definition: FFTJetVertexAdder.cc:73
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
FFTJetVertexAdder::sigmaY
const double sigmaY
Definition: FFTJetVertexAdder.cc:76
VertexFwd.h
reco::BeamSpot::x0
double x0() const
x coordinate
Definition: BeamSpot.h:61
reco::Vertex::Point
math::XYZPoint Point
point in the space
Definition: Vertex.h:40
eostools.move
def move(src, dest)
Definition: eostools.py:511
std
Definition: JetResolutionObject.h:76
HLT_FULL_cff.useBeamSpot
useBeamSpot
Definition: HLT_FULL_cff.py:26535
Vertex.h
Frameworkfwd.h
FFTJetVertexAdder::chi2
const double chi2
Definition: FFTJetVertexAdder.cc:80
Exception
Definition: hltDiff.cc:245
FFTJetVertexAdder::useBeamSpot
const bool useBeamSpot
Definition: FFTJetVertexAdder.cc:68
FFTJetVertexAdder::endJob
void endJob() override
Definition: FFTJetVertexAdder.cc:196
edm::EDProducer
Definition: EDProducer.h:35
Exception.h
cms::Exception
Definition: Exception.h:70
FFTJetVertexAdder::errY
const double errY
Definition: FFTJetVertexAdder.cc:82
View.h
ParameterSet.h
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:70
reco::BeamSpot::y0
double y0() const
y coordinate
Definition: BeamSpot.h:63
edm::Event
Definition: Event.h:73
FFTJetVertexAdder::FFTJetVertexAdder
FFTJetVertexAdder()=delete
fftjetvertexadder_cfi.errX
errX
Definition: fftjetvertexadder_cfi.py:37
edm::InputTag
Definition: InputTag.h:15
reco::Vertex
Definition: Vertex.h:35
FFTJetVertexAdder::addExistingVertices
const bool addExistingVertices
Definition: FFTJetVertexAdder.cc:69
FFTJetVertexAdder::beamSpotToken
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken
Definition: FFTJetVertexAdder.cc:63
reco::BeamSpot::BeamWidthY
double BeamWidthY() const
beam width Y
Definition: BeamSpot.h:84