CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TrackerDpgAnalysis.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: DPGAnalysis
4 // Class: TrackerDpgAnalysis
5 //
13 //
14 // Original Author: Christophe DELAERE
15 // Created: Tue Sep 23 02:11:44 CEST 2008
16 // Revised: Thu Nov 26 10:00:00 CEST 2009
17 // part of the code was inspired by http://cmssw.cvs.cern.ch/cgi-bin/cmssw.cgi/UserCode/YGao/LhcTrackAnalyzer/
18 // part of the code was inspired by
19 // other inputs from Andrea Giammanco, Gaelle Boudoul, Andrea Venturi, Steven Lowette, Gavril Giurgiu
20 // $Id: TrackerDpgAnalysis.cc,v 1.14 2013/02/27 19:49:47 wmtan Exp $
21 //
22 //
23 
24 // system include files
25 #include <memory>
26 #include <iostream>
27 #include <limits>
28 #include <utility>
29 #include <vector>
30 #include <algorithm>
31 #include <functional>
32 #include <cstring>
33 #include <sstream>
34 #include <fstream>
35 
36 // root include files
37 #include "TTree.h"
38 #include "TFile.h"
39 
40 // user include files
100 
101 // topology
103 
104 //
105 // class decleration
106 //
108 
110 public:
111  explicit TrackerDpgAnalysis(const edm::ParameterSet&);
112  ~TrackerDpgAnalysis() override;
113 
114 protected:
116  const std::vector<Trajectory>&);
117  void insertMeasurement(std::multimap<const uint32_t, std::pair<LocalPoint, double> >&,
119  double);
120  std::vector<int> onTrack(edm::Handle<edmNew::DetSetVector<SiStripCluster> >&, const reco::TrackCollection&, uint32_t);
121  void insertMeasurement(std::multimap<const uint32_t, std::pair<int, int> >&, const TrackingRecHit*, int);
122  std::map<size_t, int> inVertex(const reco::TrackCollection&, const reco::VertexCollection&, uint32_t);
123  std::vector<std::pair<double, double> > onTrackAngles(edm::Handle<edmNew::DetSetVector<SiPixelCluster> >&,
124  const std::vector<Trajectory>&);
125  void insertMeasurement(std::multimap<const uint32_t, std::pair<LocalPoint, std::pair<double, double> > >&,
127  double,
128  double);
129  std::vector<int> onTrack(edm::Handle<edmNew::DetSetVector<SiPixelCluster> >&, const reco::TrackCollection&, uint32_t);
130  void insertMeasurement(std::multimap<const uint32_t, std::pair<std::pair<float, float>, int> >&,
131  const TrackingRecHit*,
132  int);
133  std::string toStringName(uint32_t, const TrackerTopology*);
134  std::string toStringId(uint32_t);
135  double sumPtSquared(const reco::Vertex&);
136  float delay(const SiStripEventSummary&);
137  std::map<uint32_t, float> delay(const std::vector<std::string>&);
138 
139 private:
140  void beginRun(const edm::Run&, const edm::EventSetup&) override;
141  void analyze(const edm::Event&, const edm::EventSetup&) override;
142  void endJob() override;
143 
144  // ----------member data ---------------------------
145  static const int nMaxPVs_ = 50;
158  std::vector<edm::InputTag> trackLabels_;
159  std::vector<edm::EDGetTokenT<reco::TrackCollection> > trackTokens_;
160  std::vector<edm::EDGetTokenT<std::vector<Trajectory> > > trajectoryTokens_;
161  std::vector<edm::EDGetTokenT<TrajTrackAssociationCollection> > trajTrackAssoTokens_;
164  std::multimap<const uint32_t, const FedChannelConnection*> connections_;
168  TTree* clusters_;
169  TTree* pixclusters_;
170  std::vector<TTree*> tracks_;
171  std::vector<TTree*> missingHits_;
172  TTree* vertices_;
174  TTree* event_;
175  TTree* psumap_;
176  TTree* readoutmap_;
177  bool onTrack_;
178  uint32_t vertexid_;
180  uint32_t runid_;
181  uint32_t globalvertexid_;
190  float eta_, phi_, chi2_;
192  uint32_t detid_, dcuId_, type_;
195  lostHits_, ndof_;
196  uint32_t *ntracks_, *ntrajs_;
199  uint32_t nTracks_pvtx_;
207  float charge_, p_, pt_;
215  std::vector<std::string> delayFileNames_;
217  std::vector<std::string> hlNames_; // name of each HLT algorithm
218  HLTConfigProvider hltConfig_; // to get configuration for L1s/Pre
219 };
220 
221 //
222 // constructors and destructor
223 //
225  // members
226  moduleName_ = new char[256];
227  moduleId_ = new char[256];
228  PSUname_ = new char[256];
229  pset_ = iConfig;
230 
231  // enable/disable functionalities
232  functionality_offtrackClusters_ = iConfig.getUntrackedParameter<bool>("keepOfftrackClusters", true);
233  functionality_ontrackClusters_ = iConfig.getUntrackedParameter<bool>("keepOntrackClusters", true);
234  functionality_pixclusters_ = iConfig.getUntrackedParameter<bool>("keepPixelClusters", true);
235  functionality_pixvertices_ = iConfig.getUntrackedParameter<bool>("keepPixelVertices", true);
236  functionality_missingHits_ = iConfig.getUntrackedParameter<bool>("keepMissingHits", true);
237  functionality_tracks_ = iConfig.getUntrackedParameter<bool>("keepTracks", true);
238  functionality_vertices_ = iConfig.getUntrackedParameter<bool>("keepVertices", true);
239  functionality_events_ = iConfig.getUntrackedParameter<bool>("keepEvents", true);
240 
241  // parameters
242  summaryToken_ = consumes<SiStripEventSummary>(edm::InputTag("siStripDigis"));
243  clusterToken_ = consumes<edmNew::DetSetVector<SiStripCluster> >(iConfig.getParameter<edm::InputTag>("ClustersLabel"));
245  consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<edm::InputTag>("PixelClustersLabel"));
246  trackLabels_ = iConfig.getParameter<std::vector<edm::InputTag> >("TracksLabel");
248  trackLabels_, [this](edm::InputTag const& tag) { return consumes<reco::TrackCollection>(tag); });
250  trackLabels_, [this](edm::InputTag const& tag) { return consumes<std::vector<Trajectory> >(tag); });
252  trackLabels_, [this](edm::InputTag const& tag) { return consumes<TrajTrackAssociationCollection>(tag); });
253  dedx1Token_ = consumes<edm::ValueMap<reco::DeDxData> >(iConfig.getParameter<edm::InputTag>("DeDx1Label"));
254  dedx2Token_ = consumes<edm::ValueMap<reco::DeDxData> >(iConfig.getParameter<edm::InputTag>("DeDx2Label"));
255  dedx3Token_ = consumes<edm::ValueMap<reco::DeDxData> >(iConfig.getParameter<edm::InputTag>("DeDx3Label"));
256  vertexToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertexLabel"));
257  pixelVertexToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("pixelVertexLabel"));
258  bsToken_ = consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpotLabel"));
259  L1Token_ = consumes<L1GlobalTriggerReadoutRecord>(iConfig.getParameter<edm::InputTag>("L1Label"));
260  HLTTag_ = iConfig.getParameter<edm::InputTag>("HLTLabel");
261  HLTToken_ = consumes<edm::TriggerResults>(HLTTag_);
262 
263  // initialize arrays
264  size_t trackSize(trackLabels_.size());
265  ntracks_ = new uint32_t[trackSize];
266  ntrajs_ = new uint32_t[trackSize];
267  globaltrackid_ = new uint32_t[trackSize];
268  trackid_ = new uint32_t[trackSize];
269  lowPixelProbabilityFraction_ = new float[trackSize];
270  globalvertexid_ = iConfig.getParameter<uint32_t>("InitalCounter");
271  for (size_t i = 0; i < trackSize; ++i) {
272  ntracks_[i] = 0;
273  ntrajs_[i] = 0;
274  globaltrackid_[i] = iConfig.getParameter<uint32_t>("InitalCounter");
275  trackid_[i] = 0;
277  }
278 
279  // create output
281  TFileDirectory* dir = new TFileDirectory(fileService->mkdir("trackerDPG"));
282 
283  // create a TTree for clusters
284  clusters_ = dir->make<TTree>("clusters", "cluster information");
285  clusters_->Branch("eventid", &eventid_, "eventid/i");
286  clusters_->Branch("runid", &runid_, "runid/i");
287  for (size_t i = 0; i < trackSize; ++i) {
288  char buffer1[256];
289  char buffer2[256];
290  sprintf(buffer1, "trackid%lu", (unsigned long)i);
291  sprintf(buffer2, "trackid%lu/i", (unsigned long)i);
292  clusters_->Branch(buffer1, trackid_ + i, buffer2);
293  }
294  clusters_->Branch("onTrack", &onTrack_, "onTrack/O");
295  clusters_->Branch("clWidth", &clWidth_, "clWidth/F");
296  clusters_->Branch("clPosition", &clPosition_, "clPosition/F");
297  clusters_->Branch("clglobalX", &globalX_, "clglobalX/F");
298  clusters_->Branch("clglobalY", &globalY_, "clglobalY/F");
299  clusters_->Branch("clglobalZ", &globalZ_, "clglobalZ/F");
300  clusters_->Branch("angle", &angle_, "angle/F");
301  clusters_->Branch("thickness", &thickness_, "thickness/F");
302  clusters_->Branch("maxCharge", &maxCharge_, "maxCharge/F");
303  clusters_->Branch("clNormalizedCharge", &clNormalizedCharge_, "clNormalizedCharge/F");
304  clusters_->Branch("clNormalizedNoise", &clNormalizedNoise_, "clNormalizedNoise/F");
305  clusters_->Branch("clSignalOverNoise", &clSignalOverNoise_, "clSignalOverNoise/F");
306  clusters_->Branch("clCorrectedCharge", &clCorrectedCharge_, "clCorrectedCharge/F");
307  clusters_->Branch("clCorrectedSignalOverNoise", &clCorrectedSignalOverNoise_, "clCorrectedSignalOverNoise/F");
308  clusters_->Branch("clBareCharge", &clBareCharge_, "clBareCharge/F");
309  clusters_->Branch("clBareNoise", &clBareNoise_, "clBareNoise/F");
310  clusters_->Branch("stripLength", &stripLength_, "stripLength/F");
311  clusters_->Branch("detid", &detid_, "detid/i");
312  clusters_->Branch("lldChannel", &lldChannel_, "lldChannel/s");
313 
314  // create a TTree for pixel clusters
315  pixclusters_ = dir->make<TTree>("pixclusters", "pixel cluster information");
316  pixclusters_->Branch("eventid", &eventid_, "eventid/i");
317  pixclusters_->Branch("runid", &runid_, "runid/i");
318  for (size_t i = 0; i < trackSize; ++i) {
319  char buffer1[256];
320  char buffer2[256];
321  sprintf(buffer1, "trackid%lu", (unsigned long)i);
322  sprintf(buffer2, "trackid%lu/i", (unsigned long)i);
323  pixclusters_->Branch(buffer1, trackid_ + i, buffer2);
324  }
325  pixclusters_->Branch("onTrack", &onTrack_, "onTrack/O");
326  pixclusters_->Branch("clPositionX", &clPositionX_, "clPositionX/F");
327  pixclusters_->Branch("clPositionY", &clPositionY_, "clPositionY/F");
328  pixclusters_->Branch("clSize", &clSize_, "clSize/i");
329  pixclusters_->Branch("clSizeX", &clSizeX_, "clSizeX/i");
330  pixclusters_->Branch("clSizeY", &clSizeY_, "clSizeY/i");
331  pixclusters_->Branch("alpha", &alpha_, "alpha/F");
332  pixclusters_->Branch("beta", &beta_, "beta/F");
333  pixclusters_->Branch("charge", &charge_, "charge/F");
334  pixclusters_->Branch("chargeCorr", &chargeCorr_, "chargeCorr/F");
335  pixclusters_->Branch("clglobalX", &globalX_, "clglobalX/F");
336  pixclusters_->Branch("clglobalY", &globalY_, "clglobalY/F");
337  pixclusters_->Branch("clglobalZ", &globalZ_, "clglobalZ/F");
338  pixclusters_->Branch("detid", &detid_, "detid/i");
339 
340  // create a tree for tracks
341  for (size_t i = 0; i < trackSize; ++i) {
342  char buffer1[256];
343  char buffer2[256];
344  sprintf(buffer1, "tracks%lu", (unsigned long)i);
345  sprintf(buffer2, "track%lu information", (unsigned long)i);
346  TTree* thetracks_ = dir->make<TTree>(buffer1, buffer2);
347  sprintf(buffer1, "trackid%lu", (unsigned long)i);
348  sprintf(buffer2, "trackid%lu/i", (unsigned long)i);
349  thetracks_->Branch(buffer1, globaltrackid_ + i, buffer2);
350  thetracks_->Branch("eventid", &eventid_, "eventid/i");
351  thetracks_->Branch("runid", &runid_, "runid/i");
352  thetracks_->Branch("chi2", &chi2_, "chi2/F");
353  thetracks_->Branch("eta", &eta_, "eta/F");
354  thetracks_->Branch("etaerr", &etaerr_, "etaerr/F");
355  thetracks_->Branch("phi", &phi_, "phi/F");
356  thetracks_->Branch("phierr", &phierr_, "phierr/F");
357  thetracks_->Branch("dedx1", &dedx1_, "dedx1/F");
358  thetracks_->Branch("dedx2", &dedx2_, "dedx2/F");
359  thetracks_->Branch("dedx3", &dedx3_, "dedx3/F");
360  thetracks_->Branch("dedxNoM", &dedxNoM_, "dedxNoM/i");
361  thetracks_->Branch("charge", &charge_, "charge/F");
362  thetracks_->Branch("quality", &quality_, "quality/i");
363  thetracks_->Branch("foundhits", &foundhits_, "foundhits/i");
364  thetracks_->Branch("lostHits", &lostHits_, "lostHits/i");
365  thetracks_->Branch("foundhitsStrips", &foundhitsStrips_, "foundhitsStrips/i");
366  thetracks_->Branch("foundhitsPixels", &foundhitsPixels_, "foundhitsPixels/i");
367  thetracks_->Branch("losthitsStrips", &losthitsStrips_, "losthitsStrips/i");
368  thetracks_->Branch("losthitsPixels", &losthitsPixels_, "losthitsPixels/i");
369  thetracks_->Branch("p", &p_, "p/F");
370  thetracks_->Branch("pt", &pt_, "pt/F");
371  thetracks_->Branch("pterr", &pterr_, "pterr/F");
372  thetracks_->Branch("ndof", &ndof_, "ndof/i");
373  thetracks_->Branch("dz", &dz_, "dz/F");
374  thetracks_->Branch("dzerr", &dzerr_, "dzerr/F");
375  thetracks_->Branch("dzCorr", &dzCorr_, "dzCorr/F");
376  thetracks_->Branch("dxy", &dxy_, "dxy/F");
377  thetracks_->Branch("dxyerr", &dxyerr_, "dxyerr/F");
378  thetracks_->Branch("dxyCorr", &dxyCorr_, "dxyCorr/F");
379  thetracks_->Branch("qoverp", &qoverp_, "qoverp/F");
380  thetracks_->Branch("xPCA", &xPCA_, "xPCA/F");
381  thetracks_->Branch("yPCA", &yPCA_, "yPCA/F");
382  thetracks_->Branch("zPCA", &zPCA_, "zPCA/F");
383  thetracks_->Branch("nLayers", &nLayers_, "nLayers/i");
384  thetracks_->Branch("trkWeightpvtx", &trkWeightpvtx_, "trkWeightpvtx/F");
385  thetracks_->Branch("vertexid", &vertexid_, "vertexid/i");
386  tracks_.push_back(thetracks_);
387  }
388 
389  // create a tree for missing hits
390  for (size_t i = 0; i < trackSize; ++i) {
391  char buffer1[256];
392  char buffer2[256];
393  sprintf(buffer1, "misingHits%lu", (unsigned long)i);
394  sprintf(buffer2, "missing hits from track collection %lu", (unsigned long)i);
395  TTree* themissingHits_ = dir->make<TTree>(buffer1, buffer2);
396  sprintf(buffer1, "trackid%lu", (unsigned long)i);
397  sprintf(buffer2, "trackid%lu/i", (unsigned long)i);
398  themissingHits_->Branch(buffer1, globaltrackid_ + i, buffer2);
399  themissingHits_->Branch("eventid", &eventid_, "eventid/i");
400  themissingHits_->Branch("runid", &runid_, "runid/i");
401  themissingHits_->Branch("detid", &detid_, "detid/i");
402  themissingHits_->Branch("type", &type_, "type/i");
403  themissingHits_->Branch("localX", &clPositionX_, "localX/F");
404  themissingHits_->Branch("localY", &clPositionY_, "localY/F");
405  themissingHits_->Branch("globalX", &globalX_, "globalX/F");
406  themissingHits_->Branch("globalY", &globalY_, "globalY/F");
407  themissingHits_->Branch("globalZ", &globalZ_, "globalZ/F");
408  themissingHits_->Branch("measX", &measX_, "measX/F");
409  themissingHits_->Branch("measY", &measY_, "measY/F");
410  themissingHits_->Branch("errorX", &errorX_, "errorX/F");
411  themissingHits_->Branch("errorY", &errorY_, "errorY/F");
412  missingHits_.push_back(themissingHits_);
413  }
414 
415  // create a tree for the vertices
416  vertices_ = dir->make<TTree>("vertices", "vertex information");
417  vertices_->Branch("vertexid", &globalvertexid_, "vertexid/i");
418  vertices_->Branch("eventid", &eventid_, "eventid/i");
419  vertices_->Branch("runid", &runid_, "runid/i");
420  vertices_->Branch("nTracks", &nTracks_pvtx_, "nTracks/i");
421  vertices_->Branch("sumptsq", &sumptsq_pvtx_, "sumptsq/F");
422  vertices_->Branch("isValid", &isValid_pvtx_, "isValid/O");
423  vertices_->Branch("isFake", &isFake_pvtx_, "isFake/O");
424  vertices_->Branch("recx", &recx_pvtx_, "recx/F");
425  vertices_->Branch("recy", &recy_pvtx_, "recy/F");
426  vertices_->Branch("recz", &recz_pvtx_, "recz/F");
427  vertices_->Branch("recx_err", &recx_err_pvtx_, "recx_err/F");
428  vertices_->Branch("recy_err", &recy_err_pvtx_, "recy_err/F");
429  vertices_->Branch("recz_err", &recz_err_pvtx_, "recz_err/F");
430 
431  // create a tree for the vertices
432  pixelVertices_ = dir->make<TTree>("pixelVertices", "pixel vertex information");
433  pixelVertices_->Branch("eventid", &eventid_, "eventid/i");
434  pixelVertices_->Branch("runid", &runid_, "runid/i");
435  pixelVertices_->Branch("nTracks", &nTracks_pvtx_, "nTracks/i");
436  pixelVertices_->Branch("sumptsq", &sumptsq_pvtx_, "sumptsq/F");
437  pixelVertices_->Branch("isValid", &isValid_pvtx_, "isValid/O");
438  pixelVertices_->Branch("isFake", &isFake_pvtx_, "isFake/O");
439  pixelVertices_->Branch("recx", &recx_pvtx_, "recx/F");
440  pixelVertices_->Branch("recy", &recy_pvtx_, "recy/F");
441  pixelVertices_->Branch("recz", &recz_pvtx_, "recz/F");
442  pixelVertices_->Branch("recx_err", &recx_err_pvtx_, "recx_err/F");
443  pixelVertices_->Branch("recy_err", &recy_err_pvtx_, "recy_err/F");
444  pixelVertices_->Branch("recz_err", &recz_err_pvtx_, "recz_err/F");
445 
446  // create a tree for the events
447  event_ = dir->make<TTree>("events", "event information");
448  event_->Branch("eventid", &eventid_, "eventid/i");
449  event_->Branch("runid", &runid_, "runid/i");
450  event_->Branch("L1DecisionBits", L1DecisionBits_, "L1DecisionBits[192]/O");
451  event_->Branch("L1TechnicalBits", L1TechnicalBits_, "L1TechnicalBits[64]/O");
452  event_->Branch("orbit", &orbit_, "orbit/i");
453  event_->Branch("orbitL1", &orbitL1_, "orbitL1/i");
454  event_->Branch("bx", &bx_, "bx/i");
455  event_->Branch("store", &store_, "store/i");
456  event_->Branch("time", &time_, "time/i");
457  event_->Branch("delay", &delay_, "delay/F");
458  event_->Branch("lumiSegment", &lumiSegment_, "lumiSegment/s");
459  event_->Branch("physicsDeclared", &physicsDeclared_, "physicsDeclared/s");
460  event_->Branch("HLTDecisionBits", HLTDecisionBits_, "HLTDecisionBits[256]/O");
461  char buffer[256];
462  sprintf(buffer, "ntracks[%lu]/i", (unsigned long)trackSize);
463  event_->Branch("ntracks", ntracks_, buffer);
464  sprintf(buffer, "ntrajs[%lu]/i", (unsigned long)trackSize);
465  event_->Branch("ntrajs", ntrajs_, buffer);
466  sprintf(buffer, "lowPixelProbabilityFraction[%lu]/F", (unsigned long)trackSize);
467  event_->Branch("lowPixelProbabilityFraction", lowPixelProbabilityFraction_, buffer);
468  event_->Branch("nclusters", &nclusters_, "nclusters/i");
469  event_->Branch("npixClusters", &npixClusters_, "npixClusters/i");
470  event_->Branch("nclustersOntrack", &nclustersOntrack_, "nclustersOntrack/i");
471  event_->Branch("npixClustersOntrack", &npixClustersOntrack_, "npixClustersOntrack/i");
472  event_->Branch("bsX0", &bsX0_, "bsX0/F");
473  event_->Branch("bsY0", &bsY0_, "bsY0/F");
474  event_->Branch("bsZ0", &bsZ0_, "bsZ0/F");
475  event_->Branch("bsSigmaZ", &bsSigmaZ_, "bsSigmaZ/F");
476  event_->Branch("bsDxdz", &bsDxdz_, "bsDxdz/F");
477  event_->Branch("bsDydz", &bsDydz_, "bsDydz/F");
478  event_->Branch("nVertices", &nVertices_, "nVertices/i");
479  event_->Branch("thrustValue", &thrustValue_, "thrustValue/F");
480  event_->Branch("thrustX", &thrustX_, "thrustX/F");
481  event_->Branch("thrustY", &thrustY_, "thrustY/F");
482  event_->Branch("thrustZ", &thrustZ_, "thrustZ/F");
483  event_->Branch("sphericity", &sphericity_, "sphericity/F");
484  event_->Branch("planarity", &planarity_, "planarity/F");
485  event_->Branch("aplanarity", &aplanarity_, "aplanarity/F");
486  event_->Branch("MagneticField", &fBz_, "MagneticField/F");
487 
488  // cabling
489  cablingFileName_ = iConfig.getUntrackedParameter<std::string>("PSUFileName", "PSUmapping.csv");
491  iConfig.getUntrackedParameter<std::vector<std::string> >("DelayFileNames", std::vector<std::string>(0));
492  psumap_ = dir->make<TTree>("psumap", "PSU map");
493  psumap_->Branch("PSUname", PSUname_, "PSUname/C");
494  psumap_->Branch("dcuId", &dcuId_, "dcuId/i");
495  readoutmap_ = dir->make<TTree>("readoutMap", "cabling map");
496  readoutmap_->Branch("detid", &detid_, "detid/i");
497  readoutmap_->Branch("dcuId", &dcuId_, "dcuId/i");
498  readoutmap_->Branch("fecCrate", &fecCrate_, "fecCrate/s");
499  readoutmap_->Branch("fecSlot", &fecSlot_, "fecSlot/s");
500  readoutmap_->Branch("fecRing", &fecRing_, "fecRing/s");
501  readoutmap_->Branch("ccuAdd", &ccuAdd_, "ccuAdd/s");
502  readoutmap_->Branch("ccuChan", &ccuChan_, "ccuChan/s");
503  readoutmap_->Branch("lldChannel", &lldChannel_, "lldChannel/s");
504  readoutmap_->Branch("fedId", &fedId_, "fedId/s");
505  readoutmap_->Branch("fedCh", &fedCh_, "fedCh/s");
506  readoutmap_->Branch("fiberLength", &fiberLength_, "fiberLength/s");
507  readoutmap_->Branch("moduleName", moduleName_, "moduleName/C");
508  readoutmap_->Branch("moduleId", moduleId_, "moduleId/C");
509  readoutmap_->Branch("delay", &delay_, "delay/F");
510  readoutmap_->Branch("globalX", &globalX_, "globalX/F");
511  readoutmap_->Branch("globalY", &globalY_, "globalY/F");
512  readoutmap_->Branch("globalZ", &globalZ_, "globalZ/F");
513 }
514 
516  delete[] moduleName_;
517  delete[] moduleId_;
518 }
519 
520 //
521 // member functions
522 //
523 
524 // ------------ method called to for each event ------------
526  using namespace edm;
527  using namespace reco;
528  using namespace std;
529  using reco::TrackCollection;
530 
531  // load event info
532  eventid_ = iEvent.id().event();
533  runid_ = iEvent.id().run();
534  bx_ = iEvent.eventAuxiliary().bunchCrossing();
535  orbit_ = iEvent.eventAuxiliary().orbitNumber();
536  store_ = iEvent.eventAuxiliary().storeNumber();
537  time_ = iEvent.eventAuxiliary().time().value();
539 
540  // Retrieve commissioning information from "event summary", when available (for standard fine delay)
542  iEvent.getByToken(summaryToken_, summary);
543  if (summary.isValid())
544  delay_ = delay(*summary.product());
545  else
546  delay_ = 0.;
547 
548  // -- Magnetic field
550  iSetup.get<IdealMagneticFieldRecord>().get(MF);
551  const MagneticField* theMagneticField = MF.product();
552  fBz_ = fabs(theMagneticField->inTesla(GlobalPoint(0, 0, 0)).z());
553 
554  // load trigger info
556  iEvent.getByToken(L1Token_, gtrr_handle);
557  L1GlobalTriggerReadoutRecord const* gtrr = gtrr_handle.product();
558  L1GtFdlWord fdlWord = gtrr->gtFdlWord();
559  DecisionWord L1decision = fdlWord.gtDecisionWord();
560  for (int bit = 0; bit < 128; ++bit) {
561  L1DecisionBits_[bit] = L1decision[bit];
562  }
563  DecisionWordExtended L1decisionE = fdlWord.gtDecisionWordExtended();
564  for (int bit = 0; bit < 64; ++bit) {
565  L1DecisionBits_[bit + 128] = L1decisionE[bit];
566  }
567  TechnicalTriggerWord L1technical = fdlWord.gtTechnicalTriggerWord();
568  for (int bit = 0; bit < 64; ++bit) {
569  L1TechnicalBits_[bit] = L1technical[bit];
570  }
571  orbitL1_ = fdlWord.orbitNr();
572  physicsDeclared_ = fdlWord.physicsDeclared();
574  iEvent.getByToken(HLTToken_, trh);
575  size_t ntrh = trh->size();
576  for (size_t bit = 0; bit < 256; ++bit)
577  HLTDecisionBits_[bit] = bit < ntrh ? (bool)(trh->accept(bit)) : false;
578 
579  // load beamspot
580  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
581  iEvent.getByToken(bsToken_, recoBeamSpotHandle);
582  reco::BeamSpot bs = *recoBeamSpotHandle;
583  const Point beamSpot = recoBeamSpotHandle.isValid()
584  ? Point(recoBeamSpotHandle->x0(), recoBeamSpotHandle->y0(), recoBeamSpotHandle->z0())
585  : Point(0, 0, 0);
586  if (recoBeamSpotHandle.isValid()) {
587  bsX0_ = bs.x0();
588  bsY0_ = bs.y0();
589  bsZ0_ = bs.z0();
590  bsSigmaZ_ = bs.sigmaZ();
591  bsDxdz_ = bs.dxdz();
592  bsDydz_ = bs.dydz();
593  } else {
594  bsX0_ = 0.;
595  bsY0_ = 0.;
596  bsZ0_ = 0.;
597  bsSigmaZ_ = 0.;
598  bsDxdz_ = 0.;
599  bsDydz_ = 0.;
600  }
601 
602  // load primary vertex
603  static const reco::VertexCollection s_empty_vertexColl;
604  edm::Handle<reco::VertexCollection> vertexCollectionHandle;
605  iEvent.getByToken(vertexToken_, vertexCollectionHandle);
606  const reco::VertexCollection vertexColl = *(vertexCollectionHandle.product());
607  nVertices_ = 0;
608  for (reco::VertexCollection::const_iterator v = vertexColl.begin(); v != vertexColl.end(); ++v) {
609  if (v->isValid() && !v->isFake())
610  ++nVertices_;
611  }
612 
613  // load pixel vertices
614  // Pixel vertices are handled as primary vertices, but not linked to tracks.
615  edm::Handle<reco::VertexCollection> pixelVertexCollectionHandle;
616  iEvent.getByToken(pixelVertexToken_, pixelVertexCollectionHandle);
617  const reco::VertexCollection pixelVertexColl = *(pixelVertexCollectionHandle.product());
618  nPixelVertices_ = pixelVertexColl.size();
619 
620  // load the clusters
622  iEvent.getByToken(clusterToken_, clusters);
624  iEvent.getByToken(pixelclusterToken_, pixelclusters);
625 
626  // load dedx info
627  Handle<ValueMap<DeDxData> > dEdx1Handle;
628  Handle<ValueMap<DeDxData> > dEdx2Handle;
629  Handle<ValueMap<DeDxData> > dEdx3Handle;
630  try {
631  iEvent.getByToken(dedx1Token_, dEdx1Handle);
632  } catch (cms::Exception&) {
633  ;
634  }
635  try {
636  iEvent.getByToken(dedx2Token_, dEdx2Handle);
637  } catch (cms::Exception&) {
638  ;
639  }
640  try {
641  iEvent.getByToken(dedx3Token_, dEdx3Handle);
642  } catch (cms::Exception&) {
643  ;
644  }
645  const ValueMap<DeDxData> dEdxTrack1 = *dEdx1Handle.product();
646  const ValueMap<DeDxData> dEdxTrack2 = *dEdx2Handle.product();
647  const ValueMap<DeDxData> dEdxTrack3 = *dEdx3Handle.product();
648 
649  // load track collections
650  const size_t trackSize(trackLabels_.size());
651  std::vector<reco::TrackCollection> trackCollection;
652  std::vector<edm::Handle<reco::TrackCollection> > trackCollectionHandle;
653  trackCollectionHandle.resize(trackSize);
654  size_t index = 0;
655  for (std::vector<edm::EDGetTokenT<reco::TrackCollection> >::const_iterator token = trackTokens_.begin();
656  token != trackTokens_.end();
657  ++token, ++index) {
658  try {
659  iEvent.getByToken(*token, trackCollectionHandle[index]);
660  } catch (cms::Exception&) {
661  ;
662  }
663  trackCollection.push_back(*trackCollectionHandle[index].product());
664  ntracks_[index] = trackCollection[index].size();
665  }
666 
667  // load the trajectory collections
668  std::vector<std::vector<Trajectory> > trajectoryCollection;
669  std::vector<edm::Handle<std::vector<Trajectory> > > trajectoryCollectionHandle;
670  trajectoryCollectionHandle.resize(trackSize);
671  index = 0;
672  for (std::vector<edm::EDGetTokenT<std::vector<Trajectory> > >::const_iterator token = trajectoryTokens_.begin();
673  token != trajectoryTokens_.end();
674  ++token, ++index) {
675  try {
676  iEvent.getByToken(*token, trajectoryCollectionHandle[index]);
677  } catch (cms::Exception&) {
678  ;
679  }
680  trajectoryCollection.push_back(*trajectoryCollectionHandle[index].product());
681  ntrajs_[index] = trajectoryCollection[index].size();
682  }
683 
684  // load the tracks/traj association maps
685  std::vector<TrajTrackAssociationCollection> TrajToTrackMap;
686  Handle<TrajTrackAssociationCollection> trajTrackAssociationHandle;
687  for (std::vector<edm::EDGetTokenT<TrajTrackAssociationCollection> >::const_iterator token =
688  trajTrackAssoTokens_.begin();
689  token != trajTrackAssoTokens_.end();
690  ++token) {
691  try {
692  iEvent.getByToken(*token, trajTrackAssociationHandle);
693  } catch (cms::Exception&) {
694  ;
695  }
696  TrajToTrackMap.push_back(*trajTrackAssociationHandle.product());
697  }
698 
699  // sanity check
700  if (!(!trackCollection.empty() && !trajectoryCollection.empty()))
701  return;
702 
703  // build the reverse map tracks -> vertex
704  std::vector<std::map<size_t, int> > trackVertices;
705  for (size_t i = 0; i < trackSize; ++i) {
706  trackVertices.push_back(inVertex(trackCollection[0], vertexColl, globalvertexid_ + 1));
707  }
708 
709  // iterate over vertices
711  for (reco::VertexCollection::const_iterator v = vertexColl.begin(); v != vertexColl.end(); ++v) {
712  nTracks_pvtx_ = v->tracksSize();
714  isValid_pvtx_ = int(v->isValid());
715  isFake_pvtx_ = int(v->isFake());
716  recx_pvtx_ = v->x();
717  recy_pvtx_ = v->y();
718  recz_pvtx_ = v->z();
719  recx_err_pvtx_ = v->xError();
720  recy_err_pvtx_ = v->yError();
721  recz_err_pvtx_ = v->zError();
722  globalvertexid_++;
723  vertices_->Fill();
724  }
725  }
726 
727  // iterate over pixel vertices
729  for (reco::VertexCollection::const_iterator v = pixelVertexColl.begin(); v != pixelVertexColl.end(); ++v) {
730  nTracks_pvtx_ = v->tracksSize();
732  isValid_pvtx_ = int(v->isValid());
733  isFake_pvtx_ = int(v->isFake());
734  recx_pvtx_ = v->x();
735  recy_pvtx_ = v->y();
736  recz_pvtx_ = v->z();
737  recx_err_pvtx_ = v->xError();
738  recy_err_pvtx_ = v->yError();
739  recz_err_pvtx_ = v->zError();
740  pixelVertices_->Fill();
741  }
742  }
743 
744  // determine if each cluster is on a track or not, and record the local angle
745  // to do this, we use the first track/traj collection
746  std::vector<double> clusterOntrackAngles = onTrackAngles(clusters, trajectoryCollection[0]);
747  std::vector<std::pair<double, double> > pixclusterOntrackAngles =
748  onTrackAngles(pixelclusters, trajectoryCollection[0]);
749 
750  /*
751  // iterate over trajectories
752  // note: when iterating over trajectories, it might be simpler to use the tracks/trajectories association map
753  for(std::vector<Trajectory>::const_iterator traj = trajVec.begin(); traj< trajVec.end(); ++traj) {
754  }
755  // loop over all rechits from trajectories
756  //iterate over trajectories
757  for(std::vector<Trajectory>::const_iterator traj = trajVec.begin(); traj< trajVec.end(); ++traj) {
758  Trajectory::DataContainer measurements = traj->measurements();
759  // iterate over measurements
760  for(Trajectory::DataContainer::iterator meas = measurements.begin(); meas!= measurements.end(); ++meas) {
761  }
762  }
763 */
764 
765  // determine if each cluster is on a track or not, and record the trackid
766  std::vector<std::vector<int> > stripClusterOntrackIndices;
767  for (size_t i = 0; i < trackSize; ++i) {
768  stripClusterOntrackIndices.push_back(onTrack(clusters, trackCollection[i], globaltrackid_[i] + 1));
769  }
770  std::vector<std::vector<int> > pixelClusterOntrackIndices;
771  for (size_t i = 0; i < trackSize; ++i) {
772  pixelClusterOntrackIndices.push_back(onTrack(pixelclusters, trackCollection[i], globaltrackid_[i] + 1));
773  }
774  nclustersOntrack_ = count_if(
775  stripClusterOntrackIndices[0].begin(), stripClusterOntrackIndices[0].end(), [](auto c) { return c != -1; });
776  npixClustersOntrack_ = count_if(
777  pixelClusterOntrackIndices[0].begin(), pixelClusterOntrackIndices[0].end(), [](auto c) { return c != -1; });
778 
779  // iterate over tracks
780  for (size_t coll = 0; coll < trackCollection.size(); ++coll) {
781  uint32_t n_hits_barrel = 0;
782  uint32_t n_hits_lowprob = 0;
783  for (TrajTrackAssociationCollection::const_iterator it = TrajToTrackMap[coll].begin();
784  it != TrajToTrackMap[coll].end();
785  ++it) {
786  reco::TrackRef itTrack = it->val;
787  edm::Ref<std::vector<Trajectory> > traj = it->key; // bug to find type of the key
788  eta_ = itTrack->eta();
789  phi_ = itTrack->phi();
790  try { // not all track collections have the dedx info... indeed at best one.
791  dedxNoM_ = dEdxTrack1[itTrack].numberOfMeasurements();
792  dedx1_ = dEdxTrack1[itTrack].dEdx();
793  dedx2_ = dEdxTrack2[itTrack].dEdx();
794  dedx3_ = dEdxTrack3[itTrack].dEdx();
795  } catch (cms::Exception&) {
796  dedxNoM_ = 0;
797  dedx1_ = 0.;
798  dedx2_ = 0.;
799  dedx3_ = 0.;
800  }
801  charge_ = itTrack->charge();
802  quality_ = itTrack->qualityMask();
803  foundhits_ = itTrack->found();
804  lostHits_ = itTrack->lost();
805  foundhitsStrips_ = itTrack->hitPattern().numberOfValidStripHits();
806  foundhitsPixels_ = itTrack->hitPattern().numberOfValidPixelHits();
807  losthitsStrips_ = itTrack->hitPattern().numberOfLostStripHits(reco::HitPattern::TRACK_HITS);
808  losthitsPixels_ = itTrack->hitPattern().numberOfLostPixelHits(reco::HitPattern::TRACK_HITS);
809  nLayers_ = uint32_t(itTrack->hitPattern().trackerLayersWithMeasurement());
810  p_ = itTrack->p();
811  pt_ = itTrack->pt();
812  chi2_ = itTrack->chi2();
813  ndof_ = (uint32_t)itTrack->ndof();
814  dz_ = itTrack->dz();
815  dzerr_ = itTrack->dzError();
816  dzCorr_ = itTrack->dz(beamSpot);
817  dxy_ = itTrack->dxy();
818  dxyerr_ = itTrack->dxyError();
819  dxyCorr_ = itTrack->dxy(beamSpot);
820  pterr_ = itTrack->ptError();
821  etaerr_ = itTrack->etaError();
822  phierr_ = itTrack->phiError();
823  qoverp_ = itTrack->qoverp();
824  xPCA_ = itTrack->vertex().x();
825  yPCA_ = itTrack->vertex().y();
826  zPCA_ = itTrack->vertex().z();
827  try { // only one track collection (at best) is connected to the main vertex
828  if (!vertexColl.empty() && !vertexColl.begin()->isFake()) {
829  trkWeightpvtx_ = vertexColl.begin()->trackWeight(itTrack);
830  } else
831  trkWeightpvtx_ = 0.;
832  } catch (cms::Exception&) {
833  trkWeightpvtx_ = 0.;
834  }
835  globaltrackid_[coll]++;
836  std::map<size_t, int>::const_iterator theV = trackVertices[coll].find(itTrack.key());
837  vertexid_ = (theV != trackVertices[coll].end()) ? theV->second : 0;
838  // add missing hits (separate tree, common strip + pixel)
839  Trajectory::DataContainer const& measurements = traj->measurements();
841  for (Trajectory::DataContainer::const_iterator it = measurements.begin(); it != measurements.end(); ++it) {
842  TrajectoryMeasurement::ConstRecHitPointer rechit = it->recHit();
843  if (!rechit->isValid()) {
844  // detid
845  detid_ = rechit->geographicalId();
846  // status
847  type_ = rechit->getType();
848  // position
849  LocalPoint local = it->predictedState().localPosition();
850  clPositionX_ = local.x();
851  clPositionY_ = local.y();
852  // global position
853  GlobalPoint global = it->predictedState().globalPosition();
854  globalX_ = global.x();
855  globalY_ = global.y();
856  globalZ_ = global.z();
857  // position in the measurement frame
858  measX_ = 0;
859  measY_ = 0;
861  const GeomDetUnit* gdu = static_cast<const GeomDetUnit*>(tracker_->idToDetUnit(detid_));
862  if (gdu && gdu->type().isTracker()) {
863  const Topology& topo = gdu->topology();
864  MeasurementPoint meas = topo.measurementPosition(local);
865  measX_ = meas.x();
866  measY_ = meas.y();
867  }
868  }
869  // local error
870  LocalError error = it->predictedState().localError().positionError();
871  errorX_ = error.xx();
872  errorY_ = error.yy();
873  // fill
874  missingHits_[coll]->Fill();
875  }
876  }
877  }
878  // compute the fraction of low probability pixels... will be added to the event tree
879  for (trackingRecHit_iterator it = itTrack->recHitsBegin(); it != itTrack->recHitsEnd(); ++it) {
880  const TrackingRecHit* hit = &(**it);
881  const SiPixelRecHit* pixhit = dynamic_cast<const SiPixelRecHit*>(hit);
882  if (pixhit) {
883  DetId detId = pixhit->geographicalId();
884  if (detId.subdetId() == PixelSubdetector::PixelBarrel) {
885  ++n_hits_barrel;
886  double proba = pixhit->clusterProbability(0);
887  if (proba <= 0.0)
888  ++n_hits_lowprob;
889  }
890  }
891  }
892  // fill the track tree
894  tracks_[coll]->Fill();
895  }
896  lowPixelProbabilityFraction_[coll] = n_hits_barrel > 0 ? (float)n_hits_lowprob / n_hits_barrel : -1.;
897  }
898 
899  // iterate over clusters
900  nclusters_ = 0;
901  std::vector<double>::const_iterator angleIt = clusterOntrackAngles.begin();
902  uint32_t localCounter = 0;
903  for (edmNew::DetSetVector<SiStripCluster>::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
904  DSViter++) {
907  uint32_t detid = DSViter->id();
908  nclusters_ += DSViter->size();
910  for (edmNew::DetSet<SiStripCluster>::const_iterator iter = begin; iter != end;
911  ++iter, ++angleIt, ++localCounter) {
912  SiStripClusterInfo* siStripClusterInfo =
913  new SiStripClusterInfo(*iter, iSetup, detid, std::string("")); //string = quality label
914  // general quantities
915  for (size_t i = 0; i < trackSize; ++i) {
916  trackid_[i] = stripClusterOntrackIndices[i][localCounter];
917  }
918  onTrack_ = (trackid_[0] != (uint32_t)-1);
919  clWidth_ = siStripClusterInfo->width();
920  clPosition_ = siStripClusterInfo->baryStrip();
921  angle_ = *angleIt;
922  thickness_ = ((((DSViter->id() >> 25) & 0x7f) == 0xd) ||
923  ((((DSViter->id() >> 25) & 0x7f) == 0xe) && (((DSViter->id() >> 5) & 0x7) > 4)))
924  ? 500
925  : 300;
926  stripLength_ = static_cast<const StripGeomDetUnit*>(tracker_->idToDet(detid))->specificTopology().stripLength();
927  int nstrips = static_cast<const StripGeomDetUnit*>(tracker_->idToDet(detid))->specificTopology().nstrips();
928  maxCharge_ = siStripClusterInfo->maxCharge();
929  // signal and noise with gain corrections
930  clNormalizedCharge_ = siStripClusterInfo->charge();
931  clNormalizedNoise_ = siStripClusterInfo->noiseRescaledByGain();
932  clSignalOverNoise_ = siStripClusterInfo->signalOverNoise();
933  // signal and noise with gain corrections and angle corrections
934  clCorrectedCharge_ = clNormalizedCharge_ * fabs(cos(angle_)); // corrected for track angle
935  clCorrectedSignalOverNoise_ = clSignalOverNoise_ * fabs(cos(angle_)); // corrected for track angle
936  // signal and noise without gain corrections
937  clBareNoise_ = siStripClusterInfo->noise();
939  // global position
940  const StripGeomDetUnit* sgdu = static_cast<const StripGeomDetUnit*>(tracker_->idToDet(detid));
943  globalX_ = gp.x();
944  globalY_ = gp.y();
945  globalZ_ = gp.z();
946  // cabling
947  detid_ = detid;
948  lldChannel_ = 1 + (int(floor(iter->barycenter())) / 256);
949  if (lldChannel_ == 2 && nstrips == 512)
950  lldChannel_ = 3;
952  clusters_->Fill();
953  delete siStripClusterInfo;
954  }
955  }
956  }
957 
958  // iterate over pixel clusters
959  npixClusters_ = 0;
960  std::vector<std::pair<double, double> >::const_iterator pixAngleIt = pixclusterOntrackAngles.begin();
961  localCounter = 0;
962  for (edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = pixelclusters->begin();
963  DSViter != pixelclusters->end();
964  DSViter++) {
967  uint32_t detid = DSViter->id();
968  npixClusters_ += DSViter->size();
970  for (edmNew::DetSet<SiPixelCluster>::const_iterator iter = begin; iter != end;
971  ++iter, ++pixAngleIt, ++localCounter) {
972  // general quantities
973  for (size_t i = 0; i < trackSize; ++i) {
974  trackid_[i] = pixelClusterOntrackIndices[i][localCounter];
975  }
976  onTrack_ = (trackid_[0] != (uint32_t)-1);
977  clPositionX_ = iter->x();
978  clPositionY_ = iter->y();
979  clSize_ = iter->size();
980  clSizeX_ = iter->sizeX();
981  clSizeY_ = iter->sizeY();
982  alpha_ = pixAngleIt->first;
983  beta_ = pixAngleIt->second;
984  charge_ = (iter->charge()) / 1000.;
985  chargeCorr_ = charge_ * sqrt(1.0 / (1.0 / pow(tan(alpha_), 2) + 1.0 / pow(tan(beta_), 2) + 1.0)) / 1000.;
986  // global position
987  const PixelGeomDetUnit* pgdu = static_cast<const PixelGeomDetUnit*>(tracker_->idToDet(detid));
990  globalX_ = gp.x();
991  globalY_ = gp.y();
992  globalZ_ = gp.z();
993  // cabling
994  detid_ = detid;
995  // fill
996  pixclusters_->Fill();
997  }
998  }
999  }
1000 
1001  // topological quantities - uses the first track collection
1002  EventShape shape(trackCollection[0]);
1003  math::XYZTLorentzVectorF thrust = shape.thrust();
1004  thrustValue_ = thrust.t();
1005  thrustX_ = thrust.x();
1006  thrustY_ = thrust.y();
1007  thrustZ_ = thrust.z();
1008  sphericity_ = shape.sphericity();
1009  planarity_ = shape.planarity();
1010  aplanarity_ = shape.aplanarity();
1011 
1012  // fill event tree
1014  event_->Fill();
1015 }
1016 
1017 // ------------ method called once each job just before starting event loop ------------
1018 void TrackerDpgAnalysis::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
1019  //Retrieve tracker topology from geometry
1020  edm::ESHandle<TrackerTopology> tTopoHandle;
1021  iSetup.get<TrackerTopologyRcd>().get(tTopoHandle);
1022  const TrackerTopology* const tTopo = tTopoHandle.product();
1023 
1024  //geometry
1025  iSetup.get<TrackerDigiGeometryRecord>().get(tracker_);
1026 
1027  //HLT names
1028  bool changed(true);
1029  if (hltConfig_.init(iRun, iSetup, HLTTag_.process(), changed)) {
1030  if (changed) {
1032  }
1033  }
1034  int i = 0;
1035  for (std::vector<std::string>::const_iterator it = hlNames_.begin(); it < hlNames_.end(); ++it) {
1036  std::cout << (i++) << " = " << (*it) << std::endl;
1037  }
1038 
1039  // read the delay offsets for each device from input files
1040  // this is only for the so-called "random delay" run
1041  std::map<uint32_t, float> delayMap = delay(delayFileNames_);
1042  TrackerMap tmap("Delays");
1043 
1044  // cabling I (readout)
1045  iSetup.get<SiStripFedCablingRcd>().get(cabling_);
1046  auto feds = cabling_->fedIds();
1047  for (auto fedid = feds.begin(); fedid < feds.end(); ++fedid) {
1048  auto connections = cabling_->fedConnections(*fedid);
1049  for (auto conn = connections.begin(); conn < connections.end(); ++conn) {
1050  // Fill the "old" map to be used for lookup during analysis
1051  if (conn->isConnected())
1052  connections_.insert(std::make_pair(conn->detId(), new FedChannelConnection(*conn)));
1053  // Fill the standalone tree (once for all)
1054  if (conn->isConnected()) {
1055  detid_ = conn->detId();
1056  strncpy(moduleName_, toStringName(detid_, tTopo).c_str(), 256);
1057  strncpy(moduleId_, toStringId(detid_).c_str(), 256);
1058  lldChannel_ = conn->lldChannel();
1059  dcuId_ = conn->dcuId();
1060  fecCrate_ = conn->fecCrate();
1061  fecSlot_ = conn->fecSlot();
1062  fecRing_ = conn->fecRing();
1063  ccuAdd_ = conn->ccuAddr();
1064  ccuChan_ = conn->ccuChan();
1065  fedId_ = conn->fedId();
1066  fedCh_ = conn->fedCh();
1067  fiberLength_ = conn->fiberLength();
1068  delay_ = delayMap[dcuId_];
1069  const StripGeomDetUnit* sgdu = static_cast<const StripGeomDetUnit*>(tracker_->idToDet(detid_));
1071  globalX_ = gp.x();
1072  globalY_ = gp.y();
1073  globalZ_ = gp.z();
1074  readoutmap_->Fill();
1076  }
1077  }
1078  }
1079  if (!delayMap.empty())
1080  tmap.save(true, 0, 0, "delaymap.png");
1081 
1082  // cabling II (DCU map)
1083  std::ifstream cablingFile(cablingFileName_.c_str());
1084  if (cablingFile.is_open()) {
1085  char buffer[1024];
1086  cablingFile.getline(buffer, 1024);
1087  while (!cablingFile.eof()) {
1088  std::istringstream line(buffer);
1089  std::string name;
1090  // one line contains the PSU name + all dcuids connected to it.
1091  line >> name;
1092  strncpy(PSUname_, name.c_str(), 256);
1093  while (!line.eof()) {
1094  line >> dcuId_;
1095  psumap_->Fill();
1096  }
1097  cablingFile.getline(buffer, 1024);
1098  }
1099  } else {
1100  edm::LogWarning("BadConfig") << " The PSU file does not exist. The psumap tree will not be filled." << std::endl
1101  << " Looking for " << cablingFileName_.c_str() << "." << std::endl
1102  << " Please specify a valid filename through the PSUFileName untracked parameter.";
1103  }
1104 }
1105 
1106 // ------------ method called once each job just after ending the event loop ------------
1108  for (size_t i = 0; i < tracks_.size(); ++i) {
1109  char buffer[256];
1110  sprintf(buffer, "trackid%lu", (unsigned long)i);
1111  if (tracks_[i]->GetEntries())
1112  tracks_[i]->BuildIndex(buffer, "eventid");
1113  }
1114  /* not needed: missing hits is a high-level quantity
1115  for(size_t i = 0; i<missingHits_.size();++i) {
1116  char buffer[256];
1117  sprintf(buffer,"trackid%lu",(unsigned long)i);
1118  if(missingHits_[i]->GetEntries()) missingHits_[i]->BuildIndex(buffer);
1119  }
1120  */
1121  if (vertices_->GetEntries())
1122  vertices_->BuildIndex("vertexid", "eventid");
1123  if (event_->GetEntries())
1124  event_->BuildIndex("runid", "eventid");
1125  if (psumap_->GetEntries())
1126  psumap_->BuildIndex("dcuId");
1127  if (readoutmap_->GetEntries())
1128  readoutmap_->BuildIndex("detid", "lldChannel");
1129 }
1130 
1132  const std::vector<Trajectory>& trajVec) {
1133  std::vector<double> result;
1134  // first, build a list of positions and angles on trajectories
1135  std::multimap<const uint32_t, std::pair<LocalPoint, double> > onTrackPositions;
1136  for (std::vector<Trajectory>::const_iterator traj = trajVec.begin(); traj < trajVec.end(); ++traj) {
1137  Trajectory::DataContainer measurements = traj->measurements();
1138  for (Trajectory::DataContainer::iterator meas = measurements.begin(); meas != measurements.end(); ++meas) {
1139  double tla = meas->updatedState().localDirection().theta();
1140  insertMeasurement(onTrackPositions, &(*(meas->recHit())), tla);
1141  }
1142  }
1143  // then loop over the clusters to check
1144  double angle = 0.;
1145  for (edmNew::DetSetVector<SiStripCluster>::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
1146  DSViter++) {
1149  std::pair<std::multimap<uint32_t, std::pair<LocalPoint, double> >::const_iterator,
1150  std::multimap<uint32_t, std::pair<LocalPoint, double> >::const_iterator>
1151  range = onTrackPositions.equal_range(DSViter->id());
1152  const GeomDetUnit* gdu = static_cast<const GeomDetUnit*>(tracker_->idToDet(DSViter->id()));
1153  for (edmNew::DetSet<SiStripCluster>::const_iterator iter = begin; iter != end; ++iter) {
1154  angle = 0.;
1155  for (std::multimap<uint32_t, std::pair<LocalPoint, double> >::const_iterator cl = range.first; cl != range.second;
1156  ++cl) {
1157  if (fabs(gdu->topology().measurementPosition(cl->second.first).x() - iter->barycenter()) < 2) {
1158  angle = cl->second.second;
1159  }
1160  }
1161  result.push_back(angle);
1162  }
1163  }
1164  return result;
1165 }
1166 
1167 void TrackerDpgAnalysis::insertMeasurement(std::multimap<const uint32_t, std::pair<LocalPoint, double> >& collection,
1169  double tla) {
1170  if (!hit)
1171  return;
1172  const SiTrackerMultiRecHit* multihit = dynamic_cast<const SiTrackerMultiRecHit*>(hit);
1173  const SiStripRecHit2D* singlehit = dynamic_cast<const SiStripRecHit2D*>(hit);
1174  const SiStripRecHit1D* hit1d = dynamic_cast<const SiStripRecHit1D*>(hit);
1175  if (hit1d) { //...->33X
1176  collection.insert(std::make_pair(hit1d->geographicalId().rawId(), std::make_pair(hit1d->localPosition(), tla)));
1177  } else if (singlehit) { // 41X->...
1178  collection.insert(
1179  std::make_pair(singlehit->geographicalId().rawId(), std::make_pair(singlehit->localPosition(), tla)));
1180  } else if (multihit) {
1181  std::vector<const TrackingRecHit*> childs = multihit->recHits();
1182  for (std::vector<const TrackingRecHit*>::const_iterator it = childs.begin(); it != childs.end(); ++it) {
1183  insertMeasurement(collection, dynamic_cast<const TrackingRecHit*>(*it), tla);
1184  }
1185  }
1186 }
1187 
1189  const reco::TrackCollection& trackVec,
1190  uint32_t firstTrack) {
1191  std::vector<int> result;
1192  // first, build a list of positions and trackid on tracks
1193  std::multimap<const uint32_t, std::pair<int, int> > onTrackPositions;
1194  uint32_t trackid = firstTrack;
1195  for (reco::TrackCollection::const_iterator itTrack = trackVec.begin(); itTrack != trackVec.end();
1196  ++itTrack, ++trackid) {
1197  for (trackingRecHit_iterator it = itTrack->recHitsBegin(); it != itTrack->recHitsEnd(); ++it) {
1198  const TrackingRecHit* hit = &(**it);
1199  insertMeasurement(onTrackPositions, hit, trackid);
1200  }
1201  }
1202  // then loop over the clusters to check
1203  int thetrackid = -1;
1204  for (edmNew::DetSetVector<SiStripCluster>::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
1205  DSViter++) {
1208  std::pair<std::multimap<uint32_t, std::pair<int, int> >::const_iterator,
1209  std::multimap<uint32_t, std::pair<int, int> >::const_iterator>
1210  range = onTrackPositions.equal_range(DSViter->id());
1211  for (edmNew::DetSet<SiStripCluster>::const_iterator iter = begin; iter != end; ++iter) {
1212  thetrackid = -1;
1213  for (std::multimap<uint32_t, std::pair<int, int> >::const_iterator cl = range.first; cl != range.second; ++cl) {
1214  if (fabs(cl->second.first - iter->barycenter()) < 2) {
1215  thetrackid = cl->second.second;
1216  }
1217  }
1218  result.push_back(thetrackid);
1219  }
1220  }
1221  return result;
1222 }
1223 
1224 void TrackerDpgAnalysis::insertMeasurement(std::multimap<const uint32_t, std::pair<int, int> >& collection,
1225  const TrackingRecHit* hit,
1226  int trackid) {
1227  if (!hit)
1228  return;
1229  const SiTrackerMultiRecHit* multihit = dynamic_cast<const SiTrackerMultiRecHit*>(hit);
1230  const SiStripRecHit2D* singlehit = dynamic_cast<const SiStripRecHit2D*>(hit);
1231  const SiStripRecHit1D* hit1d = dynamic_cast<const SiStripRecHit1D*>(hit);
1232  if (hit1d) { // 41X->...
1233  collection.insert(
1234  std::make_pair(hit1d->geographicalId().rawId(), std::make_pair(int(hit1d->cluster()->barycenter()), trackid)));
1235  } else if (singlehit) { //...->33X
1236  collection.insert(std::make_pair(singlehit->geographicalId().rawId(),
1237  std::make_pair(int(singlehit->cluster()->barycenter()), trackid)));
1238  } else if (multihit) {
1239  std::vector<const TrackingRecHit*> childs = multihit->recHits();
1240  for (std::vector<const TrackingRecHit*>::const_iterator it = childs.begin(); it != childs.end(); ++it) {
1241  insertMeasurement(collection, *it, trackid);
1242  }
1243  }
1244 }
1245 
1248  uint32_t firstVertex) {
1249  // build reverse map track -> vertex
1250  std::map<size_t, int> output;
1251  uint32_t vertexid = firstVertex;
1252  for (reco::VertexCollection::const_iterator v = vertices.begin(); v != vertices.end(); ++v, ++vertexid) {
1253  reco::Vertex::trackRef_iterator it = v->tracks_begin();
1254  reco::Vertex::trackRef_iterator lastTrack = v->tracks_end();
1255  for (; it != lastTrack; ++it) {
1256  output[it->key()] = vertexid;
1257  }
1258  }
1259  return output;
1260 }
1261 
1262 std::vector<std::pair<double, double> > TrackerDpgAnalysis::onTrackAngles(
1263  edm::Handle<edmNew::DetSetVector<SiPixelCluster> >& clusters, const std::vector<Trajectory>& trajVec) {
1264  std::vector<std::pair<double, double> > result;
1265  // first, build a list of positions and angles on trajectories
1266  std::multimap<const uint32_t, std::pair<LocalPoint, std::pair<double, double> > > onTrackPositions;
1267  for (std::vector<Trajectory>::const_iterator traj = trajVec.begin(); traj < trajVec.end(); ++traj) {
1268  Trajectory::DataContainer measurements = traj->measurements();
1269  for (Trajectory::DataContainer::iterator meas = measurements.begin(); meas != measurements.end(); ++meas) {
1270  LocalVector localDir = meas->updatedState().localDirection();
1271  double alpha = atan2(localDir.z(), localDir.x());
1272  double beta = atan2(localDir.z(), localDir.y());
1273  insertMeasurement(onTrackPositions, &(*(meas->recHit())), alpha, beta);
1274  }
1275  }
1276  // then loop over the clusters to check
1277  double alpha = 0.;
1278  double beta = 0.;
1279  for (edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
1280  DSViter++) {
1283  for (edmNew::DetSet<SiPixelCluster>::const_iterator iter = begin; iter != end; ++iter) {
1284  alpha = 0.;
1285  beta = 0.;
1286  std::pair<std::multimap<uint32_t, std::pair<LocalPoint, std::pair<double, double> > >::const_iterator,
1287  std::multimap<uint32_t, std::pair<LocalPoint, std::pair<double, double> > >::const_iterator>
1288  range = onTrackPositions.equal_range(DSViter->id());
1289  const GeomDetUnit* gdu = static_cast<const GeomDetUnit*>(tracker_->idToDet(DSViter->id()));
1290  for (std::multimap<uint32_t, std::pair<LocalPoint, std::pair<double, double> > >::const_iterator cl = range.first;
1291  cl != range.second;
1292  ++cl) {
1293  if (fabs(gdu->topology().measurementPosition(cl->second.first).x() - iter->x()) < 2 &&
1294  fabs(gdu->topology().measurementPosition(cl->second.first).y() - iter->y()) < 2) {
1295  alpha = cl->second.second.first;
1296  beta = cl->second.second.second;
1297  }
1298  }
1299  result.push_back(std::make_pair(alpha, beta));
1300  }
1301  }
1302  return result;
1303 }
1304 
1306  std::multimap<const uint32_t, std::pair<LocalPoint, std::pair<double, double> > >& collection,
1308  double alpha,
1309  double beta) {
1310  if (!hit)
1311  return;
1312  const SiPixelRecHit* pixhit = dynamic_cast<const SiPixelRecHit*>(hit);
1313  if (pixhit) {
1314  collection.insert(std::make_pair(pixhit->geographicalId().rawId(),
1315  std::make_pair(pixhit->localPosition(), std::make_pair(alpha, beta))));
1316  }
1317 }
1318 
1320  const reco::TrackCollection& trackVec,
1321  uint32_t firstTrack) {
1322  std::vector<int> result;
1323  // first, build a list of positions and trackid on tracks
1324  std::multimap<const uint32_t, std::pair<std::pair<float, float>, int> > onTrackPositions;
1325  uint32_t trackid = firstTrack;
1326  for (reco::TrackCollection::const_iterator itTrack = trackVec.begin(); itTrack != trackVec.end();
1327  ++itTrack, ++trackid) {
1328  for (trackingRecHit_iterator it = itTrack->recHitsBegin(); it != itTrack->recHitsEnd(); ++it) {
1329  const TrackingRecHit* hit = &(**it);
1330  insertMeasurement(onTrackPositions, hit, trackid);
1331  }
1332  }
1333  // then loop over the clusters to check
1334  int thetrackid = -1;
1335  for (edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
1336  DSViter++) {
1339  for (edmNew::DetSet<SiPixelCluster>::const_iterator iter = begin; iter != end; ++iter) {
1340  thetrackid = -1;
1341  std::pair<std::multimap<uint32_t, std::pair<std::pair<float, float>, int> >::const_iterator,
1342  std::multimap<uint32_t, std::pair<std::pair<float, float>, int> >::const_iterator>
1343  range = onTrackPositions.equal_range(DSViter->id());
1344  for (std::multimap<uint32_t, std::pair<std::pair<float, float>, int> >::const_iterator cl = range.first;
1345  cl != range.second;
1346  ++cl) {
1347  if ((fabs(cl->second.first.first - iter->x()) < 2) && (fabs(cl->second.first.second - iter->y()) < 2)) {
1348  thetrackid = cl->second.second;
1349  }
1350  }
1351  result.push_back(thetrackid);
1352  }
1353  }
1354  return result;
1355 }
1356 
1358  std::multimap<const uint32_t, std::pair<std::pair<float, float>, int> >& collection,
1359  const TrackingRecHit* hit,
1360  int trackid) {
1361  if (!hit)
1362  return;
1363  const SiPixelRecHit* pixhit = dynamic_cast<const SiPixelRecHit*>(hit);
1364  if (pixhit) {
1365  collection.insert(
1366  std::make_pair(pixhit->geographicalId().rawId(),
1367  std::make_pair(std::make_pair(pixhit->cluster()->x(), pixhit->cluster()->y()), trackid)));
1368  }
1369 }
1370 
1372  SiStripDetId detid(rawid);
1373  std::string out;
1374  std::stringstream output;
1375  switch (detid.subDetector()) {
1376  case 3: {
1377  output << "TIB";
1378 
1379  output << (tTopo->tibIsZPlusSide(rawid) ? "+" : "-");
1380  output << " layer ";
1381  output << tTopo->tibLayer(rawid);
1382  output << ", string ";
1383  output << tTopo->tibString(rawid);
1384  output << (tTopo->tibIsExternalString(rawid) ? " external" : " internal");
1385  output << ", module ";
1386  output << tTopo->tibModule(rawid);
1387  if (tTopo->tibIsDoubleSide(rawid)) {
1388  output << " (double)";
1389  } else {
1390  output << (tTopo->tibIsRPhi(rawid) ? " (rphi)" : " (stereo)");
1391  }
1392  break;
1393  }
1394  case 4: {
1395  output << "TID";
1396 
1397  output << (tTopo->tidIsZPlusSide(rawid) ? "+" : "-");
1398  output << " disk ";
1399  output << tTopo->tidWheel(rawid);
1400  output << ", ring ";
1401  output << tTopo->tidRing(rawid);
1402  output << (tTopo->tidIsFrontRing(rawid) ? " front" : " back");
1403  output << ", module ";
1404  output << tTopo->tidModule(rawid);
1405  if (tTopo->tidIsDoubleSide(rawid)) {
1406  output << " (double)";
1407  } else {
1408  output << (tTopo->tidIsRPhi(rawid) ? " (rphi)" : " (stereo)");
1409  }
1410  break;
1411  }
1412  case 5: {
1413  output << "TOB";
1414 
1415  output << (tTopo->tobIsZPlusSide(rawid) ? "+" : "-");
1416  output << " layer ";
1417  output << tTopo->tobLayer(rawid);
1418  output << ", rod ";
1419  output << tTopo->tobRod(rawid);
1420  output << ", module ";
1421  output << tTopo->tobModule(rawid);
1422  if (tTopo->tobIsDoubleSide(rawid)) {
1423  output << " (double)";
1424  } else {
1425  output << (tTopo->tobIsRPhi(rawid) ? " (rphi)" : " (stereo)");
1426  }
1427  break;
1428  }
1429  case 6: {
1430  output << "TEC";
1431 
1432  output << (tTopo->tecIsZPlusSide(rawid) ? "+" : "-");
1433  output << " disk ";
1434  output << tTopo->tecWheel(rawid);
1435  output << " sector ";
1436  output << tTopo->tecPetalNumber(rawid);
1437  output << (tTopo->tecIsFrontPetal(rawid) ? " Front Petal" : " Back Petal");
1438  output << ", module ";
1439  output << tTopo->tecRing(rawid);
1440  output << tTopo->tecModule(rawid);
1441  if (tTopo->tecIsDoubleSide(rawid)) {
1442  output << " (double)";
1443  } else {
1444  output << (tTopo->tecIsRPhi(rawid) ? " (rphi)" : " (stereo)");
1445  }
1446  break;
1447  }
1448  default: {
1449  output << "UNKNOWN";
1450  }
1451  }
1452  out = output.str();
1453  return out;
1454 }
1455 
1457  std::string out;
1458  std::stringstream output;
1459  output << rawid << " (0x" << std::hex << rawid << std::dec << ")";
1460  out = output.str();
1461  return out;
1462 }
1463 
1465  double sum = 0.;
1466  double pT;
1467  for (reco::Vertex::trackRef_iterator it = v.tracks_begin(); it != v.tracks_end(); it++) {
1468  pT = (**it).pt();
1469  sum += pT * pT;
1470  }
1471  return sum;
1472 }
1473 
1475  float delay = const_cast<SiStripEventSummary&>(summary).ttcrx();
1476  uint32_t latencyCode = (const_cast<SiStripEventSummary&>(summary).layerScanned() >> 24) & 0xff;
1477  int latencyShift =
1478  latencyCode & 0x3f; // number of bunch crossings between current value and start of scan... must be positive
1479  if (latencyShift > 32)
1480  latencyShift -= 64; // allow negative values: we cover [-32,32].. should not be needed.
1481  if ((latencyCode >> 6) == 2)
1482  latencyShift -= 3; // layer in deconv, rest in peak
1483  if ((latencyCode >> 6) == 1)
1484  latencyShift += 3; // layer in peak, rest in deconv
1485  float correctedDelay =
1486  delay - (latencyShift * 25.); // shifts the delay so that 0 corresponds to the current settings.
1487  return correctedDelay;
1488 }
1489 
1490 std::map<uint32_t, float> TrackerDpgAnalysis::delay(const std::vector<std::string>& files) {
1491  // prepare output
1492  uint32_t dcuid;
1493  float delay;
1494  std::map<uint32_t, float> delayMap;
1495  //iterator over input files
1496  for (std::vector<std::string>::const_iterator file = files.begin(); file < files.end(); ++file) {
1497  // open the file
1498  std::ifstream cablingFile(file->c_str());
1499  if (cablingFile.is_open()) {
1500  char buffer[1024];
1501  // read one line
1502  cablingFile.getline(buffer, 1024);
1503  while (!cablingFile.eof()) {
1504  std::string line(buffer);
1505  size_t pos = line.find("dcuid");
1506  // one line containing dcuid
1507  if (pos != std::string::npos) {
1508  // decode dcuid
1509  std::string dcuids = line.substr(pos + 7, line.find(" ", pos) - pos - 8);
1510  std::istringstream dcuidstr(dcuids);
1511  dcuidstr >> std::hex >> dcuid;
1512  // decode delay
1513  pos = line.find("difpll");
1514  std::string diffs = line.substr(pos + 8, line.find(" ", pos) - pos - 9);
1515  std::istringstream diffstr(diffs);
1516  diffstr >> delay;
1517  // fill the map
1518  delayMap[dcuid] = delay;
1519  }
1520  // iterate
1521  cablingFile.getline(buffer, 1024);
1522  }
1523  } else {
1524  edm::LogWarning("BadConfig") << " The delay file does not exist. The delay map will not be filled properly."
1525  << std::endl
1526  << " Looking for " << file->c_str() << "." << std::endl
1527  << " Please specify valid filenames through the DelayFileNames untracked parameter.";
1528  }
1529  }
1530  return delayMap;
1531 }
1532 
1533 //define this as a plug-in
RunNumber_t run() const
Definition: EventID.h:38
ClusterRef cluster() const
uint8_t maxCharge() const
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:81
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:40
double z0() const
z coordinate
Definition: BeamSpot.h:65
T getUntrackedParameter(std::string const &, T const &) const
std::vector< TTree * > tracks_
const TechnicalTriggerWord & gtTechnicalTriggerWord() const
get/set technical trigger bits
Definition: L1GtFdlWord.h:112
float xx() const
Definition: LocalError.h:22
bool tecIsDoubleSide(const DetId &id) const
bool tobIsDoubleSide(const DetId &id) const
EventAuxiliary const & eventAuxiliary() const override
Definition: Event.h:93
void beginRun(const edm::Run &, const edm::EventSetup &) override
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
float clusterProbability(unsigned int flags=0) const
Definition: SiPixelRecHit.cc:9
T y() const
Definition: PV2DBase.h:44
bool tibIsDoubleSide(const DetId &id) const
unsigned int tibLayer(const DetId &id) const
unsigned int tibString(const DetId &id) const
unsigned int tidRing(const DetId &id) const
edm::ParameterSet pset_
virtual const GeomDetType & type() const
Definition: GeomDet.cc:69
trackRef_iterator tracks_end() const
last iterator over tracks
Definition: Vertex.cc:73
std::vector< int > onTrack(edm::Handle< edmNew::DetSetVector< SiStripCluster > > &, const reco::TrackCollection &, uint32_t)
float noise() const
double sumPtSquared(const reco::Vertex &)
static float planarity(const reco::TrackCollection &)
Definition: EventShape.cc:275
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
virtual const Topology & topology() const
Definition: GeomDet.cc:67
std::vector< const TrackingRecHit * > recHits() const override
Access to component RecHits (if any)
float noiseRescaledByGain() const
bool isTracker() const
Definition: GeomDetType.cc:21
unsigned int tecRing(const DetId &id) const
ring id
bool accept() const
Has at least one path accepted the event?
const std::vector< std::string > & triggerNames() const
names of trigger paths
edm::EDGetTokenT< edm::TriggerResults > HLTToken_
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
std::multimap< const uint32_t, const FedChannelConnection * > connections_
float baryStrip() const
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
Timestamp const & time() const
T y() const
Definition: PV3DBase.h:60
unsigned long long EventNumber_t
bool tidIsFrontRing(const DetId &id) const
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
unsigned int tidWheel(const DetId &id) const
data_type const * const_iterator
Definition: DetSetNew.h:31
static math::XYZTLorentzVectorF thrust(const reco::TrackCollection &)
Definition: EventShape.cc:124
key_type key() const
Accessor for product key.
Definition: Ref.h:250
std::string toStringId(uint32_t)
void analyze(const edm::Event &, const edm::EventSetup &) override
virtual const StripTopology & specificTopology() const
Returns a reference to the strip proxy topology.
bool tobIsRPhi(const DetId &id) const
edm::ESHandle< TrackerGeometry > tracker_
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
int bunchCrossing() const
LuminosityBlockNumber_t luminosityBlock() const
static float sphericity(const reco::TrackCollection &)
Definition: EventShape.cc:215
edm::EDGetTokenT< reco::BeamSpot > bsToken_
bool tibIsZPlusSide(const DetId &id) const
const TrackerGeomDet * idToDetUnit(DetId) const override
Return the pointer to the GeomDetUnit corresponding to a given DetId.
double dydz() const
dydz slope
Definition: BeamSpot.h:80
std::vector< bool > DecisionWordExtended
float signalOverNoise() const
Measurement2DPoint MeasurementPoint
Measurement points are two-dimensional by default.
int iEvent
Definition: GenABIO.cc:224
bool tibIsExternalString(const DetId &id) const
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
std::vector< TrajectoryMeasurement > DataContainer
Definition: Trajectory.h:40
bool tibIsRPhi(const DetId &id) const
int storeNumber() const
uint16_t charge() const
std::vector< std::string > delayFileNames_
float yy() const
Definition: LocalError.h:24
virtual LocalPoint localPosition(float strip) const =0
Class containning control, module, detector and connection information, at the level of a FED channel...
int orbitNumber() const
edm::ESHandle< SiStripFedCabling > cabling_
std::vector< bool > DecisionWord
typedefs
edm::EventNumber_t eventid_
FedsConstIterRange fedIds() const
T sqrt(T t)
Definition: SSEVec.h:19
unsigned int size() const
Get number of paths stored.
unsigned int tidModule(const DetId &id) const
void save(bool print_total=true, float minval=0., float maxval=0., std::string s="svgmap.svg", int width=1500, int height=800)
Definition: TrackerMap.cc:811
T z() const
Definition: PV3DBase.h:61
static const int nMaxPVs_
bool tobIsZPlusSide(const DetId &id) const
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
math::XYZPoint Point
virtual MeasurementPoint measurementPosition(const LocalPoint &) const =0
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
uint16_t width() const
edm::EDGetTokenT< L1GlobalTriggerReadoutRecord > L1Token_
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
static float aplanarity(const reco::TrackCollection &)
Definition: EventShape.cc:246
std::vector< bool > TechnicalTriggerWord
technical trigger bits (64 bits)
ClusterRef cluster() const
T * make(const Args &...args) const
make new ROOT object
#define end
Definition: vmac.h:39
std::vector< double > onTrackAngles(edm::Handle< edmNew::DetSetVector< SiStripCluster > > &, const std::vector< Trajectory > &)
bool isValid() const
Definition: HandleBase.h:70
bool tecIsRPhi(const DetId &id) const
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
std::vector< edm::InputTag > trackLabels_
double dxdz() const
dxdz slope
Definition: BeamSpot.h:78
HLTConfigProvider hltConfig_
unsigned int tibModule(const DetId &id) const
unsigned int tecModule(const DetId &id) const
bool tecIsFrontPetal(const DetId &id) const
void insertMeasurement(std::multimap< const uint32_t, std::pair< LocalPoint, double > > &, const TransientTrackingRecHit *, double)
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:18
edm::EDGetTokenT< SiStripEventSummary > summaryToken_
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
SubDetector subDetector() const
Definition: SiStripDetId.h:105
Definition: DetId.h:17
std::vector< edm::EDGetTokenT< reco::TrackCollection > > trackTokens_
bool tidIsRPhi(const DetId &id) const
auto vector_transform(std::vector< InputType > const &input, Function predicate) -> std::vector< typename std::remove_cv< typename std::remove_reference< decltype(predicate(input.front()))>::type >::type >
Definition: transform.h:11
edm::EDGetTokenT< edm::ValueMap< reco::DeDxData > > dedx2Token_
JetCorrectorParametersCollection coll
Definition: classes.h:10
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
Definition: TFileService.h:69
TrackerDpgAnalysis(const edm::ParameterSet &)
bool tidIsZPlusSide(const DetId &id) const
T const * product() const
Definition: Handle.h:69
edm::Service< TFileService > fileService
const L1GtFdlWord gtFdlWord(int bxInEventValue) const
get / set FDL word (record) in the GT readout record
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
edm::EDGetTokenT< reco::VertexCollection > vertexToken_
bool init(const edm::Run &iRun, const edm::EventSetup &iSetup, const std::string &processName, bool &changed)
d&#39;tor
ClusterRef cluster() const
Definition: SiPixelRecHit.h:47
double sigmaZ() const
sigma z
Definition: BeamSpot.h:76
bool tidIsDoubleSide(const DetId &id) const
const cms_uint16_t physicsDeclared() const
get/set "physics declared" bit
Definition: L1GtFdlWord.h:169
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
const DecisionWord & gtDecisionWord() const
get/set/print algorithms bits (decision word)
Definition: L1GtFdlWord.h:128
TrackingRecHit::ConstRecHitPointer ConstRecHitPointer
virtual LocalPoint localPosition(const MeasurementPoint &) const =0
std::map< size_t, int > inVertex(const reco::TrackCollection &, const reco::VertexCollection &, uint32_t)
ConnsConstIterRange fedConnections(uint16_t fed_id) const
unsigned int tobModule(const DetId &id) const
bool tecIsZPlusSide(const DetId &id) const
std::string const & process() const
Definition: InputTag.h:40
std::string toStringName(uint32_t, const TrackerTopology *)
edm::EventID id() const
Definition: EventBase.h:59
fixed size matrix
float delay(const SiStripEventSummary &)
edm::EDGetTokenT< edm::ValueMap< reco::DeDxData > > dedx1Token_
#define begin
Definition: vmac.h:32
HLT enums.
Structure Point Contains parameters of Gaussian fits to DMRs.
Definition: DMRtrends.cc:57
std::vector< TTree * > missingHits_
edm::EDGetTokenT< edmNew::DetSetVector< SiPixelCluster > > pixelclusterToken_
std::vector< edm::EDGetTokenT< TrajTrackAssociationCollection > > trajTrackAssoTokens_
T get() const
Definition: EventSetup.h:73
edm::EDGetTokenT< edm::ValueMap< reco::DeDxData > > dedx3Token_
edm::EDGetTokenT< edmNew::DetSetVector< SiStripCluster > > clusterToken_
const TrackerGeomDet * idToDet(DetId) const override
double y0() const
y coordinate
Definition: BeamSpot.h:63
LocalPoint localPosition() const final
void fill_current_val(int idmod, float current_val)
Definition: TrackerMap.cc:3259
unsigned int tecPetalNumber(const DetId &id) const
alpha
zGenParticlesMatch = cms.InputTag(""),
DetId geographicalId() const
trackRef_iterator tracks_begin() const
first iterator over tracks
Definition: Vertex.cc:71
ROOT::Math::LorentzVector< ROOT::Math::PxPyPzE4D< float > > XYZTLorentzVectorF
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:22
const DecisionWordExtended & gtDecisionWordExtended() const
get/set extended algorithms bits (extended decision word)
Definition: L1GtFdlWord.h:153
unsigned int tobRod(const DetId &id) const
T x() const
Definition: PV2DBase.h:43
T x() const
Definition: PV3DBase.h:59
unsigned int tecWheel(const DetId &id) const
T const * product() const
Definition: ESHandle.h:86
const cms_uint32_t orbitNr() const
get/set orbit number
Definition: L1GtFdlWord.h:235
TimeValue_t value() const
Definition: Timestamp.h:45
std::vector< edm::EDGetTokenT< std::vector< Trajectory > > > trajectoryTokens_
edm::EDGetTokenT< reco::VertexCollection > pixelVertexToken_
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:30
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
Definition: Vertex.h:37
std::vector< std::string > hlNames_
unsigned int tobLayer(const DetId &id) const
Definition: Run.h:45
Our base class.
Definition: SiPixelRecHit.h:23
double x0() const
x coordinate
Definition: BeamSpot.h:61
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11