CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MuonSeedCreator.cc
Go to the documentation of this file.
1 
9 #include "RecoMuon/MuonSeedGenerator/src/MuonSeedCreator.h"
10 
11 #include "RecoMuon/TransientTrackingRecHit/interface/MuonTransientTrackingRecHit.h"
12 
13 #include "RecoMuon/TrackingTools/interface/MuonPatternRecoDumper.h"
14 
15 #include "MagneticField/Engine/interface/MagneticField.h"
16 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
17 
18 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
19 #include "TrackingTools/DetLayers/interface/DetLayer.h"
20 
21 #include "DataFormats/TrajectoryState/interface/PTrajectoryStateOnDet.h"
22 #include "DataFormats/Common/interface/OwnVector.h"
23 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
24 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
25 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
26 
27 #include "FWCore/Framework/interface/EventSetup.h"
28 #include "FWCore/Framework/interface/ESHandle.h"
29 #include "FWCore/MessageLogger/interface/MessageLogger.h"
30 
31 #include "gsl/gsl_statistics.h"
32 
33 
34 /*
35  * Constructor
36  */
37 MuonSeedCreator::MuonSeedCreator(const edm::ParameterSet& pset){
38 
39  theMinMomentum = pset.getParameter<double>("minimumSeedPt");
40  theMaxMomentum = pset.getParameter<double>("maximumSeedPt");
41  defaultMomentum = pset.getParameter<double>("defaultSeedPt");
42  debug = pset.getParameter<bool>("DebugMuonSeed");
43  sysError = pset.getParameter<double>("SeedPtSystematics");
44  // load seed PT parameters
45  DT12 = pset.getParameter<std::vector<double> >("DT_12");
46  DT13 = pset.getParameter<std::vector<double> >("DT_13");
47  DT14 = pset.getParameter<std::vector<double> >("DT_14");
48  DT23 = pset.getParameter<std::vector<double> >("DT_23");
49  DT24 = pset.getParameter<std::vector<double> >("DT_24");
50  DT34 = pset.getParameter<std::vector<double> >("DT_34");
51 
52  CSC01 = pset.getParameter<std::vector<double> >("CSC_01");
53  CSC12 = pset.getParameter<std::vector<double> >("CSC_12");
54  CSC02 = pset.getParameter<std::vector<double> >("CSC_02");
55  CSC13 = pset.getParameter<std::vector<double> >("CSC_13");
56  CSC03 = pset.getParameter<std::vector<double> >("CSC_03");
57  CSC14 = pset.getParameter<std::vector<double> >("CSC_14");
58  CSC23 = pset.getParameter<std::vector<double> >("CSC_23");
59  CSC24 = pset.getParameter<std::vector<double> >("CSC_24");
60  CSC34 = pset.getParameter<std::vector<double> >("CSC_34");
61 
62  OL1213 = pset.getParameter<std::vector<double> >("OL_1213");
63  OL1222 = pset.getParameter<std::vector<double> >("OL_1222");
64  OL1232 = pset.getParameter<std::vector<double> >("OL_1232");
65  OL1213 = pset.getParameter<std::vector<double> >("OL_1213");
66  OL2222 = pset.getParameter<std::vector<double> >("OL_1222");
67 
68  SME11 = pset.getParameter<std::vector<double> >("SME_11");
69  SME12 = pset.getParameter<std::vector<double> >("SME_12");
70  SME13 = pset.getParameter<std::vector<double> >("SME_13");
71  SME21 = pset.getParameter<std::vector<double> >("SME_21");
72  SME22 = pset.getParameter<std::vector<double> >("SME_22");
73  SME31 = pset.getParameter<std::vector<double> >("SME_31");
74  SME32 = pset.getParameter<std::vector<double> >("SME_32");
75  SME41 = pset.getParameter<std::vector<double> >("SME_41");
76 
77  SMB10 = pset.getParameter<std::vector<double> >("SMB_10");
78  SMB11 = pset.getParameter<std::vector<double> >("SMB_11");
79  SMB12 = pset.getParameter<std::vector<double> >("SMB_12");
80  SMB20 = pset.getParameter<std::vector<double> >("SMB_20");
81  SMB21 = pset.getParameter<std::vector<double> >("SMB_21");
82  SMB22 = pset.getParameter<std::vector<double> >("SMB_22");
83  SMB30 = pset.getParameter<std::vector<double> >("SMB_30");
84  SMB31 = pset.getParameter<std::vector<double> >("SMB_31");
85  SMB32 = pset.getParameter<std::vector<double> >("SMB_32");
86 
87  // Load dphi scale parameters
88  CSC01_1 = pset.getParameter<std::vector<double> >("CSC_01_1_scale");
89  CSC12_1 = pset.getParameter<std::vector<double> >("CSC_12_1_scale");
90  CSC12_2 = pset.getParameter<std::vector<double> >("CSC_12_2_scale");
91  CSC12_3 = pset.getParameter<std::vector<double> >("CSC_12_3_scale");
92  CSC13_2 = pset.getParameter<std::vector<double> >("CSC_13_2_scale");
93  CSC13_3 = pset.getParameter<std::vector<double> >("CSC_13_3_scale");
94  CSC14_3 = pset.getParameter<std::vector<double> >("CSC_14_3_scale");
95  CSC23_1 = pset.getParameter<std::vector<double> >("CSC_23_1_scale");
96  CSC23_2 = pset.getParameter<std::vector<double> >("CSC_23_2_scale");
97  CSC24_1 = pset.getParameter<std::vector<double> >("CSC_24_1_scale");
98  CSC34_1 = pset.getParameter<std::vector<double> >("CSC_34_1_scale");
99 
100  DT12_1 = pset.getParameter<std::vector<double> >("DT_12_1_scale");
101  DT12_2 = pset.getParameter<std::vector<double> >("DT_12_2_scale");
102  DT13_1 = pset.getParameter<std::vector<double> >("DT_13_1_scale");
103  DT13_2 = pset.getParameter<std::vector<double> >("DT_13_2_scale");
104  DT14_1 = pset.getParameter<std::vector<double> >("DT_14_1_scale");
105  DT14_2 = pset.getParameter<std::vector<double> >("DT_14_2_scale");
106  DT23_1 = pset.getParameter<std::vector<double> >("DT_23_1_scale");
107  DT23_2 = pset.getParameter<std::vector<double> >("DT_23_2_scale");
108  DT24_1 = pset.getParameter<std::vector<double> >("DT_24_1_scale");
109  DT24_2 = pset.getParameter<std::vector<double> >("DT_24_2_scale");
110  DT34_1 = pset.getParameter<std::vector<double> >("DT_34_1_scale");
111  DT34_2 = pset.getParameter<std::vector<double> >("DT_34_2_scale");
112 
113  OL_1213 = pset.getParameter<std::vector<double> >("OL_1213_0_scale");
114  OL_1222 = pset.getParameter<std::vector<double> >("OL_1222_0_scale");
115  OL_1232 = pset.getParameter<std::vector<double> >("OL_1232_0_scale");
116  OL_2213 = pset.getParameter<std::vector<double> >("OL_2213_0_scale");
117  OL_2222 = pset.getParameter<std::vector<double> >("OL_2222_0_scale");
118 
119  SMB_10S = pset.getParameter<std::vector<double> >("SMB_10_0_scale");
120  SMB_11S = pset.getParameter<std::vector<double> >("SMB_11_0_scale");
121  SMB_12S = pset.getParameter<std::vector<double> >("SMB_12_0_scale");
122  SMB_20S = pset.getParameter<std::vector<double> >("SMB_20_0_scale");
123  SMB_21S = pset.getParameter<std::vector<double> >("SMB_21_0_scale");
124  SMB_22S = pset.getParameter<std::vector<double> >("SMB_22_0_scale");
125  SMB_30S = pset.getParameter<std::vector<double> >("SMB_30_0_scale");
126  SMB_31S = pset.getParameter<std::vector<double> >("SMB_31_0_scale");
127  SMB_32S = pset.getParameter<std::vector<double> >("SMB_32_0_scale");
128 
129  SME_11S = pset.getParameter<std::vector<double> >("SME_11_0_scale");
130  SME_12S = pset.getParameter<std::vector<double> >("SME_12_0_scale");
131  SME_13S = pset.getParameter<std::vector<double> >("SME_13_0_scale");
132  SME_21S = pset.getParameter<std::vector<double> >("SME_21_0_scale");
133  SME_22S = pset.getParameter<std::vector<double> >("SME_22_0_scale");
134 
135 }
136 
137 
138 /*
139  * Destructor
140  */
141 MuonSeedCreator::~MuonSeedCreator(){
142 
143 }
144 
145 
146 /*
147  * createSeed
148  *
149  * Note type = 1 --> CSC
150  * = 2 --> Overlap
151  * = 3 --> DT
152  */
153 TrajectorySeed MuonSeedCreator::createSeed(int type, const SegmentContainer& seg, const std::vector<int>& layers, int NShowers, int NShowerSegments ) {
154 
155  // The index of the station closest to the IP
156  int last = 0;
157 
158  double ptmean = theMinMomentum;
159  double sptmean = theMinMomentum;
160 
161  AlgebraicVector t(4);
162  AlgebraicSymMatrix mat(5,0) ;
163  LocalPoint segPos;
164 
165  // Compute the pt according to station types used;
166  if (type == 1 ) estimatePtCSC(seg, layers, ptmean, sptmean);
167  if (type == 2 ) estimatePtOverlap(seg, layers, ptmean, sptmean);
168  if (type == 3 ) estimatePtDT(seg, layers, ptmean, sptmean);
169  if (type == 4 ) estimatePtSingle(seg, layers, ptmean, sptmean);
170  // type 5 are the seeding for ME1/4
171  if (type == 5 ) estimatePtCSC(seg, layers, ptmean, sptmean);
172 
173  // for certain clear showering case, set-up the minimum value
174  if ( NShowers > 0 ) estimatePtShowering( NShowers, NShowerSegments, ptmean, sptmean );
175  //if ( NShowers > 0 ) std::cout<<" Showering happened "<<NShowers<<" times w/ "<< NShowerSegments<<std::endl; ;
176 
177 
178  // Minimal pt
179  double charge = 1.0;
180  if (ptmean < 0.) charge = -1.0;
181  if ( (charge * ptmean) < theMinMomentum ) {
182  ptmean = theMinMomentum * charge;
183  sptmean = theMinMomentum ;
184  }
185  else if ( (charge * ptmean) > theMaxMomentum ) {
186  ptmean = theMaxMomentum * charge;
187  sptmean = theMaxMomentum * 0.25 ;
188  }
189 
190  LocalTrajectoryParameters param;
191 
192  double p_err =0.0;
193  // determine the seed layer
194  int best_seg= 0;
195  double chi2_dof = 9999.0;
196  unsigned int ini_seg = 0;
197  // avoid generating seed from 1st layer(ME1/1)
198  if ( type == 5 ) ini_seg = 1;
199  for (size_t i = ini_seg ; i < seg.size(); i++) {
200  double dof = static_cast<double>(seg[i]->degreesOfFreedom());
201  if ( chi2_dof < ( seg[i]->chi2()/dof ) ) continue;
202  chi2_dof = seg[i]->chi2() / dof ;
203  best_seg = static_cast<int>(i);
204  }
205 
206 
207  if ( type==1 || type==5 || type== 4) {
208  // Fill the LocalTrajectoryParameters
210  last = best_seg;
211  // last = 0;
212  GlobalVector mom = seg[last]->globalPosition()-GlobalPoint();
213  segPos = seg[last]->localPosition();
215 
216  GlobalVector polar(GlobalVector::Spherical(mom.theta(),seg[last]->globalDirection().phi(),1.));
217 
219  polar *= fabs(ptmean)/polar.perp();
221  LocalVector segDirFromPos = seg[last]->det()->toLocal(polar);
222  int chargeI = static_cast<int>(charge);
223  LocalTrajectoryParameters param1(segPos, segDirFromPos, chargeI);
224  param = param1;
225  p_err = (sptmean*sptmean)/(polar.mag()*polar.mag()*ptmean*ptmean) ;
226  mat = seg[last]->parametersError().similarityT( seg[last]->projectionMatrix() );
227  mat[0][0]= p_err;
228  if ( type == 5 ) {
229  mat[0][0] = mat[0][0]/fabs( tan(mom.theta()) );
230  mat[1][1] = mat[1][1]/fabs( tan(mom.theta()) );
231  mat[3][3] = 2.25*mat[3][3];
232  mat[4][4] = 2.25*mat[4][4];
233  }
234  if ( type == 4 ) {
235  mat[0][0] = mat[0][0]/fabs( tan(mom.theta()) );
236  mat[1][1] = mat[1][1]/fabs( tan(mom.theta()) );
237  mat[2][2] = 2.25*mat[2][2];
238  mat[3][3] = 2.25*mat[3][3];
239  mat[4][4] = 2.25*mat[4][4];
240  }
241  double dh = fabs( seg[last]->globalPosition().eta() ) - 1.6 ;
242  if ( fabs(dh) < 0.1 && type == 1 ) {
243  mat[1][1] = 4.*mat[1][1];
244  mat[2][2] = 4.*mat[2][2];
245  mat[3][3] = 9.*mat[3][3];
246  mat[4][4] = 9.*mat[4][4];
247  }
248 
249  //if ( !highPt && type != 1 ) mat[1][1]= 2.25*mat[1][1];
250  //if ( highPt && type != 1 ) mat[3][3]= 2.25*mat[1][1];
251  //mat[2][2]= 3.*mat[2][2];
252  //mat[3][3]= 2.*mat[3][3];
253  //mat[4][4]= 2.*mat[4][4];
254  }
255  else {
256  // Fill the LocalTrajectoryParameters
258  last = 0;
259  segPos = seg[last]->localPosition();
260  GlobalVector mom = seg[last]->globalPosition()-GlobalPoint();
262  GlobalVector polar(GlobalVector::Spherical(mom.theta(),seg[last]->globalDirection().phi(),1.));
263  //GlobalVector polar(GlobalVector::Spherical(seg[last]->globalDirection().theta(),seg[last]->globalDirection().phi(),1.));
264 
266  //double ptRatio = 1.0 - (2.808/(fabs(ptmean) -1)) + (4.546/( (fabs(ptmean)-1)*(fabs(ptmean)-1)) );
267  //ptmean = ptmean*ptRatio ;
268 
270  polar *= fabs(ptmean)/polar.perp();
272  LocalVector segDirFromPos = seg[last]->det()->toLocal(polar);
273  int chargeI = static_cast<int>(charge);
274  LocalTrajectoryParameters param1(segPos, segDirFromPos, chargeI);
275  param = param1;
276  p_err = (sptmean*sptmean)/(polar.mag()*polar.mag()*ptmean*ptmean) ;
277  mat = seg[last]->parametersError().similarityT( seg[last]->projectionMatrix() );
278  //mat[0][0]= 1.44 * p_err;
279  mat[0][0]= p_err;
280  }
281 
282  if ( debug ) {
283  GlobalPoint gp = seg[last]->globalPosition();
284  float Geta = gp.eta();
285 
286  std::cout << "Type " << type << " Nsegments " << layers.size() << " ";
287  std::cout << "pt " << ptmean << " errpt " << sptmean
288  << " eta " << Geta << " charge " << charge
289  << std::endl;
290  }
291 
292  // this perform H.T() * parErr * H, which is the projection of the
293  // the measurement error (rechit rf) to the state error (TSOS rf)
294  // Legend:
295  // H => is the 4x5 projection matrix
296  // parError the 4x4 parameter error matrix of the RecHit
297 
298 
299  LocalTrajectoryError error(asSMatrix<5>(mat));
300 
301  // Create the TrajectoryStateOnSurface
302  TrajectoryStateOnSurface tsos(param, error, seg[last]->det()->surface(),&*BField);
303 
304  // Take the DetLayer on which relies the segment
305  DetId id = seg[last]->geographicalId();
306 
307  // Transform it in a TrajectoryStateOnSurface
308 
309 
310  PTrajectoryStateOnDet seedTSOS = trajectoryStateTransform::persistentState( tsos, id.rawId());
311 
312  edm::OwnVector<TrackingRecHit> container;
313  for (unsigned l=0; l<seg.size(); l++) {
314  container.push_back( seg[l]->hit()->clone() );
315  //container.push_back(seg[l]->hit());
316  }
317 
318  TrajectorySeed theSeed(seedTSOS,container,alongMomentum);
319  return theSeed;
320 }
321 
322 
323 /*
324  * estimatePtCSC
325  *
326  * Look at delta phi to determine pt as:
327  * pt = (c_1 * eta + c_2) / dphi
328  *
329  * Note that only segment pairs with one segment in the ME1 layers give sensitive results
330  *
331  */
332 void MuonSeedCreator::estimatePtCSC(const SegmentContainer& seg, const std::vector<int>& layers, double& thePt, double& theSpt ) {
333 
334  unsigned size = seg.size();
335  if (size < 2) return;
336 
337  // reverse the segment and layer container first for pure CSC case
338  //if ( layers[0] > layers[ layers.size()-1 ] ) {
339  // reverse( layers.begin(), layers.end() );
340  // reverse( seg.begin(), seg.end() );
341  //}
342 
343  std::vector<double> ptEstimate;
344  std::vector<double> sptEstimate;
345 
346  thePt = defaultMomentum;
347  theSpt = defaultMomentum;
348 
349  double pt = 0.;
350  double spt = 0.;
351  GlobalPoint segPos[2];
352 
353  int layer0 = layers[0];
354  segPos[0] = seg[0]->globalPosition();
355  float eta = fabs( segPos[0].eta() );
356  //float corr = fabs( tan(segPos[0].theta()) );
357  // use pt from vertex information
358  /*
359  if ( layer0 == 0 ) {
360  SegmentContainer seg0;
361  seg0.push_back(seg[0]);
362  std::vector<int> lyr0(1,0);
363  estimatePtSingle( seg0, lyr0, thePt, theSpt);
364  ptEstimate.push_back( thePt );
365  sptEstimate.push_back( theSpt );
366  }
367  */
368 
369  //std::cout<<" estimate CSC "<<std::endl;
370 
371  unsigned idx1 = size;
372  if (size > 1) {
373  while ( idx1 > 1 ) {
374  idx1--;
375  int layer1 = layers[idx1];
376  if (layer0 == layer1) continue;
377  segPos[1] = seg[idx1]->globalPosition();
378 
379  double dphi = segPos[0].phi() - segPos[1].phi();
380  //double temp_dphi = dphi/corr;
381  double temp_dphi = dphi;
382 
383  double sign = 1.0;
384  if (temp_dphi < 0.) {
385  temp_dphi = -1.0*temp_dphi;
386  sign = -1.0;
387  }
388 
389  // Ensure that delta phi is not too small to prevent pt from blowing up
390  if (temp_dphi < 0.0001 ) {
391  temp_dphi = 0.0001 ;
392  pt = theMaxMomentum ;
393  spt = theMaxMomentum*0.25 ;
394  ptEstimate.push_back( pt*sign );
395  sptEstimate.push_back( spt );
396  }
397  // ME1 is inner-most
398  if ( layer0 == 0 && temp_dphi >= 0.0001 ) {
399 
400  // ME1/2 is outer-most
401  if ( layer1 == 1 ) {
402  //temp_dphi = scaledPhi(temp_dphi, CSC01_1[3] );
403  pt = getPt( CSC01, eta , temp_dphi )[0];
404  spt = getPt( CSC01, eta , temp_dphi )[1];
405  }
406  // ME2 is outer-most
407  else if ( layer1 == 2 ) {
408  //temp_dphi = scaledPhi(temp_dphi, CSC12_3[3] );
409  pt = getPt( CSC02, eta , temp_dphi )[0];
410  spt = getPt( CSC02, eta , temp_dphi )[1];
411  }
412  // ME3 is outer-most
413  else if ( layer1 == 3 ) {
414  //temp_dphi = scaledPhi(temp_dphi, CSC13_3[3] );
415  pt = getPt( CSC03, eta , temp_dphi )[0];
416  spt = getPt( CSC03, eta , temp_dphi )[1];
417  }
418  // ME4 is outer-most
419  else {
420  //temp_dphi = scaledPhi(temp_dphi, CSC14_3[3]);
421  pt = getPt( CSC14, eta , temp_dphi )[0];
422  spt = getPt( CSC14, eta , temp_dphi )[1];
423  }
424  ptEstimate.push_back( pt*sign );
425  sptEstimate.push_back( spt );
426  }
427 
428  // ME1/2,ME1/3 is inner-most
429  if ( layer0 == 1 ) {
430  // ME2 is outer-most
431  if ( layer1 == 2 ) {
432 
433  //if ( eta <= 1.2 ) { temp_dphi = scaledPhi(temp_dphi, CSC12_1[3]); }
434  //if ( eta > 1.2 ) { temp_dphi = scaledPhi(temp_dphi, CSC12_2[3]); }
435  pt = getPt( CSC12, eta , temp_dphi )[0];
436  spt = getPt( CSC12, eta , temp_dphi )[1];
437  }
438  // ME3 is outer-most
439  else if ( layer1 == 3 ) {
440  temp_dphi = scaledPhi(temp_dphi, CSC13_2[3]);
441  pt = getPt( CSC13, eta , temp_dphi )[0];
442  spt = getPt( CSC13, eta , temp_dphi )[1];
443  }
444  // ME4 is outer-most
445  else {
446  temp_dphi = scaledPhi(temp_dphi, CSC14_3[3]);
447  pt = getPt( CSC14, eta , temp_dphi )[0];
448  spt = getPt( CSC14, eta , temp_dphi )[1];
449  }
450  ptEstimate.push_back( pt*sign );
451  sptEstimate.push_back( spt );
452  }
453 
454  // ME2 is inner-most
455  if ( layer0 == 2 && temp_dphi > 0.0001 ) {
456 
457  // ME4 is outer-most
458  bool ME4av =false;
459  if ( layer1 == 4 ) {
460  temp_dphi = scaledPhi(temp_dphi, CSC24_1[3]);
461  pt = getPt( CSC24, eta , temp_dphi )[0];
462  spt = getPt( CSC24, eta , temp_dphi )[1];
463  ME4av = true;
464  }
465  // ME3 is outer-most
466  else {
467  // if ME2-4 is availabe , discard ME2-3
468  if ( !ME4av ) {
469  if ( eta <= 1.7 ) { temp_dphi = scaledPhi(temp_dphi, CSC23_1[3]); }
470  if ( eta > 1.7 ) { temp_dphi = scaledPhi(temp_dphi, CSC23_2[3]); }
471  pt = getPt( CSC23, eta , temp_dphi )[0];
472  spt = getPt( CSC23, eta , temp_dphi )[1];
473  }
474  }
475  ptEstimate.push_back( pt*sign );
476  sptEstimate.push_back( spt );
477  }
478 
479  // ME3 is inner-most
480  if ( layer0 == 3 && temp_dphi > 0.0001 ) {
481 
482  temp_dphi = scaledPhi(temp_dphi, CSC34_1[3]);
483  pt = getPt( CSC34, eta , temp_dphi )[0];
484  spt = getPt( CSC34, eta , temp_dphi )[1];
485  ptEstimate.push_back( pt*sign );
486  sptEstimate.push_back( spt );
487  }
488 
489  }
490  }
491 
492  // Compute weighted average if have more than one estimator
493  if ( ptEstimate.size() > 0 ) weightedPt( ptEstimate, sptEstimate, thePt, theSpt);
494 
495 }
496 
497 
498 /*
499  * estimatePtDT
500  *
501  * Look at delta phi between segments to determine pt as:
502  * pt = (c_1 * eta + c_2) / dphi
503  */
504 void MuonSeedCreator::estimatePtDT(const SegmentContainer& seg, const std::vector<int>& layers, double& thePt, double& theSpt) {
505 
506  unsigned size = seg.size();
507  if (size < 2) return;
508 
509  std::vector<double> ptEstimate;
510  std::vector<double> sptEstimate;
511 
512  thePt = defaultMomentum;
513  theSpt = defaultMomentum;
514 
515  double pt = 0.;
516  double spt = 0.;
517  GlobalPoint segPos[2];
518 
519  int layer0 = layers[0];
520  segPos[0] = seg[0]->globalPosition();
521  float eta = fabs(segPos[0].eta());
522 
523  //std::cout<<" estimate DT "<<std::endl;
524  // inner-most layer
525  //for ( unsigned idx0 = 0; idx0 < size-1; ++idx0 ) {
526  // layer0 = layers[idx0];
527  // segPos[0] = seg[idx0]->globalPosition();
528  // outer-most layer
529  // for ( unsigned idx1 = idx0+1; idx1 < size; ++idx1 ) {
530  for ( unsigned idx1 = 1; idx1 <size ; ++idx1 ) {
531 
532  int layer1 = layers[idx1];
533  segPos[1] = seg[idx1]->globalPosition();
534 
535  //eta = fabs(segPos[1].eta());
536  //if (layer1 == -4) eta = fabs(segPos[0].eta());
537 
538  double dphi = segPos[0].phi() - segPos[1].phi();
539  double temp_dphi = dphi;
540 
541  // Ensure that delta phi is not too small to prevent pt from blowing up
542 
543  double sign = 1.0;
544  if (temp_dphi < 0.) {
545  temp_dphi = -temp_dphi;
546  sign = -1.0;
547  }
548 
549  if (temp_dphi < 0.0001 ) {
550  temp_dphi = 0.0001 ;
551  pt = theMaxMomentum ;
552  spt = theMaxMomentum*0.25 ;
553  ptEstimate.push_back( pt*sign );
554  sptEstimate.push_back( spt );
555  }
556 
557  // MB1 is inner-most
558  bool MB23av = false;
559  if (layer0 == -1 && temp_dphi > 0.0001 ) {
560  // MB2 is outer-most
561  if (layer1 == -2) {
562 
563  if ( eta <= 0.7 ) { temp_dphi = scaledPhi(temp_dphi, DT12_1[3]); }
564  if ( eta > 0.7 ) { temp_dphi = scaledPhi(temp_dphi, DT12_2[3]); }
565  pt = getPt( DT12, eta , temp_dphi )[0];
566  spt = getPt( DT12, eta , temp_dphi )[1];
567  MB23av = true;
568  }
569  // MB3 is outer-most
570  else if (layer1 == -3) {
571 
572  if ( eta <= 0.6 ) { temp_dphi = scaledPhi(temp_dphi, DT13_1[3]); }
573  if ( eta > 0.6 ) { temp_dphi = scaledPhi(temp_dphi, DT13_2[3]); }
574  pt = getPt( DT13, eta , temp_dphi )[0];
575  spt = getPt( DT13, eta , temp_dphi )[1];
576  MB23av = true;
577  }
578  // MB4 is outer-most
579  else {
580  if ( !MB23av ) {
581  if ( eta <= 0.52 ) { temp_dphi = scaledPhi(temp_dphi, DT14_1[3]); }
582  if ( eta > 0.52 ) { temp_dphi = scaledPhi(temp_dphi, DT14_2[3]); }
583  pt = getPt( DT14, eta , temp_dphi )[0];
584  spt = getPt( DT14, eta , temp_dphi )[1];
585  }
586  }
587  ptEstimate.push_back( pt*sign );
588  sptEstimate.push_back( spt );
589  }
590 
591  // MB2 is inner-most
592  if (layer0 == -2 && temp_dphi > 0.0001 ) {
593  // MB3 is outer-most
594  if ( layer1 == -3) {
595 
596  if ( eta <= 0.6 ) { temp_dphi = scaledPhi(temp_dphi, DT23_1[3]); }
597  if ( eta > 0.6 ) { temp_dphi = scaledPhi(temp_dphi, DT23_2[3]); }
598  pt = getPt( DT23, eta , temp_dphi )[0];
599  spt = getPt( DT23, eta , temp_dphi )[1];
600  }
601  // MB4 is outer-most
602  else {
603 
604  if ( eta <= 0.52 ) { temp_dphi = scaledPhi(temp_dphi, DT24_1[3]); }
605  if ( eta > 0.52 ) { temp_dphi = scaledPhi(temp_dphi, DT24_2[3]); }
606  pt = getPt( DT24, eta , temp_dphi )[0];
607  spt = getPt( DT24, eta , temp_dphi )[1];
608  }
609  ptEstimate.push_back( pt*sign );
610  sptEstimate.push_back( spt );
611  }
612 
613  // MB3 is inner-most -> only marginally useful to pick up the charge
614  if (layer0 == -3 && temp_dphi > 0.0001 ) {
615  // MB4 is outer-most
616 
617  if ( eta <= 0.51 ) { temp_dphi = scaledPhi(temp_dphi, DT34_1[3]); }
618  if ( eta > 0.51 ) { temp_dphi = scaledPhi(temp_dphi, DT34_2[3]); }
619  pt = getPt( DT34, eta , temp_dphi )[0];
620  spt = getPt( DT34, eta , temp_dphi )[1];
621  ptEstimate.push_back( pt*sign );
622  sptEstimate.push_back( spt );
623  }
624  }
625  //}
626 
627 
628  // Compute weighted average if have more than one estimator
629  if (ptEstimate.size() > 0 ) weightedPt( ptEstimate, sptEstimate, thePt, theSpt);
630 
631 }
632 
633 
634 /*
635  * estimatePtOverlap
636  *
637  */
638 void MuonSeedCreator::estimatePtOverlap(const SegmentContainer& seg, const std::vector<int>& layers, double& thePt, double& theSpt) {
639 
640  int size = layers.size();
641 
642  thePt = defaultMomentum;
643  theSpt = defaultMomentum;
644 
645  SegmentContainer segCSC;
646  std::vector<int> layersCSC;
647  SegmentContainer segDT;
648  std::vector<int> layersDT;
649 
650  // DT layers are numbered as -4 to -1, whereas CSC layers go from 0 to 4:
651  for ( unsigned j = 0; j < layers.size(); ++j ) {
652  if ( layers[j] > -1 ) {
653  segCSC.push_back(seg[j]);
654  layersCSC.push_back(layers[j]);
655  }
656  else {
657  segDT.push_back(seg[j]);
658  layersDT.push_back(layers[j]);
659  }
660  }
661 
662  std::vector<double> ptEstimate;
663  std::vector<double> sptEstimate;
664 
665  GlobalPoint segPos[2];
666  int layer0 = layers[0];
667  segPos[0] = seg[0]->globalPosition();
668  float eta = fabs(segPos[0].eta());
669  //std::cout<<" estimate OL "<<std::endl;
670 
671  if ( segDT.size() > 0 && segCSC.size() > 0 ) {
672  int layer1 = layers[size-1];
673  segPos[1] = seg[size-1]->globalPosition();
674 
675  double dphi = segPos[0].phi() - segPos[1].phi();
676  double temp_dphi = dphi;
677 
678  // Ensure that delta phi is not too small to prevent pt from blowing up
679 
680  double sign = 1.0;
681  if (temp_dphi < 0.) {
682  temp_dphi = -temp_dphi;
683  sign = -1.0;
684  }
685 
686  if (temp_dphi < 0.0001 ) {
687  temp_dphi = 0.0001 ;
688  thePt = theMaxMomentum ;
689  theSpt = theMaxMomentum*0.25 ;
690  ptEstimate.push_back( thePt*sign );
691  sptEstimate.push_back( theSpt );
692  }
693 
694  // MB1 is inner-most
695  if ( layer0 == -1 && temp_dphi > 0.0001 ) {
696  // ME1/3 is outer-most
697  if ( layer1 == 1 ) {
698  temp_dphi = scaledPhi(temp_dphi, OL_1213[3]);
699  thePt = getPt( OL1213, eta , temp_dphi )[0];
700  theSpt = getPt( OL1213, eta , temp_dphi )[1];
701  }
702  // ME2 is outer-most
703  else if ( layer1 == 2) {
704  temp_dphi = scaledPhi(temp_dphi, OL_1222[3]);
705  thePt = getPt( OL1222, eta , temp_dphi )[0];
706  theSpt = getPt( OL1222, eta , temp_dphi )[1];
707  }
708  // ME3 is outer-most
709  else {
710  temp_dphi = scaledPhi(temp_dphi, OL_1232[3]);
711  thePt = getPt( OL1232, eta , temp_dphi )[0];
712  theSpt = getPt( OL1232, eta , temp_dphi )[1];
713  }
714  ptEstimate.push_back(thePt*sign);
715  sptEstimate.push_back(theSpt);
716  }
717  // MB2 is inner-most
718  if ( layer0 == -2 && temp_dphi > 0.0001 ) {
719  // ME1/3 is outer-most
720  if ( layer1 == 1 ) {
721  temp_dphi = scaledPhi(temp_dphi, OL_2213[3]);
722  thePt = getPt( OL2213, eta , temp_dphi )[0];
723  theSpt = getPt( OL2213, eta , temp_dphi )[1];
724  ptEstimate.push_back(thePt*sign);
725  sptEstimate.push_back(theSpt);
726  }
727  // ME2 is outer-most
728  if ( layer1 == 2) {
729  temp_dphi = scaledPhi(temp_dphi, OL_2222[3]);
730  thePt = getPt( OL2222, eta , temp_dphi )[0];
731  theSpt = getPt( OL2222, eta , temp_dphi )[1];
732  }
733  }
734  }
735 
736  if ( segDT.size() > 1 ) {
737  estimatePtDT(segDT, layersDT, thePt, theSpt);
738  ptEstimate.push_back(thePt);
739  sptEstimate.push_back(theSpt);
740  }
741 
742  /*
743  // not useful ....and pt estimation is bad
744  if ( segCSC.size() > 1 ) {
745  // don't estimate pt without ME1 information
746  bool CSCLayer1=false;
747  for (unsigned i=0; i< layersCSC.size(); i++) {
748  if ( layersCSC[i]==0 || layersCSC[i]==1 ) CSCLayer1 = true;
749  }
750  if (CSCLayer1) {
751  estimatePtCSC(segCSC, layersCSC, thePt, theSpt);
752  ptEstimate.push_back(thePt);
753  sptEstimate.push_back(theSpt);
754  }
755  }
756  */
757 
758  // Compute weighted average if have more than one estimator
759  if (ptEstimate.size() > 0 ) weightedPt( ptEstimate, sptEstimate, thePt, theSpt);
760 
761 }
762 /*
763  *
764  * estimate Pt for single segment events
765  *
766  */
767 void MuonSeedCreator::estimatePtSingle(const SegmentContainer& seg, const std::vector<int>& layers, double& thePt, double& theSpt) {
768 
769  thePt = defaultMomentum;
770  theSpt = defaultMomentum;
771 
772  GlobalPoint segPos = seg[0]->globalPosition();
773  double eta = segPos.eta();
774  GlobalVector gv = seg[0]->globalDirection();
775 
776  // Psi is angle between the segment origin and segment direction
777  // Use dot product between two vectors to get Psi in global x-y plane
778  double cosDpsi = (gv.x()*segPos.x() + gv.y()*segPos.y());
779  cosDpsi /= sqrt(segPos.x()*segPos.x() + segPos.y()*segPos.y());
780  cosDpsi /= sqrt(gv.x()*gv.x() + gv.y()*gv.y());
781 
782  double axb = ( segPos.x()*gv.y() ) - ( segPos.y()*gv.x() ) ;
783  double sign = (axb < 0.) ? 1.0 : -1.0;
784 
785  double dpsi = fabs(acos(cosDpsi)) ;
786  if ( dpsi > 1.570796 ) {
787  dpsi = 3.141592 - dpsi;
788  sign = -1.*sign ;
789  }
790  if (fabs(dpsi) < 0.00005) {
791  dpsi = 0.00005;
792  }
793 
794  // the 1st layer
795  if ( layers[0] == -1 ) {
796  // MB10
797  if ( fabs(eta) < 0.3 ) {
798  dpsi = scaledPhi(dpsi, SMB_10S[3] );
799  thePt = getPt( SMB10, eta , dpsi )[0];
800  theSpt = getPt( SMB10, eta , dpsi )[1];
801  }
802  // MB11
803  if ( fabs(eta) >= 0.3 && fabs(eta) < 0.82 ) {
804  dpsi = scaledPhi(dpsi, SMB_11S[3] );
805  thePt = getPt( SMB11, eta , dpsi )[0];
806  theSpt = getPt( SMB11, eta , dpsi )[1];
807  }
808  // MB12
809  if ( fabs(eta) >= 0.82 && fabs(eta) < 1.2 ) {
810  dpsi = scaledPhi(dpsi, SMB_12S[3] );
811  thePt = getPt( SMB12, eta , dpsi )[0];
812  theSpt = getPt( SMB12, eta , dpsi )[1];
813  }
814  }
815  if ( layers[0] == 1 ) {
816  // ME13
817  if ( fabs(eta) > 0.92 && fabs(eta) < 1.16 ) {
818  dpsi = scaledPhi(dpsi, SME_13S[3] );
819  thePt = getPt( SME13, eta , dpsi )[0];
820  theSpt = getPt( SME13, eta , dpsi )[1];
821  }
822  // ME12
823  if ( fabs(eta) >= 1.16 && fabs(eta) <= 1.6 ) {
824  dpsi = scaledPhi(dpsi, SME_12S[3] );
825  thePt = getPt( SME12, eta , dpsi )[0];
826  theSpt = getPt( SME12, eta , dpsi )[1];
827  }
828  }
829  if ( layers[0] == 0 ) {
830  // ME11
831  if ( fabs(eta) > 1.6 ) {
832  dpsi = scaledPhi(dpsi, SMB_11S[3] );
833  thePt = getPt( SME11, eta , dpsi )[0];
834  theSpt = getPt( SME11, eta , dpsi )[1];
835  }
836  }
837  // the 2nd layer
838  if ( layers[0] == -2 ) {
839  // MB20
840  if ( fabs(eta) < 0.25 ) {
841  dpsi = scaledPhi(dpsi, SMB_20S[3] );
842  thePt = getPt( SMB20, eta , dpsi )[0];
843  theSpt = getPt( SMB20, eta , dpsi )[1];
844  }
845  // MB21
846  if ( fabs(eta) >= 0.25 && fabs(eta) < 0.72 ) {
847  dpsi = scaledPhi(dpsi, SMB_21S[3] );
848  thePt = getPt( SMB21, eta , dpsi )[0];
849  theSpt = getPt( SMB21, eta , dpsi )[1];
850  }
851  // MB22
852  if ( fabs(eta) >= 0.72 && fabs(eta) < 1.04 ) {
853  dpsi = scaledPhi(dpsi, SMB_22S[3] );
854  thePt = getPt( SMB22, eta , dpsi )[0];
855  theSpt = getPt( SMB22, eta , dpsi )[1];
856  }
857  }
858  if ( layers[0] == 2 ) {
859  // ME22
860  if ( fabs(eta) > 0.95 && fabs(eta) <= 1.6 ) {
861  dpsi = scaledPhi(dpsi, SME_22S[3] );
862  thePt = getPt( SME22, eta , dpsi )[0];
863  theSpt = getPt( SME22, eta , dpsi )[1];
864  }
865  // ME21
866  if ( fabs(eta) > 1.6 && fabs(eta) < 2.45 ) {
867  dpsi = scaledPhi(dpsi, SME_21S[3] );
868  thePt = getPt( SME21, eta , dpsi )[0];
869  theSpt = getPt( SME21, eta , dpsi )[1];
870  }
871  }
872 
873  // the 3rd layer
874  if ( layers[0] == -3 ) {
875  // MB30
876  if ( fabs(eta) <= 0.22 ) {
877  dpsi = scaledPhi(dpsi, SMB_30S[3] );
878  thePt = getPt( SMB30, eta , dpsi )[0];
879  theSpt = getPt( SMB30, eta , dpsi )[1];
880  }
881  // MB31
882  if ( fabs(eta) > 0.22 && fabs(eta) <= 0.6 ) {
883  dpsi = scaledPhi(dpsi, SMB_31S[3] );
884  thePt = getPt( SMB31, eta , dpsi )[0];
885  theSpt = getPt( SMB31, eta , dpsi )[1];
886  }
887  // MB32
888  if ( fabs(eta) > 0.6 && fabs(eta) < 0.95 ) {
889  dpsi = scaledPhi(dpsi, SMB_32S[3] );
890  thePt = getPt( SMB32, eta , dpsi )[0];
891  theSpt = getPt( SMB32, eta , dpsi )[1];
892  }
893  }
894  thePt = fabs(thePt)*sign;
895  theSpt = fabs(theSpt);
896 
897  return;
898 }
899 
900 // setup the minimum value for obvious showering cases
901 void MuonSeedCreator::estimatePtShowering(int& NShowers, int& NShowerSegments, double& thePt, double& theSpt) {
902 
903  if ( NShowers > 2 && thePt < 300. ) {
904  thePt = 800. ;
905  theSpt = 200. ;
906  }
907  if ( NShowers == 2 && NShowerSegments > 11 && thePt < 150. ) {
908  thePt = 280. ;
909  theSpt = 70. ;
910  }
911  if ( NShowers == 2 && NShowerSegments <= 11 && thePt < 50. ) {
912  thePt = 80.;
913  theSpt = 40. ;
914  }
915  if ( NShowers == 1 && NShowerSegments <= 5 && thePt < 10. ) {
916  thePt = 16. ;
917  theSpt = 8. ;
918  }
919 
920 }
921 
922 /*
923  * weightedPt
924  *
925  * Look at delta phi between segments to determine pt as:
926  * pt = (c_1 * eta + c_2) / dphi
927  */
928 void MuonSeedCreator::weightedPt(const std::vector<double>& ptEstimate, const std::vector<double>& sptEstimate, double& thePt, double& theSpt) {
929 
930 
931  int size = ptEstimate.size();
932 
933  // If only one element, by-pass computation below
934  if (size < 2) {
935  thePt = ptEstimate[0];
936  theSpt = sptEstimate[0];
937  return;
938  }
939 
940  double charge = 0.;
941  // If have more than one pt estimator, first look if any estimator is biased
942  // by comparing the charge of each estimator
943 
944  for ( unsigned j = 0; j < ptEstimate.size(); j++ ) {
945  //std::cout<<" weighting pt: "<< ptEstimate[j] <<std::endl;
946  if ( ptEstimate[j] < 0. ) {
947  // To prevent from blowing up, add 0.1
948  charge -= 1. * (ptEstimate[j]*ptEstimate[j]) / (sptEstimate[j]*sptEstimate[j] ); // weight by relative error on pt
949  } else {
950  charge += 1. * (ptEstimate[j]*ptEstimate[j]) / (sptEstimate[j]*sptEstimate[j] ); // weight by relative error on pt
951  }
952  }
953 
954  // No need to normalize as we want to know only sign ( + or - )
955  if (charge < 0.) {
956  charge = -1.;
957  } else {
958  charge = 1.;
959  }
960 
961  //int n = 0;
962  double weightPtSum = 0.;
963  double sigmaSqr_sum = 0.;
964 
965  // Now, we want to compute average Pt using estimators with "correct" charge
966  // This is to remove biases
967  for ( unsigned j = 0; j < ptEstimate.size(); ++j ) {
968  //if ( (minpt_ratio < 0.5) && (fabs(ptEstimate[j]) < 5.0) ) continue;
969  //if ( ptEstimate[j] * charge > 0. ) {
970  //n++;
971  sigmaSqr_sum += 1.0 / (sptEstimate[j]*sptEstimate[j]);
972  weightPtSum += fabs(ptEstimate[j])/(sptEstimate[j]*sptEstimate[j]);
973  //}
974  }
975  /*
976  if (n < 1) {
977  thePt = defaultMomentum*charge;
978  theSpt = defaultMomentum;
979  return;
980  }
981  */
982  // Compute weighted mean and error
983 
984  thePt = (charge*weightPtSum) / sigmaSqr_sum;
985  theSpt = sqrt( 1.0 / sigmaSqr_sum ) ;
986 
987  //std::cout<<" final weighting : "<< thePt <<" ~ "<< fabs( theSpt/thePt ) <<std::endl;
988 
989  return;
990 }
991 
992 std::vector<double> MuonSeedCreator::getPt(const std::vector<double>& vPara, double eta, double dPhi ) {
993 
994  double h = fabs(eta);
995  double estPt = ( vPara[0] + vPara[1]*h + vPara[2]*h*h ) / dPhi;
996  double estSPt = ( vPara[3] + vPara[4]*h + vPara[5]*h*h ) * estPt;
997  std::vector<double> paraPt ;
998  paraPt.push_back( estPt );
999  paraPt.push_back( estSPt ) ;
1000 
1001  //std::cout<<" pt:"<<estPt<<" +/-"<< estSPt<<" h:"<<eta<<" df:"<<dPhi<<std::endl;
1002  return paraPt ;
1003 }
1004 
1005 double MuonSeedCreator::scaledPhi( double dphi, double t1) {
1006 
1007  if (dphi != 0. ) {
1008 
1009  double oPhi = 1./dphi ;
1010  dphi = dphi /( 1. + t1/( oPhi + 10. ) ) ;
1011  return dphi ;
1012 
1013  } else {
1014  return dphi ;
1015  }
1016 
1017 }
1018 
MuonSeedCreator::OL2213
std::vector< double > OL2213
Definition: MuonSeedCreator.h:113
type
type
Definition: HCALResponse.h:21
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
MuonSeedCreator::DT13_2
std::vector< double > DT13_2
Definition: MuonSeedCreator.h:152
i
int i
Definition: DBlmapReader.cc:9
MuonSeedCreator::SMB32
std::vector< double > SMB32
Definition: MuonSeedCreator.h:133
MuonSeedCreator::CSC23_1
std::vector< double > CSC23_1
Definition: MuonSeedCreator.h:144
MuonSeedCreator::SMB_30S
std::vector< double > SMB_30S
Definition: MuonSeedCreator.h:174
MuonSeedCreator::SMB_21S
std::vector< double > SMB_21S
Definition: MuonSeedCreator.h:172
MuonSeedCreator::DT14_2
std::vector< double > DT14_2
Definition: MuonSeedCreator.h:154
LayerTriplets::layers
std::vector< LayerSetAndLayers > layers(const SeedingLayerSetsHits &sets)
Definition: LayerTriplets.cc:4
MuonSeedCreator::DT23_2
std::vector< double > DT23_2
Definition: MuonSeedCreator.h:156
MuonSeedCreator::DT34_1
std::vector< double > DT34_1
Definition: MuonSeedCreator.h:159
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
MuonSeedCreator::DT14_1
std::vector< double > DT14_1
Definition: MuonSeedCreator.h:153
MuonSeedCreator::SME_22S
std::vector< double > SME_22S
Definition: MuonSeedCreator.h:182
MuonSeedCreator::SMB22
std::vector< double > SMB22
Definition: MuonSeedCreator.h:130
MuonSeedCreator::OL1232
std::vector< double > OL1232
Definition: MuonSeedCreator.h:112
MuonSeedCreator::DT34_2
std::vector< double > DT34_2
Definition: MuonSeedCreator.h:160
MuonSeedCreator::CSC24
std::vector< double > CSC24
Definition: MuonSeedCreator.h:107
MuonSeedCreator::estimatePtDT
void estimatePtDT(const SegmentContainer &seg, const std::vector< int > &layers, double &pt, double &spt)
Estimate transverse momentum of track from DT measurements.
Definition: MuonSeedCreator.cc:504
MuonSeedCreator::weightedPt
void weightedPt(const std::vector< double > &ptEstimate, const std::vector< double > &sptEstimate, double &ptAvg, double &sptAvg)
Compute weighted mean pt from different pt estimators.
Definition: MuonSeedCreator.cc:928
MuonSeedCreator::CSC23_2
std::vector< double > CSC23_2
Definition: MuonSeedCreator.h:145
MuonSeedCreator::SMB20
std::vector< double > SMB20
Definition: MuonSeedCreator.h:128
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
GlobalPoint
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
MuonSeedCreator::CSC14
std::vector< double > CSC14
Definition: MuonSeedCreator.h:105
trajectoryStateTransform::persistentState
PTrajectoryStateOnDet persistentState(const TrajectoryStateOnSurface &ts, unsigned int detid)
Definition: TrajectoryStateTransform.cc:16
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
MuonSeedCreator::DT34
std::vector< double > DT34
Definition: MuonSeedCreator.h:98
MuonSeedCreator::CSC12_2
std::vector< double > CSC12_2
Definition: MuonSeedCreator.h:139
MuonSeedCreator::CSC13_2
std::vector< double > CSC13_2
Definition: MuonSeedCreator.h:141
MuonSeedCreator::CSC12_1
std::vector< double > CSC12_1
Definition: MuonSeedCreator.h:138
MuonSeedCreator::SME12
std::vector< double > SME12
Definition: MuonSeedCreator.h:117
MuonSeedCreator::SME22
std::vector< double > SME22
Definition: MuonSeedCreator.h:120
MuonSeedCreator::createSeed
TrajectorySeed createSeed(int type, const SegmentContainer &seg, const std::vector< int > &layers, int NShower, int NShowerSeg)
Create a seed from set of segments.
Definition: MuonSeedCreator.cc:153
MuonSeedCreator::SegmentContainer
MuonTransientTrackingRecHit::MuonRecHitContainer SegmentContainer
Definition: MuonSeedCreator.h:34
eta
T eta() const
Definition: Basic3DVectorLD.h:177
MuonSeedCreator::SME21
std::vector< double > SME21
Definition: MuonSeedCreator.h:119
DeDxDiscriminatorTools::charge
double charge(const std::vector< uint8_t > &Ampls)
Definition: DeDxDiscriminatorTools.cc:43
MuonSeedCreator::estimatePtOverlap
void estimatePtOverlap(const SegmentContainer &seg, const std::vector< int > &layers, double &pt, double &spt)
Estimate transverse momentum of track from CSC + DT measurements.
Definition: MuonSeedCreator.cc:638
MuonSeedCreator::CSC13_3
std::vector< double > CSC13_3
Definition: MuonSeedCreator.h:142
MuonSeedCreator::MuonSeedCreator
MuonSeedCreator(const edm::ParameterSet &pset)
Constructor.
Definition: MuonSeedCreator.cc:37
MuonSeedCreator::BField
const MagneticField * BField
Definition: MuonSeedCreator.h:90
MuonSeedCreator::estimatePtSingle
void estimatePtSingle(const SegmentContainer &seg, const std::vector< int > &layers, double &pt, double &spt)
Estimate transverse momentum of track from single segment.
Definition: MuonSeedCreator.cc:767
MuonSeedCreator::CSC34
std::vector< double > CSC34
Definition: MuonSeedCreator.h:108
MuonSeedCreator::estimatePtCSC
void estimatePtCSC(const SegmentContainer &seg, const std::vector< int > &layers, double &pt, double &spt)
Estimate transverse momentum of track from CSC measurements.
Definition: MuonSeedCreator.cc:332
MuonSeedCreator::SME13
std::vector< double > SME13
Definition: MuonSeedCreator.h:118
PV3DBase::theta
Geom::Theta< T > theta() const
Definition: PV3DBase.h:75
MuonSeedCreator::DT13
std::vector< double > DT13
Definition: MuonSeedCreator.h:94
MuonSeedCreator::CSC12_3
std::vector< double > CSC12_3
Definition: MuonSeedCreator.h:140
edm::OwnVector::push_back
void push_back(D *&d)
Definition: OwnVector.h:274
MuonSeedCreator::SMB10
std::vector< double > SMB10
Definition: MuonSeedCreator.h:125
MuonSeedCreator::CSC12
std::vector< double > CSC12
Definition: MuonSeedCreator.h:101
MuonSeedCreator::DT12_2
std::vector< double > DT12_2
Definition: MuonSeedCreator.h:150
MuonSeedCreator::OL_1222
std::vector< double > OL_1222
Definition: MuonSeedCreator.h:163
MuonSeedCreator::OL1222
std::vector< double > OL1222
Definition: MuonSeedCreator.h:111
dPhi
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
MuonSeedCreator::estimatePtShowering
void estimatePtShowering(int &NShowers, int &NShowerSeg, double &pt, double &spt)
Estimate transverse momentum of track from showering segment.
Definition: MuonSeedCreator.cc:901
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
MuonSeedCreator::CSC23
std::vector< double > CSC23
Definition: MuonSeedCreator.h:106
MuonSeedCreator::CSC13
std::vector< double > CSC13
Definition: MuonSeedCreator.h:103
MuonSeedCreator::DT23_1
std::vector< double > DT23_1
Definition: MuonSeedCreator.h:155
funct::tan
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
MuonSeedCreator::SMB30
std::vector< double > SMB30
Definition: MuonSeedCreator.h:131
j
int j
Definition: DBlmapReader.cc:9
MuonSeedCreator::SMB_10S
std::vector< double > SMB_10S
Definition: MuonSeedCreator.h:168
h
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
MuonSeedCreator::DT24_2
std::vector< double > DT24_2
Definition: MuonSeedCreator.h:158
MuonSeedCreator::SMB_31S
std::vector< double > SMB_31S
Definition: MuonSeedCreator.h:175
MuonSeedCreator::SMB_32S
std::vector< double > SMB_32S
Definition: MuonSeedCreator.h:176
MuonSeedCreator::DT12_1
std::vector< double > DT12_1
Definition: MuonSeedCreator.h:149
MuonSeedCreator::CSC14_3
std::vector< double > CSC14_3
Definition: MuonSeedCreator.h:143
MuonSeedCreator::DT12
std::vector< double > DT12
Definition: MuonSeedCreator.h:93
MuonSeedCreator::SMB_20S
std::vector< double > SMB_20S
Definition: MuonSeedCreator.h:171
MuonSeedCreator::SMB_11S
std::vector< double > SMB_11S
Definition: MuonSeedCreator.h:169
MuonSeedCreator::OL1213
std::vector< double > OL1213
Definition: MuonSeedCreator.h:110
MuonSeedCreator::getPt
std::vector< double > getPt(const std::vector< double > &vParameters, double eta, double dPhi)
Compute pt from parameters.
Definition: MuonSeedCreator.cc:992
DetId
Definition: DetId.h:18
AlgebraicVector
CLHEP::HepVector AlgebraicVector
Definition: AlgebraicObjects.h:12
MuonSeedCreator::SMB12
std::vector< double > SMB12
Definition: MuonSeedCreator.h:127
MuonSeedCreator::SME_21S
std::vector< double > SME_21S
Definition: MuonSeedCreator.h:181
MuonSeedCreator::SME_12S
std::vector< double > SME_12S
Definition: MuonSeedCreator.h:179
MuonSeedCreator::OL2222
std::vector< double > OL2222
Definition: MuonSeedCreator.h:114
clone
TEveGeoShape * clone(const TEveElement *element, TEveElement *parent)
Definition: eve_macros.cc:135
MuonSeedCreator::SME_11S
std::vector< double > SME_11S
Definition: MuonSeedCreator.h:178
MuonSeedCreator::DT14
std::vector< double > DT14
Definition: MuonSeedCreator.h:95
PV3DBase::eta
T eta() const
Definition: PV3DBase.h:76
MuonSeedCreator::DT13_1
std::vector< double > DT13_1
Definition: MuonSeedCreator.h:151
MuonSeedCreator::SME31
std::vector< double > SME31
Definition: MuonSeedCreator.h:121
MuonSeedCreator::DT24
std::vector< double > DT24
Definition: MuonSeedCreator.h:97
MuonSeedCreator::OL_2213
std::vector< double > OL_2213
Definition: MuonSeedCreator.h:165
MuonSeedCreator::DT23
std::vector< double > DT23
Definition: MuonSeedCreator.h:96
MuonSeedCreator::SME32
std::vector< double > SME32
Definition: MuonSeedCreator.h:122
MuonSeedCreator::CSC24_1
std::vector< double > CSC24_1
Definition: MuonSeedCreator.h:146
MuonSeedCreator::SMB21
std::vector< double > SMB21
Definition: MuonSeedCreator.h:129
AlgebraicSymMatrix
CLHEP::HepSymMatrix AlgebraicSymMatrix
Definition: AlgebraicObjects.h:14
MuonSeedCreator::CSC01_1
std::vector< double > CSC01_1
Definition: MuonSeedCreator.h:137
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
MuonSeedCreator::scaledPhi
double scaledPhi(double dphi, double t1)
Scale the dPhi from segment position.
Definition: MuonSeedCreator.cc:1005
MuonSeedCreator::~MuonSeedCreator
~MuonSeedCreator()
Destructor.
Definition: MuonSeedCreator.cc:141
MuonSeedCreator::SME_13S
std::vector< double > SME_13S
Definition: MuonSeedCreator.h:180
MuonSeedCreator::OL_2222
std::vector< double > OL_2222
Definition: MuonSeedCreator.h:166
MuonSeedCreator::SMB31
std::vector< double > SMB31
Definition: MuonSeedCreator.h:132
MuonSeedCreator::CSC03
std::vector< double > CSC03
Definition: MuonSeedCreator.h:104
MuonSeedCreator::CSC02
std::vector< double > CSC02
Definition: MuonSeedCreator.h:102
MuonSeedCreator::DT24_1
std::vector< double > DT24_1
Definition: MuonSeedCreator.h:157
MuonSeedCreator::CSC34_1
std::vector< double > CSC34_1
Definition: MuonSeedCreator.h:147
MuonSeedCreator::SME11
std::vector< double > SME11
Definition: MuonSeedCreator.h:116
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
MuonSeedCreator::SMB11
std::vector< double > SMB11
Definition: MuonSeedCreator.h:126
MuonSeedCreator::OL_1232
std::vector< double > OL_1232
Definition: MuonSeedCreator.h:164
MuonSeedCreator::SMB_12S
std::vector< double > SMB_12S
Definition: MuonSeedCreator.h:170
MuonSeedCreator::SME41
std::vector< double > SME41
Definition: MuonSeedCreator.h:123
findQualityFiles.size
tuple size
Write out results.
Definition: findQualityFiles.py:442
CSC01
Definition: QTest.h:566
cuy.dh
tuple dh
Definition: cuy.py:353
MuonSeedCreator::OL_1213
std::vector< double > OL_1213
Definition: MuonSeedCreator.h:162
MuonSeedCreator::SMB_22S
std::vector< double > SMB_22S
Definition: MuonSeedCreator.h:173
Generated for CMSSW Reference Manual by doxygen 1.8.5