Project

General

Profile

Bug #1185 » PatternRecognitionTest.cc

Dobbs, Adam, 12 December 2012 15:05

 
1
/* This file is part of MAUS: http://micewww.pp.rl.ac.uk:8080/projects/maus
2
 *
3
 * MAUS is free software: you can redistribute it and/or modify
4
 * it under the terms of the GNU General Public License as published by
5
 * the Free Software Foundation, either version 3 of the License, or
6
 * (at your option) any later version.
7
 *
8
 * MAUS is distributed in the hope that it will be useful,
9
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
11
 * GNU General Public License for more details.
12
 *
13
 * You should have received a copy of the GNU General Public License
14
 * along with MAUS.  If not, see <http://www.gnu.org/licenses/>.
15
 *
16
 */
17

    
18
#include <cmath>
19

    
20
#include "gtest/gtest.h"
21

    
22
#include "src/common_cpp/Recon/SciFi/PatternRecognition.hh"
23
#include "src/common_cpp/DataStructure/SciFiSpacePoint.hh"
24
#include "src/common_cpp/DataStructure/SciFiStraightPRTrack.hh"
25
#include "src/common_cpp/DataStructure/SciFiEvent.hh"
26
#include "src/common_cpp/DataStructure/ThreeVector.hh"
27

    
28
namespace MAUS {
29

    
30
class PatternRecognitionTest : public ::testing::Test {
31
 protected:
32
  PatternRecognitionTest()  {}
33
  virtual ~PatternRecognitionTest() {}
34
  virtual void SetUp()    {}
35
  virtual void TearDown() {}
36

    
37
  void set_up_spacepoints() {
38
    _sp1 = new SciFiSpacePoint();
39
    _sp2 = new SciFiSpacePoint();
40
    _sp3 = new SciFiSpacePoint();
41
    _sp4 = new SciFiSpacePoint();
42
    _sp5 = new SciFiSpacePoint();
43

    
44
    ThreeVector pos(-68.24883333333334, -57.810948479361, -0.652299999999741);
45
    _sp1->set_position(pos);
46
    _sp1->set_tracker(0);
47
    _sp1->set_station(1);
48
    _sp1->set_type("triplet");
49
    _sp1->set_used(false);
50

    
51
    pos.set(-62.84173333333334, -67.17694825239995, -200.6168999999991);
52
    _sp2->set_position(pos);
53
    _sp2->set_tracker(0);
54
    _sp2->set_station(2);
55
    _sp2->set_type("triplet");
56
    _sp2->set_used(false);
57

    
58
    pos.set(-56.99676666666667, -76.0964980027428, -450.4798999999994);
59
    _sp3->set_position(pos);
60
    _sp3->set_tracker(0);
61
    _sp3->set_station(3);
62
    _sp3->set_type("triplet");
63
    _sp3->set_used(false);
64

    
65
    pos.set(-47.89523333333333, -87.75184770769343, -750.4801999999991);
66
    _sp4->set_position(pos);
67
    _sp4->set_tracker(0);
68
    _sp4->set_station(4);
69
    _sp4->set_type("triplet");
70
    _sp4->set_used(false);
71

    
72
    pos.set(-35.86799999999999, -99.22774738994798, -1100.410099999999);
73
    _sp5->set_position(pos);
74
    _sp5->set_tracker(0);
75
    _sp5->set_station(5);
76
    _sp5->set_type("triplet");
77
    _sp5->set_used(false);
78
  }
79

    
80
  SciFiSpacePoint *_sp1;
81
  SciFiSpacePoint *_sp2;
82
  SciFiSpacePoint *_sp3;
83
  SciFiSpacePoint *_sp4;
84
  SciFiSpacePoint *_sp5;
85
};
86

    
87
TEST_F(PatternRecognitionTest, test_process_good) {
88

    
89
  PatternRecognition pr;
90

    
91
  // Set up the spacepoints vector
92
  set_up_spacepoints();
93
  std::vector<SciFiSpacePoint*> spnts;
94
  spnts.push_back(_sp5);
95
  spnts.push_back(_sp2);
96
  spnts.push_back(_sp3);
97
  spnts.push_back(_sp1);
98
  spnts.push_back(_sp4);
99

    
100
  // For a straight fit
101
  // ------------------
102
  SciFiEvent evt1;
103
  evt1.set_spacepoints(spnts);
104

    
105
  pr.process(false, true, evt1); // Helical off, Straight on
106

    
107
  std::vector<SciFiStraightPRTrack> strks = evt1.straightprtracks();
108
  std::vector<SciFiHelicalPRTrack> htrks = evt1.helicalprtracks();
109

    
110
  // The track parameters that should be reconstructed from the spacepoints
111
  int num_points = 5;
112

    
113
  double line_y0 = -58.85201389;
114
  double line_x0 = -68.94108927;
115
  double line_my = 0.03755825;
116
  double line_mx = -0.02902014;
117
  double line_x_chisq = 22.87148204;
118
  double line_y_chisq = 20.99052559;
119

    
120
  // Check it matches to within a tolerance epsilon
121
  double epsilon = 0.001;
122
  ASSERT_EQ(static_cast<unsigned int>(1), strks.size());
123
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
124
  EXPECT_NEAR(line_x0, strks[0].get_x0(), epsilon);
125
  EXPECT_NEAR(line_mx, strks[0].get_mx(), epsilon);
126
  EXPECT_NEAR(line_x_chisq, strks[0].get_x_chisq(), epsilon);
127
  EXPECT_NEAR(line_y0, strks[0].get_y0(), epsilon);
128
  EXPECT_NEAR(line_my, strks[0].get_my(), epsilon);
129
  EXPECT_NEAR(line_y_chisq, strks[0].get_y_chisq(), epsilon);
130
  EXPECT_EQ(num_points, strks[0].get_num_points());
131

    
132
  // For a helical fit
133
  //------------------
134
  _sp1->set_used(false);
135
  _sp2->set_used(false);
136
  _sp3->set_used(false);
137
  _sp4->set_used(false);
138
  _sp5->set_used(false);
139

    
140
  pr.process(true, false, evt1); // Helical on, Straight off
141

    
142
  strks = evt1.straightprtracks();
143
  htrks = evt1.helicalprtracks();
144

    
145
  double helix_x0 = -68.2488;
146
  double helix_y0 = -57.8109;
147
  double helix_R = 136.335;
148
  double helix_dsdz = -0.0470962; // Need to check this value is physical
149

    
150
  ASSERT_EQ(static_cast<unsigned int>(1), htrks.size());
151
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
152
  EXPECT_NEAR(helix_x0, htrks[0].get_x0(), epsilon);
153
  EXPECT_NEAR(helix_y0, htrks[0].get_y0(), epsilon);
154
  EXPECT_NEAR(helix_R, htrks[0].get_R(), epsilon);
155
  EXPECT_NEAR(helix_dsdz, htrks[0].get_dsdz(), epsilon);
156
  EXPECT_EQ(num_points, htrks[0].get_num_points());
157
}
158

    
159
TEST_F(PatternRecognitionTest, test_make_tracks) {
160

    
161
  // Set up the spacepoints vector
162
  set_up_spacepoints();
163
  std::vector<SciFiSpacePoint*> spnts;
164
  spnts.push_back(_sp5);
165
  spnts.push_back(_sp2);
166
  spnts.push_back(_sp1);
167

    
168
  PatternRecognition pr;
169
  int n_stations = 5;
170

    
171
  // Set up the spacepoints by station 2D vector
172
  std::vector< std::vector<SciFiSpacePoint*> > spnts_by_station(n_stations);
173
  pr.sort_by_station(spnts, spnts_by_station);
174

    
175
  SciFiEvent evt;
176
  bool track_type = 0; // Straight tracks
177
  int tracker_num = 0;
178

    
179
  // The track parameters that should be reconstructed from the spacepoints
180
  double x_chisq = 22.87148204;
181
  double y_chisq = 20.99052559;
182
  double y0 = -58.85201389;
183
  double x0 = -68.94108927;
184
  double my = 0.03755825;
185
  double mx = -0.02902014;
186

    
187
  // Make a 3 point track
188
  // ---------------------
189
  pr.make_all_tracks(track_type, tracker_num, spnts_by_station, evt);
190
  std::vector<SciFiStraightPRTrack> strks = evt.straightprtracks();
191
  std::vector<SciFiHelicalPRTrack> htrks = evt.helicalprtracks();
192

    
193
  // Check it matches to within a tolerance
194
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
195
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
196
  EXPECT_EQ(3, strks[0].get_num_points());
197
  EXPECT_NEAR(x0, strks[0].get_x0(), 1);
198
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.001);
199
  EXPECT_NEAR(0.9, strks[0].get_x_chisq(), 0.1);
200
  EXPECT_NEAR(y0, strks[0].get_y0(), 1);
201
  EXPECT_NEAR(my, strks[0].get_my(), 0.001);
202
  EXPECT_NEAR(13.3, strks[0].get_y_chisq(), 0.1);
203

    
204
  // Make a 4 point track
205
  // ---------------------
206
  spnts.push_back(_sp3);
207
  _sp1->set_used(false);
208
  _sp2->set_used(false);
209
  _sp3->set_used(false);
210
  _sp5->set_used(false);
211

    
212
  spnts_by_station.clear();
213
  spnts_by_station.resize(0);
214
  spnts_by_station.resize(n_stations);
215
  pr.sort_by_station(spnts, spnts_by_station);
216
  strks.resize(0);
217
  evt.set_straightprtrack(strks);
218
  pr.make_all_tracks(track_type, tracker_num, spnts_by_station, evt);
219
  strks = evt.straightprtracks();
220
  htrks = evt.helicalprtracks();
221

    
222
  // Check it matches to within a tolerance
223
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
224
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
225
  EXPECT_EQ(4, strks[0].get_num_points());
226
  EXPECT_NEAR(x0, strks[0].get_x0(), 1);
227
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.001);
228
  EXPECT_NEAR(17.5, strks[0].get_x_chisq(), 0.1);
229
  EXPECT_NEAR(y0, strks[0].get_y0(), 1);
230
  EXPECT_NEAR(my, strks[0].get_my(), 0.001);
231
  EXPECT_NEAR(16.0, strks[0].get_y_chisq(), 0.1);
232

    
233
  // Make a 5 point track
234
  // ---------------------
235
  spnts.push_back(_sp4);
236
  _sp1->set_used(false);
237
  _sp2->set_used(false);
238
  _sp3->set_used(false);
239
  _sp4->set_used(false);
240
  _sp5->set_used(false);
241

    
242
  spnts_by_station.clear();
243
  spnts_by_station.resize(0);
244
  spnts_by_station.resize(n_stations);
245
  pr.sort_by_station(spnts, spnts_by_station);
246
  strks.resize(0);
247
  evt.set_straightprtrack(strks);
248
  pr.make_all_tracks(track_type, tracker_num, spnts_by_station, evt);
249
  strks = evt.straightprtracks();
250
  htrks = evt.helicalprtracks();
251

    
252
  // Check it matches to within a tolerance
253
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
254
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
255
  EXPECT_EQ(5, strks[0].get_num_points());
256
  EXPECT_NEAR(x0, strks[0].get_x0(), 1);
257
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.001);
258
  EXPECT_NEAR(x_chisq, strks[0].get_x_chisq(), 0.1);
259
  EXPECT_NEAR(y0, strks[0].get_y0(), 1);
260
  EXPECT_NEAR(my, strks[0].get_my(), 0.001);
261
  EXPECT_NEAR(y_chisq, strks[0].get_y_chisq(), 0.1);
262
}
263

    
264
TEST_F(PatternRecognitionTest, test_make_4pt_tracks) {
265

    
266
  // Set up the spacepoints vector
267
  set_up_spacepoints();
268
  std::vector<SciFiSpacePoint*> spnts;
269
  spnts.push_back(_sp5);
270
  spnts.push_back(_sp3);
271
  spnts.push_back(_sp2);
272
  spnts.push_back(_sp1);
273

    
274
  PatternRecognition pr;
275
  int n_stations = 5;
276

    
277
  // Set up the spacepoints by station 2D vector
278
  std::vector< std::vector<SciFiSpacePoint*> > spnts_by_station(n_stations);
279
  pr.sort_by_station(spnts, spnts_by_station);
280

    
281
  SciFiEvent evt;
282
  bool track_type = 0; // Straight tracks
283

    
284
  // The track parameters that should be reconstructed from the spacepoints in 5 pt track case
285
  int num_points = 4;
286
  double y0 = -58.85201389;
287
  double x0 = -68.94108927;
288
  double my = 0.03755825;
289
  double mx = -0.02902014;
290
  std::vector<SciFiHelicalPRTrack> htrks;
291
  std::vector<SciFiStraightPRTrack> strks;
292

    
293
  // Make a 4 point track with 4 spacepoints
294
  // ---------------------------------------
295
  _sp1->set_used(false);
296
  _sp2->set_used(false);
297
  _sp3->set_used(false);
298
  _sp5->set_used(false);
299

    
300
  spnts_by_station.clear();
301
  spnts_by_station.resize(0);
302
  spnts_by_station.resize(n_stations);
303
  pr.sort_by_station(spnts, spnts_by_station);
304
  strks.resize(0);
305

    
306
  pr.make_4tracks(track_type, spnts_by_station, strks, htrks);
307

    
308
  // Check it matches to within a tolerance
309
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
310
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
311
  EXPECT_EQ(num_points, strks[0].get_num_points());
312
  EXPECT_NEAR(x0, strks[0].get_x0(), 3);
313
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.002);
314
  EXPECT_NEAR(17.5, strks[0].get_x_chisq(), 0.1);
315
  EXPECT_NEAR(y0, strks[0].get_y0(), 3);
316
  EXPECT_NEAR(my, strks[0].get_my(), 0.005);
317
  EXPECT_NEAR(15.9, strks[0].get_y_chisq(), 0.1);
318

    
319
    // Make a 4 point track with 5 spacepoints
320
  // ---------------------------------------
321
  spnts.push_back(_sp4);
322
  _sp1->set_used(false);
323
  _sp2->set_used(false);
324
  _sp3->set_used(false);
325
  _sp4->set_used(false);
326
  _sp5->set_used(false);
327

    
328
  spnts_by_station.clear();
329
  spnts_by_station.resize(0);
330
  spnts_by_station.resize(n_stations);
331
  pr.sort_by_station(spnts, spnts_by_station);
332
  strks.resize(0);
333

    
334
  pr.make_4tracks(track_type, spnts_by_station, strks, htrks);
335

    
336
  // Check it matches to within a tolerance
337
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
338
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
339
  EXPECT_EQ(num_points, strks[0].get_num_points());
340
  EXPECT_NEAR(x0, strks[0].get_x0(), 3);
341
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.01); // Needed a wider tolerance than the others
342
  EXPECT_NEAR(16.3, strks[0].get_x_chisq(), 0.1);
343
  EXPECT_NEAR(y0, strks[0].get_y0(), 3);
344
  EXPECT_NEAR(my, strks[0].get_my(), 0.005);
345
  EXPECT_NEAR(4.5, strks[0].get_y_chisq(), 0.1);
346
}
347

    
348
TEST_F(PatternRecognitionTest, test_make_3pt_tracks) {
349

    
350
  // Set up the spacepoints vector
351
  set_up_spacepoints();
352
  std::vector<SciFiSpacePoint*> spnts;
353
  spnts.push_back(_sp5);
354
  spnts.push_back(_sp2);
355
  spnts.push_back(_sp1);
356

    
357
  PatternRecognition pr;
358
  int n_stations = 5;
359

    
360
  // Set up the spacepoints by station 2D vector
361
  std::vector< std::vector<SciFiSpacePoint*> > spnts_by_station(n_stations);
362
  pr.sort_by_station(spnts, spnts_by_station);
363

    
364
  SciFiEvent evt;
365
  bool track_type = 0; // Straight tracks
366

    
367
  // The track parameters that should be reconstructed from the spacepoints in 5 pt track case
368
  int num_points = 3;
369
  double y0 = -58.85201389;
370
  double x0 = -68.94108927;
371
  double my = 0.03755825;
372
  double mx = -0.02902014;
373
  std::vector<SciFiHelicalPRTrack> htrks;
374
  std::vector<SciFiStraightPRTrack> strks;
375

    
376
  // Make a 3 point track with 3 spacepoints
377
  // ---------------------------------------
378
  pr.make_3tracks(track_type, spnts_by_station, strks, htrks);
379

    
380
  // Check it matches to within a tolerance
381
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
382
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
383
  EXPECT_EQ(num_points, strks[0].get_num_points());
384
  EXPECT_NEAR(x0, strks[0].get_x0(), 3);
385
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.002);
386
  EXPECT_NEAR(0.9, strks[0].get_x_chisq(), 0.1);
387
  EXPECT_NEAR(y0, strks[0].get_y0(), 3);
388
  EXPECT_NEAR(my, strks[0].get_my(), 0.005);
389
  EXPECT_NEAR(13.3, strks[0].get_y_chisq(), 0.1);
390

    
391
  // Make a 3 point track with 4 spacepoints
392
  // ---------------------------------------
393
  spnts.push_back(_sp3);
394
  _sp1->set_used(false);
395
  _sp2->set_used(false);
396
  _sp3->set_used(false);
397
  _sp5->set_used(false);
398

    
399
  spnts_by_station.clear();
400
  spnts_by_station.resize(0);
401
  spnts_by_station.resize(n_stations);
402
  pr.sort_by_station(spnts, spnts_by_station);
403
  strks.resize(0);
404

    
405
  pr.make_3tracks(track_type, spnts_by_station, strks, htrks);
406

    
407
  // Check it matches to within a tolerance
408
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
409
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
410
  EXPECT_EQ(num_points, strks[0].get_num_points());
411
  EXPECT_NEAR(x0, strks[0].get_x0(), 3);
412
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.002);
413
  EXPECT_NEAR(11.3, strks[0].get_x_chisq(), 0.1);
414
  EXPECT_NEAR(y0, strks[0].get_y0(), 3);
415
  EXPECT_NEAR(my, strks[0].get_my(), 0.005);
416
  EXPECT_NEAR(0.0015, strks[0].get_y_chisq(), 0.001);
417

    
418
    // Make a 3 point track with 5 spacepoints
419
  // ---------------------------------------
420
  spnts.push_back(_sp4);
421
  _sp1->set_used(false);
422
  _sp2->set_used(false);
423
  _sp3->set_used(false);
424
  _sp4->set_used(false);
425
  _sp5->set_used(false);
426

    
427
  spnts_by_station.clear();
428
  spnts_by_station.resize(0);
429
  spnts_by_station.resize(n_stations);
430
  pr.sort_by_station(spnts, spnts_by_station);
431
  strks.resize(0);
432

    
433
  pr.make_3tracks(track_type, spnts_by_station, strks, htrks);
434

    
435
  // Check it matches to within a tolerance
436
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
437
  EXPECT_EQ(static_cast<unsigned int>(0), htrks.size());
438
  EXPECT_EQ(num_points, strks[0].get_num_points());
439
  EXPECT_NEAR(x0, strks[0].get_x0(), 3);
440
  EXPECT_NEAR(mx, strks[0].get_mx(), 0.01); // Needed a wider tolerance than the others
441
  EXPECT_NEAR(1.8, strks[0].get_x_chisq(), 0.1);
442
  EXPECT_NEAR(y0, strks[0].get_y0(), 3);
443
  EXPECT_NEAR(my, strks[0].get_my(), 0.005);
444
  EXPECT_NEAR(4.2, strks[0].get_y_chisq(), 0.1);
445
}
446

    
447
TEST_F(PatternRecognitionTest, test_make_straight_tracks) {
448

    
449
  int n_stations = 5;
450
  PatternRecognition pr;
451

    
452
  // Set up the spacepoints vector
453
  set_up_spacepoints();
454
  std::vector<SciFiSpacePoint*> spnts;
455
  spnts.push_back(_sp5);
456
  spnts.push_back(_sp2);
457
  spnts.push_back(_sp3);
458
  spnts.push_back(_sp1);
459
  spnts.push_back(_sp4);
460

    
461
  // Set up the spacepoints by station 2D vector
462
  std::vector< std::vector<SciFiSpacePoint*> > spnts_by_station(n_stations);
463
  pr.sort_by_station(spnts, spnts_by_station);
464

    
465
  // Check the spacepoints have setup correctly
466
  EXPECT_EQ(_sp1, spnts_by_station[0][0]);
467
  EXPECT_EQ(_sp2, spnts_by_station[1][0]);
468
  EXPECT_EQ(_sp3, spnts_by_station[2][0]);
469
  EXPECT_EQ(_sp4, spnts_by_station[3][0]);
470
  EXPECT_EQ(_sp5, spnts_by_station[4][0]);
471
  EXPECT_EQ(-68.24883333333334, spnts_by_station[0][0]->get_position().x());
472

    
473
  std::vector<int> ignore_stations;
474
  std::vector<SciFiStraightPRTrack> strks;
475

    
476
  // The track parameters that should be reconstructed from the spacepoints
477
  int num_points = 5;
478
  double x_chisq = 22.87148204;
479
  double y_chisq = 20.99052559;
480
  double y0 = -58.85201389;
481
  double x0 = -68.94108927;
482
  double my = 0.03755825;
483
  double mx = -0.02902014;
484

    
485
  // Make the track from the spacepoints
486
  pr.make_straight_tracks(num_points, ignore_stations, spnts_by_station, strks);
487

    
488
  // Check it matches to within a tolerance epsilon
489
  double epsilon = 0.000001;
490
  EXPECT_EQ(static_cast<unsigned int>(1), strks.size());
491
  EXPECT_NEAR(x0, strks[0].get_x0(), epsilon);
492
  EXPECT_NEAR(mx, strks[0].get_mx(), epsilon);
493
  EXPECT_NEAR(x_chisq, strks[0].get_x_chisq(), epsilon);
494
  EXPECT_NEAR(y0, strks[0].get_y0(), epsilon);
495
  EXPECT_NEAR(my, strks[0].get_my(), epsilon);
496
  EXPECT_NEAR(y_chisq, strks[0].get_y_chisq(), epsilon);
497
}
498

    
499
TEST_F(PatternRecognitionTest, test_set_ignore_stations) {
500

    
501
  PatternRecognition pr;
502
  std::vector<int> ignore_stations(0);
503
  int is1, is2;
504

    
505
  EXPECT_TRUE(pr.set_ignore_stations(ignore_stations, is1, is2));
506
  EXPECT_EQ(-1, is1);
507
  EXPECT_EQ(-1, is2);
508

    
509
  ignore_stations.push_back(4);
510
  EXPECT_TRUE(pr.set_ignore_stations(ignore_stations, is1, is2));
511
  EXPECT_EQ(4, is1);
512
  EXPECT_EQ(-1, is2);
513

    
514
  ignore_stations.push_back(1);
515
  EXPECT_TRUE(pr.set_ignore_stations(ignore_stations, is1, is2));
516
  EXPECT_EQ(4, is1);
517
  EXPECT_EQ(1, is2);
518

    
519
  ignore_stations.push_back(2);
520
  EXPECT_FALSE(pr.set_ignore_stations(ignore_stations, is1, is2));
521
  EXPECT_EQ(-1, is1);
522
  EXPECT_EQ(-1, is2);
523
}
524

    
525
TEST_F(PatternRecognitionTest, test_set_seed_stations) {
526

    
527
  PatternRecognition pr;
528
  std::vector<int> ignore_stations(0);
529
  int outer_st_num, inner_st_num, mid_st_num;
530

    
531
  // 5 pt track case
532
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
533
  EXPECT_EQ(4, outer_st_num);
534
  EXPECT_EQ(0, inner_st_num);
535
  EXPECT_EQ(2, mid_st_num);
536

    
537
  // 4 pt track case
538
  ignore_stations.push_back(4);
539
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
540
  EXPECT_EQ(3, outer_st_num);
541
  EXPECT_EQ(0, inner_st_num);
542
  EXPECT_EQ(2, mid_st_num);
543

    
544
  ignore_stations[0] = 3;
545
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
546
  EXPECT_EQ(4, outer_st_num);
547
  EXPECT_EQ(0, inner_st_num);
548
  EXPECT_EQ(2, mid_st_num);
549

    
550
  ignore_stations[0] = 2;
551
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
552
  EXPECT_EQ(4, outer_st_num);
553
  EXPECT_EQ(0, inner_st_num);
554
  EXPECT_EQ(1, mid_st_num);
555

    
556
  ignore_stations[0] = 1;
557
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
558
  EXPECT_EQ(4, outer_st_num);
559
  EXPECT_EQ(0, inner_st_num);
560
  EXPECT_EQ(2, mid_st_num);
561

    
562
  ignore_stations[0] = 0;
563
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
564
  EXPECT_EQ(4, outer_st_num);
565
  EXPECT_EQ(1, inner_st_num);
566
  EXPECT_EQ(2, mid_st_num);
567

    
568
  // 3 pt track case (not testing all possible combinations)
569
  ignore_stations[0] = 4;
570
  ignore_stations.push_back(3);
571
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
572
  EXPECT_EQ(2, outer_st_num);
573
  EXPECT_EQ(0, inner_st_num);
574
  EXPECT_EQ(1, mid_st_num);
575

    
576
  ignore_stations[0] = 4;
577
  ignore_stations[1] = 2;
578
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
579
  EXPECT_EQ(3, outer_st_num);
580
  EXPECT_EQ(0, inner_st_num);
581
  EXPECT_EQ(1, mid_st_num);
582

    
583
  ignore_stations[0] = 4;
584
  ignore_stations[1] = 1;
585
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
586
  EXPECT_EQ(3, outer_st_num);
587
  EXPECT_EQ(0, inner_st_num);
588
  EXPECT_EQ(2, mid_st_num);
589

    
590
  ignore_stations[0] = 4;
591
  ignore_stations[1] = 0;
592
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
593
  EXPECT_EQ(3, outer_st_num);
594
  EXPECT_EQ(1, inner_st_num);
595
  EXPECT_EQ(2, mid_st_num);
596

    
597
  ignore_stations[0] = 3;
598
  ignore_stations[1] = 4;
599
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
600
  EXPECT_EQ(2, outer_st_num);
601
  EXPECT_EQ(0, inner_st_num);
602
  EXPECT_EQ(1, mid_st_num);
603

    
604
  ignore_stations[0] = 3;
605
  ignore_stations[1] = 2;
606
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
607
  EXPECT_EQ(4, outer_st_num);
608
  EXPECT_EQ(0, inner_st_num);
609
  EXPECT_EQ(1, mid_st_num);
610

    
611
  ignore_stations[0] = 3;
612
  ignore_stations[1] = 1;
613
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
614
  EXPECT_EQ(4, outer_st_num);
615
  EXPECT_EQ(0, inner_st_num);
616
  EXPECT_EQ(2, mid_st_num);
617

    
618
  ignore_stations[0] = 3;
619
  ignore_stations[1] = 0;
620
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
621
  EXPECT_EQ(4, outer_st_num);
622
  EXPECT_EQ(1, inner_st_num);
623
  EXPECT_EQ(2, mid_st_num);
624

    
625
  ignore_stations[0] = 2;
626
  ignore_stations[1] = 1;
627
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
628
  EXPECT_EQ(4, outer_st_num);
629
  EXPECT_EQ(0, inner_st_num);
630
  EXPECT_EQ(3, mid_st_num);
631

    
632
  ignore_stations[0] = 2;
633
  ignore_stations[1] = 0;
634
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
635
  EXPECT_EQ(4, outer_st_num);
636
  EXPECT_EQ(1, inner_st_num);
637
  EXPECT_EQ(3, mid_st_num);
638

    
639
  ignore_stations[0] = 1;
640
  ignore_stations[1] = 0;
641
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
642
  EXPECT_EQ(4, outer_st_num);
643
  EXPECT_EQ(2, inner_st_num);
644
  EXPECT_EQ(3, mid_st_num);
645

    
646
  ignore_stations[0] = 0;
647
  ignore_stations[1] = 1;
648
  EXPECT_TRUE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
649
  EXPECT_EQ(4, outer_st_num);
650
  EXPECT_EQ(2, inner_st_num);
651
  EXPECT_EQ(3, mid_st_num);
652

    
653
  // Test input error cases
654

    
655
  // Same two stations numbers entered
656
  ignore_stations[0] = 4;
657
  ignore_stations.push_back(4);
658
  EXPECT_FALSE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
659

    
660
  // Out of bounds station numbers
661
  ignore_stations.resize(1);
662
  ignore_stations[0] = 12;
663
  EXPECT_FALSE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
664

    
665
  ignore_stations[0] = -1;
666
  EXPECT_FALSE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
667

    
668
  // Too large ignore_stations vector
669
  ignore_stations.resize(3);
670
  EXPECT_FALSE(pr.set_seed_stations(ignore_stations, outer_st_num, inner_st_num, mid_st_num));
671
}
672

    
673

    
674
TEST_F(PatternRecognitionTest, test_set_end_stations) {
675

    
676
  PatternRecognition pr;
677
  std::vector<int> ignore_stations(0);
678
  int outer_st_num, inner_st_num;
679

    
680
  // 5 pt track case
681
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
682
  EXPECT_EQ(4, outer_st_num);
683
  EXPECT_EQ(0, inner_st_num);
684

    
685
  // 4 pt track case
686
  ignore_stations.push_back(4);
687
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
688
  EXPECT_EQ(3, outer_st_num);
689
  EXPECT_EQ(0, inner_st_num);
690

    
691
  ignore_stations[0] = 3;
692
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
693
  EXPECT_EQ(4, outer_st_num);
694
  EXPECT_EQ(0, inner_st_num);
695

    
696
  ignore_stations[0] = 2;
697
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
698
  EXPECT_EQ(4, outer_st_num);
699
  EXPECT_EQ(0, inner_st_num);
700

    
701
  ignore_stations[0] = 1;
702
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
703
  EXPECT_EQ(4, outer_st_num);
704
  EXPECT_EQ(0, inner_st_num);
705

    
706
  ignore_stations[0] = 0;
707
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
708
  EXPECT_EQ(4, outer_st_num);
709
  EXPECT_EQ(1, inner_st_num);
710

    
711
  // 3 pt track case (not testing all possible combinations)
712
  ignore_stations[0] = 4;
713
  ignore_stations.push_back(3);
714
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
715
  EXPECT_EQ(2, outer_st_num);
716
  EXPECT_EQ(0, inner_st_num);
717

    
718
  ignore_stations[0] = 4;
719
  ignore_stations[1] = 2;
720
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
721
  EXPECT_EQ(3, outer_st_num);
722
  EXPECT_EQ(0, inner_st_num);
723

    
724
  ignore_stations[0] = 4;
725
  ignore_stations[1] = 1;
726
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
727
  EXPECT_EQ(3, outer_st_num);
728
  EXPECT_EQ(0, inner_st_num);
729

    
730
  ignore_stations[0] = 4;
731
  ignore_stations[1] = 0;
732
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
733
  EXPECT_EQ(3, outer_st_num);
734
  EXPECT_EQ(1, inner_st_num);
735

    
736
  ignore_stations[0] = 3;
737
  ignore_stations[1] = 4;
738
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
739
  EXPECT_EQ(2, outer_st_num);
740
  EXPECT_EQ(0, inner_st_num);
741

    
742
  ignore_stations[0] = 3;
743
  ignore_stations[1] = 2;
744
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
745
  EXPECT_EQ(4, outer_st_num);
746
  EXPECT_EQ(0, inner_st_num);
747

    
748
  ignore_stations[0] = 3;
749
  ignore_stations[1] = 1;
750
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
751
  EXPECT_EQ(4, outer_st_num);
752
  EXPECT_EQ(0, inner_st_num);
753

    
754
  ignore_stations[0] = 3;
755
  ignore_stations[1] = 0;
756
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
757
  EXPECT_EQ(4, outer_st_num);
758
  EXPECT_EQ(1, inner_st_num);
759

    
760
  ignore_stations[0] = 2;
761
  ignore_stations[1] = 1;
762
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
763
  EXPECT_EQ(4, outer_st_num);
764
  EXPECT_EQ(0, inner_st_num);
765

    
766
  ignore_stations[0] = 2;
767
  ignore_stations[1] = 0;
768
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
769
  EXPECT_EQ(4, outer_st_num);
770
  EXPECT_EQ(1, inner_st_num);
771

    
772
  ignore_stations[0] = 1;
773
  ignore_stations[1] = 0;
774
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
775
  EXPECT_EQ(4, outer_st_num);
776
  EXPECT_EQ(2, inner_st_num);
777

    
778
  ignore_stations[0] = 0;
779
  ignore_stations[1] = 1;
780
  EXPECT_TRUE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
781
  EXPECT_EQ(4, outer_st_num);
782
  EXPECT_EQ(2, inner_st_num);
783

    
784
  // Test input error cases
785

    
786
  // Same two stations numbers entered
787
  ignore_stations[0] = 4;
788
  ignore_stations.push_back(4);
789
  EXPECT_FALSE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
790

    
791
  // Out of bounds station numbers
792
  ignore_stations.resize(1);
793
  ignore_stations[0] = 12;
794
  EXPECT_FALSE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
795

    
796
  ignore_stations[0] = -1;
797
  EXPECT_FALSE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
798

    
799
  // Too large ignore_stations vector
800
  ignore_stations.resize(3);
801
  EXPECT_FALSE(pr.set_end_stations(ignore_stations, outer_st_num, inner_st_num));
802
}
803

    
804
TEST_F(PatternRecognitionTest, test_sort_by_station) {
805

    
806
  PatternRecognition pr;
807

    
808
  SciFiSpacePoint *sp1 = new SciFiSpacePoint();
809
  SciFiSpacePoint *sp2 = new SciFiSpacePoint();
810
  SciFiSpacePoint *sp3 = new SciFiSpacePoint();
811
  SciFiSpacePoint *sp4 = new SciFiSpacePoint();
812
  SciFiSpacePoint *sp5 = new SciFiSpacePoint();
813

    
814
  sp1->set_station(1);
815
  sp2->set_station(2);
816
  sp3->set_station(3);
817
  sp4->set_station(4);
818
  sp5->set_station(5);
819

    
820
  std::vector<SciFiSpacePoint*> spnts;
821
  spnts.push_back(sp5);
822
  spnts.push_back(sp2);
823
  spnts.push_back(sp3);
824
  spnts.push_back(sp1);
825
  spnts.push_back(sp4);
826

    
827
  std::vector< std::vector<SciFiSpacePoint*> > spnts_by_station(5);
828
  pr.sort_by_station(spnts, spnts_by_station);
829
  EXPECT_EQ(sp1, spnts_by_station[0][0]);
830
  EXPECT_EQ(sp2, spnts_by_station[1][0]);
831
  EXPECT_EQ(sp3, spnts_by_station[2][0]);
832
  EXPECT_EQ(sp4, spnts_by_station[3][0]);
833
  EXPECT_EQ(sp5, spnts_by_station[4][0]);
834
}
835

    
836
TEST_F(PatternRecognitionTest, test_stations_with_unused_sp) {
837

    
838
  // Set up spacepoints, leaving station 3 empty to check function copes with an empty station
839
  SciFiSpacePoint *sp1 = new SciFiSpacePoint();
840
  SciFiSpacePoint *sp2 = new SciFiSpacePoint();
841
  // SciFiSpacePoint *sp3 = new SciFiSpacePoint();
842
  SciFiSpacePoint *sp4 = new SciFiSpacePoint();
843
  SciFiSpacePoint *sp4_1 = new SciFiSpacePoint();
844
  SciFiSpacePoint *sp5 = new SciFiSpacePoint();
845

    
846
  sp1->set_station(1);
847
  sp2->set_station(2);
848
  // sp3->set_station(3);
849
  sp4->set_station(4);
850
  sp4_1->set_station(4);
851
  sp5->set_station(5);
852

    
853
  sp1->set_used(true);
854
  sp2->set_used(true);
855
  // sp3->set_used(true);
856
  sp4->set_used(false);
857
  sp4_1->set_used(true);
858
  sp5->set_used(false);
859

    
860
  std::vector<SciFiSpacePoint*> spnts;
861
  spnts.push_back(sp5);
862
  spnts.push_back(sp2);
863
  // spnts.push_back(sp3);
864
  spnts.push_back(sp1);
865
  spnts.push_back(sp4);
866
  spnts.push_back(sp4_1);
867

    
868
  SpacePoint2dPArray spnts_by_station(5);
869

    
870
  PatternRecognition pr;
871
  pr.sort_by_station(spnts, spnts_by_station);
872
  ASSERT_EQ(static_cast<unsigned int>(5), spnts_by_station.size());
873
  ASSERT_EQ(static_cast<unsigned int>(1), spnts_by_station[0].size());
874
  ASSERT_EQ(static_cast<unsigned int>(1), spnts_by_station[1].size());
875
  ASSERT_EQ(static_cast<unsigned int>(2), spnts_by_station[3].size());
876
  ASSERT_EQ(static_cast<unsigned int>(1), spnts_by_station[4].size());
877

    
878
  std::vector<int> stations_hit, stations_not_hit;
879
  pr.stations_with_unused_spnts(spnts_by_station, stations_hit, stations_not_hit);
880

    
881
  ASSERT_EQ(static_cast<unsigned int>(2), stations_hit.size());
882
  ASSERT_EQ(static_cast<unsigned int>(3), stations_not_hit.size());
883
  EXPECT_EQ(3, stations_hit[0]);
884
  EXPECT_EQ(4, stations_hit[1]);
885
  EXPECT_EQ(0, stations_not_hit[0]);
886
  EXPECT_EQ(1, stations_not_hit[1]);
887
  EXPECT_EQ(2, stations_not_hit[2]);
888

    
889
  int stats_with_unused = pr.num_stations_with_unused_spnts(spnts_by_station);
890
  EXPECT_EQ(2, stats_with_unused);
891
}
892

    
893
TEST_F(PatternRecognitionTest, test_circle_fit) {
894

    
895
  PatternRecognition pr;
896

    
897
  // Set up spacepoints from an MC helical track
898
  SciFiSpacePoint *sp1 = new SciFiSpacePoint();
899
  SciFiSpacePoint *sp2 = new SciFiSpacePoint();
900
  SciFiSpacePoint *sp3 = new SciFiSpacePoint();
901
  SciFiSpacePoint *sp4 = new SciFiSpacePoint();
902
  SciFiSpacePoint *sp5 = new SciFiSpacePoint();
903

    
904
  sp1->set_station(1);
905
  sp2->set_station(2);
906
  sp3->set_station(3);
907
  sp4->set_station(4);
908
  sp5->set_station(5);
909

    
910
  ThreeVector pos(14.1978, 9.05992, 0.6523);
911
  sp1->set_position(pos);
912
  pos.set(-7.97067, 10.3542, 200.652);
913
  sp2->set_position(pos);
914
  pos.set(-11.4578, -16.3941, 450.652);
915
  sp3->set_position(pos);
916
  pos.set(19.9267, -12.0799, 750.652);
917
  sp4->set_position(pos);
918
  pos.set(-5.47983, 12.9427, 1100.65);
919
  sp5->set_position(pos);
920

    
921
  std::vector<SciFiSpacePoint*> spnts;
922
  spnts.push_back(sp5);
923
  spnts.push_back(sp2);
924
  spnts.push_back(sp3);
925
  spnts.push_back(sp1);
926
  spnts.push_back(sp4);
927

    
928
  SimpleCircle circle;
929
  bool good_radius = pr.circle_fit(spnts, circle);
930

    
931
  double epsilon = 0.01;
932

    
933
  ASSERT_TRUE(good_radius);
934
  EXPECT_NEAR(circle.get_x0(), 2.56, epsilon);
935
  EXPECT_NEAR(circle.get_y0(), -4.62, epsilon);
936
  EXPECT_NEAR(circle.get_R(), 18.56, epsilon);
937
  EXPECT_NEAR(circle.get_chisq(), 0.0994, epsilon);
938
}
939

    
940
TEST_F(PatternRecognitionTest, test_linear_fit) {
941

    
942
  PatternRecognition pr;
943

    
944
  // Test with a simple line, c = 2, m = 1, with three points, small errors
945
  std::vector<double> x, y, y_err;
946
  x.push_back(1);
947
  x.push_back(2);
948
  x.push_back(3);
949
  y.push_back(3);
950
  y.push_back(4);
951
  y.push_back(5);
952
  y_err.push_back(0.05);
953
  y_err.push_back(0.05);
954
  y_err.push_back(0.05);
955

    
956
  SimpleLine line;
957
  pr.linear_fit(x, y, y_err, line);
958

    
959
  double epsilon = 0.00001;
960

    
961
  EXPECT_NEAR(2.0, line.get_c(), epsilon);
962
  EXPECT_NEAR(1.0, line.get_m(), epsilon);
963
}
964

    
965
TEST_F(PatternRecognitionTest, test_calculate_dipangle) {
966

    
967
  PatternRecognition pr;
968

    
969
  // Set up spacepoints from an MC helical track
970
  SciFiSpacePoint *sp1 = new SciFiSpacePoint();
971
  SciFiSpacePoint *sp2 = new SciFiSpacePoint();
972
  SciFiSpacePoint *sp3 = new SciFiSpacePoint();
973
  SciFiSpacePoint *sp4 = new SciFiSpacePoint();
974
  SciFiSpacePoint *sp5 = new SciFiSpacePoint();
975

    
976
  sp1->set_station(1);
977
  sp2->set_station(2);
978
  sp3->set_station(3);
979
  sp4->set_station(4);
980
  sp5->set_station(5);
981

    
982
  ThreeVector pos(14.1978, 9.05992, 0.6523);
983
  sp1->set_position(pos);
984
  pos.set(-7.97067, 10.3542, 200.652);
985
  sp2->set_position(pos);
986
  pos.set(-11.4578, -16.3941, 450.652);
987
  sp3->set_position(pos);
988
  pos.set(19.9267, -12.0799, 750.652);
989
  sp4->set_position(pos);
990
  pos.set(-5.47983, 12.9427, 1100.65);
991
  sp5->set_position(pos);
992

    
993
  std::vector<SciFiSpacePoint*> spnts;
994
  spnts.push_back(sp1);
995
  spnts.push_back(sp2);
996
  spnts.push_back(sp3);
997
  spnts.push_back(sp4);
998
  spnts.push_back(sp5);
999

    
1000
  SimpleCircle circle;
1001
  bool good_radius = pr.circle_fit(spnts, circle);
1002

    
1003
  double epsilon = 0.01;
1004

    
1005
  ASSERT_TRUE(good_radius);
1006
  EXPECT_NEAR(circle.get_x0(), 2.56, epsilon);
1007
  EXPECT_NEAR(circle.get_y0(), -4.62, epsilon);
1008
  EXPECT_NEAR(circle.get_R(), 18.56, epsilon);
1009
  EXPECT_NEAR(circle.get_chisq(), 0.0994, epsilon);
1010

    
1011
  SimpleLine line_sz;
1012
  std::vector<double> dphi;
1013
  double phi_0;
1014

    
1015
  pr.calculate_dipangle(spnts, circle, dphi, line_sz, phi_0);
1016

    
1017
  EXPECT_NEAR(line_sz.get_c(), -1.09, epsilon);
1018
  EXPECT_NEAR(line_sz.get_m(), 0.126, epsilon);
1019
  EXPECT_NEAR(line_sz.get_chisq(), 0.440, epsilon);
1020
}
1021

    
1022
TEST_F(PatternRecognitionTest, test_AB_ratio) {
1023

    
1024
  PatternRecognition pr;
1025
  double phi_i = 1.0;
1026
  double phi_j = 0.5;
1027
  double z_i = 200.0;
1028
  double z_j = 450.0;
1029

    
1030
  double epsilon = 0.01;
1031

    
1032
  bool result = pr.AB_ratio(phi_i, phi_j, z_i, z_j);
1033
  ASSERT_TRUE(result);
1034
  EXPECT_NEAR(phi_i, 7.28319, epsilon);
1035
  EXPECT_NEAR(phi_j, 6.783, epsilon);
1036
  EXPECT_EQ(z_i, 200.0);
1037
  EXPECT_EQ(z_j, 450.0);
1038
}
1039

    
1040
/*
1041
TEST_F(PatternRecognitionTest, test_process_bad) {
1042

    
1043
  PatternRecognition pr;
1044

    
1045
  // Now give it random spacepoints and check it produces no tracks
1046
  SciFiSpacePoint *sp1 = new SciFiSpacePoint();
1047
  SciFiSpacePoint *sp2 = new SciFiSpacePoint();
1048
  SciFiSpacePoint *sp3 = new SciFiSpacePoint();
1049
  SciFiSpacePoint *sp4 = new SciFiSpacePoint();
1050
  SciFiSpacePoint *sp5 = new SciFiSpacePoint();
1051

    
1052
  ThreeVector pos(-50.0, -50.0, -0.652299999999741);
1053
  sp1->set_position(pos);
1054
  sp1->set_tracker(0);
1055
  sp1->set_station(1);
1056
  sp1->set_type("triplet");
1057

    
1058
  pos.set(-20, -70, -200.6168999999991);
1059
  sp2->set_position(pos);
1060
  sp2->set_tracker(0);
1061
  sp2->set_station(2);
1062
  sp2->set_type("triplet");
1063

    
1064
  pos.set(40.0, -80, -450.4798999999994);
1065
  sp3->set_position(pos);
1066
  sp3->set_tracker(0);
1067
  sp3->set_station(3);
1068
  sp3->set_type("triplet");
1069

    
1070
  pos.set(-10.0, 90, -750.4801999999991);
1071
  sp4->set_position(pos);
1072
  sp4->set_tracker(0);
1073
  sp4->set_station(4);
1074
  sp4->set_type("triplet");
1075

    
1076
  pos.set(1.0, 20.0, -1100.410099999999);
1077
  sp5->set_position(pos);
1078
  sp5->set_tracker(0);
1079
  sp5->set_station(5);
1080
  sp5->set_type("triplet");
1081

    
1082
  // Set up the spacepoints vector
1083
  std::vector<SciFiSpacePoint*> spnts;
1084
  spnts.push_back(sp5);
1085
  spnts.push_back(sp2);
1086
  spnts.push_back(sp3);
1087
  spnts.push_back(sp1);
1088
  spnts.push_back(sp4);
1089

    
1090
  SciFiEvent evt;
1091
  evt.set_spacepoints(spnts);
1092

    
1093
  pr.process(evt);
1094

    
1095
  std::vector<SciFiStraightPRTrack> strks = evt.straightprtracks();
1096
  std::vector<SciFiHelicalPRTrack> htrks = evt.helicalprtracks();
1097

    
1098
  SciFiHelicalPRTrack htrk = htrks[0];
1099
  std::cerr << " x0 is " << htrk.get_x0() << std::endl;
1100
  std::cerr << " y0 is " << htrk.get_y0() << std::endl;
1101
  std::cerr << " z0 is " << htrk.get_z0() << std::endl;
1102
  std::cerr << " phi0 is " << htrk.get_phi0() << std::endl;
1103
  std::cerr << " psi0 is " << htrk.get_psi0() << std::endl;
1104
  std::cerr << " ds/dz is " << htrk.get_dsdz() << std::endl;
1105
  std::cerr << " R is " << htrk.get_R() << std::endl;
1106
  std::cerr << " line_sz_chi2 is " << htrk.get_line_sz_chisq() << std::endl;
1107
  std::cerr << " circle_x0 is " << htrk.get_circle_x0() << std::endl;
1108
  std::cerr << " circle_y0 is " << htrk.get_circle_y0() << std::endl;
1109
  std::cerr << " circle_chi2 is " << htrk.get_circle_chisq() << std::endl;
1110
  std::cerr << " chi2 is " << htrk.get_chisq() << std::endl;
1111
  std::cerr << " chi2_dof is " << htrk.get_chisq_dof() << std::endl;
1112
  std::cerr << " num_points is " << htrk.get_num_points() << std::endl;
1113

    
1114
  EXPECT_EQ(0, strks.size());
1115
  EXPECT_EQ(0, htrks.size());
1116
}
1117
*/
1118

    
1119

    
1120
} // ~namespace MAUS
(2-2/2)