CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFTrackTransformer.cc
Go to the documentation of this file.
1 //
2 // -*- C++ -*-
3 // Package: PFTracking
4 // Class: PFTrackTransformer
5 //
6 // Original Author: Michele Pioppi
7 // Other Author: Daniele Benedetti
8 
9 #include "RecoParticleFlow/PFTracking/interface/PFTrackTransformer.h"
10 
11 #include "FWCore/MessageLogger/interface/MessageLogger.h"
12 #include "FWCore/Framework/interface/ESHandle.h"
13 
14 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
15 #include "DataFormats/ParticleFlowReco/interface/PFRecTrackFwd.h"
16 #include "DataFormats/TrackReco/interface/Track.h"
17 
18 #include "TrackingTools/PatternTools/interface/Trajectory.h"
19 #include "FastSimulation/BaseParticlePropagator/interface/BaseParticlePropagator.h"
20 
21 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
22 // Add by Daniele
23 #include "TrackingTools/GsfTools/interface/MultiGaussianStateTransform.h"
24 #include "TrackingTools/GsfTools/interface/MultiGaussianState1D.h"
25 #include "TrackingTools/GsfTools/interface/GaussianSumUtilities1D.h"
26 #include "RecoParticleFlow/PFTracking/interface/PFGsfHelper.h"
27 
28 using namespace std;
29 using namespace reco;
30 using namespace edm;
31 
32 
33 
34 PFTrackTransformer::PFTrackTransformer(const math::XYZVector& B):B_(B){
35  LogInfo("PFTrackTransformer")<<"PFTrackTransformer built";
36 
37  onlyprop_=false;
38 }
39 
40 PFTrackTransformer::~PFTrackTransformer(){
41 
42 }
43 
44 
45 bool
46 PFTrackTransformer::addPoints( reco::PFRecTrack& pftrack,
47  const reco::Track& track,
48  const Trajectory& traj,
49  bool msgwarning) const {
50 
51  LogDebug("PFTrackTransformer")<<"Trajectory propagation started";
52  using namespace reco;
53  using namespace std;
54 
55  float PT= track.pt();
56  float pfmass= (pftrack.algoType()==reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
57  float pfenergy=sqrt((pfmass*pfmass)+(track.p()*track.p()));
58  // closest approach
59  BaseParticlePropagator theParticle =
60  BaseParticlePropagator(
61  RawParticle(XYZTLorentzVector(track.px(),
62  track.py(),
63  track.pz(),
64  pfenergy),
65  XYZTLorentzVector(track.vertex().x(),
66  track.vertex().y(),
67  track.vertex().z(),
68  0.)),
69  0.,0.,B_.z());
70 
71  theParticle.setCharge(track.charge());
72  float pfoutenergy=sqrt((pfmass*pfmass)+track.outerMomentum().Mag2());
73  BaseParticlePropagator theOutParticle =
74  BaseParticlePropagator(
75  RawParticle(XYZTLorentzVector(track.outerMomentum().x(),
76  track.outerMomentum().y(),
77  track.outerMomentum().z(),
78  pfoutenergy),
79  XYZTLorentzVector(track.outerPosition().x(),
80  track.outerPosition().y(),
81  track.outerPosition().z(),
82  0.)),
83  0.,0.,B_.z());
84  theOutParticle.setCharge(track.charge());
85 
86 
87  math::XYZTLorentzVector momClosest
88  = math::XYZTLorentzVector(track.px(), track.py(),
89  track.pz(), track.p());
90  math::XYZPoint posClosest = track.vertex();
91 
92  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ClosestApproach,
93  posClosest,momClosest));
94 
95 
96  //BEAMPIPE
97  theParticle.setPropagationConditions(pfGeometry_.outerRadius(PFGeometry::BeamPipe),
98  pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
99  theParticle.propagate();
100  if(theParticle.getSuccess()!=0)
101  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::BeamPipeOrEndVertex,
102  math::XYZPoint(theParticle.vertex()),
103  math::XYZTLorentzVector(theParticle.momentum())));
104  else {
105  PFTrajectoryPoint dummyMaxSh;
106  pftrack.addPoint(dummyMaxSh);
107  }
108 
109 
110 
111  //trajectory points
112 
113  if (!onlyprop_){
114  bool direction =(traj.direction() == alongMomentum);
115  vector<TrajectoryMeasurement> measurements =traj.measurements();
116  int iTrajFirst = (direction) ? 0 : measurements.size() - 1;
117  int increment = (direction) ? +1 : -1;
118  int iTrajLast = (direction) ? int(measurements.size()) : -1;
119 
120 
121  for (int iTraj = iTrajFirst; iTraj != iTrajLast; iTraj += increment) {
122  GlobalPoint v=measurements[iTraj].updatedState().globalPosition();
123  GlobalVector p=measurements[iTraj].updatedState().globalMomentum();
124  unsigned int iid=measurements[iTraj].recHit()->det()->geographicalId().rawId();
125  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
126  math::XYZPoint(v.x(), v.y(), v.z()),
127  math::XYZTLorentzVector(p.x(),p.y(),p.z(),p.mag())));
128  }
129  }
130 
131  bool isBelowPS=false;
132  theOutParticle.propagateToPreshowerLayer1(false);
133  if(theOutParticle.getSuccess()!=0)
134  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
135  math::XYZPoint(theOutParticle.vertex()),
136  math::XYZTLorentzVector(theOutParticle.momentum())));
137  else {
138  PFTrajectoryPoint dummyPS1;
139  pftrack.addPoint(dummyPS1);
140  }
141 
142 
143  theOutParticle.propagateToPreshowerLayer2(false);
144  if(theOutParticle.getSuccess()!=0){
145  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
146  math::XYZPoint(theOutParticle.vertex()),
147  math::XYZTLorentzVector(theOutParticle.momentum())));
148  isBelowPS=true;
149  } else {
150  PFTrajectoryPoint dummyPS2;
151  pftrack.addPoint(dummyPS2);
152  }
153 
154  theOutParticle.propagateToEcalEntrance(false);
155 
156  if(theOutParticle.getSuccess()!=0){
157  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
158  math::XYZPoint(theOutParticle.vertex()),
159  math::XYZTLorentzVector(theOutParticle.momentum())));
160  double ecalShowerDepth
161  = PFCluster::getDepthCorrection(theOutParticle.momentum().E(),
162  isBelowPS,
163  false);
164 
165  math::XYZPoint meanShower=math::XYZPoint(theOutParticle.vertex())+
166  math::XYZTLorentzVector(theOutParticle.momentum()).Vect().Unit()*ecalShowerDepth;
167 
168  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
169  meanShower,
170  math::XYZTLorentzVector(theOutParticle.momentum())));}
171  else {
172  if (PT>5. && msgwarning)
173  LogWarning("PFTrackTransformer")<<"KF TRACK "<<pftrack<< " PROPAGATION TO THE ECAL HAS FAILED";
174  PFTrajectoryPoint dummyECAL;
175  pftrack.addPoint(dummyECAL);
176  PFTrajectoryPoint dummyMaxSh;
177  pftrack.addPoint(dummyMaxSh);
178  }
179 
180 
181 
182  //HCAL entrance
183  theOutParticle.propagateToHcalEntrance(false);
184  if(theOutParticle.getSuccess()!=0)
185  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
186  math::XYZPoint(theOutParticle.vertex()),
187  math::XYZTLorentzVector(theOutParticle.momentum())));
188  else{
189  if (PT>5.&& msgwarning)
190  LogWarning("PFTrackTransformer")<<"KF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
191  PFTrajectoryPoint dummyHCALentrance;
192  pftrack.addPoint(dummyHCALentrance);
193  }
194 
195  //HCAL exit
196  // theOutParticle.setMagneticField(0); //Show we propagate as straight line inside HCAL ?
197  theOutParticle.propagateToHcalExit(false);
198  if(theOutParticle.getSuccess()!=0)
199  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
200  math::XYZPoint(theOutParticle.vertex()),
201  math::XYZTLorentzVector(theOutParticle.momentum())));
202  else{
203  if (PT>5.&& msgwarning)
204  LogWarning("PFTrackTransformer")<<"KF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL EXIT HAS FAILED";
205  PFTrajectoryPoint dummyHCALexit;
206  pftrack.addPoint(dummyHCALexit);
207  }
208 
209 
210  //HO layer0
211  // if (abs(theOutParticle.vertex().z())<550) {
212  if ( PT>3.0) { //Same value is used in PFBlockAlgo::link( case PFBlockLink::TRACKandHO:
213  theOutParticle.setMagneticField(0);
214  theOutParticle.setCharge(0);
215  theOutParticle.propagateToHOLayer(false);
216  if(theOutParticle.getSuccess()!=0) {
217  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
218  math::XYZPoint(theOutParticle.vertex()),
219  math::XYZTLorentzVector(theOutParticle.momentum())));
220  } else {
221  if (PT>5. && abs(theOutParticle.Z()) < 700.25 && msgwarning)
222  LogWarning("PFTrackTransformer")<<"KF TRACK "<<pftrack<< " PROPAGATION TO THE HO HAS FAILED";
223  PFTrajectoryPoint dummyHOLayer;
224  pftrack.addPoint(dummyHOLayer);
225  }
226  }
227 
228  return true;
229 }
230 bool
231 PFTrackTransformer::addPointsAndBrems( reco::GsfPFRecTrack& pftrack,
232  const reco::Track& track,
233  const Trajectory& traj,
234  const bool& GetMode) const {
235 
236  float PT= track.pt();
237  // Trajectory for each trajectory point
238 
239  bool direction =(traj.direction() == alongMomentum);
240  vector<TrajectoryMeasurement> measurements =traj.measurements();
241  int iTrajFirst = (direction) ? 0 : measurements.size() - 1;
242  int increment = (direction) ? +1 : -1;
243  int iTrajLast = (direction) ? int(measurements.size()) : -1;
244 
245 
246  unsigned int iTrajPos = 0;
247  for (int iTraj = iTrajFirst; iTraj != iTrajLast; iTraj += increment) {
248 
249  GlobalPoint v=measurements[iTraj].updatedState().globalPosition();
250  PFGsfHelper* PFGsf = new PFGsfHelper(measurements[iTraj]);
251  //if (PFGsf->isValid()){
252  bool ComputeMODE = GetMode;
253  GlobalVector p = PFGsf->computeP(ComputeMODE);
254  double DP = PFGsf->fittedDP();
255  double SigmaDP = PFGsf->sigmafittedDP();
256  unsigned int iid=measurements[iTraj].recHit()->det()->geographicalId().rawId();
257  delete PFGsf;
258 
259  // -------------------------- Fill GSF Track -------------------------------------
260 
261 
262  // float pfmass= (pftrack.algoType()==reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
263  float ptot = sqrt((p.x()*p.x())+(p.y()*p.y())+(p.z()*p.z()));
264  float pfenergy= ptot;
265 
266  if (iTraj == iTrajFirst) {
267 
268  math::XYZTLorentzVector momClosest
269  = math::XYZTLorentzVector(p.x(), p.y(),
270  p.z(), ptot);
271  math::XYZPoint posClosest = track.vertex();
272  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ClosestApproach,
273  posClosest,momClosest));
274 
275  BaseParticlePropagator theInnerParticle =
276  BaseParticlePropagator(
277  RawParticle(XYZTLorentzVector(p.x(),
278  p.y(),
279  p.z(),
280  pfenergy),
281  XYZTLorentzVector(track.vertex().x(),
282  track.vertex().y(),
283  track.vertex().z(),
284  0.)), //DANIELE Same thing v.x(),v.y(),v.()?
285  0.,0.,B_.z());
286  theInnerParticle.setCharge(track.charge());
287 
288  //BEAMPIPE
289  theInnerParticle.setPropagationConditions(pfGeometry_.outerRadius(PFGeometry::BeamPipe),
290  pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
291  theInnerParticle.propagate();
292  if(theInnerParticle.getSuccess()!=0)
293  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::BeamPipeOrEndVertex,
294  math::XYZPoint(theInnerParticle.vertex()),
295  math::XYZTLorentzVector(theInnerParticle.momentum())));
296  else {
297  PFTrajectoryPoint dummyMaxSh;
298  pftrack.addPoint(dummyMaxSh);
299  }
300 
301  // First Point for the trajectory == Vertex ??
302  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
303  math::XYZPoint(v.x(), v.y(), v.z()),
304  math::XYZTLorentzVector(p.x(),p.y(),p.z(),p.mag())));
305 
306 
307  }
308  if (iTraj != iTrajFirst && iTraj != (abs(iTrajLast)-1)) {
309  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
310  math::XYZPoint(v.x(), v.y(), v.z()),
311  math::XYZTLorentzVector(p.x(),p.y(),p.z(),p.mag())));
312 
313 
314  }
315  if (iTraj == (abs(iTrajLast)-1)) {
316 
317  // Last Trajectory Meas
318  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
319  math::XYZPoint(v.x(), v.y(), v.z()),
320  math::XYZTLorentzVector(p.x(),p.y(),p.z(),p.mag())));
321 
322 
323 
324 
325  BaseParticlePropagator theOutParticle =
326  BaseParticlePropagator(
327  RawParticle(XYZTLorentzVector(p.x(),
328  p.y(),
329  p.z(),
330  pfenergy),
331  XYZTLorentzVector(v.x(),
332  v.y(),
333  v.z(),
334  0.)),
335  0.,0.,B_.z());
336  theOutParticle.setCharge(track.charge());
337  bool isBelowPS=false;
338  theOutParticle.propagateToPreshowerLayer1(false);
339  if(theOutParticle.getSuccess()!=0)
340  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
341  math::XYZPoint(theOutParticle.vertex()),
342  math::XYZTLorentzVector(theOutParticle.momentum())));
343  else {
344  PFTrajectoryPoint dummyPS1;
345  pftrack.addPoint(dummyPS1);
346  }
347 
348 
349  theOutParticle.propagateToPreshowerLayer2(false);
350  if(theOutParticle.getSuccess()!=0){
351  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
352  math::XYZPoint(theOutParticle.vertex()),
353  math::XYZTLorentzVector(theOutParticle.momentum())));
354  isBelowPS=true;
355  } else {
356  PFTrajectoryPoint dummyPS2;
357  pftrack.addPoint(dummyPS2);
358  }
359 
360  theOutParticle.propagateToEcalEntrance(false);
361 
362  if(theOutParticle.getSuccess()!=0){
363  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
364  math::XYZPoint(theOutParticle.vertex()),
365  math::XYZTLorentzVector(theOutParticle.momentum())));
366  double ecalShowerDepth
367  = PFCluster::getDepthCorrection(theOutParticle.momentum().E(),
368  isBelowPS,
369  false);
370 
371  math::XYZPoint meanShower=math::XYZPoint(theOutParticle.vertex())+
372  math::XYZTLorentzVector(theOutParticle.momentum()).Vect().Unit()*ecalShowerDepth;
373 
374  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
375  meanShower,
376  math::XYZTLorentzVector(theOutParticle.momentum())));}
377  else {
378  if (PT>5.)
379  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE ECAL HAS FAILED";
380  PFTrajectoryPoint dummyECAL;
381  pftrack.addPoint(dummyECAL);
382  PFTrajectoryPoint dummyMaxSh;
383  pftrack.addPoint(dummyMaxSh);
384  }
385 
386 
387 
388  //HCAL entrance
389  theOutParticle.propagateToHcalEntrance(false);
390  if(theOutParticle.getSuccess()!=0)
391  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
392  math::XYZPoint(theOutParticle.vertex()),
393  math::XYZTLorentzVector(theOutParticle.momentum())));
394  else{
395  if (PT>5.)
396  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
397  PFTrajectoryPoint dummyHCALentrance;
398  pftrack.addPoint(dummyHCALentrance);
399  }
400  //HCAL exit
401  theOutParticle.propagateToHcalExit(false);
402  if(theOutParticle.getSuccess()!=0)
403  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
404  math::XYZPoint(theOutParticle.vertex()),
405  math::XYZTLorentzVector(theOutParticle.momentum())));
406  else{
407  if (PT>5.)
408  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL EXIT HAS FAILED";
409  PFTrajectoryPoint dummyHCALexit;
410  pftrack.addPoint(dummyHCALexit);
411  }
412 
413  //HO Layer0
414  if ( abs(theOutParticle.vertex().z())<550) {
415  theOutParticle.setMagneticField(0);
416  theOutParticle.propagateToHOLayer(false);
417  if(theOutParticle.getSuccess()!=0)
418  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
419  math::XYZPoint(theOutParticle.vertex()),
420  math::XYZTLorentzVector(theOutParticle.momentum())));
421  else{
422  if (PT>5. && abs(theOutParticle.Z()) < 700.25 )
423  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE HO HAS FAILED";
424  PFTrajectoryPoint dummyHOLayer;
425  pftrack.addPoint(dummyHOLayer);
426  }
427  }
428  }
429 
430  // -------------------------- END GSF Track -------------------------------------
431 
432  // -------------------------- Fill Brem "Track" ---------------------------------
433  // Fill the brem for each traj point
434 
435  //check that the vertex of the brem is in the tracker volume
436  if ((v.perp()>110) ||(fabs(v.z())>280)) continue;
437  unsigned int iTrajPoint = iTrajPos + 2;
438  if(iid%2 == 1) iTrajPoint = 99;
439 
440  PFBrem brem(DP,SigmaDP,iTrajPoint);
441 
442 
443  GlobalVector p_gamma= p*(fabs(DP)/p.mag()); // Direction from the electron (tangent), DP without any sign!;
444  float e_gamma = fabs(DP); // DP = pout-pin so could be negative
445  BaseParticlePropagator theBremParticle =
446  BaseParticlePropagator(
447  RawParticle(XYZTLorentzVector(p_gamma.x(),
448  p_gamma.y(),
449  p_gamma.z(),
450  e_gamma),
451  XYZTLorentzVector(v.x(),
452  v.y(),
453  v.z(),
454  0.)),
455  0.,0.,B_.z());
456  int gamma_charge = 0;
457  theBremParticle.setCharge(gamma_charge);
458 
459 
460  // add TrajectoryPoint for Brem, PS, ECAL, ECALShowMax, HCAL
461  // Brem Entrance PS Layer1
462 
463  PFTrajectoryPoint dummyClosest; // Added just to have the right number order in PFTrack.cc
464  brem.addPoint(dummyClosest);
465 
466 
467  PFTrajectoryPoint dummyBeamPipe; // Added just to have the right number order in PFTrack.cc
468  brem.addPoint(dummyBeamPipe);
469 
470 
471 
472  bool isBelowPS=false;
473  theBremParticle.propagateToPreshowerLayer1(false);
474  if(theBremParticle.getSuccess()!=0)
475  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
476  math::XYZPoint(theBremParticle.vertex()),
477  math::XYZTLorentzVector(theBremParticle.momentum())));
478  else {
479  PFTrajectoryPoint dummyPS1;
480  brem.addPoint(dummyPS1);
481  }
482 
483  // Brem Entrance PS Layer 2
484 
485  theBremParticle.propagateToPreshowerLayer2(false);
486  if(theBremParticle.getSuccess()!=0){
487  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
488  math::XYZPoint(theBremParticle.vertex()),
489  math::XYZTLorentzVector(theBremParticle.momentum())));
490  isBelowPS=true;
491  } else {
492  PFTrajectoryPoint dummyPS2;
493  brem.addPoint(dummyPS2);
494  }
495 
496  theBremParticle.propagateToEcalEntrance(false);
497 
498  if(theBremParticle.getSuccess()!=0){
499  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
500  math::XYZPoint(theBremParticle.vertex()),
501  math::XYZTLorentzVector(theBremParticle.momentum())));
502  double ecalShowerDepth
503  = PFCluster::getDepthCorrection(theBremParticle.momentum().E(),
504  isBelowPS,
505  false);
506 
507  math::XYZPoint meanShower=math::XYZPoint(theBremParticle.vertex())+
508  math::XYZTLorentzVector(theBremParticle.momentum()).Vect().Unit()*ecalShowerDepth;
509 
510  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
511  meanShower,
512  math::XYZTLorentzVector(theBremParticle.momentum())));}
513  else {
514  if ((DP>5.) && ((DP/SigmaDP)>3))
515  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE ECAL HAS FAILED";
516  PFTrajectoryPoint dummyECAL;
517  brem.addPoint(dummyECAL);
518  PFTrajectoryPoint dummyMaxSh;
519  brem.addPoint(dummyMaxSh);
520  }
521 
522 
523 
524  //HCAL entrance
525  theBremParticle.propagateToHcalEntrance(false);
526  if(theBremParticle.getSuccess()!=0)
527  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
528  math::XYZPoint(theBremParticle.vertex()),
529  math::XYZTLorentzVector(theBremParticle.momentum())));
530  else{
531  if ((DP>5.) && ((DP/SigmaDP)>3))
532  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
533  PFTrajectoryPoint dummyHCALentrance;
534  brem.addPoint(dummyHCALentrance);
535  }
536 
537  //HCAL exit
538  theBremParticle.propagateToHcalExit(false);
539  if(theBremParticle.getSuccess()!=0)
540  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
541  math::XYZPoint(theBremParticle.vertex()),
542  math::XYZTLorentzVector(theBremParticle.momentum())));
543  else{
544  if ((DP>5.) && ((DP/SigmaDP)>3))
545  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL EXIT HAS FAILED";
546  PFTrajectoryPoint dummyHCALexit;
547  brem.addPoint(dummyHCALexit);
548  }
549 
550  //HO Layer0
551  if ( abs(theBremParticle.vertex().z())<550.0) {
552  theBremParticle.setMagneticField(0);
553  theBremParticle.propagateToHOLayer(false);
554  if(theBremParticle.getSuccess()!=0)
555  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
556  math::XYZPoint(theBremParticle.vertex()),
557  math::XYZTLorentzVector(theBremParticle.momentum())));
558  else {
559  if ((DP>5.) && ((DP/SigmaDP)>3) && abs(theBremParticle.Z()) < 700.25 )
560  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE H0 HAS FAILED";
561  PFTrajectoryPoint dummyHOLayer;
562  brem.addPoint(dummyHOLayer);
563  }
564  }
565  brem.calculatePositionREP();
566  pftrack.addBrem(brem);
567  iTrajPos++;
568  }
569  return true;
570 }
571 
572 
573 
574 bool
575 PFTrackTransformer::addPointsAndBrems( reco::GsfPFRecTrack& pftrack,
576  const reco::GsfTrack& track,
577  const MultiTrajectoryStateTransform& mtjstate) const {
578 
579  // float PT= track.pt();
580  unsigned int iTrajPos = 0;
581  unsigned int iid = 0; // not anymore saved
582 
583 
584  // ***************************** INNER State *************************************
585  TrajectoryStateOnSurface inTSOS = mtjstate.innerStateOnSurface((track));
586  TrajectoryStateOnSurface outTSOS = mtjstate.outerStateOnSurface((track));
587 
588  if(!inTSOS.isValid() || !outTSOS.isValid()) {
589  if(!inTSOS.isValid())
590  LogWarning("PFTrackTransformer")<<" INNER TSOS NOT VALID ";
591  if(!outTSOS.isValid())
592  LogWarning("PFTrackTransformer")<<" OUTER TSOS NOT VALID ";
593  return false;
594  }
595 
596  GlobalVector InMom;
597  GlobalPoint InPos;
598  if(inTSOS.isValid()) {
599  mtsMode_->momentumFromModeCartesian(inTSOS,InMom);
600  mtsMode_->positionFromModeCartesian(inTSOS,InPos);
601  }
602  else {
603  InMom = GlobalVector(track.pxMode(),track.pyMode(),track.pzMode());
604  InPos = GlobalPoint(0.,0.,0.);
605  }
606 
607  // float pfmass= (pftrack.algoType()==reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
608  float ptot = sqrt((InMom.x()*InMom.x())+(InMom.y()*InMom.y())+(InMom.z()*InMom.z()));
609  float pfenergy= ptot;
610 
611  math::XYZTLorentzVector momClosest
612  = math::XYZTLorentzVector(InMom.x(), InMom.y(),
613  InMom.z(), ptot);
614  math::XYZPoint posClosest = track.vertex();
615  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ClosestApproach,
616  posClosest,momClosest));
617 
618  BaseParticlePropagator theInnerParticle =
619  BaseParticlePropagator( RawParticle(XYZTLorentzVector(InMom.x(),
620  InMom.y(),
621  InMom.z(),
622  pfenergy),
623  XYZTLorentzVector(track.vertex().x(),
624  track.vertex().y(),
625  track.vertex().z(),
626  0.)), //DANIELE Same thing v.x(),v.y(),v.()?
627  0.,0.,B_.z());
628  theInnerParticle.setCharge(track.charge()); // Use the chargeMode ??
629  //BEAMPIPE
630  theInnerParticle.setPropagationConditions(pfGeometry_.outerRadius(PFGeometry::BeamPipe),
631  pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
632  theInnerParticle.propagate();
633  if(theInnerParticle.getSuccess()!=0)
634  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::BeamPipeOrEndVertex,
635  math::XYZPoint(theInnerParticle.vertex()),
636  math::XYZTLorentzVector(theInnerParticle.momentum())));
637  else {
638  PFTrajectoryPoint dummyBeam;
639  pftrack.addPoint(dummyBeam);
640  }
641 
642 
643  // first tjpoint
644  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
645  math::XYZPoint(InPos.x(),InPos.y(), InPos.z()),
646  math::XYZTLorentzVector(InMom.x(),InMom.y(),InMom.z(),InMom.mag())));
647 
648 
649  //######### Photon at INNER State ##########
650 
651 
652  unsigned int iTrajPoint = iTrajPos + 2;
653  double dp_tang = ptot;
654  double sdp_tang = track.ptModeError()*(track.pMode()/track.ptMode());
655  PFBrem brem(dp_tang,sdp_tang,iTrajPoint);
656  BaseParticlePropagator theBremParticle =
657  BaseParticlePropagator(
658  RawParticle(XYZTLorentzVector(InMom.x(),
659  InMom.y(),
660  InMom.z(),
661  dp_tang),
662  XYZTLorentzVector(InPos.x(),
663  InPos.y(),
664  InPos.z(),
665  0.)),
666  0.,0.,B_.z());
667  int gamma_charge = 0;
668  theBremParticle.setCharge(gamma_charge);
669  // add TrajectoryPoint for Brem, PS, ECAL, ECALShowMax, HCAL
670  // Brem Entrance PS Layer1
671  PFTrajectoryPoint dummyClosest; // Added just to have the right number order in PFTrack.cc
672  brem.addPoint(dummyClosest);
673 
674 
675  PFTrajectoryPoint dummyBeamPipe; // Added just to have the right number order in PFTrack.cc
676  brem.addPoint(dummyBeamPipe);
677 
678 
679 
680  bool isBelowPS=false;
681  theBremParticle.propagateToPreshowerLayer1(false);
682  if(theBremParticle.getSuccess()!=0)
683  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
684  math::XYZPoint(theBremParticle.vertex()),
685  math::XYZTLorentzVector(theBremParticle.momentum())));
686  else {
687  PFTrajectoryPoint dummyPS1;
688  brem.addPoint(dummyPS1);
689  }
690 
691  // Brem Entrance PS Layer 2
692 
693  theBremParticle.propagateToPreshowerLayer2(false);
694  if(theBremParticle.getSuccess()!=0){
695  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
696  math::XYZPoint(theBremParticle.vertex()),
697  math::XYZTLorentzVector(theBremParticle.momentum())));
698  isBelowPS=true;
699  } else {
700  PFTrajectoryPoint dummyPS2;
701  brem.addPoint(dummyPS2);
702  }
703 
704  theBremParticle.propagateToEcalEntrance(false);
705 
706  if(theBremParticle.getSuccess()!=0){
707  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
708  math::XYZPoint(theBremParticle.vertex()),
709  math::XYZTLorentzVector(theBremParticle.momentum())));
710 
711  // for the first brem give a low default DP of 100 MeV.
712  double EDepthCorr = 0.01;
713  double ecalShowerDepth
714  = PFCluster::getDepthCorrection(EDepthCorr,
715  isBelowPS,
716  false);
717 
718  math::XYZPoint meanShower=math::XYZPoint(theBremParticle.vertex())+
719  math::XYZTLorentzVector(theBremParticle.momentum()).Vect().Unit()*ecalShowerDepth;
720 
721  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
722  meanShower,
723  math::XYZTLorentzVector(theBremParticle.momentum())));}
724  else {
725  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
726  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE ECAL HAS FAILED";
727  PFTrajectoryPoint dummyECAL;
728  brem.addPoint(dummyECAL);
729  PFTrajectoryPoint dummyMaxSh;
730  brem.addPoint(dummyMaxSh);
731  }
732 
733 
734 
735  //HCAL entrance
736  theBremParticle.propagateToHcalEntrance(false);
737  if(theBremParticle.getSuccess()!=0)
738  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
739  math::XYZPoint(theBremParticle.vertex()),
740  math::XYZTLorentzVector(theBremParticle.momentum())));
741  else{
742  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
743  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
744  PFTrajectoryPoint dummyHCALentrance;
745  brem.addPoint(dummyHCALentrance);
746  }
747 
748  //HCAL exit
749  theBremParticle.propagateToHcalExit(false);
750  if(theBremParticle.getSuccess()!=0)
751  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
752  math::XYZPoint(theBremParticle.vertex()),
753  math::XYZTLorentzVector(theBremParticle.momentum())));
754  else{
755  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
756  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL EXIT HAS FAILED";
757  PFTrajectoryPoint dummyHCALexit;
758  brem.addPoint(dummyHCALexit);
759  }
760 
761  //HO Layer0
762  if ( abs(theBremParticle.vertex().z())<550) {
763  theBremParticle.setMagneticField(0);
764  theBremParticle.propagateToHOLayer(false);
765  if(theBremParticle.getSuccess()!=0)
766  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
767  math::XYZPoint(theBremParticle.vertex()),
768  math::XYZTLorentzVector(theBremParticle.momentum())));
769  else{
770  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3) && abs(theBremParticle.Z()) < 700.25 )
771  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE H0 HAS FAILED";
772  PFTrajectoryPoint dummyHOLayer;
773  brem.addPoint(dummyHOLayer);
774  }
775  }
776 
777  brem.calculatePositionREP();
778  pftrack.addBrem(brem);
779  iTrajPos++;
780 
781 
782 
783 
784  // ***************************** INTERMIDIATE State *************************************
785  //From the new Wolfgang code
786 
787  // To think if the cout should be removed.
788  if(track.gsfExtra()->tangentsSize() == 0)
789  LogError("PFTrackTransformer")
790  <<"BE CAREFUL: Gsf Tangents not stored in the event. You need to re-reco the particle-flow with RecoToDisplay_cfg.py and not RecoToDisplay_NoTracking_cfg.py ";
791 
792 
793  vector<GsfTangent> gsftang = track.gsfExtra()->tangents();
794  for(unsigned int iTang = 0; iTang < track.gsfExtra()->tangentsSize(); iTang++) {
795 
796  dp_tang = gsftang[iTang].deltaP().value();
797  sdp_tang = gsftang[iTang].deltaP().error();
798 
799  //check that the vertex of the brem is in the tracker volume
800  if ((sqrt(gsftang[iTang].position().x()*gsftang[iTang].position().x()
801  + gsftang[iTang].position().y()*gsftang[iTang].position().y())>110)
802  ||(fabs(gsftang[iTang].position().z())>280)) continue;
803 
804  iTrajPoint = iTrajPos + 2;
805  PFBrem brem(dp_tang,sdp_tang,iTrajPoint);
806 
807 
808 
809  GlobalVector p_tang= GlobalVector(gsftang[iTang].momentum().x(),
810  gsftang[iTang].momentum().y(),
811  gsftang[iTang].momentum().z());
812 
813 
814  // ###### track tj points
815  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
816  math::XYZPoint(gsftang[iTang].position().x(),gsftang[iTang].position().y(),gsftang[iTang].position().z()),
817  math::XYZTLorentzVector(p_tang.x(),p_tang.y(),p_tang.z(),p_tang.mag())));
818 
819 
820  //rescale
821  GlobalVector p_gamma = p_tang *(fabs(dp_tang)/p_tang.mag());
822 
823  // GlobalVector
824 
825 
826  double e_gamma = fabs(dp_tang); // DP = pout-pin so could be negative
827  theBremParticle = BaseParticlePropagator(
828  RawParticle(XYZTLorentzVector(p_gamma.x(),
829  p_gamma.y(),
830  p_gamma.z(),
831  e_gamma),
832  XYZTLorentzVector(gsftang[iTang].position().x(),
833  gsftang[iTang].position().y(),
834  gsftang[iTang].position().z(),
835  0.)),
836  0.,0.,B_.z());
837 
838  theBremParticle.setCharge(gamma_charge);
839 
840 
841  PFTrajectoryPoint dummyClosest; // Added just to have the right number order in PFTrack.cc
842  brem.addPoint(dummyClosest);
843 
844 
845  PFTrajectoryPoint dummyBeamPipe; // Added just to have the right number order in PFTrack.cc
846  brem.addPoint(dummyBeamPipe);
847 
848 
849 
850  isBelowPS=false;
851  theBremParticle.propagateToPreshowerLayer1(false);
852  if(theBremParticle.getSuccess()!=0)
853  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
854  math::XYZPoint(theBremParticle.vertex()),
855  math::XYZTLorentzVector(theBremParticle.momentum())));
856  else {
857  PFTrajectoryPoint dummyPS1;
858  brem.addPoint(dummyPS1);
859  }
860 
861  // Brem Entrance PS Layer 2
862 
863  theBremParticle.propagateToPreshowerLayer2(false);
864  if(theBremParticle.getSuccess()!=0){
865  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
866  math::XYZPoint(theBremParticle.vertex()),
867  math::XYZTLorentzVector(theBremParticle.momentum())));
868  isBelowPS=true;
869  } else {
870  PFTrajectoryPoint dummyPS2;
871  brem.addPoint(dummyPS2);
872  }
873 
874  theBremParticle.propagateToEcalEntrance(false);
875 
876  if(theBremParticle.getSuccess()!=0){
877  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
878  math::XYZPoint(theBremParticle.vertex()),
879  math::XYZTLorentzVector(theBremParticle.momentum())));
880 
881  double ecalShowerDepth
882  = PFCluster::getDepthCorrection(theBremParticle.momentum().E(),
883  isBelowPS,
884  false);
885 
886  math::XYZPoint meanShower=math::XYZPoint(theBremParticle.vertex())+
887  math::XYZTLorentzVector(theBremParticle.momentum()).Vect().Unit()*ecalShowerDepth;
888 
889  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
890  meanShower,
891  math::XYZTLorentzVector(theBremParticle.momentum())));}
892  else {
893  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
894  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE ECAL HAS FAILED";
895  PFTrajectoryPoint dummyECAL;
896  brem.addPoint(dummyECAL);
897  PFTrajectoryPoint dummyMaxSh;
898  brem.addPoint(dummyMaxSh);
899  }
900 
901 
902 
903  //HCAL entrance
904  theBremParticle.propagateToHcalEntrance(false);
905  if(theBremParticle.getSuccess()!=0)
906  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
907  math::XYZPoint(theBremParticle.vertex()),
908  math::XYZTLorentzVector(theBremParticle.momentum())));
909  else{
910  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
911  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
912  PFTrajectoryPoint dummyHCALentrance;
913  brem.addPoint(dummyHCALentrance);
914  }
915 
916  //HCAL exit
917  theBremParticle.propagateToHcalExit(false);
918  if(theBremParticle.getSuccess()!=0)
919  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
920  math::XYZPoint(theBremParticle.vertex()),
921  math::XYZTLorentzVector(theBremParticle.momentum())));
922  else{
923  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
924  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL EXIT HAS FAILED";
925  PFTrajectoryPoint dummyHCALexit;
926  brem.addPoint(dummyHCALexit);
927  }
928 
929  //HO Layer0
930  if ( abs(theBremParticle.vertex().z())<550) {
931  theBremParticle.setMagneticField(0);
932  theBremParticle.propagateToHOLayer(false);
933  if(theBremParticle.getSuccess()!=0)
934  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
935  math::XYZPoint(theBremParticle.vertex()),
936  math::XYZTLorentzVector(theBremParticle.momentum())));
937  else{
938  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3) && abs(theBremParticle.Z()) < 700.25 )
939  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE H0 HAS FAILED";
940  PFTrajectoryPoint dummyHOLayer;
941  brem.addPoint(dummyHOLayer);
942  }
943  }
944 
945  brem.calculatePositionREP();
946  pftrack.addBrem(brem);
947  iTrajPos++;
948  }
949 
950 
951 
952 
953  // ***************************** OUTER State *************************************
954 
955  if(outTSOS.isValid()) {
956  GlobalVector OutMom;
957  GlobalPoint OutPos;
958 
959  // DANIELE ????? if the out is not valid maybe take the last tangent?
960  // From Wolfgang. It should be always valid
961 
962  mtsMode_->momentumFromModeCartesian(outTSOS,OutMom);
963  mtsMode_->positionFromModeCartesian(outTSOS,OutPos);
964 
965 
966 
967  // last tjpoint
968  pftrack.addPoint(PFTrajectoryPoint(iid,-1,
969  math::XYZPoint(OutPos.x(),OutPos.y(), OutPos.z()),
970  math::XYZTLorentzVector(OutMom.x(),OutMom.y(),OutMom.z(),OutMom.mag())));
971 
972 
973  float ptot_out = sqrt((OutMom.x()*OutMom.x())+(OutMom.y()*OutMom.y())+(OutMom.z()*OutMom.z()));
974  float pTtot_out = sqrt((OutMom.x()*OutMom.x())+(OutMom.y()*OutMom.y()));
975  float pfenergy_out = ptot_out;
976  BaseParticlePropagator theOutParticle =
977  BaseParticlePropagator( RawParticle(XYZTLorentzVector(OutMom.x(),
978  OutMom.y(),
979  OutMom.z(),
980  pfenergy_out),
981  XYZTLorentzVector(OutPos.x(),
982  OutPos.y(),
983  OutPos.z(),
984  0.)),
985  0.,0.,B_.z());
986  theOutParticle.setCharge(track.charge());
987  isBelowPS=false;
988  theOutParticle.propagateToPreshowerLayer1(false);
989  if(theOutParticle.getSuccess()!=0)
990  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
991  math::XYZPoint(theOutParticle.vertex()),
992  math::XYZTLorentzVector(theOutParticle.momentum())));
993  else {
994  PFTrajectoryPoint dummyPS1;
995  pftrack.addPoint(dummyPS1);
996  }
997 
998 
999  theOutParticle.propagateToPreshowerLayer2(false);
1000  if(theOutParticle.getSuccess()!=0){
1001  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
1002  math::XYZPoint(theOutParticle.vertex()),
1003  math::XYZTLorentzVector(theOutParticle.momentum())));
1004  isBelowPS=true;
1005  } else {
1006  PFTrajectoryPoint dummyPS2;
1007  pftrack.addPoint(dummyPS2);
1008  }
1009 
1010  theOutParticle.propagateToEcalEntrance(false);
1011 
1012  if(theOutParticle.getSuccess()!=0){
1013  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
1014  math::XYZPoint(theOutParticle.vertex()),
1015  math::XYZTLorentzVector(theOutParticle.momentum())));
1016  double EDepthCorr = 0.01;
1017  double ecalShowerDepth
1018  = PFCluster::getDepthCorrection(EDepthCorr,
1019  isBelowPS,
1020  false);
1021 
1022  math::XYZPoint meanShower=math::XYZPoint(theOutParticle.vertex())+
1023  math::XYZTLorentzVector(theOutParticle.momentum()).Vect().Unit()*ecalShowerDepth;
1024 
1025  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
1026  meanShower,
1027  math::XYZTLorentzVector(theOutParticle.momentum())));}
1028  else {
1029  if (pTtot_out>5.)
1030  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE ECAL HAS FAILED";
1031  PFTrajectoryPoint dummyECAL;
1032  pftrack.addPoint(dummyECAL);
1033  PFTrajectoryPoint dummyMaxSh;
1034  pftrack.addPoint(dummyMaxSh);
1035  }
1036 
1037 
1038 
1039  //HCAL entrance
1040  theOutParticle.propagateToHcalEntrance(false);
1041  if(theOutParticle.getSuccess()!=0)
1042  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
1043  math::XYZPoint(theOutParticle.vertex()),
1044  math::XYZTLorentzVector(theOutParticle.momentum())));
1045  else{
1046  if (pTtot_out>5.)
1047  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
1048  PFTrajectoryPoint dummyHCALentrance;
1049  pftrack.addPoint(dummyHCALentrance);
1050  }
1051  //HCAL exit
1052  theOutParticle.propagateToHcalExit(false);
1053  if(theOutParticle.getSuccess()!=0)
1054  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
1055  math::XYZPoint(theOutParticle.vertex()),
1056  math::XYZTLorentzVector(theOutParticle.momentum())));
1057  else{
1058  if (pTtot_out>5.)
1059  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<< " PROPAGATION TO THE HCAL EXIT HAS FAILED";
1060  PFTrajectoryPoint dummyHCALexit;
1061  pftrack.addPoint(dummyHCALexit);
1062  }
1063 
1064  //HO Layer0
1065  if ( abs(theOutParticle.vertex().z())<550) {
1066  theOutParticle.setMagneticField(0);
1067  theOutParticle.propagateToHOLayer(false);
1068  if(theOutParticle.getSuccess()!=0)
1069  pftrack.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
1070  math::XYZPoint(theOutParticle.vertex()),
1071  math::XYZTLorentzVector(theOutParticle.momentum())));
1072  else{
1073  if ( pTtot_out > 5. && abs(theOutParticle.Z()) < 700.25 )
1074  LogWarning("PFTrackTransformer")<<"GSF TRACK "<<pftrack<<" PROPAGATION TO THE HO HAS FAILED";
1075  PFTrajectoryPoint dummyHOLayer;
1076  pftrack.addPoint(dummyHOLayer);
1077  }
1078  }
1079  //######## Photon at the OUTER State ##########
1080 
1081  dp_tang = OutMom.mag();
1082  // for the moment same inner error just for semplicity
1083  sdp_tang = track.ptModeError()*(track.pMode()/track.ptMode());
1084  iTrajPoint = iTrajPos + 2;
1085  PFBrem brem(dp_tang,sdp_tang,iTrajPoint);
1086 
1087  theBremParticle =
1088  BaseParticlePropagator( RawParticle(XYZTLorentzVector(OutMom.x(),
1089  OutMom.y(),
1090  OutMom.z(),
1091  dp_tang),
1092  XYZTLorentzVector(OutPos.x(),
1093  OutPos.y(),
1094  OutPos.z(),
1095  0.)),
1096  0.,0.,B_.z());
1097  theBremParticle.setCharge(gamma_charge);
1098 
1099 
1100  PFTrajectoryPoint dummyClosest; // Added just to have the right number order in PFTrack.cc
1101  brem.addPoint(dummyClosest);
1102 
1103 
1104  PFTrajectoryPoint dummyBeamPipe; // Added just to have the right number order in PFTrack.cc
1105  brem.addPoint(dummyBeamPipe);
1106 
1107 
1108 
1109  isBelowPS=false;
1110  theBremParticle.propagateToPreshowerLayer1(false);
1111  if(theBremParticle.getSuccess()!=0)
1112  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS1,
1113  math::XYZPoint(theBremParticle.vertex()),
1114  math::XYZTLorentzVector(theBremParticle.momentum())));
1115  else {
1116  PFTrajectoryPoint dummyPS1;
1117  brem.addPoint(dummyPS1);
1118  }
1119 
1120  // Brem Entrance PS Layer 2
1121 
1122  theBremParticle.propagateToPreshowerLayer2(false);
1123  if(theBremParticle.getSuccess()!=0){
1124  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::PS2,
1125  math::XYZPoint(theBremParticle.vertex()),
1126  math::XYZTLorentzVector(theBremParticle.momentum())));
1127  isBelowPS=true;
1128  } else {
1129  PFTrajectoryPoint dummyPS2;
1130  brem.addPoint(dummyPS2);
1131  }
1132 
1133  theBremParticle.propagateToEcalEntrance(false);
1134 
1135  if(theBremParticle.getSuccess()!=0){
1136  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALEntrance,
1137  math::XYZPoint(theBremParticle.vertex()),
1138  math::XYZTLorentzVector(theBremParticle.momentum())));
1139  double ecalShowerDepth
1140  = PFCluster::getDepthCorrection(theBremParticle.momentum().E(),
1141  isBelowPS,
1142  false);
1143 
1144  math::XYZPoint meanShower=math::XYZPoint(theBremParticle.vertex())+
1145  math::XYZTLorentzVector(theBremParticle.momentum()).Vect().Unit()*ecalShowerDepth;
1146 
1147  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::ECALShowerMax,
1148  meanShower,
1149  math::XYZTLorentzVector(theBremParticle.momentum())));}
1150  else {
1151  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
1152  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE ECAL HAS FAILED";
1153  PFTrajectoryPoint dummyECAL;
1154  brem.addPoint(dummyECAL);
1155  PFTrajectoryPoint dummyMaxSh;
1156  brem.addPoint(dummyMaxSh);
1157  }
1158 
1159 
1160 
1161  //HCAL entrance
1162  theBremParticle.propagateToHcalEntrance(false);
1163  if(theBremParticle.getSuccess()!=0)
1164  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALEntrance,
1165  math::XYZPoint(theBremParticle.vertex()),
1166  math::XYZTLorentzVector(theBremParticle.momentum())));
1167  else{
1168  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
1169  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
1170  PFTrajectoryPoint dummyHCALentrance;
1171  brem.addPoint(dummyHCALentrance);
1172  }
1173  //HCAL exit
1174  theBremParticle.propagateToHcalExit(false);
1175  if(theBremParticle.getSuccess()!=0)
1176  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HCALExit,
1177  math::XYZPoint(theBremParticle.vertex()),
1178  math::XYZTLorentzVector(theBremParticle.momentum())));
1179  else{
1180  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3))
1181  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE HCAL EXIT HAS FAILED";
1182  PFTrajectoryPoint dummyHCALexit;
1183  brem.addPoint(dummyHCALexit);
1184  }
1185 
1186  //HO Layer0
1187  if ( abs(theBremParticle.vertex().z())<550) {
1188  theBremParticle.setMagneticField(0);
1189  theBremParticle.propagateToHOLayer(false);
1190  if(theBremParticle.getSuccess()!=0)
1191  brem.addPoint(PFTrajectoryPoint(-1,PFTrajectoryPoint::HOLayer,
1192  math::XYZPoint(theBremParticle.vertex()),
1193  math::XYZTLorentzVector(theBremParticle.momentum())));
1194  else{
1195  if ((dp_tang>5.) && ((dp_tang/sdp_tang)>3) && abs(theBremParticle.Z()) < 700.25 )
1196  LogWarning("PFTrackTransformer")<<"BREM "<<brem<<" PROPAGATION TO THE H0 HAS FAILED";
1197  PFTrajectoryPoint dummyHOLayer;
1198  brem.addPoint(dummyHOLayer);
1199  }
1200  }
1201  brem.calculatePositionREP();
1202  pftrack.addBrem(brem);
1203  iTrajPos++;
1204  }
1205 
1206  return true;
1207 }
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:501
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:578
RawParticle::setCharge
void setCharge(float q)
set the MEASURED charge
Definition: RawParticle.cc:139
reco::PFRecTrack
reconstructed track used as an input to particle flow
Definition: PFRecTrack.h:22
reco::GsfTrack::pyMode
double pyMode() const
y coordinate of momentum vector from mode
Definition: GsfTrack.h:53
MultiTrajectoryStateMode::positionFromModeCartesian
bool positionFromModeCartesian(const TrajectoryStateOnSurface tsos, GlobalPoint &position) const
Definition: MultiTrajectoryStateMode.cc:60
BaseParticlePropagator::propagateToPreshowerLayer1
bool propagateToPreshowerLayer1(bool first=true)
Definition: BaseParticlePropagator.cc:446
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:72
reco::GsfTrack::pMode
double pMode() const
momentum vector magnitude from mode
Definition: GsfTrack.h:47
PT
#define PT
Definition: PFJetBenchmark.cc:34
PFTrackTransformer::B_
math::XYZVector B_
B field.
Definition: PFTrackTransformer.h:68
PFTrackTransformer::addPointsAndBrems
bool addPointsAndBrems(reco::GsfPFRecTrack &pftrack, const reco::Track &track, const Trajectory &traj, const bool &GetMode) const
Definition: PFTrackTransformer.cc:231
PFTrackTransformer::PFTrackTransformer
PFTrackTransformer(const math::XYZVector &)
Definition: PFTrackTransformer.cc:34
PFTrackTransformer::addPoints
bool addPoints(reco::PFRecTrack &pftrack, const reco::Track &track, const Trajectory &traj, bool msgwarning=true) const
Add points to a PFTrack. return false if a TSOS is invalid.
Definition: PFTrackTransformer.cc:46
GlobalPoint
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
BaseParticlePropagator::getSuccess
int getSuccess() const
Has propagation been performed and was barrel or endcap reached ?
Definition: BaseParticlePropagator.h:297
PFGeometry::outerZ
float outerZ(PFGeometry::Layers_t layer) const
return outer position along z axis of a given layer
Definition: PFGeometry.h:68
BaseParticlePropagator::setPropagationConditions
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
Definition: BaseParticlePropagator.cc:350
BaseParticlePropagator::setMagneticField
void setMagneticField(double b)
Set the magnetic field.
Definition: BaseParticlePropagator.h:300
reco::TrackBase::px
double px() const
x coordinate of momentum vector
Definition: TrackBase.h:590
MultiTrajectoryStateMode::momentumFromModeCartesian
bool momentumFromModeCartesian(const TrajectoryStateOnSurface tsos, GlobalVector &momentum) const
Definition: MultiTrajectoryStateMode.cc:13
dt_dqm_sourceclient_common_cff.reco
tuple reco
Definition: dt_dqm_sourceclient_common_cff.py:106
reco::Track::outerPosition
const math::XYZPoint & outerPosition() const
position of the outermost hit
Definition: Track.h:65
PFGsfHelper::sigmafittedDP
double sigmafittedDP() const
Definition: PFGsfHelper.cc:140
MultiTrajectoryStateTransform::outerStateOnSurface
TrajectoryStateOnSurface outerStateOnSurface(const reco::GsfTrack &tk) const
Definition: MultiTrajectoryStateTransform.cc:20
PFTrackTransformer::mtsMode_
const MultiTrajectoryStateMode * mtsMode_
Definition: PFTrackTransformer.h:69
reco::GsfTrack::ptModeError
double ptModeError() const
error on Pt (set to 1000 TeV if charge==0 for safety) from mode
Definition: GsfTrack.h:81
double_binary
Definition: DDStreamer.cc:227
Trajectory::direction
PropagationDirection const & direction() const
Definition: Trajectory.cc:118
reco::GsfPFRecTrack::addBrem
void addBrem(const reco::PFBrem &brem)
add a Bremsstrahlung photon
Definition: GsfPFRecTrack.cc:22
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
Trajectory::measurements
DataContainer const & measurements() const
Definition: Trajectory.h:203
reco::TrackBase::vertex
const Point & vertex() const
reference point on the track. This method is DEPRECATED, please use referencePoint() instead ...
Definition: TrackBase.h:650
RawParticle::momentum
const XYZTLorentzVector & momentum() const
the momentum fourvector
Definition: RawParticle.h:286
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:67
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
reco::TrackBase::pt
double pt() const
track transverse momentum
Definition: TrackBase.h:584
BaseParticlePropagator::propagateToHcalExit
bool propagateToHcalExit(bool first=true)
Definition: BaseParticlePropagator.cc:566
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
RawParticle::Z
double Z() const
z of vertex
Definition: RawParticle.h:276
reco::GsfTrack::pzMode
double pzMode() const
z coordinate of momentum vector from mode
Definition: GsfTrack.h:55
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
PFGeometry::outerRadius
float outerRadius(PFGeometry::Layers_t layer) const
return outer radius of a given layer
Definition: PFGeometry.h:60
reco::PFTrack::addPoint
void addPoint(const reco::PFTrajectoryPoint &trajPt)
Definition: PFTrack.cc:42
BaseParticlePropagator::propagateToEcalEntrance
bool propagateToEcalEntrance(bool first=true)
Definition: BaseParticlePropagator.cc:487
reco::TrackBase::pz
double pz() const
z coordinate of momentum vector
Definition: TrackBase.h:602
RawParticle::vertex
const XYZTLorentzVector & vertex() const
the vertex fourvector
Definition: RawParticle.h:285
reco::GsfTrack::gsfExtra
const GsfTrackExtraRef & gsfExtra() const
reference to &quot;extra&quot; object
Definition: GsfTrack.h:32
reco::Track::outerMomentum
const math::XYZVector & outerMomentum() const
momentum vector at the outermost hit position
Definition: Track.h:70
math::XYZVector
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
reco::PFTrack::calculatePositionREP
void calculatePositionREP()
Definition: PFTrack.cc:68
reco::GsfTrack::pxMode
double pxMode() const
x coordinate of momentum vector from mode
Definition: GsfTrack.h:51
PFGsfHelper::fittedDP
double fittedDP() const
Definition: PFGsfHelper.cc:136
AlCaHLTBitMon_ParallelJobs.p
tuple p
Definition: AlCaHLTBitMon_ParallelJobs.py:152
MultiTrajectoryStateTransform::innerStateOnSurface
TrajectoryStateOnSurface innerStateOnSurface(const reco::GsfTrack &tk) const
Definition: MultiTrajectoryStateTransform.cc:26
BaseParticlePropagator::propagateToHcalEntrance
bool propagateToHcalEntrance(bool first=true)
Definition: BaseParticlePropagator.cc:513
reco::PFRecTrack::algoType
unsigned int algoType() const
Definition: PFRecTrack.h:47
position
static int position[264][3]
Definition: ReadPGInfo.cc:509
reco::TrackBase::charge
int charge() const
track electric charge
Definition: TrackBase.h:530
reco::PFTrajectoryPoint
A PFTrack holds several trajectory points, which basically contain the position and momentum of a tra...
Definition: PFTrajectoryPoint.h:26
PFGsfHelper::computeP
GlobalVector computeP(bool ComputeMode) const
Definition: PFGsfHelper.cc:130
BaseParticlePropagator::propagateToHOLayer
bool propagateToHOLayer(bool first=true)
Definition: BaseParticlePropagator.cc:584
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
reco::GsfTrack::ptMode
double ptMode() const
track transverse momentum from mode
Definition: GsfTrack.h:49
BaseParticlePropagator::propagateToPreshowerLayer2
bool propagateToPreshowerLayer2(bool first=true)
Definition: BaseParticlePropagator.cc:466
reco::TrackBase::py
double py() const
y coordinate of momentum vector
Definition: TrackBase.h:596
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10
XYZTLorentzVector
math::XYZTLorentzVector XYZTLorentzVector
Definition: RawParticle.h:15
Generated for CMSSW Reference Manual by doxygen 1.8.5