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