CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
VZeroFinder.cc
Go to the documentation of this file.
2 
4 
6 
8 
9 #include <utility>
10 using namespace std;
11 
12 /*****************************************************************************/
14  (const edm::EventSetup& es,
15  const edm::ParameterSet& pset)
16 {
17  // Get track-pair level cuts
18  maxDca = pset.getParameter<double>("maxDca");
19 
20  // Get mother cuts
21  minCrossingRadius = pset.getParameter<double>("minCrossingRadius");
22  maxCrossingRadius = pset.getParameter<double>("maxCrossingRadius");
23  maxImpactMother = pset.getParameter<double>("maxImpactMother");
24 
25  // Get magnetic field
27  es.get<IdealMagneticFieldRecord>().get(magField);
28  theMagField = magField.product();
29 }
30 
31 /*****************************************************************************/
33 {
34 }
35 
36 /*****************************************************************************/
38  (const reco::Track& track)
39 {
40  GlobalPoint position(track.vertex().x(),
41  track.vertex().y(),
42  track.vertex().z());
43 
44  GlobalVector momentum(track.momentum().x(),
45  track.momentum().y(),
46  track.momentum().z());
47 
49  track.charge(),theMagField);
50 
51  return gtp;
52 }
53 
54 /*****************************************************************************/
56 {
57  double pt = p.perp();
58  return GlobalVector( -pt*cos(a), -pt*sin(a), p.z());
59 }
60 
61 /*****************************************************************************/
63  const reco::Track& negTrack,
66 {
67 /*
68  LogTrace("MinBiasTracking") << " [VZeroFinder] tracks"
69  << " +" << posTrack.algoName() << " " << posTrack.d0()
70  << " -" << negTrack.algoName() << " " << negTrack.d0();
71 */
72 
73  // Get trajectories
74  GlobalTrajectoryParameters posGtp = getGlobalTrajectoryParameters(posTrack);
75  GlobalTrajectoryParameters negGtp = getGlobalTrajectoryParameters(negTrack);
76 
77  // Two track minimum distance
79 
80  // Closest approach
81  theMinimum.calculate(posGtp,negGtp);
82 
83  // Closest points
84  pair<GlobalPoint, GlobalPoint> points = theMinimum.points();
85 
86  // Momenta at closest point
87  pair<GlobalVector,GlobalVector> momenta;
88  momenta.first = rotate(posGtp.momentum(), theMinimum.firstAngle() );
89  momenta.second = rotate(negGtp.momentum(), theMinimum.secondAngle());
90 
91  // dcaR
92  float dca = theMinimum.distance();
93  GlobalPoint crossing = theMinimum.crossingPoint();
94 
95 /*
96  LogTrace("MinBiasTracking") << fixed << setprecision(2)
97  << " [VZeroFinder] dca = "<<dca<<" (<"<<maxDca<<")";
98  LogTrace("MinBiasTracking") << fixed << setprecision(2)
99  << " [VZeroFinder] crossR = "<<crossing.perp()<<" ("
100  <<minCrossingRadius<<"< <"<<maxCrossingRadius<<")";
101 */
102 
103  if(dca < maxDca &&
104  crossing.perp() > minCrossingRadius &&
105  crossing.perp() < maxCrossingRadius)
106  {
107  // Momentum of the mother
108  GlobalVector momentum = momenta.first + momenta.second;
109  float impact = -1.;
110 
111  if(vertices->size() > 0)
112  {
113  // Impact parameter of the mother wrt vertices, choose smallest
114  for(reco::VertexCollection::const_iterator
115  vertex = vertices->begin(); vertex!= vertices->end(); vertex++)
116  {
117  GlobalVector r(crossing.x(),
118  crossing.y(),
119  crossing.z() - vertex->position().z());
120  GlobalVector p(momentum.x(),momentum.y(),momentum.z());
121 
122  GlobalVector b = r - (r*p)*p / p.mag2();
123 
124  float im = b.mag();
125  if(im < impact || vertex == vertices->begin())
126  impact = im;
127  }
128  }
129  else
130  {
131  // Impact parameter of the mother in the plane
132  GlobalVector r_(crossing.x(),crossing.y(),0);
133  GlobalVector p_(momentum.x(),momentum.y(),0);
134 
135  GlobalVector b_ = r_ - (r_*p_)*p_ / p_.mag2();
136  impact = b_.mag();
137  }
138 
139 /*
140  LogTrace("MinBiasTracking") << fixed << setprecision(2)
141  << " [VZeroFinder] impact = "<<impact<<" (<"<<maxImpactMother<<")";
142 */
143 
144  if(impact < maxImpactMother)
145  {
146  // Armenteros variables
147  float pt =
148  (momenta.first.cross(momenta.second)).mag()/momentum.mag();
149  float alpha =
150  (momenta.first.mag2() - momenta.second.mag2())/momentum.mag2();
151 
152  // Fill data
153  data.dca = dca;
154  data.crossingPoint = crossing;
155  data.impactMother = impact;
156  data.armenterosPt = pt;
157  data.armenterosAlpha = alpha;
158  data.momenta = momenta;
159 
160 /*
161  LogTrace("MinBiasTracking") << fixed << setprecision(2)
162  << " [VZeroFinder] success {alpha = "<<alpha<<", pt = "<<pt<<"}";
163 */
164 
165  return true;
166  }
167  }
168 
169  return false;
170 }
171 
virtual float distance() const
T getParameter(std::string const &) const
math::GlobalPoint crossingPoint
Definition: VZero.h:18
virtual bool calculate(const TrajectoryStateOnSurface &sta, const TrajectoryStateOnSurface &stb)
float alpha
Definition: AMPTWrapper.h:95
T mag2() const
Definition: PV3DBase.h:60
T perp() const
Definition: PV3DBase.h:66
const Vector & momentum() const
track momentum vector
Definition: TrackBase.h:149
virtual GlobalPoint crossingPoint() const
GlobalTrajectoryParameters getGlobalTrajectoryParameters(const reco::Track &track)
Definition: VZeroFinder.cc:38
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
T y() const
Definition: PV3DBase.h:57
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
float armenterosPt
Definition: VZero.h:17
std::pair< GlobalVector, GlobalVector > momenta
Definition: VZero.h:19
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
T mag() const
Definition: PV3DBase.h:61
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
float dca
Definition: VZero.h:17
T z() const
Definition: PV3DBase.h:58
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
float armenterosAlpha
Definition: VZero.h:17
tuple pset
Definition: CrabTask.py:85
const Point & vertex() const
reference point on the track. This method is DEPRECATED, please use referencePoint() instead ...
Definition: TrackBase.h:155
GlobalVector rotate(const GlobalVector &p, double a)
Definition: VZeroFinder.cc:55
float impactMother
Definition: VZero.h:17
const T & get() const
Definition: EventSetup.h:55
double b
Definition: hdecay.h:120
def rotate
Definition: svgfig.py:704
#define begin
Definition: vmac.h:31
double a
Definition: hdecay.h:121
VZeroFinder(const edm::EventSetup &es, const edm::ParameterSet &pset)
Definition: VZeroFinder.cc:14
int charge() const
track electric charge
Definition: TrackBase.h:112
bool checkTrackPair(const reco::Track &posTrack, const reco::Track &negTrack, const reco::VertexCollection *vertices, reco::VZeroData &data)
Definition: VZeroFinder.cc:62
T x() const
Definition: PV3DBase.h:56
virtual std::pair< GlobalPoint, GlobalPoint > points() const
Global3DVector GlobalVector
Definition: GlobalVector.h:10