ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
FullBilloirVertexFitterTests.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file FullBilloirVertexFitterTests.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2019 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #include <boost/test/data/test_case.hpp>
10 #include <boost/test/tools/output_test_stream.hpp>
11 #include <boost/test/unit_test.hpp>
12 
20 #include "Acts/Utilities/Units.hpp"
24 
25 namespace bdata = boost::unit_test::data;
26 using namespace Acts::UnitLiterals;
27 
28 namespace Acts {
29 namespace Test {
30 
32 using Linearizer =
33  HelicalTrackLinearizer<Propagator<EigenStepper<ConstantBField>>>;
34 
35 // Create a test context
38 
41 BOOST_AUTO_TEST_CASE(billoir_vertex_fitter_empty_input_test) {
42  // Set up constant B-Field
43  ConstantBField bField(0.0, 0.0, 1_T);
44 
45  // Set up Eigenstepper
47 
48  // Set up propagator with void navigator
49  auto propagator =
50  std::make_shared<Propagator<EigenStepper<ConstantBField>>>(stepper);
51 
52  Linearizer::Config ltConfig(bField, propagator);
53  Linearizer linearizer(ltConfig);
54 
55  // Set up Billoir Vertex Fitter
58  vertexFitterCfg);
59 
60  // Constraint for vertex fit
61  Vertex<BoundParameters> myConstraint;
62  // Some abitrary values
63  SpacePointSymMatrix myCovMat = SpacePointSymMatrix::Zero();
64  myCovMat(0, 0) = 30.;
65  myCovMat(1, 1) = 30.;
66  myCovMat(2, 2) = 30.;
67  myCovMat(3, 3) = 30.;
68  myConstraint.setFullCovariance(std::move(myCovMat));
69  myConstraint.setFullPosition(SpacePointVector(0, 0, 0, 0));
70 
71  const std::vector<const BoundParameters*> emptyVector;
72 
74  myConstraint);
75 
76  Vertex<BoundParameters> fittedVertex =
77  billoirFitter.fit(emptyVector, linearizer, vfOptions).value();
78 
79  Vector3D origin(0., 0., 0.);
80  BOOST_CHECK_EQUAL(fittedVertex.position(), origin);
81 
82  SpacePointSymMatrix zeroMat = SpacePointSymMatrix::Zero();
83  BOOST_CHECK_EQUAL(fittedVertex.fullCovariance(), zeroMat);
84 
85  fittedVertex = billoirFitter.fit(emptyVector, linearizer, vfOptions).value();
86 
87  BOOST_CHECK_EQUAL(fittedVertex.position(), origin);
88  BOOST_CHECK_EQUAL(fittedVertex.fullCovariance(), zeroMat);
89 }
90 
91 // Vertex x/y position distribution
92 std::uniform_real_distribution<> vXYDist(-0.1_mm, 0.1_mm);
93 // Vertex z position distribution
94 std::uniform_real_distribution<> vZDist(-20_mm, 20_mm);
95 // Track d0 distribution
96 std::uniform_real_distribution<> d0Dist(-0.01_mm, 0.01_mm);
97 // Track z0 distribution
98 std::uniform_real_distribution<> z0Dist(-0.2_mm, 0.2_mm);
99 // Track pT distribution
100 std::uniform_real_distribution<> pTDist(0.4_GeV, 10_GeV);
101 // Track phi distribution
102 std::uniform_real_distribution<> phiDist(-M_PI, M_PI);
103 // Track theta distribution
104 std::uniform_real_distribution<> thetaDist(1.0, M_PI - 1.0);
105 // Track charge helper distribution
106 std::uniform_real_distribution<> qDist(-1, 1);
107 // Track IP resolution distribution
108 std::uniform_real_distribution<> resIPDist(0., 100_um);
109 // Track angular distribution
110 std::uniform_real_distribution<> resAngDist(0., 0.1);
111 // Track q/p resolution distribution
112 std::uniform_real_distribution<> resQoPDist(-0.1, 0.1);
113 // Number of tracks distritbution
114 std::uniform_int_distribution<> nTracksDist(3, 10);
115 
120 BOOST_AUTO_TEST_CASE(billoir_vertex_fitter_defaulttrack_test) {
121  bool debugMode = false;
122  // Set up RNG
123  int mySeed = 31415;
124  std::mt19937 gen(mySeed);
125  // Set up constant B-Field
126  ConstantBField bField(0.0, 0.0, 1_T);
127 
128  // Set up Eigenstepper
130  // Set up propagator with void navigator
131  auto propagator =
132  std::make_shared<Propagator<EigenStepper<ConstantBField>>>(stepper);
133 
134  Linearizer::Config ltConfig(bField, propagator);
135  Linearizer linearizer(ltConfig);
136 
137  // Number of events
138  const int nEvents = 10;
139  for (int eventIdx = 0; eventIdx < nEvents; ++eventIdx) {
140  // Number of tracks
141  unsigned int nTracks = nTracksDist(gen);
142 
143  // Set up Billoir Vertex Fitter
145  vertexFitterCfg;
147  vertexFitterCfg);
148  // Constraint for vertex fit
149  Vertex<BoundParameters> myConstraint;
150  // Some abitrary values
151  SpacePointSymMatrix myCovMat = SpacePointSymMatrix::Zero();
152  myCovMat(0, 0) = 30.;
153  myCovMat(1, 1) = 30.;
154  myCovMat(2, 2) = 30.;
155  myCovMat(3, 3) = 30.;
156  myConstraint.setFullCovariance(std::move(myCovMat));
157  myConstraint.setFullPosition(SpacePointVector(0, 0, 0, 0));
159 
160  VertexingOptions<BoundParameters> vfOptionsConstr(
161  geoContext, magFieldContext, myConstraint);
162  // Create position of vertex and perigee surface
163  double x = vXYDist(gen);
164  double y = vXYDist(gen);
165  double z = vZDist(gen);
166 
167  Vector3D vertexPosition(x, y, z);
168  std::shared_ptr<PerigeeSurface> perigeeSurface =
169  Surface::makeShared<PerigeeSurface>(Vector3D(0., 0., 0.));
170  // Calculate d0 and z0 corresponding to vertex position
171  double d0V = sqrt(x * x + y * y);
172  double z0V = z;
173 
174  // Start constructing nTracks tracks in the following
175  std::vector<BoundParameters> tracks;
176 
177  // Construct random track emerging from vicinity of vertex position
178  // Vector to store track objects used for vertex fit
179  for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
180  // Construct positive or negative charge randomly
181  double q = qDist(gen) < 0 ? -1. : 1.;
182 
183  // Construct random track parameters
184  BoundVector paramVec;
185  paramVec << d0V + d0Dist(gen), z0V + z0Dist(gen), phiDist(gen),
186  thetaDist(gen), q / pTDist(gen), 0.;
187 
188  // Resolutions
189  double resD0 = resIPDist(gen);
190  double resZ0 = resIPDist(gen);
191  double resPh = resAngDist(gen);
192  double resTh = resAngDist(gen);
193  double resQp = resQoPDist(gen);
194 
195  // Fill vector of track objects with simple covariance matrix
196  Covariance covMat;
197 
198  covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
199  0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
200  0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
201  tracks.push_back(BoundParameters(geoContext, std::move(covMat), paramVec,
202  perigeeSurface));
203  }
204 
205  std::vector<const BoundParameters*> tracksPtr;
206  for (const auto& trk : tracks) {
207  tracksPtr.push_back(&trk);
208  }
209 
210  // Do the actual fit with 4 tracks without constraint
211  Vertex<BoundParameters> fittedVertex =
212  billoirFitter.fit(tracksPtr, linearizer, vfOptions).value();
213  if (fittedVertex.tracks().size() > 0) {
214  CHECK_CLOSE_ABS(fittedVertex.position(), vertexPosition, 1_mm);
215  }
216  // Do the fit with a constraint
217  Vertex<BoundParameters> fittedVertexConstraint =
218  billoirFitter.fit(tracksPtr, linearizer, vfOptionsConstr).value();
219  if (fittedVertexConstraint.tracks().size() > 0) {
220  CHECK_CLOSE_ABS(fittedVertexConstraint.position(), vertexPosition, 1_mm);
221  }
222 
223  if (debugMode) {
224  std::cout << "Fitting nTracks: " << nTracks << std::endl;
225  std::cout << "True Vertex: " << x << ", " << y << ", " << z << std::endl;
226  std::cout << "Fitted Vertex: " << fittedVertex.position() << std::endl;
227  std::cout << "Fitted constraint Vertex: "
228  << fittedVertexConstraint.position() << std::endl;
229  }
230  }
231 }
232 
233 // Dummy user-defined InputTrack type
234 struct InputTrack {
235  InputTrack(const BoundParameters& params) : m_parameters(params) {}
236 
237  const BoundParameters& parameters() const { return m_parameters; }
238 
239  // store e.g. link to original objects here
240 
241  private:
242  BoundParameters m_parameters;
243 };
244 
249 BOOST_AUTO_TEST_CASE(billoir_vertex_fitter_usertrack_test) {
250  bool debugMode = false;
251 
252  // Set up RNG
253  int mySeed = 31415;
254  std::mt19937 gen(mySeed);
255 
256  // Set up constant B-Field
257  ConstantBField bField(0.0, 0.0, 1_T);
258 
259  // Set up Eigenstepper
261 
262  // Set up propagator with void navigator
263  auto propagator =
264  std::make_shared<Propagator<EigenStepper<ConstantBField>>>(stepper);
265 
266  Linearizer::Config ltConfig(bField, propagator);
267  Linearizer linearizer(ltConfig);
268 
269  const int nEvents = 10;
270 
271  for (int eventIdx = 0; eventIdx < nEvents; ++eventIdx) {
272  unsigned int nTracks = nTracksDist(gen);
273 
274  // Create a custom std::function to extract BoundParameters from
275  // user-defined InputTrack
276  std::function<BoundParameters(InputTrack)> extractParameters =
277  [](InputTrack params) { return params.parameters(); };
278 
279  // Set up Billoir Vertex Fitter
282  vertexFitterCfg, extractParameters);
283 
284  // Constraint for vertex fit
285  Vertex<InputTrack> myConstraint;
286  // Some abitrary values
287  SpacePointSymMatrix myCovMat = SpacePointSymMatrix::Zero();
288  myCovMat(0, 0) = 30.;
289  myCovMat(1, 1) = 30.;
290  myCovMat(2, 2) = 30.;
291  myCovMat(3, 3) = 30.;
292  myConstraint.setFullCovariance(std::move(myCovMat));
293  myConstraint.setFullPosition(SpacePointVector(0, 0, 0, 0));
294 
296 
298  myConstraint);
299 
300  // Create position of vertex and perigee surface
301  double x = vXYDist(gen);
302  double y = vXYDist(gen);
303  double z = vZDist(gen);
304 
305  Vector3D vertexPosition(x, y, z);
306  std::shared_ptr<PerigeeSurface> perigeeSurface =
307  Surface::makeShared<PerigeeSurface>(Vector3D(0., 0., 0.));
308 
309  // Calculate d0 and z0 corresponding to vertex position
310  double d0V = sqrt(x * x + y * y);
311  double z0V = z;
312 
313  // Start constructing nTracks tracks in the following
314  std::vector<InputTrack> tracks;
315 
316  // Construct random track emerging from vicinity of vertex position
317  // Vector to store track objects used for vertex fit
318  for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
319  // Construct positive or negative charge randomly
320  double q = qDist(gen) < 0 ? -1. : 1.;
321 
322  // Construct random track parameters
323  BoundVector paramVec;
324  paramVec << d0V + d0Dist(gen), z0V + z0Dist(gen), phiDist(gen),
325  thetaDist(gen), q / pTDist(gen), 0.;
326 
327  // Resolutions
328  double resD0 = resIPDist(gen);
329  double resZ0 = resIPDist(gen);
330  double resPh = resAngDist(gen);
331  double resTh = resAngDist(gen);
332  double resQp = resQoPDist(gen);
333 
334  // Fill vector of track objects with simple covariance matrix
335  Covariance covMat;
336 
337  covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
338  0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
339  0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
340  tracks.push_back(InputTrack(BoundParameters(geoContext, std::move(covMat),
341  paramVec, perigeeSurface)));
342  }
343 
344  std::vector<const InputTrack*> tracksPtr;
345  for (const auto& trk : tracks) {
346  tracksPtr.push_back(&trk);
347  }
348 
349  // Do the actual fit with 4 tracks without constraint
350  Vertex<InputTrack> fittedVertex =
351  billoirFitter.fit(tracksPtr, linearizer, vfOptions).value();
352  if (fittedVertex.tracks().size() > 0) {
353  CHECK_CLOSE_ABS(fittedVertex.position(), vertexPosition, 1_mm);
354  }
355  // Do the fit with a constraint
356  Vertex<InputTrack> fittedVertexConstraint =
357  billoirFitter.fit(tracksPtr, linearizer, vfOptionsConstr).value();
358  if (fittedVertexConstraint.tracks().size() > 0) {
359  CHECK_CLOSE_ABS(fittedVertexConstraint.position(), vertexPosition, 1_mm);
360  }
361 
362  if (debugMode) {
363  std::cout << "Fitting nTracks: " << nTracks << std::endl;
364  std::cout << "True Vertex: " << x << ", " << y << ", " << z << std::endl;
365  std::cout << "Fitted Vertex: " << fittedVertex.position() << std::endl;
366  std::cout << "Fitted constraint Vertex: "
367  << fittedVertexConstraint.position() << std::endl;
368  }
369  }
370 }
371 
372 } // namespace Test
373 } // namespace Acts