CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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
 
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
 
double value () const
 
void write (std::ofstream &output)
 
virtual ~CSCPairResidualsConstraint ()
 
- 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
Enumerator
kModePhiy 
kModePhiPos 
kModePhiz 
kModeRadius 

Definition at line 37 of file CSCPairResidualsConstraint.h.

37  {
38  kModePhiy,
39  kModePhiPos,
40  kModePhiz,
41  kModeRadius
42  };

Constructor & Destructor Documentation

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), m_id_i(id_i), m_id_j(id_j)
32  , m_sum1(0.), m_sumx(0.), m_sumy(0.), m_sumxx(0.), m_sumyy(0.), m_sumxy(0.), m_sumN(0)
33  , m_Zplane(1000.), m_iZ1(1000.), m_iZ6(1000.), m_jZ1(1000.), m_jZ6(1000.), m_cscGeometry(NULL), m_propagator(NULL)
34  {};
CSCPairConstraint::CSCPairConstraint
CSCPairConstraint(int i, int j, double value, double error)
Definition: CSCPairConstraint.h:16
NULL
#define NULL
Definition: scimark2.h:8
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:21
virtual CSCPairResidualsConstraint::~CSCPairResidualsConstraint ( )
inlinevirtual

Definition at line 35 of file CSCPairResidualsConstraint.h.

35 {};

Member Function Documentation

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

Definition at line 93 of file CSCPairResidualsConstraint.cc.

References assert(), calculatePhi(), CSCDetId::chamber(), MuonSubdetId::CSC, CSCDetId, dphidzFromTrack(), CSCDetId::endcap(), TrackingRecHit::geographicalId(), 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(), alignCSCRings::r, CSCDetId::ring(), mathSSE::sqrt(), CSCDetId::station(), puppiForMET_cff::weight, and z.

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

Definition at line 535 of file CSCPairResidualsConstraint.cc.

References angle(), assert(), funct::cos(), CSCDetId, CSCDetId::endcap(), TrackingRecHit::geographicalId(), CSCLayer::geometry(), TrackingRecHit::globalPosition(), CSCGeometry::layer(), TrackingRecHit::localPosition(), TrackingRecHit::localPositionError(), m_cscGeometry, M_PI, CSCLayerGeometry::nearestStrip(), PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), dttmaxenums::R, alignCSCRings::r, CSCDetId::ring(), funct::sin(), CSCDetId::station(), CSCLayerGeometry::stripAngle(), PV3DBase< T, PVType, FrameType >::x(), LocalError::xx(), LocalError::xy(), create_public_lumi_plots::xy, PV3DBase< T, PVType, FrameType >::y(), and LocalError::yy().

Referenced by addTrack(), and isFiducial().

535  {
536  align::LocalPoint pos = hit->localPosition();
537  DetId id = hit->geographicalId();
538  CSCDetId cscid = CSCDetId(id.rawId());
539 
540  double r = 0.;
541  if (globalPhi) {
542  phi = hit->globalPosition().phi();
543  r = hit->globalPosition().perp();
544 
545 // double sinAngle = sin(phi);
546 // double cosAngle = cos(phi);
547 // double xx = hit->globalPositionError().cxx();
548 // double xy = hit->globalPositionError().cyx();
549 // double yy = hit->globalPositionError().cyy();
550 // phierr2 = (xx*cosAngle*cosAngle + 2.*xy*sinAngle*cosAngle + yy*sinAngle*sinAngle) / (r*r);
551  }
552  else {
553  // these constants are related to the way CSC chambers are built--- really constant!
554  const double R_ME11 = 181.5;
555  const double R_ME12 = 369.7;
556  const double R_ME21 = 242.7;
557  const double R_ME31 = 252.7;
558  const double R_ME41 = 262.65;
559  const double R_MEx2 = 526.5;
560 
561  double R = 0.;
562  if (cscid.station() == 1 && (cscid.ring() == 1 || cscid.ring() == 4)) R = R_ME11;
563  else if (cscid.station() == 1 && cscid.ring() == 2) R = R_ME12;
564  else if (cscid.station() == 2 && cscid.ring() == 1) R = R_ME21;
565  else if (cscid.station() == 3 && cscid.ring() == 1) R = R_ME31;
566  else if (cscid.station() == 4 && cscid.ring() == 1) R = R_ME41;
567  else if (cscid.station() > 1 && cscid.ring() == 2) R = R_MEx2;
568  else assert(false);
569  r = (pos.y() + R);
570 
571  phi = atan2(pos.x(), r);
572 
573  if (cscid.endcap() == 1 && cscid.station() >= 3) phi *= -1;
574  else if (cscid.endcap() == 2 && cscid.station() <= 2) phi *= -1;
575  }
576 
577  int strip = m_cscGeometry->layer(id)->geometry()->nearestStrip(pos);
578  double angle = m_cscGeometry->layer(id)->geometry()->stripAngle(strip) - M_PI/2.;
579  double sinAngle = sin(angle);
580  double cosAngle = cos(angle);
581  double xx = hit->localPositionError().xx();
582  double xy = hit->localPositionError().xy();
583  double yy = hit->localPositionError().yy();
584  phierr2 = (xx*cosAngle*cosAngle + 2.*xy*sinAngle*cosAngle + yy*sinAngle*sinAngle) / (r*r);
585 
586  if (doRphi) {
587  phi *= r;
588  phierr2 *= r*r;
589  }
590 }
LocalError::xx
float xx() const
Definition: LocalError.h:24
CSCDetId
CSCDetId
Definition: RecoLocalMuon_RECO.doi:36
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
assert
assert(m_qm.get())
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
create_public_lumi_plots.xy
tuple xy
Definition: create_public_lumi_plots.py:1811
TrackingRecHit::globalPosition
virtual GlobalPoint globalPosition() const
Definition: TrackingRecHit.cc:47
CSCDetId::endcap
int endcap() const
Definition: CSCDetId.h:93
LocalError::xy
float xy() const
Definition: LocalError.h:25
LocalError::yy
float yy() const
Definition: LocalError.h:26
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:51
CSCLayerGeometry::nearestStrip
int nearestStrip(const LocalPoint &lp) const
Definition: CSCLayerGeometry.h:94
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:75
DetId
Definition: DetId.h:18
TrackingRecHit::localPositionError
virtual LocalError localPositionError() const =0
CSCGeometry::layer
const CSCLayer * layer(CSCDetId id) const
Return the layer corresponding to given DetId.
Definition: CSCGeometry.cc:124
CSCDetId::station
int station() const
Definition: CSCDetId.h:86
TrackingRecHit::geographicalId
DetId geographicalId() const
Definition: TrackingRecHit.h:116
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
TrackingRecHit::localPosition
virtual LocalPoint localPosition() const =0
CSCLayerGeometry::stripAngle
float stripAngle(int strip) const
Definition: CSCLayerGeometry.cc:140
CSCLayer::geometry
const CSCLayerGeometry * geometry() const
Definition: CSCLayer.h:47
angle
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11
void CSCPairResidualsConstraint::configure ( CSCOverlapsAlignmentAlgorithm *  parent)

Definition at line 43 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(), mergeVDriftHistosByStation::name, NULL, and indexGen::title.

43  {
44  m_parent = parent;
45 
46  if (m_parent->m_makeHistograms) {
47  edm::Service<TFileService> tFileService;
48 
49  std::stringstream name, name2, name3, title;
50  title << "i =" << m_id_i << " j =" << m_id_j;
51 
52  name << "slopeResiduals_" << m_identifier;
53  m_slopeResiduals = tFileService->make<TH1F>(name.str().c_str(), title.str().c_str(), 300, -30., 30.);
54 
55  name2 << "offsetResiduals_" << m_identifier;
56  m_offsetResiduals = tFileService->make<TH1F>(name2.str().c_str(), title.str().c_str(), 300, -30., 30.);
57 
58  name3 << "radial_" << m_identifier;
59  m_radial = tFileService->make<TH1F>(name3.str().c_str(), title.str().c_str(), 700, 0., 700.);
60  }
61  else {
62  m_slopeResiduals = NULL;
63  m_offsetResiduals = NULL;
64  m_radial = NULL;
65  }
66 }
NULL
#define NULL
Definition: scimark2.h:8
TFileService::make
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
bool CSCPairResidualsConstraint::dphidzFromTrack ( const std::vector< TrajectoryMeasurement > &  measurements,
const reco::TransientTrack &  track,
const TrackTransformer *  trackTransformer,
double &  drphidz 
)
protected

Definition at line 385 of file CSCPairResidualsConstraint.cc.

References angle(), CSCDetId::chamber(), TrajectoryStateCombiner::combine(), funct::cos(), MuonSubdetId::CSC, CSCDetId, 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(), mergeVDriftHistosByStation::stations, pileupDistInMC::total, HLT_FULL_cff::trajectories, TrackTransformer::transform(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by addTrack().

385  {
386  // make a list of hits on all chambers *other* than the ones associated with this constraint
387  std::map<int,int> stations;
388  int total = 0;
389  TransientTrackingRecHit::ConstRecHitContainer cscHits;
390  for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements.begin(); measurement != measurements.end(); ++measurement) {
391  DetId id = measurement->recHit()->geographicalId();
392  if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
393  CSCDetId cscid(id.rawId());
394  CSCDetId chamberId(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber(), 0);
395  if (m_parent->m_combineME11 && cscid.station() == 1 && cscid.ring() == 4) chamberId = CSCDetId(cscid.endcap(), 1, 1, cscid.chamber(), 0);
396 
397  if (chamberId != m_id_i && chamberId != m_id_j) {
398  int station = (cscid.endcap() == 1 ? 1 : -1) * cscid.station();
399  if (stations.find(station) == stations.end()) {
400  stations[station] = 0;
401  }
402  stations[station]++;
403  total++;
404 
405  cscHits.push_back(measurement->recHit());
406  }
407  }
408  }
409 
410  // for the fit to be reliable, it needs to cross multiple stations
411  int numStations = 0;
412  for (std::map<int,int>::const_iterator station = stations.begin(); station != stations.end(); ++station) {
413  if (station->second >= m_parent->m_minHitsPerChamber) {
414  numStations++;
415  }
416  }
417 
418  if (numStations >= m_parent->m_minStationsInTrackRefits) {
419  // refit the track with these hits
420  std::vector<Trajectory> trajectories = trackTransformer->transform(track, cscHits);
421 
422  if (!trajectories.empty()) {
423  const std::vector<TrajectoryMeasurement> &measurements2 = trajectories.begin()->measurements();
424 
425  // find the closest TSOS to the Z plane (on both sides)
426  bool found_plus = false;
427  bool found_minus = false;
428  TrajectoryStateOnSurface tsos_plus, tsos_minus;
429  for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements2.begin(); measurement != measurements2.end(); ++measurement) {
430  double z = measurement->recHit()->globalPosition().z();
431  if (z > m_Zplane) {
432  if (!found_plus || fabs(z - m_Zplane) < fabs(tsos_plus.globalPosition().z() - m_Zplane)) {
433  tsos_plus = TrajectoryStateCombiner().combine(measurement->forwardPredictedState(), measurement->backwardPredictedState());
434  }
435  if (tsos_plus.isValid()) found_plus = true;
436  }
437  else {
438  if (!found_minus || fabs(z - m_Zplane) < fabs(tsos_minus.globalPosition().z() - m_Zplane)) {
439  tsos_minus = TrajectoryStateCombiner().combine(measurement->forwardPredictedState(), measurement->backwardPredictedState());
440  }
441  if (tsos_minus.isValid()) found_minus = true;
442  }
443  }
444 
445  // propagate from the closest TSOS to the Z plane (from both sides, if possible)
446  TrajectoryStateOnSurface from_plus, from_minus;
447  if (found_plus) {
448  from_plus = m_propagator->propagate(tsos_plus, *m_Zsurface);
449  }
450  if (found_minus) {
451  from_minus = m_propagator->propagate(tsos_minus, *m_Zsurface);
452  }
453 
454  // if you have two sides, merge them
455  TrajectoryStateOnSurface merged;
456  if (found_plus && from_plus.isValid() && found_minus && from_minus.isValid()) {
457  merged = TrajectoryStateCombiner().combine(from_plus, from_minus);
458  }
459  else if (found_plus && from_plus.isValid()) {
460  merged = from_plus;
461  }
462  else if (found_minus && from_minus.isValid()) {
463  merged = from_minus;
464  }
465  else return false;
466 
467  // if, after all that, we have a good fit-and-propagation, report the direction
468  if (merged.isValid()) {
469  double angle = merged.globalPosition().phi() + M_PI/2.;
470  GlobalVector direction = merged.globalDirection();
471  double dxdz = direction.x() / direction.z();
472  double dydz = direction.y() / direction.z();
473  dphidz = (dxdz*cos(angle) + dydz*sin(angle)) / merged.globalPosition().perp();
474  return true;
475  }
476 
477  } // end if refit successful
478  } // end if enough hits
479  return false;
480 }
CSCDetId
CSCDetId
Definition: RecoLocalMuon_RECO.doi:36
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
TrajectoryStateCombiner::combine
TSOS combine(const TSOS &pTsos1, const TSOS &pTsos2) const
Definition: TrajectoryStateCombiner.cc:7
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:83
MuonSubdetId::CSC
static const int CSC
Definition: MuonSubdetId.h:13
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
HLT_FULL_cff.trajectories
tuple trajectories
Definition: HLT_FULL_cff.py:7275
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:51
TrackingRecHit::ConstRecHitContainer
std::vector< ConstRecHitPointer > ConstRecHitContainer
Definition: TrackingRecHit.h:35
DetId
Definition: DetId.h:18
Propagator::propagate
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:53
TrackTransformer::transform
virtual std::vector< Trajectory > transform(const reco::Track &) const
Convert a reco::Track into Trajectory.
Definition: TrackTransformer.cc:165
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
TrajectoryStateOnSurface::globalDirection
GlobalVector globalDirection() const
Definition: TrajectoryStateOnSurface.h:89
angle
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11
double CSCPairResidualsConstraint::error ( ) const
virtual

Reimplemented from CSCPairConstraint.

Definition at line 21 of file CSCPairResidualsConstraint.cc.

References assert(), 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 argparse.ArgumentParser::_get_option_tuples(), python.rootplot.argparse.ArgumentParser::_get_option_tuples(), argparse.ArgumentParser::_parse_known_args(), python.rootplot.argparse.ArgumentParser::_parse_known_args(), argparse.ArgumentParser::_parse_optional(), python.rootplot.argparse.ArgumentParser::_parse_optional(), argparse.ArgumentParser::_read_args_from_files(), python.rootplot.argparse.ArgumentParser::_read_args_from_files(), argparse.ArgumentParser::add_subparsers(), python.rootplot.argparse.ArgumentParser::add_subparsers(), argparse.ArgumentParser::parse_args(), python.rootplot.argparse.ArgumentParser::parse_args(), argparse.ArgumentParser::parse_known_args(), and python.rootplot.argparse.ArgumentParser::parse_known_args().

21  {
22  if (m_parent->m_errorFromRMS) {
23  assert(m_sum1 > 0.);
24  return sqrt((m_sumyy/m_sum1) - pow(m_sumy/m_sum1, 2))/sqrt(m_sumN);
25  }
26  else {
27  double delta = (m_sum1*m_sumxx) - (m_sumx*m_sumx);
28  assert(delta > 0.);
29  if (m_parent->m_mode == kModePhiy || m_parent->m_mode == kModePhiPos || m_parent->m_mode == kModeRadius) {
30  return sqrt(m_sumxx/delta);
31  }
32  else if (m_parent->m_mode == kModePhiz) {
33  return sqrt(m_sum1/delta);
34  }
35  else assert(false);
36  }
37 }
delta
dbl * delta
Definition: mlp_gen.cc:36
assert
assert(m_qm.get())
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
CSCDetId CSCPairResidualsConstraint::id_i ( ) const
inline

Definition at line 46 of file CSCPairResidualsConstraint.h.

References m_id_i.

Referenced by CSCChamberFitter::radiusCorrection().

46 { return m_id_i; };
CSCDetId CSCPairResidualsConstraint::id_j ( ) const
inline

Definition at line 47 of file CSCPairResidualsConstraint.h.

References m_id_j.

47 { return m_id_j; };
bool CSCPairResidualsConstraint::isFiducial ( std::vector< const TransientTrackingRecHit * > &  hits,
bool  is_i 
)
protected

Definition at line 592 of file CSCPairResidualsConstraint.cc.

References assert(), calculatePhi(), delta, 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(), puppiForMET_cff::weight, and z.

Referenced by addTrack().

592  {
593  // these constants are related to the way CSC chambers are built--- really constant!
594  const double cut_ME11 = 0.086;
595  const double cut_ME12 = 0.090;
596  const double cut_MEx1 = 0.180;
597  const double cut_MEx2 = 0.090;
598 
599  double sum1 = 0.;
600  double sumx = 0.;
601  double sumy = 0.;
602  double sumxx = 0.;
603  double sumxy = 0.;
604  for (std::vector<const TransientTrackingRecHit*>::const_iterator hit = hits.begin(); hit != hits.end(); ++hit) {
605  double phi, phierr2;
606  calculatePhi(*hit, phi, phierr2);
607  double z = (is_i ? m_cscGeometry->idToDet(m_id_i)->surface() : m_cscGeometry->idToDet(m_id_j)->surface()).toLocal((*hit)->globalPosition()).z();
608 
609  if (m_parent->m_makeHistograms) {
610  if (m_id_i.station() == 1 && (m_id_i.ring() == 1 || m_id_i.ring() == 4)) {
611  m_parent->m_fiducial_ME11->Fill(fabs(phi), sqrt(phierr2));
612  }
613  else if (m_id_i.station() == 1 && m_id_i.ring() == 2) {
614  m_parent->m_fiducial_ME12->Fill(fabs(phi), sqrt(phierr2));
615  }
616  else if (m_id_i.station() > 1 && m_id_i.ring() == 1) {
617  m_parent->m_fiducial_MEx1->Fill(fabs(phi), sqrt(phierr2));
618  }
619  else if (m_id_i.station() > 1 && m_id_i.ring() == 2) {
620  m_parent->m_fiducial_MEx2->Fill(fabs(phi), sqrt(phierr2));
621  }
622  }
623 
624  double weight = 1.;
625  if (m_parent->m_useHitWeights) weight = 1./phierr2;
626  sum1 += weight;
627  sumx += weight * z;
628  sumy += weight * phi;
629  sumxx += weight * z*z;
630  sumxy += weight * z*phi;
631  }
632  double delta = (sum1*sumxx) - (sumx*sumx);
633  if (delta == 0.) return false;
634  double intercept = ((sumxx*sumy) - (sumx*sumxy))/delta;
635  double slope = ((sum1*sumxy) - (sumx*sumy))/delta;
636 
637  double phi1 = intercept + slope*(is_i ? m_iZ1 : m_jZ1);
638  double phi6 = intercept + slope*(is_i ? m_iZ6 : m_jZ6);
639 
640  if (m_id_i.station() == 1 && (m_id_i.ring() == 1 || m_id_i.ring() == 4)) {
641  return (fabs(phi1) < cut_ME11 && fabs(phi6) < cut_ME11);
642  }
643  else if (m_id_i.station() == 1 && m_id_i.ring() == 2) {
644  return (fabs(phi1) < cut_ME12 && fabs(phi6) < cut_ME12);
645  }
646  else if (m_id_i.station() > 1 && m_id_i.ring() == 1) {
647  return (fabs(phi1) < cut_MEx1 && fabs(phi6) < cut_MEx1);
648  }
649  else if (m_id_i.station() > 1 && m_id_i.ring() == 2) {
650  return (fabs(phi1) < cut_MEx2 && fabs(phi6) < cut_MEx2);
651  }
652  else assert(false);
653 }
delta
dbl * delta
Definition: mlp_gen.cc:36
slope
static const double slope[3]
Definition: CastorTimeSlew.cc:6
assert
assert(m_qm.get())
puppiForMET_cff.weight
tuple weight
Definition: puppiForMET_cff.py:35
toLocal
LocalVector toLocal(const reco::Track::Vector &v, const Surface &s)
Definition: ConversionProducer.h:207
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:40
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
CSCPairResidualsConstraint::calculatePhi
void calculatePhi(const TransientTrackingRecHit *hit, double &phi, double &phierr2, bool doRphi=false, bool globalPhi=false)
Definition: CSCPairResidualsConstraint.cc:535
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:75
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
CSCDetId::station
int station() const
Definition: CSCDetId.h:86
CSCGeometry::idToDet
virtual const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:99
double CSCPairResidualsConstraint::radius ( bool  is_i) const
inline

Definition at line 49 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().

49 { return m_cscGeometry->idToDet((is_i ? m_id_i : m_id_j))->surface().position().perp(); };
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:40
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:42
CSCGeometry::idToDet
virtual const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:99
void CSCPairResidualsConstraint::read ( std::vector< std::ifstream * > &  input,
std::vector< std::string > &  filenames 
)

Definition at line 488 of file CSCPairResidualsConstraint.cc.

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

Referenced by edmIntegrityCheck.PublishToFileSystem::get(), Vispa.Plugins.EdmBrowser.EdmDataAccessor.EdmDataAccessor::goto(), and Vispa.Plugins.EdmBrowser.EdmDataAccessor.EdmDataAccessor::setFilterBranches().

488  {
489  m_sumN = 0;
490  m_sum1 = 0.;
491  m_sumx = 0.;
492  m_sumy = 0.;
493  m_sumxx = 0.;
494  m_sumyy = 0.;
495  m_sumxy = 0.;
496 
497  std::vector<std::ifstream*>::const_iterator inputiter = input.begin();
498  std::vector<std::string>::const_iterator filename = filenames.begin();
499  for (; inputiter != input.end(); ++inputiter, ++filename) {
500  int linenumber = 0;
501  bool touched = false;
502  while (!(*inputiter)->eof()) {
503  linenumber++;
504  std::string name, eoln;
505  unsigned int identifier;
506  int i, j;
507  int sumN;
508  double sum1, sumx, sumy, sumxx, sumyy, sumxy;
509 
510  (**inputiter) >> name >> identifier >> i >> j >> sumN >> sum1 >> sumx >> sumy >> sumxx >> sumyy >> sumxy >> eoln;
511 
512  if (!(*inputiter)->eof() && (name != "CSCPairResidualsConstraint" || eoln != "EOLN")) throw cms::Exception("CorruptTempFile") << "Temporary file " << *filename << " is incorrectly formatted on line " << linenumber << std::endl;
513 
514  if (identifier == m_identifier) {
515  if (i != m_i || j != m_j) throw cms::Exception("CorruptTempFile") << "Wrong (i,j) for CSCPairResidualsConstraint " << m_identifier << " (" << m_i << "," << m_j << ") in file " << *filename << " on line " << linenumber << std::endl;
516  touched = true;
517 
518  m_sumN += sumN;
519  m_sum1 += sum1;
520  m_sumx += sumx;
521  m_sumy += sumy;
522  m_sumxx += sumxx;
523  m_sumyy += sumyy;
524  m_sumxy += sumxy;
525  }
526  }
527 
528  (*inputiter)->clear();
529  (*inputiter)->seekg(0, std::ios::beg);
530 
531  if (!touched) throw cms::Exception("CorruptTempFile") << "CSCPairResidualsConstraint " << m_identifier << " is missing from file " << *filename << std::endl;
532  }
533 }
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
input
static std::string const input
Definition: EdmProvDump.cc:44
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:21
lut2db_cfg.filename
tuple filename
Definition: lut2db_cfg.py:20
void CSCPairResidualsConstraint::setPropagator ( const Propagator *  propagator)

Definition at line 89 of file CSCPairResidualsConstraint.cc.

References m_propagator, and HLT_FULL_cff::propagator.

89  {
90  m_propagator = propagator;
91 }
HLT_FULL_cff.propagator
tuple propagator
Definition: HLT_FULL_cff.py:700
void CSCPairResidualsConstraint::setZplane ( const CSCGeometry *  cscGeometry)

Definition at line 68 of file CSCPairResidualsConstraint.cc.

References Plane::build(), CSCDetId::chamber(), CSCDetId::endcap(), 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().

68  {
69  m_cscGeometry = cscGeometry;
70 
71  m_Zplane = (m_cscGeometry->idToDet(m_id_i)->surface().position().z() + m_cscGeometry->idToDet(m_id_j)->surface().position().z()) / 2.;
72  m_averageRadius = (m_cscGeometry->idToDet(m_id_i)->surface().position().perp() + m_cscGeometry->idToDet(m_id_j)->surface().position().perp()) / 2.;
73 
74  m_iZ = m_cscGeometry->idToDet(m_id_i)->surface().position().z();
75  m_jZ = m_cscGeometry->idToDet(m_id_j)->surface().position().z();
76 
77  CSCDetId i1(m_id_i.endcap(), m_id_i.station(), m_id_i.ring(), m_id_i.chamber(), 1);
78  CSCDetId i6(m_id_i.endcap(), m_id_i.station(), m_id_i.ring(), m_id_i.chamber(), 6);
79  CSCDetId j1(m_id_j.endcap(), m_id_j.station(), m_id_j.ring(), m_id_j.chamber(), 1);
80  CSCDetId j6(m_id_j.endcap(), m_id_j.station(), m_id_j.ring(), m_id_j.chamber(), 6);
81  m_iZ1 = m_cscGeometry->idToDet(m_id_i)->surface().toLocal(m_cscGeometry->idToDet(i1)->surface().position()).z();
82  m_iZ6 = m_cscGeometry->idToDet(m_id_i)->surface().toLocal(m_cscGeometry->idToDet(i6)->surface().position()).z();
83  m_jZ1 = m_cscGeometry->idToDet(m_id_j)->surface().toLocal(m_cscGeometry->idToDet(j1)->surface().position()).z();
84  m_jZ6 = m_cscGeometry->idToDet(m_id_j)->surface().toLocal(m_cscGeometry->idToDet(j6)->surface().position()).z();
85 
86  m_Zsurface = Plane::build(Plane::PositionType(0., 0., m_Zplane), Plane::RotationType());
87 }
CSCDetId::chamber
int chamber() const
Definition: CSCDetId.h:68
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
GeomDet::surface
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:40
CSCDetId::endcap
int endcap() const
Definition: CSCDetId.h:93
GloballyPositioned::toLocal
LocalPoint toLocal(const GlobalPoint &gp) const
Definition: GloballyPositioned.h:125
Plane::build
static PlanePointer build(Args &&...args)
Definition: Plane.h:33
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:75
CSCDetId::station
int station() const
Definition: CSCDetId.h:86
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:42
CSCGeometry::idToDet
virtual const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:99
bool CSCPairResidualsConstraint::valid ( ) const
virtual

Reimplemented from CSCPairConstraint.

Definition at line 39 of file CSCPairResidualsConstraint.cc.

References CSCOverlapsAlignmentAlgorithm::m_minTracksPerOverlap, m_parent, and m_sumN.

Referenced by CSCChamberFitter::radiusCorrection().

39  {
40  return (m_sumN >= m_parent->m_minTracksPerOverlap);
41 }
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
double CSCPairResidualsConstraint::value ( ) const
virtual

Reimplemented from CSCPairConstraint.

Definition at line 9 of file CSCPairResidualsConstraint.cc.

References assert(), 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(), Vispa.Views.PropertyView.FileProperty::buttonClicked(), Types.string::configValue(), Types.FileInPath::configValue(), Mixins.UsingBlock::dumpPython(), Types.int32::insertInto(), Types.uint32::insertInto(), Types.int64::insertInto(), Types.uint64::insertInto(), Mixins.UsingBlock::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(), Vispa.Views.PropertyView.FileProperty::labelDoubleClicked(), 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  }
15  else if (m_parent->m_mode == kModePhiz) {
16  return ((m_sum1*m_sumxy) - (m_sumx*m_sumy))/delta;
17  }
18  else assert(false);
19 }
delta
dbl * delta
Definition: mlp_gen.cc:36
assert
assert(m_qm.get())
CSCPairResidualsConstraint::m_parent
CSCOverlapsAlignmentAlgorithm * m_parent
Definition: CSCPairResidualsConstraint.h:69
void CSCPairResidualsConstraint::write ( std::ofstream &  output)

Definition at line 482 of file CSCPairResidualsConstraint.cc.

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

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

482  {
483  output << std::setprecision(14) << std::fixed;
484  output << "CSCPairResidualsConstraint " << m_identifier << " " << i() << " " << j() << " "
485  << m_sumN << " " << m_sum1 << " " << m_sumx << " " << m_sumy << " " << m_sumxx << " " << m_sumyy << " " << m_sumxy << " EOLN" << std::endl;
486 }
CSCPairConstraint::i
virtual int i() const
Definition: CSCPairConstraint.h:20
CSCPairConstraint::j
virtual int j() const
Definition: CSCPairConstraint.h:21
convertSQLitetoXML_cfg.output
tuple output
Definition: convertSQLitetoXML_cfg.py:32

Member Data Documentation

double CSCPairResidualsConstraint::m_averageRadius
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

const CSCGeometry* CSCPairResidualsConstraint::m_cscGeometry
protected

Definition at line 72 of file CSCPairResidualsConstraint.h.

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

CSCDetId CSCPairResidualsConstraint::m_id_i
protected

Definition at line 65 of file CSCPairResidualsConstraint.h.

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

CSCDetId CSCPairResidualsConstraint::m_id_j
protected

Definition at line 65 of file CSCPairResidualsConstraint.h.

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

unsigned int CSCPairResidualsConstraint::m_identifier
protected

Definition at line 64 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_iZ
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

double CSCPairResidualsConstraint::m_iZ1
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

double CSCPairResidualsConstraint::m_iZ6
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

double CSCPairResidualsConstraint::m_jZ
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and setZplane().

double CSCPairResidualsConstraint::m_jZ1
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

double CSCPairResidualsConstraint::m_jZ6
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

Referenced by isFiducial(), and setZplane().

TH1F* CSCPairResidualsConstraint::m_offsetResiduals
protected

Definition at line 77 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

CSCOverlapsAlignmentAlgorithm* CSCPairResidualsConstraint::m_parent
protected

Definition at line 69 of file CSCPairResidualsConstraint.h.

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

const Propagator* CSCPairResidualsConstraint::m_propagator
protected

Definition at line 73 of file CSCPairResidualsConstraint.h.

Referenced by dphidzFromTrack(), and setPropagator().

TH1F* CSCPairResidualsConstraint::m_radial
protected

Definition at line 78 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

TH1F* CSCPairResidualsConstraint::m_slopeResiduals
protected

Definition at line 76 of file CSCPairResidualsConstraint.h.

Referenced by addTrack(), and configure().

double CSCPairResidualsConstraint::m_sum1
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

int CSCPairResidualsConstraint::m_sumN
protected

Definition at line 67 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_sumx
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_sumxx
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_sumxy
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_sumy
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_sumyy
protected

Definition at line 66 of file CSCPairResidualsConstraint.h.

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

double CSCPairResidualsConstraint::m_Zplane
protected

Definition at line 71 of file CSCPairResidualsConstraint.h.

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

Plane::PlanePointer CSCPairResidualsConstraint::m_Zsurface
protected

Definition at line 74 of file CSCPairResidualsConstraint.h.

Referenced by dphidzFromTrack(), and setZplane().

Generated for CMSSW Reference Manual by doxygen 1.8.5