CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
AlignmentMonitorGeneric.cc
Go to the documentation of this file.
3 
7 #include "TObject.h"
8 
9 #include <TString.h>
10 
12  : AlignmentMonitorBase(cfg, iC, "AlignmentMonitorGeneric") {}
13 
15  std::vector<std::string> residNames; // names of residual histograms
16 
17  residNames.push_back("x hit residuals pos track");
18  residNames.push_back("x hit residuals neg track");
19  residNames.push_back("y hit residuals pos track");
20  residNames.push_back("y hit residuals neg track");
21 
22  auto alignableObjectId = AlignableObjectId::commonObjectIdProvider(pTracker(), pMuon());
23 
24  const auto& alignables = pStore()->alignables();
25 
26  unsigned int nAlignable = alignables.size();
27  unsigned int nResidName = residNames.size();
28 
29  for (unsigned int i = 0; i < nAlignable; ++i) {
30  const Alignable* ali = alignables[i];
31 
32  Hist1Ds& hists = m_resHists[ali];
33 
34  hists.resize(nResidName, nullptr);
35 
36  align::ID id = ali->id();
38 
39  for (unsigned int n = 0; n < nResidName; ++n) {
40  const std::string& name = residNames[n];
41 
42  TString histName(name.c_str());
43  histName += Form("_%s_%d", alignableObjectId.idToString(type), id);
44  histName.ReplaceAll(" ", "");
45 
46  TString histTitle(name.c_str());
47  histTitle += Form(" for %s with ID %d (subdet %d)", alignableObjectId.idToString(type), id, DetId(id).subdetId());
48 
49  hists[n] = book1D(std::string("/iterN/") + std::string(name) + std::string("/"),
50  std::string(histName.Data()),
51  std::string(histTitle.Data()),
52  nBin_,
53  -5.,
54  5.);
55  }
56  }
57 
58  m_trkHists.resize(6, nullptr);
59 
60  m_trkHists[0] = book1D("/iterN/", "pt", "track p_{t} (GeV)", nBin_, 0.0, 100.0);
61  m_trkHists[1] = book1D("/iterN/", "eta", "track #eta", nBin_, -3.0, 3.0);
62  m_trkHists[2] = book1D("/iterN/", "phi", "track #phi", nBin_, -M_PI, M_PI);
63  m_trkHists[3] = book1D("/iterN/", "d0", "track d0 (cm)", nBin_, -0.02, 0.02);
64  m_trkHists[4] = book1D("/iterN/", "dz", "track dz (cm)", nBin_, -20.0, 20.0);
65  m_trkHists[5] = book1D("/iterN/", "chi2", "track #chi^{2}/dof", nBin_, 0.0, 20.0);
66 }
67 
69  const edm::EventSetup&,
71  TrajectoryStateCombiner tsoscomb;
72 
73  for (unsigned int t = 0; t < tracks.size(); ++t) {
74  const reco::Track* track = tracks[t].second;
75 
76  float charge = tracks[t].second->charge();
77 
78  const std::vector<TrajectoryMeasurement>& meass = tracks[t].first->measurements();
79 
80  for (unsigned int m = 0; m < meass.size(); ++m) {
81  const TrajectoryMeasurement& meas = meass[m];
82  const TransientTrackingRecHit& hit = *meas.recHit();
83 
84  if (hit.isValid()) {
86 
87  while (ali) {
88  std::map<const Alignable*, Hist1Ds>::iterator h = m_resHists.find(ali);
89  if (h != m_resHists.end()) {
91 
92  align::LocalVector res = tsos.localPosition() - hit.localPosition();
93  LocalError err1 = tsos.localError().positionError();
94  LocalError err2 = hit.localPositionError(); // CPE+APE
95 
96  // subtract APEs from err2 (if existing) from covariance matrix
97  auto det = static_cast<const TrackerGeomDet*>(hit.det());
98  const auto localAPE = det->localAlignmentError();
99  if (localAPE.valid()) {
100  err2 = LocalError(err2.xx() - localAPE.xx(), err2.xy() - localAPE.xy(), err2.yy() - localAPE.yy());
101  }
102 
103  float errX = std::sqrt(err1.xx() + err2.xx());
104  float errY = std::sqrt(err1.yy() + err2.yy());
105 
106  h->second[charge > 0 ? 0 : 1]->Fill(res.x() / errX);
107  h->second[charge > 0 ? 2 : 3]->Fill(res.y() / errY);
108  }
109  ali = ali->mother();
110  }
111  }
112  }
113 
114  m_trkHists[0]->Fill(track->pt());
115  m_trkHists[1]->Fill(track->eta());
116  m_trkHists[2]->Fill(track->phi());
117  m_trkHists[3]->Fill(track->d0());
118  m_trkHists[4]->Fill(track->dz());
119  m_trkHists[5]->Fill(track->normalizedChi2());
120  }
121 }
122 
align::ID id() const
Return the ID of Alignable, i.e. DetId of &#39;first&#39; component GeomDet(Unit).
Definition: Alignable.h:180
float xx() const
Definition: LocalError.h:22
AlignableNavigator * pNavigator()
AlignmentParameterStore * pStore()
tuple cfg
Definition: looper.py:296
ConstRecHitPointer const & recHit() const
double d0() const
dxy parameter in perigee convention (d0 = -dxy)
Definition: TrackBase.h:611
uint16_t *__restrict__ id
std::map< const Alignable *, Hist1Ds > m_resHists
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:593
static AlignableObjectId commonObjectIdProvider(const AlignableObjectId &, const AlignableObjectId &)
uint32_t ID
Definition: Definitions.h:24
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:649
AlignmentMonitorGeneric(const edm::ParameterSet &, edm::ConsumesCollector iC)
AlignableTracker * pTracker()
auto const & tracks
cannot be loose
void event(const edm::Event &, const edm::EventSetup &, const ConstTrajTrackPairCollection &) override
Called for each event (by &quot;run()&quot;): may be reimplemented.
LocalError positionError() const
void book() override
Book or retrieve histograms; MUST be reimplemented.
std::vector< TH1F * > Hist1Ds
float xy() const
Definition: LocalError.h:23
int iEvent
Definition: GenABIO.cc:224
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:652
float yy() const
Definition: LocalError.h:24
const GeomDet * det() const
T sqrt(T t)
Definition: SSEVec.h:19
double pt() const
track transverse momentum
Definition: TrackBase.h:637
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
virtual StructureType alignableObjectId() const =0
Return the alignable type identifier.
const LocalTrajectoryError & localError() const
std::vector< ConstTrajTrackPair > ConstTrajTrackPairCollection
TH1F * book1D(std::string dir, std::string name, std::string title, int nchX, double lowX, double highX)
TrajectoryStateOnSurface const & forwardPredictedState() const
Access to forward predicted state (from fitter or builder)
#define M_PI
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:622
Definition: DetId.h:17
virtual LocalError localPositionError() const =0
bool isValid() const
list hists
Definition: compare.py:318
#define DEFINE_EDM_PLUGIN(factory, type, name)
DetId geographicalId() const
static const unsigned int nBin_
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
virtual LocalPoint localPosition() const =0
AlignableDetOrUnitPtr alignableFromDetId(const DetId &detid)
Returns AlignableDetOrUnitPtr corresponding to given DetId.
Alignable * mother() const
Return pointer to container alignable (if any)
Definition: Alignable.h:91
LocalError const & localAlignmentError() const
Return local alligment error.
const align::Alignables & alignables(void) const
get all alignables
TrajectoryStateOnSurface const & backwardPredictedState() const
Access to backward predicted state (from smoother)