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