CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
TrackingRecoMaterialAnalyser Class Reference
Inheritance diagram for TrackingRecoMaterialAnalyser:
DQMEDAnalyzer edm::stream::EDAnalyzer< edm::RunSummaryCache< dqmDetails::NoCache >, edm::LuminosityBlockSummaryCache< dqmDetails::NoCache > > edm::stream::EDAnalyzerBase edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &, const edm::EventSetup &) override
 
virtual void bookHistograms (DQMStore::IBooker &i, edm::Run const &, edm::EventSetup const &) override
 
 TrackingRecoMaterialAnalyser (const edm::ParameterSet &)
 
virtual ~TrackingRecoMaterialAnalyser ()
 
- Public Member Functions inherited from DQMEDAnalyzer
virtual void beginRun (edm::Run const &, edm::EventSetup const &) final
 
virtual void beginStream (edm::StreamID id) final
 
virtual void dqmBeginRun (edm::Run const &, edm::EventSetup const &)
 
 DQMEDAnalyzer (void)
 
virtual void endLuminosityBlockSummary (edm::LuminosityBlock const &, edm::EventSetup const &, dqmDetails::NoCache *) const final
 
virtual void endRunSummary (edm::Run const &, edm::EventSetup const &, dqmDetails::NoCache *) const final
 
uint32_t streamId () const
 
- Public Member Functions inherited from edm::stream::EDAnalyzer< edm::RunSummaryCache< dqmDetails::NoCache >, edm::LuminosityBlockSummaryCache< dqmDetails::NoCache > >
 EDAnalyzer ()=default
 
- Public Member Functions inherited from edm::stream::EDAnalyzerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzerBase ()
 
ModuleDescription const & moduleDescription () const
 
virtual ~EDAnalyzerBase ()
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector
< ProductResolverIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

bool isDoubleSided (SiStripDetId strip_id)
 

Private Attributes

const edm::EDGetTokenT
< reco::BeamSpot
beamspotToken_
 
MonitorElementdeltaP_in_out_vs_eta_
 
MonitorElementdeltaP_in_out_vs_eta_2d_
 
MonitorElementdeltaP_in_out_vs_eta_vs_phi_2d_
 
MonitorElementdeltaP_in_out_vs_z_
 
MonitorElementdeltaP_in_out_vs_z_2d_
 
MonitorElementdeltaPl_in_out_vs_eta_
 
MonitorElementdeltaPl_in_out_vs_z_
 
MonitorElementdeltaPt_in_out_2d_
 
MonitorElementdeltaPt_in_out_vs_eta_
 
MonitorElementdeltaPt_in_out_vs_z_
 
std::string folder_
 
MonitorElementhisto_RZ_
 
MonitorElementhisto_RZ_Ori_
 
std::unordered_map
< std::string, MonitorElement * > 
histosEta_
 
std::unordered_map
< std::string, MonitorElement * > 
histosOriEta_
 
MonitorElementP_vs_eta_2d_
 
TrackTransformer refitter_
 
const edm::EDGetTokenT
< reco::TrackCollection
tracksToken_
 
bool usePV_
 
const edm::EDGetTokenT
< reco::VertexCollection
verticesToken_
 

Additional Inherited Members

- Public Types inherited from edm::stream::EDAnalyzer< edm::RunSummaryCache< dqmDetails::NoCache >, edm::LuminosityBlockSummaryCache< dqmDetails::NoCache > >
typedef CacheContexts< T...> CacheTypes
 
typedef CacheTypes::GlobalCache GlobalCache
 
typedef AbilityChecker< T...> HasAbility
 
typedef
CacheTypes::LuminosityBlockCache 
LuminosityBlockCache
 
typedef
LuminosityBlockContextT
< LuminosityBlockCache,
RunCache, GlobalCache
LuminosityBlockContext
 
typedef
CacheTypes::LuminosityBlockSummaryCache 
LuminosityBlockSummaryCache
 
typedef CacheTypes::RunCache RunCache
 
typedef RunContextT< RunCache,
GlobalCache
RunContext
 
typedef CacheTypes::RunSummaryCache RunSummaryCache
 
- Public Types inherited from edm::stream::EDAnalyzerBase
typedef EDAnalyzerAdaptorBase ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from DQMEDAnalyzer
static std::shared_ptr
< dqmDetails::NoCache
globalBeginLuminosityBlockSummary (edm::LuminosityBlock const &, edm::EventSetup const &, LuminosityBlockContext const *)
 
static std::shared_ptr
< dqmDetails::NoCache
globalBeginRunSummary (edm::Run const &, edm::EventSetup const &, RunContext const *)
 
static void globalEndLuminosityBlockSummary (edm::LuminosityBlock const &, edm::EventSetup const &, LuminosityBlockContext const *, dqmDetails::NoCache *)
 
static void globalEndRunSummary (edm::Run const &, edm::EventSetup const &, RunContext const *, dqmDetails::NoCache *)
 
- Static Public Member Functions inherited from edm::stream::EDAnalyzerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

Definition at line 36 of file TrackingRecoMaterialAnalyser.cc.

Constructor & Destructor Documentation

TrackingRecoMaterialAnalyser::TrackingRecoMaterialAnalyser ( const edm::ParameterSet iPSet)
explicit

Definition at line 72 of file TrackingRecoMaterialAnalyser.cc.

72  :
73  refitter_(iPSet),
74  tracksToken_(consumes<reco::TrackCollection>(iPSet.getParameter<edm::InputTag>("tracks"))),
75  beamspotToken_(consumes<reco::BeamSpot>(iPSet.getParameter<edm::InputTag>("beamspot"))),
76  verticesToken_(mayConsume<reco::VertexCollection>(iPSet.getParameter<edm::InputTag>("vertices"))),
77  usePV_(iPSet.getParameter<bool>("usePV")),
78  folder_(iPSet.getParameter<std::string>("folder")),
79  histo_RZ_(0),
80  histo_RZ_Ori_(0),
91  P_vs_eta_2d_(0)
92 {
93 }
T getParameter(std::string const &) const
const edm::EDGetTokenT< reco::VertexCollection > verticesToken_
const edm::EDGetTokenT< reco::TrackCollection > tracksToken_
const edm::EDGetTokenT< reco::BeamSpot > beamspotToken_
TrackingRecoMaterialAnalyser::~TrackingRecoMaterialAnalyser ( void  )
virtual

Definition at line 96 of file TrackingRecoMaterialAnalyser.cc.

97 {
98 }

Member Function Documentation

void TrackingRecoMaterialAnalyser::analyze ( const edm::Event event,
const edm::EventSetup setup 
)
overridevirtual

Implements edm::stream::EDAnalyzerBase.

Definition at line 191 of file TrackingRecoMaterialAnalyser.cc.

References funct::abs(), assert(), SiPixelRawToDigiRegional_cfi::beamSpot, beamspotToken_, deltaP_in_out_vs_eta_, deltaP_in_out_vs_eta_2d_, deltaP_in_out_vs_eta_vs_phi_2d_, deltaP_in_out_vs_z_, deltaP_in_out_vs_z_2d_, deltaPl_in_out_vs_eta_, deltaPl_in_out_vs_z_, deltaPt_in_out_2d_, deltaPt_in_out_vs_eta_, deltaPt_in_out_vs_z_, MonitorElement::Fill(), edm::EventSetup::get(), histo_RZ_, histo_RZ_Ori_, histosEta_, histosOriEta_, SurfaceOrientation::inner, isDoubleSided(), reco::Vertex::isFake(), MediumProperties::isValid(), LogTrace, Surface::mediumProperties(), TrackingRecHit::missing, reco::Vertex::ndof(), SurfaceOrientation::outer, p2, P_vs_eta_2d_, reco::Vertex::position(), MetAnalyzer::pv(), MediumProperties::radLen(), dt_dqm_sourceclient_common_cff::reco, refitter_, sDETS, TrackTransformer::setServices(), sLAYS, mathSSE::sqrt(), testEve_cfg::tracks, tracksToken_, TrackTransformer::transform(), usePV_, findQualityFiles::v, TrackingRecHit::valid, HLT_FULL_cff::vertices, verticesToken_, MediumProperties::xi(), and reco::Vertex::z().

193 {
194  using namespace edm;
195  using namespace reco;
196  using namespace std;
197 
198  refitter_.setServices(setup);
199 
202  ESHandle<TrackerTopology> trk_topology;
203 
204  // Get the TrackerTopology
205  setup.get<TrackerTopologyRcd>().get(trk_topology);
206 
207  // Get Tracks
208  event.getByToken(tracksToken_, tracks);
209  if (!tracks.isValid() || tracks->size() == 0) {
210  LogInfo("TrackingRecoMaterialAnalyser") << "Invalid or empty track collection" << endl;
211  return;
212  }
213 
214  // Track Selector
215  auto selector = [](const Track &track, const reco::Vertex::Point &pv) -> bool {
216  return (track.quality(track.qualityByName("highPurity"))
217  && track.dxy(pv) < 0.01
218  && track.hitPattern().numberOfLostTrackerHits(HitPattern::MISSING_OUTER_HITS) == 0);
219  };
220 
221  // Get BeamSpot
223  event.getByToken(beamspotToken_, beamSpot);
224  // Bail out if missing
225  if (!beamSpot.isValid())
226  return;
227 
228  reco::Vertex::Point pv(beamSpot->position());
229  if (usePV_) {
230  event.getByToken(verticesToken_, vertices);
231  if (vertices.isValid() && vertices->size() != 0) {
232  // Since we need to use eta and Z information from the tracks, in
233  // order not to have the reco material distribution washed out due
234  // to geometrical effects, we need to confine the PV to some small
235  // region.
236  const Vertex &v = (*vertices)[0];
237  if (!v.isFake() && v.ndof() > 4
238  && std::fabs(v.z()) < 24
239  && std::fabs(v.position().rho()) < 2)
240  pv = v.position();
241  }
242  }
243 
244  // Main idea:
245  // * select first good tracks in input, according to reasonable criteria
246  // * refit the tracks so that we have access to the TrajectoryMeasurements
247  // that internally have all the information about the TSOS on all
248  // crossed layers. We need the refit track and not the original one so
249  // that we are able to correctly compute the path travelled by the track
250  // within the detector, using its updated TSOS. The material description
251  // can in principle be derived also directly from the rechits, via the
252  // det()[(->)GeomDet *]->mediumProperties chain, but that would simply give the
253  // face values, not the "real" ones used while propagating the track.
254  // * Loop on all measurements, extract the information about the TSOS,
255  // its surface and its mediumProperties
256  // * Make plots for the untouched material properties, but also for the
257  // ones corrected by the track direction, since the material properties,
258  // according to the documentation, should refer to normal incidence of
259  // the track, which is seldom the case, according to the current direction
260  TrajectoryStateOnSurface current_tsos;
261  DetId current_det;
262  for (auto const track : *tracks) {
263  if (!selector(track, pv))
264  continue;
265  auto const inner = track.innerMomentum();
266  auto const outer = track.outerMomentum();
267  deltaP_in_out_vs_eta_->Fill(inner.eta(), inner.R() - outer.R());
268  deltaP_in_out_vs_z_->Fill(track.outerZ(), inner.R() - outer.R());
269  deltaP_in_out_vs_eta_2d_->Fill(inner.eta(), inner.R() - outer.R());
270  deltaP_in_out_vs_eta_vs_phi_2d_->Fill(inner.eta(), inner.phi(), inner.R() - outer.R());
271  deltaP_in_out_vs_z_2d_->Fill(track.outerZ(), inner.R() - outer.R());
272  deltaPt_in_out_vs_eta_->Fill(inner.eta(), inner.rho() - outer.rho());
273  deltaPt_in_out_vs_z_->Fill(track.outerZ(), inner.rho() - outer.rho());
274  deltaPl_in_out_vs_eta_->Fill(inner.eta(), inner.z() - outer.z());
275  deltaPl_in_out_vs_z_->Fill(track.outerZ(), inner.z() - outer.z());
276  deltaPt_in_out_2d_->Fill(track.outerZ(), track.outerPosition().rho(), inner.rho() - outer.rho());
277  P_vs_eta_2d_->Fill(track.eta(), track.p());
278  vector<Trajectory> traj = refitter_.transform(track);
279  if (traj.size() > 1 || traj.size() == 0)
280  continue;
281  for (auto const &tm : traj.front().measurements()) {
282  if (tm.recHit().get() &&
283  (tm.recHitR().type() == TrackingRecHit::valid ||
284  tm.recHitR().type() == TrackingRecHit::missing)) {
285  current_tsos = tm.updatedState().isValid() ? tm.updatedState() : tm.forwardPredictedState();
286  auto const & localP = current_tsos.localMomentum();
287  current_det = tm.recHit()->geographicalId();
288  const Surface& surface = current_tsos.surface();
289  assert(tm.recHit()->surface() == &surface);
290  if (!surface.mediumProperties().isValid()) {
291  LogError("TrackingRecoMaterialAnalyser")
292  << "Medium properties for material linked to detector"
293  << " are invalid at: "
294  << current_tsos.globalPosition() << " "
295  << (SiStripDetId)current_det << endl;
296  assert(0);
297  continue;
298  }
299  float p2 = localP.mag2();
300  float xf = std::abs(std::sqrt(p2)/localP.z());
301  float ori_xi = surface.mediumProperties().xi();
302  float ori_radLen = surface.mediumProperties().radLen();
303  float xi = ori_xi*xf;
304  float radLen = ori_radLen*xf;
305 
306  // NOTA BENE: THIS ASSUMES THAT THE TRACKS HAVE BEEN PRODUCED
307  // WITH SPLITTING, AS IS THE CASE FOR THE generalTracks
308  // collection.
309 
310  // Since there are double-sided (glued) modules all over the
311  // tracker, the material budget has been internally
312  // partitioned in two equal components, so that each single
313  // layer will receive half of the correct radLen. For this
314  // reason, only for the double-sided components, we rescale
315  // the obtained radLen by 2.
316 
317  // In particular see code here:
318  // http://cmslxr.fnal.gov/dxr/CMSSW_8_0_5/source/Geometry/TrackerGeometryBuilder/src/TrackerGeomBuilderFromGeometricDet.cc#213
319  // where, in the SiStrip Tracker, if the module has a partner
320  // (i.e. it's a glued detector) the plane is built with a
321  // scaling of 0.5. The actual plane is built few lines below:
322  // http://cmslxr.fnal.gov/dxr/CMSSW_8_0_5/source/Geometry/TrackerGeometryBuilder/src/TrackerGeomBuilderFromGeometricDet.cc#287
323 
324  if (isDoubleSided(current_det)) {
325  LogTrace("TrackingRecoMaterialAnalyser") << "Eta: " << track.eta() << " "
326  << sDETS[current_det.subdetId()]+sLAYS[trk_topology->layer(current_det)]
327  << " has ori_radLen: " << ori_radLen << " and ori_xi: " << xi
328  << " and has radLen: " << radLen << " and xi: " << xi << endl;
329  ori_radLen *= 2.;
330  radLen *= 2.;
331  }
332 
333  histosOriEta_[sDETS[current_det.subdetId()]+sLAYS[trk_topology->layer(current_det)]]->Fill(current_tsos.globalPosition().eta(), ori_radLen);
334  histosEta_[sDETS[current_det.subdetId()]+sLAYS[trk_topology->layer(current_det)]]->Fill(current_tsos.globalPosition().eta(), radLen);
335  histo_RZ_Ori_->Fill(current_tsos.globalPosition().z(), current_tsos.globalPosition().perp(), ori_radLen);
336  histo_RZ_->Fill(current_tsos.globalPosition().z(), current_tsos.globalPosition().perp(), radLen);
337  LogInfo("TrackingRecoMaterialAnalyser") << "Eta: " << track.eta() << " "
338  << sDETS[current_det.subdetId()]+sLAYS[trk_topology->layer(current_det)]
339  << " has ori_radLen: " << ori_radLen << " and ori_xi: " << xi
340  << " and has radLen: " << radLen << " and xi: " << xi << endl;
341  }
342  }
343  }
344 }
const edm::EDGetTokenT< reco::VertexCollection > verticesToken_
float radLen() const
assert(m_qm.get())
std::unordered_map< std::string, MonitorElement * > histosEta_
const edm::EDGetTokenT< reco::TrackCollection > tracksToken_
const edm::EDGetTokenT< reco::BeamSpot > beamspotToken_
const Point & position() const
position
Definition: Vertex.h:109
static const std::vector< std::string > sDETS
void Fill(long long x)
std::unordered_map< std::string, MonitorElement * > histosOriEta_
T sqrt(T t)
Definition: SSEVec.h:18
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double z() const
z coordinate
Definition: Vertex.h:115
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
virtual void setServices(const edm::EventSetup &)
set the services needed by the TrackTransformer
double p2[4]
Definition: TauolaWrapper.h:90
#define LogTrace(id)
double ndof() const
Definition: Vertex.h:105
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:17
Definition: DetId.h:18
bool isFake() const
Definition: Vertex.h:72
tuple tracks
Definition: testEve_cfg.py:39
float xi() const
const T & get() const
Definition: EventSetup.h:56
virtual std::vector< Trajectory > transform(const reco::Track &) const
Convert a reco::Track into Trajectory.
bool isDoubleSided(SiStripDetId strip_id)
static const std::vector< std::string > sLAYS
bool isValid() const
const MediumProperties & mediumProperties() const
Definition: Surface.h:112
void TrackingRecoMaterialAnalyser::bookHistograms ( DQMStore::IBooker i,
edm::Run const &  ,
edm::EventSetup const &  setup 
)
overridevirtual

Implements DQMEDAnalyzer.

Definition at line 101 of file TrackingRecoMaterialAnalyser.cc.

References DQMStore::IBooker::book2D(), DQMStore::IBooker::bookProfile(), DQMStore::IBooker::bookProfile2D(), data, deltaP_in_out_vs_eta_, deltaP_in_out_vs_eta_2d_, deltaP_in_out_vs_eta_vs_phi_2d_, deltaP_in_out_vs_z_, deltaP_in_out_vs_z_2d_, deltaPl_in_out_vs_eta_, deltaPl_in_out_vs_z_, deltaPt_in_out_2d_, deltaPt_in_out_vs_eta_, deltaPt_in_out_vs_z_, folder_, edm::EventSetup::get(), histo_RZ_, histo_RZ_Ori_, histosEta_, histosOriEta_, relval_steps::key, P_vs_eta_2d_, sDETS, DQMStore::IBooker::setCurrentFolder(), and SiPixelPhase1TrackClusters_cfi::title.

103  {
104  using namespace std;
105  edm::ESHandle<TrackerGeometry> trackerGeometry;
106  setup.get<TrackerDigiGeometryRecord>().get(trackerGeometry);
107 
108  ibook.setCurrentFolder(folder_);
109 
110  // Histogram to store the radiation length map, in the R-Z plane,
111  // gathering numbers directly from the trackerRecoMaterial.xml
112  // file. The numbers are not corrected for the track angle.
113  histo_RZ_Ori_ = ibook.bookProfile2D("OriRadLen", "Original_RadLen",
114  600, -300., 300, 120, 0., 120., 0., 1.);
115 
116  // Histogram to store the radiation length map, as before, but
117  // correcting the numbers with the track angle. This represents the
118  // real material seen by the track.
119  histo_RZ_ = ibook.bookProfile2D("RadLen", "RadLen",
120  600, -300., 300, 120, 0., 120., 0., 1.);
121 
122  // Histogram to show the deltaP Out-In in the eta-phi plane.
123  deltaP_in_out_vs_eta_vs_phi_2d_ = ibook.bookProfile2D("DeltaP_in_out_vs_eta_vs_phi_2d",
124  "DeltaP_in_out_vs_eta_vs_phi_2d",
125  100, -3.0, 3.0,
126  100, -3.15, 3.15,
127  0., 100.);
128 
129  // Histogram to show the deltaP Out-In vs eta.
130  deltaP_in_out_vs_eta_2d_ = ibook.book2D("DeltaP_in_out_vs_eta_2d", "DeltaP_in_out_vs_eta_2d",
131  100, -3.0, 3.0, 100, 0., 1);
132 
133  // Histogram to show the deltaP Out-In vs Z. The Z coordinate is the
134  // one computed at the outermost hit.
135  deltaP_in_out_vs_z_2d_ = ibook.book2D("DeltaP_in_out_vs_z_2d", "DeltaP_in_out_vs_z_2d",
136  600, -300, 300, 200., -1, 1.);
137 
138  // Histogram to show the deltaP Out-In vs eta. The eta is the one
139  // computed at the innermost hit.
140  deltaP_in_out_vs_eta_ = ibook.bookProfile("DeltaP_in_out_vs_eta", "DeltaP_in_out_vs_eta",
141  100, -3.0, 3.0, -100., 100.);
142  deltaP_in_out_vs_z_ = ibook.bookProfile("DeltaP_in_out_vs_z", "DeltaP_in_out_vs_z",
143  600, -300, 300, -100., 100.);
144 
145  // Histogram to show the delta_Pt Out-In vs eta. The eta is the one
146  // computed at the innermost hit.
147  deltaPt_in_out_vs_eta_ = ibook.bookProfile("DeltaPt_in_out_vs_eta", "DeltaPt_in_out_vs_eta",
148  100, -3.0, 3.0, -100., 100.);
149 
150  // Histogram to show the delta_Pt Out-In vs Z. The Z is the one
151  // computed at the outermost hit.
152  deltaPt_in_out_vs_z_ = ibook.bookProfile("DeltaPt_in_out_vs_z", "DeltaPt_in_out_vs_z",
153  600, -300, 300, -100., 100);
154 
155  // Histogram to show the delta_Pl Out-In vs eta. The eta is the one
156  // computed at the innermost hit.
157  deltaPl_in_out_vs_eta_ = ibook.bookProfile("DeltaPz_in_out_vs_eta", "DeltaPz_in_out_vs_eta",
158  100, -3.0, 3.0, -100., 100.);
159 
160  // Histogram to show the delta_Pl Out-In vs Z. The Z is the one
161  // computed at the outermost hit.
162  deltaPl_in_out_vs_z_ = ibook.bookProfile("DeltaPz_in_out_vs_z", "DeltaPz_in_out_vs_z",
163  600, -300, 300, -100., 100.);
164 
165  // Histogram to show the delta_Pt Out-In in the Z-R plane. Z and R
166  // are related to the outermost hit.
167  deltaPt_in_out_2d_ = ibook.bookProfile2D("DeltaPt 2D", "DeltaPt 2D",
168  600, -300., 300, 120, 0., 120., -100., 100.);
169 
170  // Histogram to show the distribution of p vs eta for all tracks.
171  P_vs_eta_2d_ = ibook.book2D("P_vs_eta_2d", "P_vs_eta_2d",
172  100, -3.0, 3.0, 100., 0., 5.);
173  char title[50];
174  char key[20];
175  for (unsigned int det = 1; det < sDETS.size(); ++det ) {
176  for (unsigned int sub_det = 1;
177  sub_det <= trackerGeometry->numberOfLayers(det); ++sub_det) {
178  memset(title, 0, sizeof(title));
179  snprintf(title, sizeof(title), "Original_RadLen_vs_Eta_%s%d", sDETS[det].data(), sub_det);
180  snprintf(key, sizeof(key), "%s%d", sDETS[det].data(), sub_det);
181  histosOriEta_.insert(make_pair<string, MonitorElement*>(key,
182  ibook.bookProfile(title, title, 250, -5.0, 5.0, 0., 1.)));
183  snprintf(title, sizeof(title), "RadLen_vs_Eta_%s%d", sDETS[det].data(), sub_det);
184  histosEta_.insert(make_pair<string, MonitorElement*>(key,
185  ibook.bookProfile(title, title, 250, -5.0, 5.0, 0., 1.)));
186  }
187  }
188 }
std::unordered_map< std::string, MonitorElement * > histosEta_
static const std::vector< std::string > sDETS
std::unordered_map< std::string, MonitorElement * > histosOriEta_
string key
FastSim: produces sample of signal events, overlayed with premixed minbias events.
char data[epos_bytes_allocation]
Definition: EPOS_Wrapper.h:82
bool TrackingRecoMaterialAnalyser::isDoubleSided ( SiStripDetId  strip_id)
inlineprivate

Definition at line 45 of file TrackingRecoMaterialAnalyser.cc.

References SiStripDetId::glued(), SiStripDetId::subDetector(), and SiStripDetId::UNKNOWN.

Referenced by analyze().

45  {
46  return ((strip_id.subDetector() != SiStripDetId::UNKNOWN) && strip_id.glued());
47  }
SubDetector subDetector() const
Definition: SiStripDetId.h:114
uint32_t glued() const
Definition: SiStripDetId.h:155

Member Data Documentation

const edm::EDGetTokenT<reco::BeamSpot> TrackingRecoMaterialAnalyser::beamspotToken_
private

Definition at line 50 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze().

MonitorElement* TrackingRecoMaterialAnalyser::deltaP_in_out_vs_eta_
private

Definition at line 59 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaP_in_out_vs_eta_2d_
private

Definition at line 61 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaP_in_out_vs_eta_vs_phi_2d_
private

Definition at line 62 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaP_in_out_vs_z_
private

Definition at line 60 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaP_in_out_vs_z_2d_
private

Definition at line 63 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaPl_in_out_vs_eta_
private

Definition at line 66 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaPl_in_out_vs_z_
private

Definition at line 67 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaPt_in_out_2d_
private

Definition at line 58 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaPt_in_out_vs_eta_
private

Definition at line 64 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::deltaPt_in_out_vs_z_
private

Definition at line 65 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

std::string TrackingRecoMaterialAnalyser::folder_
private

Definition at line 53 of file TrackingRecoMaterialAnalyser.cc.

Referenced by bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::histo_RZ_
private

Definition at line 56 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::histo_RZ_Ori_
private

Definition at line 57 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

std::unordered_map<std::string, MonitorElement *> TrackingRecoMaterialAnalyser::histosEta_
private

Definition at line 55 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

std::unordered_map<std::string, MonitorElement *> TrackingRecoMaterialAnalyser::histosOriEta_
private

Definition at line 54 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* TrackingRecoMaterialAnalyser::P_vs_eta_2d_
private

Definition at line 68 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze(), and bookHistograms().

TrackTransformer TrackingRecoMaterialAnalyser::refitter_
private

Definition at line 48 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze().

const edm::EDGetTokenT<reco::TrackCollection> TrackingRecoMaterialAnalyser::tracksToken_
private

Definition at line 49 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze().

bool TrackingRecoMaterialAnalyser::usePV_
private

Definition at line 52 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze().

const edm::EDGetTokenT<reco::VertexCollection> TrackingRecoMaterialAnalyser::verticesToken_
private

Definition at line 51 of file TrackingRecoMaterialAnalyser.cc.

Referenced by analyze().