CMS 3D CMS Logo
List of all members | Public Types | Public Member Functions | Protected Member Functions | Protected Attributes
CSCPairResidualsConstraint Class Reference

#include <CSCPairResidualsConstraint.h>

Inheritance diagram for CSCPairResidualsConstraint:
CSCPairConstraint

Public Types

enum  { kModePhiy, kModePhiPos, kModePhiz, kModeRadius }
 

Public Member Functions

bool addTrack (const std::vector< TrajectoryMeasurement > &measurements, const reco::TransientTrack &track, const TrackTransformer *trackTransformer)
 
void configure (CSCOverlapsAlignmentAlgorithm *parent)
 
 CSCPairResidualsConstraint (unsigned int identifier, int i, int j, CSCDetId id_i, CSCDetId id_j)
 
double error () const override
 
CSCDetId id_i () const
 
CSCDetId id_j () const
 
double radius (bool is_i) const
 
void read (std::vector< std::ifstream *> &input, std::vector< std::string > &filenames)
 
void setPropagator (const Propagator *propagator)
 
void setZplane (const CSCGeometry *cscGeometry)
 
bool valid () const override
 
double value () const override
 
void write (std::ofstream &output)
 
 ~CSCPairResidualsConstraint () override
 
- Public Member Functions inherited from CSCPairConstraint
 CSCPairConstraint (int i, int j, double value, double error)
 
virtual int i () const
 
virtual int j () const
 
virtual ~CSCPairConstraint ()
 

Protected Member Functions

void calculatePhi (const TransientTrackingRecHit *hit, double &phi, double &phierr2, bool doRphi=false, bool globalPhi=false)
 
bool dphidzFromTrack (const std::vector< TrajectoryMeasurement > &measurements, const reco::TransientTrack &track, const TrackTransformer *trackTransformer, double &drphidz)
 
bool isFiducial (std::vector< const TransientTrackingRecHit *> &hits, bool is_i)
 

Protected Attributes

double m_averageRadius
 
const CSCGeometry * m_cscGeometry
 
CSCDetId m_id_i
 
CSCDetId m_id_j
 
unsigned int m_identifier
 
double m_iZ
 
double m_iZ1
 
double m_iZ6
 
double m_jZ
 
double m_jZ1
 
double m_jZ6
 
TH1F * m_offsetResiduals
 
CSCOverlapsAlignmentAlgorithm * m_parent
 
const Propagator * m_propagator
 
TH1F * m_radial
 
TH1F * m_slopeResiduals
 
double m_sum1
 
int m_sumN
 
double m_sumx
 
double m_sumxx
 
double m_sumxy
 
double m_sumy
 
double m_sumyy
 
double m_Zplane
 
Plane::PlanePointer m_Zsurface
 
- Protected Attributes inherited from CSCPairConstraint
double m_error
 
int m_i
 
int m_j
 
double m_value
 

Detailed Description

Date
2010/05/06 05:57:33
Revision
1.4
Author
J. Pivarski - Texas A&M University pivar.nosp@m.ski@.nosp@m.physi.nosp@m.cs.t.nosp@m.amu.e.nosp@m.du

Definition at line 27 of file CSCPairResidualsConstraint.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
Enumerator
kModePhiy 
kModePhiPos 
kModePhiz 
kModeRadius 

Definition at line 50 of file CSCPairResidualsConstraint.h.

50 { kModePhiy, kModePhiPos, kModePhiz, kModeRadius };

Constructor & Destructor Documentation

◆ CSCPairResidualsConstraint()

CSCPairResidualsConstraint::CSCPairResidualsConstraint ( unsigned int  identifier,
int  i,
int  j,
CSCDetId  id_i,
CSCDetId  id_j 
)
inline

Definition at line 29 of file CSCPairResidualsConstraint.h.

30  : CSCPairConstraint(i, j, 0., 0.),
31  m_identifier(identifier),
32  m_id_i(id_i),
33  m_id_j(id_j),
34  m_sum1(0.),
35  m_sumx(0.),
36  m_sumy(0.),
37  m_sumxx(0.),
38  m_sumyy(0.),
39  m_sumxy(0.),
40  m_sumN(0),
41  m_Zplane(1000.),
42  m_iZ1(1000.),
43  m_iZ6(1000.),
44  m_jZ1(1000.),
45  m_jZ6(1000.),
46  m_cscGeometry(nullptr),
47  m_propagator(nullptr){};
CSCPairConstraint::CSCPairConstraint
CSCPairConstraint(int i, int j, double value, double error)
Definition: CSCPairConstraint.h:16
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:19

◆ ~CSCPairResidualsConstraint()

CSCPairResidualsConstraint::~CSCPairResidualsConstraint ( )
inlineoverride

Definition at line 48 of file CSCPairResidualsConstraint.h.

48 {};

Member Function Documentation

◆ addTrack()

bool CSCPairResidualsConstraint::addTrack ( const std::vector< TrajectoryMeasurement > &  measurements,
const reco::TransientTrack &  track,
const TrackTransformer *  trackTransformer 
)

Definition at line 89 of file CSCPairResidualsConstraint.cc.

References cms::cuda::assert(), calculatePhi(), CSCDetId::chamber(), MuonSubdetId::CSC, dphidzFromTrack(), CSCDetId::endcap(), isFiducial(), kModePhiPos, kModePhiy, kModePhiz, kModeRadius, m_averageRadius, CSCOverlapsAlignmentAlgorithm::m_combineME11, CSCOverlapsAlignmentAlgorithm::m_drdz, CSCOverlapsAlignmentAlgorithm::m_fiducial, CSCOverlapsAlignmentAlgorithm::m_hitsPerChamber, m_id_i, m_id_j, m_iZ, m_jZ, CSCOverlapsAlignmentAlgorithm::m_makeHistograms, CSCOverlapsAlignmentAlgorithm::m_maxdrdz, CSCOverlapsAlignmentAlgorithm::m_minHitsPerChamber, CSCOverlapsAlignmentAlgorithm::m_mode, CSCOverlapsAlignmentAlgorithm::m_occupancy, m_offsetResiduals, CSCOverlapsAlignmentAlgorithm::m_offsetResiduals, CSCOverlapsAlignmentAlgorithm::m_offsetResiduals_normalized, CSCOverlapsAlignmentAlgorithm::m_offsetResiduals_weighted, m_parent, m_radial, CSCOverlapsAlignmentAlgorithm::m_slope, CSCOverlapsAlignmentAlgorithm::m_slope_MEm1, CSCOverlapsAlignmentAlgorithm::m_slope_MEm2, CSCOverlapsAlignmentAlgorithm::m_slope_MEm3, CSCOverlapsAlignmentAlgorithm::m_slope_MEm4, CSCOverlapsAlignmentAlgorithm::m_slope_MEp1, CSCOverlapsAlignmentAlgorithm::m_slope_MEp2, CSCOverlapsAlignmentAlgorithm::m_slope_MEp3, CSCOverlapsAlignmentAlgorithm::m_slope_MEp4, CSCOverlapsAlignmentAlgorithm::m_slopeFromTrackRefit, m_slopeResiduals, CSCOverlapsAlignmentAlgorithm::m_slopeResiduals, CSCOverlapsAlignmentAlgorithm::m_slopeResiduals_normalized, CSCOverlapsAlignmentAlgorithm::m_slopeResiduals_weighted, m_sum1, m_sumN, m_sumx, m_sumxx, m_sumxy, m_sumy, m_sumyy, CSCOverlapsAlignmentAlgorithm::m_truncateOffsetResid, CSCOverlapsAlignmentAlgorithm::m_truncateSlopeResid, CSCOverlapsAlignmentAlgorithm::m_useHitWeights, CSCOverlapsAlignmentAlgorithm::m_useTrackWeights, m_Zplane, DetId::Muon, phi, funct::pow(), MjjMonitor_cfi::quantity, alignCSCRings::r, CSCDetId::ring(), mathSSE::sqrt(), CSCDetId::station(), HLT_2022v15_cff::track, mps_merge::weight, and z.

91  {
92  std::vector<const TransientTrackingRecHit *> hits_i;
93  std::vector<const TransientTrackingRecHit *> hits_j;
94 
95  for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements.begin();
96  measurement != measurements.end();
97  ++measurement) {
98  const TransientTrackingRecHit *hit = &*(measurement->recHit());
99 
100  DetId id = hit->geographicalId();
101  if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
102  CSCDetId cscid(id.rawId());
103  CSCDetId chamberId(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber(), 0);
104  if (m_parent->m_combineME11 && cscid.station() == 1 && cscid.ring() == 4)
105  chamberId = CSCDetId(cscid.endcap(), 1, 1, cscid.chamber(), 0);
106 
107  if (chamberId == m_id_i)
108  hits_i.push_back(hit);
109  if (chamberId == m_id_j)
110  hits_j.push_back(hit);
111  }
112  }
113 
114  if (m_parent->m_makeHistograms) {
115  m_parent->m_hitsPerChamber->Fill(hits_i.size());
116  m_parent->m_hitsPerChamber->Fill(hits_j.size());
117  }
118 
119  // require minimum number of hits (if the requirement is too low (~2), some NANs might result...)
120  if (int(hits_i.size()) < m_parent->m_minHitsPerChamber || int(hits_j.size()) < m_parent->m_minHitsPerChamber)
121  return false;
122 
123  // maybe require segments to be fiducial
124  if (m_parent->m_fiducial && !(isFiducial(hits_i, true) && isFiducial(hits_j, false)))
125  return false;
126 
127  double intercept_i = 0.;
128  double interceptError2_i = 0.;
129  double slope_i = 0.;
130  double slopeError2_i = 0.;
131  double intercept_j = 0.;
132  double interceptError2_j = 0.;
133  double slope_j = 0.;
134  double slopeError2_j = 0.;
135 
136  // if slopeFromTrackRefit, then you'll need to refit the whole track without this station's hits;
137  // need at least two other stations for that to be reliable
138  if (m_parent->m_slopeFromTrackRefit) {
139  double dphidz;
140  if (dphidzFromTrack(measurements, track, trackTransformer, dphidz)) {
141  double sum1_i = 0.;
142  double sumy_i = 0.;
143  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_i.begin(); hit != hits_i.end();
144  ++hit) {
145  double phi, phierr2;
146  calculatePhi(*hit, phi, phierr2, false, true);
147  double z = (*hit)->globalPosition().z() - m_Zplane;
148 
149  double weight = 1.;
150  if (m_parent->m_useHitWeights)
151  weight = 1. / phierr2;
152  sum1_i += weight;
153  sumy_i += weight * (phi - z * dphidz);
154  }
155 
156  double sum1_j = 0.;
157  double sumy_j = 0.;
158  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_j.begin(); hit != hits_j.end();
159  ++hit) {
160  double phi, phierr2;
161  calculatePhi(*hit, phi, phierr2, false, true);
162  double z = (*hit)->globalPosition().z() - m_Zplane;
163 
164  double weight = 1.;
165  if (m_parent->m_useHitWeights)
166  weight = 1. / phierr2;
167  sum1_j += weight;
168  sumy_j += weight * (phi - z * dphidz);
169  }
170 
171  if (sum1_i != 0. && sum1_j != 0.) {
172  slope_i = slope_j = dphidz;
173 
174  intercept_i = sumy_i / sum1_i;
175  interceptError2_i = 1. / sum1_i;
176 
177  intercept_j = sumy_j / sum1_j;
178  interceptError2_j = 1. / sum1_j;
179  } else
180  return false;
181  }
182  }
183 
184  else { // not slopeFromTrackRefit
185  double sum1_i = 0.;
186  double sumx_i = 0.;
187  double sumy_i = 0.;
188  double sumxx_i = 0.;
189  double sumxy_i = 0.;
190  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_i.begin(); hit != hits_i.end();
191  ++hit) {
192  double phi, phierr2;
193  calculatePhi(*hit, phi, phierr2, false, true);
194  double z = (*hit)->globalPosition().z() - m_Zplane;
195 
196  double weight = 1.;
197  if (m_parent->m_useHitWeights)
198  weight = 1. / phierr2;
199  sum1_i += weight;
200  sumx_i += weight * z;
201  sumy_i += weight * phi;
202  sumxx_i += weight * z * z;
203  sumxy_i += weight * z * phi;
204  }
205 
206  double sum1_j = 0.;
207  double sumx_j = 0.;
208  double sumy_j = 0.;
209  double sumxx_j = 0.;
210  double sumxy_j = 0.;
211  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_j.begin(); hit != hits_j.end();
212  ++hit) {
213  double phi, phierr2;
214  calculatePhi(*hit, phi, phierr2, false, true);
215  double z = (*hit)->globalPosition().z() - m_Zplane;
216 
217  double weight = 1.;
218  if (m_parent->m_useHitWeights)
219  weight = 1. / phierr2;
220  sum1_j += weight;
221  sumx_j += weight * z;
222  sumy_j += weight * phi;
223  sumxx_j += weight * z * z;
224  sumxy_j += weight * z * phi;
225  }
226 
227  double delta_i = (sum1_i * sumxx_i) - (sumx_i * sumx_i);
228  double delta_j = (sum1_j * sumxx_j) - (sumx_j * sumx_j);
229  if (delta_i != 0. && delta_j != 0.) {
230  intercept_i = ((sumxx_i * sumy_i) - (sumx_i * sumxy_i)) / delta_i;
231  interceptError2_i = sumxx_i / delta_i;
232  slope_i = ((sum1_i * sumxy_i) - (sumx_i * sumy_i)) / delta_i;
233  slopeError2_i = sum1_i / delta_i;
234 
235  intercept_j = ((sumxx_j * sumy_j) - (sumx_j * sumxy_j)) / delta_j;
236  interceptError2_j = sumxx_j / delta_j;
237  slope_j = ((sum1_j * sumxy_j) - (sumx_j * sumy_j)) / delta_j;
238  slopeError2_j = sum1_j / delta_j;
239  } else
240  return false;
241  }
242 
243  // from hits on the two chambers, determine radial_intercepts separately and radial_slope together
244  double sum1_ri = 0.;
245  double sumx_ri = 0.;
246  double sumy_ri = 0.;
247  double sumxx_ri = 0.;
248  double sumxy_ri = 0.;
249  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_i.begin(); hit != hits_i.end(); ++hit) {
250  double r = (*hit)->globalPosition().perp();
251  double z = (*hit)->globalPosition().z() - m_Zplane;
252  sum1_ri += 1.;
253  sumx_ri += z;
254  sumy_ri += r;
255  sumxx_ri += z * z;
256  sumxy_ri += z * r;
257  }
258  double radial_delta_i = (sum1_ri * sumxx_ri) - (sumx_ri * sumx_ri);
259  if (radial_delta_i == 0.)
260  return false;
261  double radial_slope_i = ((sum1_ri * sumxy_ri) - (sumx_ri * sumy_ri)) / radial_delta_i;
262  double radial_intercept_i =
263  ((sumxx_ri * sumy_ri) - (sumx_ri * sumxy_ri)) / radial_delta_i + radial_slope_i * (m_iZ - m_Zplane);
264 
265  double sum1_rj = 0.;
266  double sumx_rj = 0.;
267  double sumy_rj = 0.;
268  double sumxx_rj = 0.;
269  double sumxy_rj = 0.;
270  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits_j.begin(); hit != hits_j.end(); ++hit) {
271  double r = (*hit)->globalPosition().perp();
272  double z = (*hit)->globalPosition().z() - m_Zplane;
273  sum1_rj += 1.;
274  sumx_rj += z;
275  sumy_rj += r;
276  sumxx_rj += z * z;
277  sumxy_rj += z * r;
278  }
279  double radial_delta_j = (sum1_rj * sumxx_rj) - (sumx_rj * sumx_rj);
280  if (radial_delta_j == 0.)
281  return false;
282  double radial_slope_j = ((sum1_rj * sumxy_rj) - (sumx_rj * sumy_rj)) / radial_delta_j;
283  double radial_intercept_j =
284  ((sumxx_rj * sumy_rj) - (sumx_rj * sumxy_rj)) / radial_delta_j + radial_slope_j * (m_jZ - m_Zplane);
285 
286  double radial_delta = ((sum1_ri + sum1_rj) * (sumxx_ri + sumxx_rj)) - ((sumx_ri + sumx_rj) * (sumx_ri + sumx_rj));
287  if (radial_delta == 0.)
288  return false;
289  double radial_intercept =
290  (((sumxx_ri + sumxx_rj) * (sumy_ri + sumy_rj)) - ((sumx_ri + sumx_rj) * (sumxy_ri + sumxy_rj))) / radial_delta;
291  double radial_slope =
292  (((sum1_ri + sum1_rj) * (sumxy_ri + sumxy_rj)) - ((sumx_ri + sumx_rj) * (sumy_ri + sumy_rj))) / radial_delta;
293 
294  if (m_parent->m_makeHistograms) {
295  m_parent->m_drdz->Fill(radial_slope);
296  }
297  if (m_parent->m_maxdrdz > 0. && fabs(radial_slope) > m_parent->m_maxdrdz)
298  return false;
299 
300  double quantity = 0.;
301  double quantityError2 = 0.;
302  if (m_parent->m_mode == kModePhiy) { // phiy comes from track d(rphi)/dz
303  quantity = (slope_i * radial_intercept_i) - (slope_j * radial_intercept_j);
304  quantityError2 = (slopeError2_i)*pow(radial_intercept_i, 2) + (slopeError2_j)*pow(radial_intercept_j, 2);
305  } else if (m_parent->m_mode == kModePhiPos || m_parent->m_mode == kModeRadius) { // phipos comes from phi intercepts
306  quantity = intercept_i - intercept_j;
307  quantityError2 = interceptError2_i + interceptError2_j;
308  } else if (m_parent->m_mode == kModePhiz) { // phiz comes from the slope of rphi intercepts
309  quantity = (intercept_i - intercept_j) * radial_intercept;
310  quantityError2 = (interceptError2_i + interceptError2_j) * pow(radial_intercept, 2);
311  } else
312  assert(false);
313 
314  if (quantityError2 == 0.)
315  return false;
316 
317  double slopeResid = ((slope_i * radial_intercept_i) - (slope_j * radial_intercept_j)) * 1000.;
318  double slopeResidError2 =
319  ((slopeError2_i)*pow(radial_intercept_i, 2) + (slopeError2_j)*pow(radial_intercept_j, 2)) * 1000. * 1000.;
320  double offsetResid = (intercept_i - intercept_j) * radial_intercept * 10.;
321  double offsetResidError2 = (interceptError2_i + interceptError2_j) * pow(radial_intercept, 2) * 10. * 10.;
322 
323  if (m_parent->m_truncateSlopeResid > 0. && fabs(slopeResid) > m_parent->m_truncateSlopeResid)
324  return false;
325  if (m_parent->m_truncateOffsetResid > 0. && fabs(offsetResid) > m_parent->m_truncateOffsetResid)
326  return false;
327 
328  double weight = 1.;
329  if (m_parent->m_useTrackWeights)
330  weight = 1. / quantityError2;
331 
332  // fill the running sums for this CSCPairResidualsConstraint
333  m_sumN += 1;
334  m_sum1 += weight;
335  m_sumx += weight * (radial_intercept - m_averageRadius);
336  m_sumy += weight * quantity;
337  m_sumxx += weight * pow(radial_intercept - m_averageRadius, 2);
338  m_sumyy += weight * quantity * quantity;
339  m_sumxy += weight * (radial_intercept - m_averageRadius) * quantity;
340 
341  if (m_parent->m_makeHistograms) {
342  double rphi_slope_i = slope_i * radial_intercept_i;
343  double rphi_slope_j = slope_j * radial_intercept_j;
344 
345  if (m_parent->m_slopeFromTrackRefit) {
346  m_parent->m_slope->Fill(rphi_slope_i); // == rphi_slope_j
347 
348  if (m_id_i.endcap() == 1 && m_id_i.station() == 4)
349  m_parent->m_slope_MEp4->Fill(rphi_slope_i);
350  if (m_id_i.endcap() == 1 && m_id_i.station() == 3)
351  m_parent->m_slope_MEp3->Fill(rphi_slope_i);
352  if (m_id_i.endcap() == 1 && m_id_i.station() == 2)
353  m_parent->m_slope_MEp2->Fill(rphi_slope_i);
354  if (m_id_i.endcap() == 1 && m_id_i.station() == 1)
355  m_parent->m_slope_MEp1->Fill(rphi_slope_i);
356  if (m_id_i.endcap() == 2 && m_id_i.station() == 1)
357  m_parent->m_slope_MEm1->Fill(rphi_slope_i);
358  if (m_id_i.endcap() == 2 && m_id_i.station() == 2)
359  m_parent->m_slope_MEm2->Fill(rphi_slope_i);
360  if (m_id_i.endcap() == 2 && m_id_i.station() == 3)
361  m_parent->m_slope_MEm3->Fill(rphi_slope_i);
362  if (m_id_i.endcap() == 2 && m_id_i.station() == 4)
363  m_parent->m_slope_MEm4->Fill(rphi_slope_i);
364  } else {
365  m_parent->m_slope->Fill(rphi_slope_i);
366  m_parent->m_slope->Fill(rphi_slope_j);
367 
368  if (m_id_i.endcap() == 1 && m_id_i.station() == 4) {
369  m_parent->m_slope_MEp4->Fill(rphi_slope_i);
370  m_parent->m_slope_MEp4->Fill(rphi_slope_j);
371  }
372  if (m_id_i.endcap() == 1 && m_id_i.station() == 3) {
373  m_parent->m_slope_MEp3->Fill(rphi_slope_i);
374  m_parent->m_slope_MEp3->Fill(rphi_slope_j);
375  }
376  if (m_id_i.endcap() == 1 && m_id_i.station() == 2) {
377  m_parent->m_slope_MEp2->Fill(rphi_slope_i);
378  m_parent->m_slope_MEp2->Fill(rphi_slope_j);
379  }
380  if (m_id_i.endcap() == 1 && m_id_i.station() == 1) {
381  m_parent->m_slope_MEp1->Fill(rphi_slope_i);
382  m_parent->m_slope_MEp1->Fill(rphi_slope_j);
383  }
384  if (m_id_i.endcap() == 2 && m_id_i.station() == 1) {
385  m_parent->m_slope_MEm1->Fill(rphi_slope_i);
386  m_parent->m_slope_MEm1->Fill(rphi_slope_j);
387  }
388  if (m_id_i.endcap() == 2 && m_id_i.station() == 2) {
389  m_parent->m_slope_MEm2->Fill(rphi_slope_i);
390  m_parent->m_slope_MEm2->Fill(rphi_slope_j);
391  }
392  if (m_id_i.endcap() == 2 && m_id_i.station() == 3) {
393  m_parent->m_slope_MEm3->Fill(rphi_slope_i);
394  m_parent->m_slope_MEm3->Fill(rphi_slope_j);
395  }
396  if (m_id_i.endcap() == 2 && m_id_i.station() == 4) {
397  m_parent->m_slope_MEm4->Fill(rphi_slope_i);
398  m_parent->m_slope_MEm4->Fill(rphi_slope_j);
399  }
400  }
401 
402  m_slopeResiduals->Fill(slopeResid);
403  m_offsetResiduals->Fill(offsetResid);
404  m_radial->Fill(radial_intercept);
405 
406  m_parent->m_slopeResiduals->Fill(slopeResid);
407  m_parent->m_slopeResiduals_weighted->Fill(slopeResid, 1. / slopeResidError2);
408  m_parent->m_slopeResiduals_normalized->Fill(slopeResid / sqrt(slopeResidError2));
409 
410  m_parent->m_offsetResiduals->Fill(offsetResid);
411  m_parent->m_offsetResiduals_weighted->Fill(offsetResid, 1. / offsetResidError2);
412  m_parent->m_offsetResiduals_normalized->Fill(offsetResid / sqrt(offsetResidError2));
413 
414  double ringbin = 0;
415  if (m_id_i.endcap() == 2 && m_id_i.station() == 4 && m_id_i.ring() == 2)
416  ringbin = 1.5;
417  else if (m_id_i.endcap() == 2 && m_id_i.station() == 4 && m_id_i.ring() == 1)
418  ringbin = 2.5;
419  else if (m_id_i.endcap() == 2 && m_id_i.station() == 3 && m_id_i.ring() == 2)
420  ringbin = 3.5;
421  else if (m_id_i.endcap() == 2 && m_id_i.station() == 3 && m_id_i.ring() == 1)
422  ringbin = 4.5;
423  else if (m_id_i.endcap() == 2 && m_id_i.station() == 2 && m_id_i.ring() == 2)
424  ringbin = 5.5;
425  else if (m_id_i.endcap() == 2 && m_id_i.station() == 2 && m_id_i.ring() == 1)
426  ringbin = 6.5;
427  else if (m_id_i.endcap() == 2 && m_id_i.station() == 1 && m_id_i.ring() == 3)
428  ringbin = 7.5;
429  else if (m_id_i.endcap() == 2 && m_id_i.station() == 1 && m_id_i.ring() == 2)
430  ringbin = 8.5;
431  else if (m_id_i.endcap() == 2 && m_id_i.station() == 1 && m_id_i.ring() == 1)
432  ringbin = 9.5;
433  else if (m_id_i.endcap() == 2 && m_id_i.station() == 1 && m_id_i.ring() == 4)
434  ringbin = 10.5;
435  else if (m_id_i.endcap() == 1 && m_id_i.station() == 1 && m_id_i.ring() == 4)
436  ringbin = 11.5;
437  else if (m_id_i.endcap() == 1 && m_id_i.station() == 1 && m_id_i.ring() == 1)
438  ringbin = 12.5;
439  else if (m_id_i.endcap() == 1 && m_id_i.station() == 1 && m_id_i.ring() == 2)
440  ringbin = 13.5;
441  else if (m_id_i.endcap() == 1 && m_id_i.station() == 1 && m_id_i.ring() == 3)
442  ringbin = 14.5;
443  else if (m_id_i.endcap() == 1 && m_id_i.station() == 2 && m_id_i.ring() == 1)
444  ringbin = 15.5;
445  else if (m_id_i.endcap() == 1 && m_id_i.station() == 2 && m_id_i.ring() == 2)
446  ringbin = 16.5;
447  else if (m_id_i.endcap() == 1 && m_id_i.station() == 3 && m_id_i.ring() == 1)
448  ringbin = 17.5;
449  else if (m_id_i.endcap() == 1 && m_id_i.station() == 3 && m_id_i.ring() == 2)
450  ringbin = 18.5;
451  else if (m_id_i.endcap() == 1 && m_id_i.station() == 4 && m_id_i.ring() == 1)
452  ringbin = 19.5;
453  else if (m_id_i.endcap() == 1 && m_id_i.station() == 4 && m_id_i.ring() == 2)
454  ringbin = 20.5;
455  m_parent->m_occupancy->Fill(m_id_i.chamber() + 0.5, ringbin);
456  }
457 
458  return true;
459 }
CSCPairResidualsConstraint::isFiducial
bool isFiducial(std::vector< const TransientTrackingRecHit *> &hits, bool is_i)
Definition: CSCPairResidualsConstraint.cc:694
weight
Definition: weight.py:1
cms::cuda::assert
assert(be >=bs)
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
CSCPairResidualsConstraint::calculatePhi
void calculatePhi(const TransientTrackingRecHit *hit, double &phi, double &phierr2, bool doRphi=false, bool globalPhi=false)
Definition: CSCPairResidualsConstraint.cc:628
CSCDetId::chamber
int chamber() const
Definition: CSCDetId.h:62
DetId
Definition: DetId.h:17
CSCDetId::station
int station() const
Definition: CSCDetId.h:79
CSCDetId::endcap
int endcap() const
Definition: CSCDetId.h:85
MjjMonitor_cfi.quantity
quantity
Definition: MjjMonitor_cfi.py:10
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85
HLT_2022v15_cff.track
track
Definition: HLT_2022v15_cff.py:9650
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
CSCPairResidualsConstraint::dphidzFromTrack
bool dphidzFromTrack(const std::vector< TrajectoryMeasurement > &measurements, const reco::TransientTrack &track, const TrackTransformer *trackTransformer, double &drphidz)
Definition: CSCPairResidualsConstraint.cc:461
mps_merge.weight
weight
Definition: mps_merge.py:88
MuonSubdetId::CSC
static constexpr int CSC
Definition: MuonSubdetId.h:12
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ calculatePhi()

void CSCPairResidualsConstraint::calculatePhi ( const TransientTrackingRecHit *  hit,
double &  phi,
double &  phierr2,
bool  doRphi = false,
bool  globalPhi = false 
)
protected

Definition at line 628 of file CSCPairResidualsConstraint.cc.

References angle(), cms::cuda::assert(), funct::cos(), CSCDetId::endcap(), CSCLayer::geometry(), CSCGeometry::layer(), m_cscGeometry, M_PI, CSCLayerGeometry::nearestStrip(), phi, dttmaxenums::R, alignCSCRings::r, CSCDetId::ring(), funct::sin(), CSCDetId::station(), digitizers_cfi::strip, CSCLayerGeometry::stripAngle(), geometryCSVtoXML::xx, geometryCSVtoXML::xy, and geometryCSVtoXML::yy.

Referenced by addTrack(), and isFiducial().

629  {
630  align::LocalPoint pos = hit->localPosition();
631  DetId id = hit->geographicalId();
632  CSCDetId cscid = CSCDetId(id.rawId());
633 
634  double r = 0.;
635  if (globalPhi) {
636  phi = hit->globalPosition().phi();
637  r = hit->globalPosition().perp();
638 
639  // double sinAngle = sin(phi);
640  // double cosAngle = cos(phi);
641  // double xx = hit->globalPositionError().cxx();
642  // double xy = hit->globalPositionError().cyx();
643  // double yy = hit->globalPositionError().cyy();
644  // phierr2 = (xx*cosAngle*cosAngle + 2.*xy*sinAngle*cosAngle + yy*sinAngle*sinAngle) / (r*r);
645  } else {
646  // these constants are related to the way CSC chambers are built--- really constant!
647  const double R_ME11 = 181.5;
648  const double R_ME12 = 369.7;
649  const double R_ME21 = 242.7;
650  const double R_ME31 = 252.7;
651  const double R_ME41 = 262.65;
652  const double R_MEx2 = 526.5;
653 
654  double R = 0.;
655  if (cscid.station() == 1 && (cscid.ring() == 1 || cscid.ring() == 4))
656  R = R_ME11;
657  else if (cscid.station() == 1 && cscid.ring() == 2)
658  R = R_ME12;
659  else if (cscid.station() == 2 && cscid.ring() == 1)
660  R = R_ME21;
661  else if (cscid.station() == 3 && cscid.ring() == 1)
662  R = R_ME31;
663  else if (cscid.station() == 4 && cscid.ring() == 1)
664  R = R_ME41;
665  else if (cscid.station() > 1 && cscid.ring() == 2)
666  R = R_MEx2;
667  else
668  assert(false);
669  r = (pos.y() + R);
670 
671  phi = atan2(pos.x(), r);
672 
673  if (cscid.endcap() == 1 && cscid.station() >= 3)
674  phi *= -1;
675  else if (cscid.endcap() == 2 && cscid.station() <= 2)
676  phi *= -1;
677  }
678 
679  int strip = m_cscGeometry->layer(id)->geometry()->nearestStrip(pos);
680  double angle = m_cscGeometry->layer(id)->geometry()->stripAngle(strip) - M_PI / 2.;
681  double sinAngle = sin(angle);
682  double cosAngle = cos(angle);
683  double xx = hit->localPositionError().xx();
684  double xy = hit->localPositionError().xy();
685  double yy = hit->localPositionError().yy();
686  phierr2 = (xx * cosAngle * cosAngle + 2. * xy * sinAngle * cosAngle + yy * sinAngle * sinAngle) / (r * r);
687 
688  if (doRphi) {
689  phi *= r;
690  phierr2 *= r * r;
691  }
692 }
CSCLayerGeometry::nearestStrip
int nearestStrip(const LocalPoint &lp) const
Definition: CSCLayerGeometry.h:96
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
cms::cuda::assert
assert(be >=bs)
CSCLayerGeometry::stripAngle
float stripAngle(int strip) const
Definition: CSCLayerGeometry.cc:147
CSCLayer::geometry
const CSCLayerGeometry * geometry() const
Definition: CSCLayer.h:44
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:50
DetId
Definition: DetId.h:17
CSCDetId::station
int station() const
Definition: CSCDetId.h:79
digitizers_cfi.strip
strip
Definition: digitizers_cfi.py:19
CSCDetId::endcap
int endcap() const
Definition: CSCDetId.h:85
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
CSCGeometry::layer
const CSCLayer * layer(CSCDetId id) const
Return the layer corresponding to given DetId.
Definition: CSCGeometry.cc:105
angle
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11

◆ configure()

void CSCPairResidualsConstraint::configure ( CSCOverlapsAlignmentAlgorithm *  parent)

Definition at line 38 of file CSCPairResidualsConstraint.cc.

References m_id_i, m_id_j, m_identifier, CSCOverlapsAlignmentAlgorithm::m_makeHistograms, m_offsetResiduals, m_parent, m_radial, m_slopeResiduals, TFileService::make(), Skims_PA_cff::name, hfnoseParametersInitialization_cfi::name2, class-composition::parent, and runGCPTkAlMap::title.

38  {
39  m_parent = parent;
40 
41  if (m_parent->m_makeHistograms) {
42  edm::Service<TFileService> tFileService;
43 
44  std::stringstream name, name2, name3, title;
45  title << "i =" << m_id_i << " j =" << m_id_j;
46 
47  name << "slopeResiduals_" << m_identifier;
48  m_slopeResiduals = tFileService->make<TH1F>(name.str().c_str(), title.str().c_str(), 300, -30., 30.);
49 
50  name2 << "offsetResiduals_" << m_identifier;
51  m_offsetResiduals = tFileService->make<TH1F>(name2.str().c_str(), title.str().c_str(), 300, -30., 30.);
52 
53  name3 << "radial_" << m_identifier;
54  m_radial = tFileService->make<TH1F>(name3.str().c_str(), title.str().c_str(), 700, 0., 700.);
55  } else {
56  m_slopeResiduals = nullptr;
57  m_offsetResiduals = nullptr;
58  m_radial = nullptr;
59  }
60 }
hfnoseParametersInitialization_cfi.name2
name2
Definition: hfnoseParametersInitialization_cfi.py:7
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85
TFileService::make
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
Skims_PA_cff.name
name
Definition: Skims_PA_cff.py:17

◆ dphidzFromTrack()

bool CSCPairResidualsConstraint::dphidzFromTrack ( const std::vector< TrajectoryMeasurement > &  measurements,
const reco::TransientTrack &  track,
const TrackTransformer *  trackTransformer,
double &  drphidz 
)
protected

Definition at line 461 of file CSCPairResidualsConstraint.cc.

References angle(), CSCDetId::chamber(), TrajectoryStateCombiner::combine(), funct::cos(), MuonSubdetId::CSC, dxdz, dydz, CSCDetId::endcap(), TrajectoryStateOnSurface::globalDirection(), TrajectoryStateOnSurface::globalPosition(), TrajectoryStateOnSurface::isValid(), CSCOverlapsAlignmentAlgorithm::m_combineME11, m_id_i, m_id_j, CSCOverlapsAlignmentAlgorithm::m_minHitsPerChamber, CSCOverlapsAlignmentAlgorithm::m_minStationsInTrackRefits, m_parent, M_PI, m_propagator, m_Zplane, m_Zsurface, DetId::Muon, PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), Propagator::propagate(), CSCDetId::ring(), funct::sin(), relativeConstraints::station, CSCDetId::station(), dqmMemoryStats::total, HLT_2022v15_cff::track, FastTrackerRecHitMaskProducer_cfi::trajectories, TrackTransformer::transform(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by addTrack().

464  {
465  // make a list of hits on all chambers *other* than the ones associated with this constraint
466  std::map<int, int> stations;
467  int total = 0;
468  TransientTrackingRecHit::ConstRecHitContainer cscHits;
469  for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements.begin();
470  measurement != measurements.end();
471  ++measurement) {
472  DetId id = measurement->recHit()->geographicalId();
473  if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
474  CSCDetId cscid(id.rawId());
475  CSCDetId chamberId(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber(), 0);
476  if (m_parent->m_combineME11 && cscid.station() == 1 && cscid.ring() == 4)
477  chamberId = CSCDetId(cscid.endcap(), 1, 1, cscid.chamber(), 0);
478 
479  if (chamberId != m_id_i && chamberId != m_id_j) {
480  int station = (cscid.endcap() == 1 ? 1 : -1) * cscid.station();
481  if (stations.find(station) == stations.end()) {
482  stations[station] = 0;
483  }
484  stations[station]++;
485  total++;
486 
487  cscHits.push_back(measurement->recHit());
488  }
489  }
490  }
491 
492  // for the fit to be reliable, it needs to cross multiple stations
493  int numStations = 0;
494  for (std::map<int, int>::const_iterator station = stations.begin(); station != stations.end(); ++station) {
495  if (station->second >= m_parent->m_minHitsPerChamber) {
496  numStations++;
497  }
498  }
499 
500  if (numStations >= m_parent->m_minStationsInTrackRefits) {
501  // refit the track with these hits
502  std::vector<Trajectory> trajectories = trackTransformer->transform(track, cscHits);
503 
504  if (!trajectories.empty()) {
505  const std::vector<TrajectoryMeasurement> &measurements2 = trajectories.begin()->measurements();
506 
507  // find the closest TSOS to the Z plane (on both sides)
508  bool found_plus = false;
509  bool found_minus = false;
510  TrajectoryStateOnSurface tsos_plus, tsos_minus;
511  for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements2.begin();
512  measurement != measurements2.end();
513  ++measurement) {
514  double z = measurement->recHit()->globalPosition().z();
515  if (z > m_Zplane) {
516  if (!found_plus || fabs(z - m_Zplane) < fabs(tsos_plus.globalPosition().z() - m_Zplane)) {
517  tsos_plus = TrajectoryStateCombiner().combine(measurement->forwardPredictedState(),
518  measurement->backwardPredictedState());
519  }
520  if (tsos_plus.isValid())
521  found_plus = true;
522  } else {
523  if (!found_minus || fabs(z - m_Zplane) < fabs(tsos_minus.globalPosition().z() - m_Zplane)) {
524  tsos_minus = TrajectoryStateCombiner().combine(measurement->forwardPredictedState(),
525  measurement->backwardPredictedState());
526  }
527  if (tsos_minus.isValid())
528  found_minus = true;
529  }
530  }
531 
532  // propagate from the closest TSOS to the Z plane (from both sides, if possible)
533  TrajectoryStateOnSurface from_plus, from_minus;
534  if (found_plus) {
535  from_plus = m_propagator->propagate(tsos_plus, *m_Zsurface);
536  }
537  if (found_minus) {
538  from_minus = m_propagator->propagate(tsos_minus, *m_Zsurface);
539  }
540 
541  // if you have two sides, merge them
542  TrajectoryStateOnSurface merged;
543  if (found_plus && from_plus.isValid() && found_minus && from_minus.isValid()) {
544  merged = TrajectoryStateCombiner().combine(from_plus, from_minus);
545  } else if (found_plus && from_plus.isValid()) {
546  merged = from_plus;
547  } else if (found_minus && from_minus.isValid()) {
548  merged = from_minus;
549  } else
550  return false;
551 
552  // if, after all that, we have a good fit-and-propagation, report the direction
553  if (merged.isValid()) {
554  double angle = merged.globalPosition().phi() + M_PI / 2.;
555  GlobalVector direction = merged.globalDirection();
556  double dxdz = direction.x() / direction.z();
557  double dydz = direction.y() / direction.z();
558  dphidz = (dxdz * cos(angle) + dydz * sin(angle)) / merged.globalPosition().perp();
559  return true;
560  }
561 
562  } // end if refit successful
563  } // end if enough hits
564  return false;
565 }
dydz
float dydz
Definition: beamSpotDipStandalone.cc:59
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
dxdz
float dxdz
Definition: beamSpotDipStandalone.cc:58
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Propagator::propagate
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:50
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:65
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:50
TrackingRecHit::ConstRecHitContainer
std::vector< ConstRecHitPointer > ConstRecHitContainer
Definition: TrackingRecHit.h:32
FastTrackerRecHitMaskProducer_cfi.trajectories
trajectories
Definition: FastTrackerRecHitMaskProducer_cfi.py:7
DetId
Definition: DetId.h:17
TrajectoryStateOnSurface::globalDirection
GlobalVector globalDirection() const
Definition: TrajectoryStateOnSurface.h:67
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85
HLT_2022v15_cff.track
track
Definition: HLT_2022v15_cff.py:9650
MuonSubdetId::CSC
static constexpr int CSC
Definition: MuonSubdetId.h:12
TrackTransformer::transform
std::vector< Trajectory > transform(const reco::Track &) const override
Convert a reco::Track into Trajectory.
Definition: TrackTransformer.cc:154
TrajectoryStateCombiner::combine
TSOS combine(const TSOS &pTsos1, const TSOS &pTsos2) const
Definition: TrajectoryStateCombiner.cc:6
angle
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11

◆ error()

double CSCPairResidualsConstraint::error ( ) const
overridevirtual

Reimplemented from CSCPairConstraint.

Definition at line 20 of file CSCPairResidualsConstraint.cc.

References cms::cuda::assert(), dumpMFGeometry_cfg::delta, kModePhiPos, kModePhiy, kModePhiz, kModeRadius, CSCOverlapsAlignmentAlgorithm::m_errorFromRMS, CSCOverlapsAlignmentAlgorithm::m_mode, m_parent, m_sum1, m_sumN, m_sumx, m_sumxx, m_sumy, m_sumyy, funct::pow(), and mathSSE::sqrt().

Referenced by Page1Parser.Page1Parser::check_for_whole_start_tag().

20  {
21  if (m_parent->m_errorFromRMS) {
22  assert(m_sum1 > 0.);
23  return sqrt((m_sumyy / m_sum1) - pow(m_sumy / m_sum1, 2)) / sqrt(m_sumN);
24  } else {
25  double delta = (m_sum1 * m_sumxx) - (m_sumx * m_sumx);
26  assert(delta > 0.);
27  if (m_parent->m_mode == kModePhiy || m_parent->m_mode == kModePhiPos || m_parent->m_mode == kModeRadius) {
28  return sqrt(m_sumxx / delta);
29  } else if (m_parent->m_mode == kModePhiz) {
30  return sqrt(m_sum1 / delta);
31  } else
32  assert(false);
33  }
34 }
cms::cuda::assert
assert(be >=bs)
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ id_i()

CSCDetId CSCPairResidualsConstraint::id_i ( ) const
inline

Definition at line 54 of file CSCPairResidualsConstraint.h.

References m_id_i.

Referenced by CSCChamberFitter::radiusCorrection().

54 { return m_id_i; };

◆ id_j()

CSCDetId CSCPairResidualsConstraint::id_j ( ) const
inline

Definition at line 55 of file CSCPairResidualsConstraint.h.

References m_id_j.

55 { return m_id_j; };

◆ isFiducial()

bool CSCPairResidualsConstraint::isFiducial ( std::vector< const TransientTrackingRecHit *> &  hits,
bool  is_i 
)
protected

Definition at line 694 of file CSCPairResidualsConstraint.cc.

References cms::cuda::assert(), calculatePhi(), dumpMFGeometry_cfg::delta, hfClusterShapes_cfi::hits, CSCGeometry::idToDet(), m_cscGeometry, CSCOverlapsAlignmentAlgorithm::m_fiducial_ME11, CSCOverlapsAlignmentAlgorithm::m_fiducial_ME12, CSCOverlapsAlignmentAlgorithm::m_fiducial_MEx1, CSCOverlapsAlignmentAlgorithm::m_fiducial_MEx2, m_id_i, m_id_j, m_iZ1, m_iZ6, m_jZ1, m_jZ6, CSCOverlapsAlignmentAlgorithm::m_makeHistograms, m_parent, CSCOverlapsAlignmentAlgorithm::m_useHitWeights, phi, CSCDetId::ring(), slope, mathSSE::sqrt(), CSCDetId::station(), GeomDet::surface(), toLocal(), mps_merge::weight, and z.

Referenced by addTrack().

694  {
695  // these constants are related to the way CSC chambers are built--- really constant!
696  const double cut_ME11 = 0.086;
697  const double cut_ME12 = 0.090;
698  const double cut_MEx1 = 0.180;
699  const double cut_MEx2 = 0.090;
700 
701  double sum1 = 0.;
702  double sumx = 0.;
703  double sumy = 0.;
704  double sumxx = 0.;
705  double sumxy = 0.;
706  for (std::vector<const TransientTrackingRecHit *>::const_iterator hit = hits.begin(); hit != hits.end(); ++hit) {
707  double phi, phierr2;
708  calculatePhi(*hit, phi, phierr2);
709  double z = (is_i ? m_cscGeometry->idToDet(m_id_i)->surface() : m_cscGeometry->idToDet(m_id_j)->surface())
710  .toLocal((*hit)->globalPosition())
711  .z();
712 
713  if (m_parent->m_makeHistograms) {
714  if (m_id_i.station() == 1 && (m_id_i.ring() == 1 || m_id_i.ring() == 4)) {
715  m_parent->m_fiducial_ME11->Fill(fabs(phi), sqrt(phierr2));
716  } else if (m_id_i.station() == 1 && m_id_i.ring() == 2) {
717  m_parent->m_fiducial_ME12->Fill(fabs(phi), sqrt(phierr2));
718  } else if (m_id_i.station() > 1 && m_id_i.ring() == 1) {
719  m_parent->m_fiducial_MEx1->Fill(fabs(phi), sqrt(phierr2));
720  } else if (m_id_i.station() > 1 && m_id_i.ring() == 2) {
721  m_parent->m_fiducial_MEx2->Fill(fabs(phi), sqrt(phierr2));
722  }
723  }
724 
725  double weight = 1.;
726  if (m_parent->m_useHitWeights)
727  weight = 1. / phierr2;
728  sum1 += weight;
729  sumx += weight * z;
730  sumy += weight * phi;
731  sumxx += weight * z * z;
732  sumxy += weight * z * phi;
733  }
734  double delta = (sum1 * sumxx) - (sumx * sumx);
735  if (delta == 0.)
736  return false;
737  double intercept = ((sumxx * sumy) - (sumx * sumxy)) / delta;
738  double slope = ((sum1 * sumxy) - (sumx * sumy)) / delta;
739 
740  double phi1 = intercept + slope * (is_i ? m_iZ1 : m_jZ1);
741  double phi6 = intercept + slope * (is_i ? m_iZ6 : m_jZ6);
742 
743  if (m_id_i.station() == 1 && (m_id_i.ring() == 1 || m_id_i.ring() == 4)) {
744  return (fabs(phi1) < cut_ME11 && fabs(phi6) < cut_ME11);
745  } else if (m_id_i.station() == 1 && m_id_i.ring() == 2) {
746  return (fabs(phi1) < cut_ME12 && fabs(phi6) < cut_ME12);
747  } else if (m_id_i.station() > 1 && m_id_i.ring() == 1) {
748  return (fabs(phi1) < cut_MEx1 && fabs(phi6) < cut_MEx1);
749  } else if (m_id_i.station() > 1 && m_id_i.ring() == 2) {
750  return (fabs(phi1) < cut_MEx2 && fabs(phi6) < cut_MEx2);
751  } else
752  assert(false);
753 }
slope
static const double slope[3]
Definition: CastorTimeSlew.cc:6
hfClusterShapes_cfi.hits
hits
Definition: hfClusterShapes_cfi.py:5
weight
Definition: weight.py:1
cms::cuda::assert
assert(be >=bs)
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
CSCPairResidualsConstraint::calculatePhi
void calculatePhi(const TransientTrackingRecHit *hit, double &phi, double &phierr2, bool doRphi=false, bool globalPhi=false)
Definition: CSCPairResidualsConstraint.cc:628
toLocal
LocalVector toLocal(const reco::Track::Vector &v, const Surface &s)
Definition: ConversionProducer.cc:202
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
CSCDetId::station
int station() const
Definition: CSCDetId.h:79
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
mps_merge.weight
weight
Definition: mps_merge.py:88
CSCGeometry::idToDet
const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:91

◆ radius()

double CSCPairResidualsConstraint::radius ( bool  is_i) const
inline

Definition at line 57 of file CSCPairResidualsConstraint.h.

References CSCGeometry::idToDet(), m_cscGeometry, m_id_i, m_id_j, PV3DBase< T, PVType, FrameType >::perp(), GloballyPositioned< T >::position(), and GeomDet::surface().

Referenced by CSCChamberFitter::radiusCorrection().

57  {
58  return m_cscGeometry->idToDet((is_i ? m_id_i : m_id_j))->surface().position().perp();
59  };
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:36
CSCGeometry::idToDet
const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:91

◆ read()

void CSCPairResidualsConstraint::read ( std::vector< std::ifstream *> &  input,
std::vector< std::string > &  filenames 
)

Definition at line 574 of file CSCPairResidualsConstraint.cc.

References Exception, corrVsCorr::filename, CSCPairConstraint::i(), input, CSCPairConstraint::j(), CSCPairConstraint::m_i, m_identifier, CSCPairConstraint::m_j, m_sum1, m_sumN, m_sumx, m_sumxx, m_sumxy, m_sumy, m_sumyy, Skims_PA_cff::name, and AlCaHLTBitMon_QueryRunRegistry::string.

Referenced by edmIntegrityCheck.PublishToFileSystem::get().

574  {
575  m_sumN = 0;
576  m_sum1 = 0.;
577  m_sumx = 0.;
578  m_sumy = 0.;
579  m_sumxx = 0.;
580  m_sumyy = 0.;
581  m_sumxy = 0.;
582 
583  std::vector<std::ifstream *>::const_iterator inputiter = input.begin();
584  std::vector<std::string>::const_iterator filename = filenames.begin();
585  for (; inputiter != input.end(); ++inputiter, ++filename) {
586  int linenumber = 0;
587  bool touched = false;
588  while (!(*inputiter)->eof()) {
589  linenumber++;
590  std::string name, eoln;
591  unsigned int identifier;
592  int i, j;
593  int sumN;
594  double sum1, sumx, sumy, sumxx, sumyy, sumxy;
595 
596  (**inputiter) >> name >> identifier >> i >> j >> sumN >> sum1 >> sumx >> sumy >> sumxx >> sumyy >> sumxy >> eoln;
597 
598  if (!(*inputiter)->eof() && (name != "CSCPairResidualsConstraint" || eoln != "EOLN"))
599  throw cms::Exception("CorruptTempFile")
600  << "Temporary file " << *filename << " is incorrectly formatted on line " << linenumber << std::endl;
601 
602  if (identifier == m_identifier) {
603  if (i != m_i || j != m_j)
604  throw cms::Exception("CorruptTempFile")
605  << "Wrong (i,j) for CSCPairResidualsConstraint " << m_identifier << " (" << m_i << "," << m_j
606  << ") in file " << *filename << " on line " << linenumber << std::endl;
607  touched = true;
608 
609  m_sumN += sumN;
610  m_sum1 += sum1;
611  m_sumx += sumx;
612  m_sumy += sumy;
613  m_sumxx += sumxx;
614  m_sumyy += sumyy;
615  m_sumxy += sumxy;
616  }
617  }
618 
619  (*inputiter)->clear();
620  (*inputiter)->seekg(0, std::ios::beg);
621 
622  if (!touched)
623  throw cms::Exception("CorruptTempFile")
624  << "CSCPairResidualsConstraint " << m_identifier << " is missing from file " << *filename << std::endl;
625  }
626 }
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
input
static std::string const input
Definition: EdmProvDump.cc:50
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:19
Skims_PA_cff.name
name
Definition: Skims_PA_cff.py:17

◆ setPropagator()

void CSCPairResidualsConstraint::setPropagator ( const Propagator *  propagator)

Definition at line 87 of file CSCPairResidualsConstraint.cc.

References m_propagator, and TrackCandidateProducer_cfi::propagator.

87 { m_propagator = propagator; }
TrackCandidateProducer_cfi.propagator
propagator
Definition: TrackCandidateProducer_cfi.py:17

◆ setZplane()

void CSCPairResidualsConstraint::setZplane ( const CSCGeometry *  cscGeometry)

Definition at line 62 of file CSCPairResidualsConstraint.cc.

References Plane::build(), CSCDetId::chamber(), CSCDetId::endcap(), testProducerWithPsetDescEmpty_cfi::i1, CSCGeometry::idToDet(), m_averageRadius, m_cscGeometry, m_id_i, m_id_j, m_iZ, m_iZ1, m_iZ6, m_jZ, m_jZ1, m_jZ6, m_Zplane, m_Zsurface, PV3DBase< T, PVType, FrameType >::perp(), GloballyPositioned< T >::position(), CSCDetId::ring(), CSCDetId::station(), GeomDet::surface(), GloballyPositioned< T >::toLocal(), z, and PV3DBase< T, PVType, FrameType >::z().

62  {
63  m_cscGeometry = cscGeometry;
64 
65  m_Zplane = (m_cscGeometry->idToDet(m_id_i)->surface().position().z() +
66  m_cscGeometry->idToDet(m_id_j)->surface().position().z()) /
67  2.;
68  m_averageRadius = (m_cscGeometry->idToDet(m_id_i)->surface().position().perp() +
69  m_cscGeometry->idToDet(m_id_j)->surface().position().perp()) /
70  2.;
71 
72  m_iZ = m_cscGeometry->idToDet(m_id_i)->surface().position().z();
73  m_jZ = m_cscGeometry->idToDet(m_id_j)->surface().position().z();
74 
75  CSCDetId i1(m_id_i.endcap(), m_id_i.station(), m_id_i.ring(), m_id_i.chamber(), 1);
76  CSCDetId i6(m_id_i.endcap(), m_id_i.station(), m_id_i.ring(), m_id_i.chamber(), 6);
77  CSCDetId j1(m_id_j.endcap(), m_id_j.station(), m_id_j.ring(), m_id_j.chamber(), 1);
78  CSCDetId j6(m_id_j.endcap(), m_id_j.station(), m_id_j.ring(), m_id_j.chamber(), 6);
79  m_iZ1 = m_cscGeometry->idToDet(m_id_i)->surface().toLocal(m_cscGeometry->idToDet(i1)->surface().position()).z();
80  m_iZ6 = m_cscGeometry->idToDet(m_id_i)->surface().toLocal(m_cscGeometry->idToDet(i6)->surface().position()).z();
81  m_jZ1 = m_cscGeometry->idToDet(m_id_j)->surface().toLocal(m_cscGeometry->idToDet(j1)->surface().position()).z();
82  m_jZ6 = m_cscGeometry->idToDet(m_id_j)->surface().toLocal(m_cscGeometry->idToDet(j6)->surface().position()).z();
83 
84  m_Zsurface = Plane::build(Plane::PositionType(0., 0., m_Zplane), Plane::RotationType());
85 }
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
GloballyPositioned::toLocal
LocalPoint toLocal(const GlobalPoint &gp) const
Definition: GloballyPositioned.h:98
Plane::build
static PlanePointer build(Args &&... args)
Definition: Plane.h:33
CSCDetId::chamber
int chamber() const
Definition: CSCDetId.h:62
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:36
CSCDetId::station
int station() const
Definition: CSCDetId.h:79
CSCDetId::endcap
int endcap() const
Definition: CSCDetId.h:85
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
CSCGeometry::idToDet
const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:91

◆ valid()

bool CSCPairResidualsConstraint::valid ( ) const
overridevirtual

Reimplemented from CSCPairConstraint.

Definition at line 36 of file CSCPairResidualsConstraint.cc.

References CSCOverlapsAlignmentAlgorithm::m_minTracksPerOverlap, m_parent, and m_sumN.

Referenced by CSCChamberFitter::radiusCorrection().

36 { return (m_sumN >= m_parent->m_minTracksPerOverlap); }
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85

◆ value()

double CSCPairResidualsConstraint::value ( ) const
overridevirtual

Reimplemented from CSCPairConstraint.

Definition at line 9 of file CSCPairResidualsConstraint.cc.

References cms::cuda::assert(), dumpMFGeometry_cfg::delta, kModePhiPos, kModePhiy, kModePhiz, kModeRadius, CSCOverlapsAlignmentAlgorithm::m_mode, m_parent, m_sum1, m_sumx, m_sumxx, m_sumxy, and m_sumy.

Referenced by Types.int32::__nonzero__(), Types.uint32::__nonzero__(), Types.int64::__nonzero__(), Types.uint64::__nonzero__(), Types.double::__nonzero__(), Types.bool::__nonzero__(), Types.string::__nonzero__(), average.Average::average(), Types.string::configValue(), Types.FileInPath::configValue(), Mixins.UsingBlock::dumpPython(), Mixins.UsingBlock::insertInto(), Types.int32::insertInto(), Types.uint32::insertInto(), Types.int64::insertInto(), Types.uint64::insertInto(), Types.double::insertInto(), Types.bool::insertInto(), Types.string::insertInto(), Types.FileInPath::insertInto(), Types.vint32::insertInto(), Types.vuint32::insertInto(), Types.vint64::insertInto(), Types.vuint64::insertInto(), Types.vdouble::insertInto(), Types.vbool::insertInto(), Types.vstring::insertInto(), and CSCChamberFitter::radiusCorrection().

9  {
10  double delta = (m_sum1 * m_sumxx) - (m_sumx * m_sumx);
11  assert(delta > 0.);
12  if (m_parent->m_mode == kModePhiy || m_parent->m_mode == kModePhiPos || m_parent->m_mode == kModeRadius) {
13  return ((m_sumxx * m_sumy) - (m_sumx * m_sumxy)) / delta;
14  } else if (m_parent->m_mode == kModePhiz) {
15  return ((m_sum1 * m_sumxy) - (m_sumx * m_sumy)) / delta;
16  } else
17  assert(false);
18 }
cms::cuda::assert
assert(be >=bs)
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:85

◆ write()

void CSCPairResidualsConstraint::write ( std::ofstream &  output)

Definition at line 567 of file CSCPairResidualsConstraint.cc.

References alignBH_cfg::fixed, CSCPairConstraint::i(), CSCPairConstraint::j(), m_identifier, m_sum1, m_sumN, m_sumx, m_sumxx, m_sumxy, m_sumy, m_sumyy, and convertSQLitetoXML_cfg::output.

Referenced by pkg.AbstractPkg::generate(), and querying.connection::write_and_commit().

567  {
568  output << std::setprecision(14) << std::fixed;
569  output << "CSCPairResidualsConstraint " << m_identifier << " " << i() << " " << j() << " " << m_sumN << " " << m_sum1
570  << " " << m_sumx << " " << m_sumy << " " << m_sumxx << " " << m_sumyy << " " << m_sumxy << " EOLN"
571  << std::endl;
572 }
alignBH_cfg.fixed
fixed
Definition: alignBH_cfg.py:54
convertSQLitetoXML_cfg.output
output
Definition: convertSQLitetoXML_cfg.py:72
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:19

Member Data Documentation

◆ m_averageRadius

double CSCPairResidualsConstraint::m_averageRadius
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

◆ m_cscGeometry

const CSCGeometry* CSCPairResidualsConstraint::m_cscGeometry
protected

Definition at line 88 of file CSCPairResidualsConstraint.h.

Referenced by calculatePhi(), isFiducial(), radius(), and setZplane().

◆ m_id_i

CSCDetId CSCPairResidualsConstraint::m_id_i
protected

Definition at line 81 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), configure(), dphidzFromTrack(), id_i(), isFiducial(), radius(), and setZplane().

◆ m_id_j

CSCDetId CSCPairResidualsConstraint::m_id_j
protected

Definition at line 81 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), configure(), dphidzFromTrack(), id_j(), isFiducial(), radius(), and setZplane().

◆ m_identifier

unsigned int CSCPairResidualsConstraint::m_identifier
protected

Definition at line 80 of file CSCPairResidualsConstraint.h.

Referenced by configure(), read(), and write().

◆ m_iZ

double CSCPairResidualsConstraint::m_iZ
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

◆ m_iZ1

double CSCPairResidualsConstraint::m_iZ1
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

◆ m_iZ6

double CSCPairResidualsConstraint::m_iZ6
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

◆ m_jZ

double CSCPairResidualsConstraint::m_jZ
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

◆ m_jZ1

double CSCPairResidualsConstraint::m_jZ1
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

◆ m_jZ6

double CSCPairResidualsConstraint::m_jZ6
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

◆ m_offsetResiduals

TH1F* CSCPairResidualsConstraint::m_offsetResiduals
protected

Definition at line 93 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

◆ m_parent

CSCOverlapsAlignmentAlgorithm* CSCPairResidualsConstraint::m_parent
protected

Definition at line 85 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), configure(), dphidzFromTrack(), error(), isFiducial(), valid(), and value().

◆ m_propagator

const Propagator* CSCPairResidualsConstraint::m_propagator
protected

Definition at line 89 of file CSCPairResidualsConstraint.h.

Referenced by dphidzFromTrack(), and setPropagator().

◆ m_radial

TH1F* CSCPairResidualsConstraint::m_radial
protected

Definition at line 94 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

◆ m_slopeResiduals

TH1F* CSCPairResidualsConstraint::m_slopeResiduals
protected

Definition at line 92 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

◆ m_sum1

double CSCPairResidualsConstraint::m_sum1
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), value(), and write().

◆ m_sumN

int CSCPairResidualsConstraint::m_sumN
protected

Definition at line 83 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), valid(), and write().

◆ m_sumx

double CSCPairResidualsConstraint::m_sumx
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), value(), and write().

◆ m_sumxx

double CSCPairResidualsConstraint::m_sumxx
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), value(), and write().

◆ m_sumxy

double CSCPairResidualsConstraint::m_sumxy
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), read(), value(), and write().

◆ m_sumy

double CSCPairResidualsConstraint::m_sumy
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), value(), and write().

◆ m_sumyy

double CSCPairResidualsConstraint::m_sumyy
protected

Definition at line 82 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), error(), read(), and write().

◆ m_Zplane

double CSCPairResidualsConstraint::m_Zplane
protected

Definition at line 87 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), dphidzFromTrack(), and setZplane().

◆ m_Zsurface

Plane::PlanePointer CSCPairResidualsConstraint::m_Zsurface
protected

Definition at line 90 of file CSCPairResidualsConstraint.h.

Referenced by dphidzFromTrack(), and setZplane().

Generated for CMSSW Reference Manual by doxygen 1.8.14