← Back to branch summary

~chris-rogers/maus/emr_mc_digitization

« back to all changes in this revision

Viewing changes to tests/cpp_unit/Recon/Global/DataStructureHelperTest.cc

  • Committer: Chris Rogers
  • Date: 2014-04-16 11:48:45 UTC
  • mfrom: (707 merge)
  • mto: This revision was merged to the branch mainline in revision 711.
  • Revision ID: chris.rogers@stfc.ac.uk-20140416114845-h3u3q7pdcxkxvovs
Update to trunk

Show diffs side-by-side

added added

removed removed

Lines of Context:
tests/cpp_unit/Recon/Global/DataStructureHelperTest.cc
 
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
/* Author: Peter Lane
 
18
 */
 
19
 
 
20
#include <cmath>
 
21
#include <iostream>
 
22
#include <sstream>
 
23
#include <string>
 
24
 
 
25
#include <cstdlib>
 
26
 
 
27
#include "CLHEP/Vector/ThreeVector.h"
 
28
#include "gtest/gtest.h"
 
29
#include "json/reader.h"
 
30
#include "json/value.h"
 
31
 
 
32
#include "Config/MiceModule.hh"
 
33
#include "DataStructure/Global/ReconEnums.hh"
 
34
#include "DataStructure/Global/TrackPoint.hh"
 
35
#include "Globals/GlobalsManager.hh"
 
36
#include "Optics/PhaseSpaceVector.hh"
 
37
#include "Recon/Global/DataStructureHelper.hh"
 
38
#include "Recon/Global/Detector.hh"
 
39
#include "Utils/Exception.hh"
 
40
#include "Utils/Globals.hh"
 
41
 
 
42
using MAUS::recon::global::DataStructureHelper;
 
43
namespace Recon = MAUS::recon::global;
 
44
namespace GlobalDS = MAUS::DataStructure::Global;
 
45
class DataStructureHelperTest : public testing::Test {
 
46
 public:
 
47
  DataStructureHelperTest() {
 
48
  }
 
49
 
 
50
  ~DataStructureHelperTest() {
 
51
  }
 
52
 
 
53
 protected:
 
54
};
 
55
 
 
56
// *************************************************
 
57
// DataStructureHelperTest static const initializations
 
58
// *************************************************
 
59
 
 
60
// ***********
 
61
// test cases
 
62
// ***********
 
63
 
 
64
TEST_F(DataStructureHelperTest, FindModulesByName) {
 
65
  MiceModule root(NULL, "root");
 
66
  MiceModule * orwell = new MiceModule(&root, "Orwell");
 
67
  MiceModule * george_orwell = new MiceModule(orwell, "George");
 
68
  MiceModule * curie = new MiceModule(&root, "Curie");
 
69
  MiceModule * marie_curie = new MiceModule(curie, "Marie");
 
70
  MiceModule * washington = new MiceModule(&root, "Washington");
 
71
  MiceModule * george_washington = new MiceModule(washington, "George");
 
72
  MiceModule * nother = new MiceModule(&root, "Nother");
 
73
  MiceModule * emmy_nother = new MiceModule(curie, "Nother");
 
74
  MiceModule * scott = new MiceModule(&root, "Scott");
 
75
  MiceModule * c_scott = new MiceModule(scott, "C");
 
76
  MiceModule * george_c_scott = new MiceModule(c_scott, "George");
 
77
  root.addDaughter(orwell);
 
78
  orwell->addDaughter(george_orwell);
 
79
  root.addDaughter(curie);
 
80
  curie->addDaughter(marie_curie);
 
81
  root.addDaughter(washington);
 
82
  washington->addDaughter(george_washington);
 
83
  root.addDaughter(nother);
 
84
  nother->addDaughter(emmy_nother);
 
85
  root.addDaughter(scott);
 
86
  scott->addDaughter(c_scott);
 
87
  c_scott->addDaughter(george_c_scott);
 
88
 
 
89
  const std::vector<const MiceModule *> georges
 
90
    = DataStructureHelper::GetInstance().FindModulesByName(&root, "George");
 
91
  const size_t number_of_georges = georges.size();
 
92
  ASSERT_EQ(number_of_georges, static_cast<size_t>(3));
 
93
 
 
94
  const std::string names[3] = {
 
95
    "root/Orwell/George",
 
96
    "root/Washington/George",
 
97
    "root/Scott/C/George"
 
98
  };
 
99
  for (size_t index = 0; index < 3; ++index) {
 
100
    std::string actual_name = georges[index]->fullName();
 
101
    std::string expected_name = names[index];
 
102
    EXPECT_EQ(expected_name, actual_name);
 
103
  }
 
104
}
 
105
 
 
106
TEST_F(DataStructureHelperTest, GetDetectorZPosition) {
 
107
  const GlobalDS::DetectorPoint detector_ids[14] = {
 
108
    GlobalDS::kTOF0, GlobalDS::kTOF1,
 
109
    GlobalDS::kTracker0_1, GlobalDS::kTracker0_2, GlobalDS::kTracker0_3,
 
110
    GlobalDS::kTracker0_4, GlobalDS::kTracker0_5,
 
111
    GlobalDS::kTracker1_1, GlobalDS::kTracker1_2, GlobalDS::kTracker1_3,
 
112
    GlobalDS::kTracker1_4, GlobalDS::kTracker1_5,
 
113
    GlobalDS::kTOF2,
 
114
    GlobalDS::kEMR  // used for testing unsupported detector failure
 
115
  };
 
116
 
 
117
  MAUS::Globals * old_globals = MAUS::Globals::GetInstance();
 
118
  // Deep copy since GlobalsManager deletes Globals' old modules
 
119
  // on call to SetReconstructionMiceModules()
 
120
  MiceModule * old_modules
 
121
    = MiceModule::deepCopy(*old_globals->GetReconstructionMiceModules());
 
122
  MiceModule * recon_modules = new MiceModule(NULL, "root");
 
123
 
 
124
  // *** Test missing module exception ***
 
125
  MAUS::GlobalsManager::SetReconstructionMiceModules(recon_modules);
 
126
  for (int index = 0; index < 14; ++index) {
 
127
    bool success = false;
 
128
    try {
 
129
      DataStructureHelper::GetInstance().GetDetectorZPosition(
 
130
        detector_ids[index]);
 
131
    } catch (MAUS::Exception& exc) {
 
132
      success = true;
 
133
    }
 
134
    EXPECT_TRUE(success);
 
135
  }
 
136
 
 
137
  // *** Test duplicate module exception ***
 
138
  recon_modules = new MiceModule(NULL, "root");
 
139
  for (int index = 0; index < 3; ++index) {
 
140
    std::stringstream detector_name;
 
141
    detector_name << "TOF" << index << ".dat";
 
142
    recon_modules->addDaughter(
 
143
      new MiceModule(recon_modules, detector_name.str()));
 
144
    recon_modules->addDaughter(
 
145
      new MiceModule(recon_modules, detector_name.str()));
 
146
  }
 
147
  for (int index = 0; index < 2; ++index) {
 
148
    for (int station = 0; station < 5; ++station) {
 
149
      std::stringstream detector_name;
 
150
      detector_name << "Tracker";
 
151
      if (index > 0) detector_name << index;
 
152
      detector_name << "Station" << (station+1) << ".dat";
 
153
      recon_modules->addDaughter(
 
154
        new MiceModule(recon_modules, detector_name.str()));
 
155
      recon_modules->addDaughter(
 
156
        new MiceModule(recon_modules, detector_name.str()));
 
157
    }
 
158
  }
 
159
  MAUS::GlobalsManager::SetReconstructionMiceModules(recon_modules);
 
160
  for (int index = 0; index < 13; ++index) {
 
161
    bool success = false;
 
162
    try {
 
163
      DataStructureHelper::GetInstance().GetDetectorZPosition(
 
164
        detector_ids[index]);
 
165
    } catch (MAUS::Exception& exc) {
 
166
      success = true;
 
167
    }
 
168
    EXPECT_TRUE(success);
 
169
  }
 
170
 
 
171
  // create a new root module since the next SetReconstructionMiceModules() call
 
172
  // will delete the original used during the above missing module tests
 
173
  recon_modules = new MiceModule(NULL, "root");
 
174
 
 
175
  // Setup a minimal MiceModule tree manually
 
176
  MiceModule * tof0 = new MiceModule(recon_modules, "TOF0.dat");
 
177
  tof0->addPropertyHep3Vector("Position", Hep3Vector(0., 0., -2000.));
 
178
  recon_modules->addDaughter(tof0);
 
179
 
 
180
  MiceModule * tof1 = new MiceModule(recon_modules, "TOF1.dat");
 
181
  tof1->addPropertyHep3Vector("Position", Hep3Vector(0., 0., -900.));
 
182
  recon_modules->addDaughter(tof1);
 
183
  MiceModule * tof1_detector = new MiceModule(tof1, "TOF1Detector.dat");
 
184
  tof1_detector->addPropertyHep3Vector("Position", Hep3Vector(0., 0., -100.));
 
185
  tof1->addDaughter(tof1_detector);
 
186
 
 
187
  MiceModule * tracker0 = new MiceModule(recon_modules, "Tracker0.dat");
 
188
  tracker0->addPropertyHep3Vector("Position", Hep3Vector(0., 0., 1000.));
 
189
  recon_modules->addDaughter(tracker0);
 
190
  MiceModule * tracker0_stations[5];
 
191
  for (int index = 0; index < 5; ++index) {
 
192
    std::stringstream detector_name;
 
193
    detector_name << "TrackerStation" << (index+1) << ".dat";
 
194
    tracker0_stations[index] = new MiceModule(tracker0,
 
195
                                              detector_name.str());
 
196
    tracker0_stations[index]->addPropertyHep3Vector("Position",
 
197
      Hep3Vector(0., 0., 10*(index-2)));
 
198
    tracker0->addDaughter(tracker0_stations[index]);
 
199
 
 
200
    MiceModule * tracker_view_w = NULL;
 
201
    if (index == 4) {  // Station 5
 
202
      tracker_view_w = new MiceModule(tracker0_stations[index],
 
203
                                      "Tracker0Station5ViewW.dat");
 
204
    } else {
 
205
      tracker_view_w = new MiceModule(tracker0_stations[index],
 
206
                                      "TrackerViewW.dat");
 
207
    }
 
208
    tracker_view_w->addPropertyHep3Vector("Position", Hep3Vector(0., 0., 0.));
 
209
    tracker0_stations[index]->addDaughter(tracker_view_w);
 
210
  }
 
211
 
 
212
  MiceModule * tracker1 = new MiceModule(recon_modules, "Tracker1.dat");
 
213
  tracker1->addPropertyHep3Vector("Position", Hep3Vector(0., 0., 2000.));
 
214
  recon_modules->addDaughter(tracker1);
 
215
  MiceModule * tracker1_stations[5];
 
216
  for (int index = 0; index < 5; ++index) {
 
217
    std::stringstream detector_name;
 
218
    detector_name << "Tracker1Station" << (index+1) << ".dat";
 
219
    tracker1_stations[index] = new MiceModule(tracker1,
 
220
                                              detector_name.str());
 
221
    tracker1_stations[index]->addPropertyHep3Vector("Position",
 
222
      Hep3Vector(0., 0., 10*(index-2)));
 
223
    tracker1->addDaughter(tracker1_stations[index]);
 
224
 
 
225
    MiceModule * tracker_view_w = new MiceModule(tracker1_stations[index],
 
226
                                                 "TrackerViewW.dat");
 
227
    tracker_view_w->addPropertyHep3Vector("Position", Hep3Vector(0., 0., 0.));
 
228
    tracker1_stations[index]->addDaughter(tracker_view_w);
 
229
  }
 
230
 
 
231
  MiceModule * tof2 = new MiceModule(recon_modules, "TOF2.dat");
 
232
  tof2->addPropertyHep3Vector("Position", Hep3Vector(0., 0., 3100.));
 
233
  recon_modules->addDaughter(tof2);
 
234
  MiceModule * tof2_detector = new MiceModule(tof2, "TOF2Detector.dat");
 
235
  tof2_detector->addPropertyHep3Vector("Position", Hep3Vector(0., 0., -100.));
 
236
  tof2->addDaughter(tof2_detector);
 
237
 
 
238
  MAUS::GlobalsManager::SetReconstructionMiceModules(recon_modules);
 
239
 
 
240
  // *** Test getting z positions for all TOF and Tracker stations
 
241
  const double z_positions[13] = {
 
242
    -2000, -1000,                 // TOF0, TOF1
 
243
    980, 990, 1000, 1010, 1020,   // Tracker0, Station1-5
 
244
    1980, 1990, 2000, 2010, 2020, // Tracker1, Station1-5
 
245
    3000                          // TOF2
 
246
  };
 
247
  for (int index = 0; index < 13; ++index) {
 
248
    const double detector_z = DataStructureHelper::GetInstance()
 
249
                            .GetDetectorZPosition(detector_ids[index]);
 
250
    EXPECT_EQ(z_positions[index], detector_z);
 
251
  }
 
252
 
 
253
  // Restore the original modules
 
254
  MAUS::GlobalsManager::SetReconstructionMiceModules(old_modules);
 
255
}
 
256
 
 
257
TEST_F(DataStructureHelperTest, GetDetectorAttributes) {
 
258
  std::stringstream configuration_string;
 
259
  configuration_string << "{"
 
260
    << "\"global_recon_detector_attributes\":["
 
261
      << "{"
 
262
        << "\"id\":2,"
 
263
        << "\"uncertainties\":["
 
264
          << "[4.0, 0.0, 0.0, 0.0, 0.0, 0.0],"
 
265
          << "[0.0, 1.44e7, 0.0, 0.0, 0.0, 0.0],"
 
266
          << "[0.0, 0.0, 300.0, 0.0, 0.0, 0.0],"
 
267
          << "[0.0, 0.0, 0.0, 1.0e7, 0.0, 0.0],"
 
268
          << "[0.0, 0.0, 0.0, 0.0, 300.0, 0.0],"
 
269
          << "[0.0, 0.0, 0.0, 0.0, 0.0, 1.0e7]"
 
270
        << "]"
 
271
      << "},"
 
272
      << "{"
 
273
        << "\"id\":6,"
 
274
        << "\"uncertainties\":["
 
275
          << "[4.0, 0.0, 0.0, 0.0, 0.0, 0.0],"
 
276
          << "[0.0, 1.44e7, 0.0, 0.0, 0.0, 0.0],"
 
277
          << "[0.0, 0.0, 300.0, 0.0, 0.0, 0.0],"
 
278
          << "[0.0, 0.0, 0.0, 1.0e7, 0.0, 0.0],"
 
279
          << "[0.0, 0.0, 0.0, 0.0, 300.0, 0.0],"
 
280
          << "[0.0, 0.0, 0.0, 0.0, 0.0, 1.0e7]"
 
281
        << "]"
 
282
      << "}]}";
 
283
  MAUS::recon::global::DetectorMap detectors;
 
284
  bool success = true;
 
285
  try {
 
286
    const Json::Value configuration
 
287
      = JsonWrapper::StringToJson(configuration_string.str());
 
288
 
 
289
    DataStructureHelper::GetInstance().GetDetectorAttributes(
 
290
        configuration, detectors);
 
291
  } catch (MAUS::Exception& exception) {
 
292
    success = false;
 
293
  } catch (std::exception& exc) {
 
294
    success = false;
 
295
  }
 
296
  ASSERT_TRUE(success);
 
297
 
 
298
  Recon::DetectorMap::iterator detector = detectors.begin();
 
299
  const GlobalDS::DetectorPoint ids[2] = {GlobalDS::kTOF0, GlobalDS::kTOF1};
 
300
  const double uncertainty_data[36] = {
 
301
    4.0, 0.0, 0.0, 0.0, 0.0, 0.0,
 
302
    0.0, 1.44e7, 0.0, 0.0, 0.0, 0.0,
 
303
    0.0, 0.0, 300.0, 0.0, 0.0, 0.0,
 
304
    0.0, 0.0, 0.0, 1.0e7, 0.0, 0.0,
 
305
    0.0, 0.0, 0.0, 0.0, 300.0, 0.0,
 
306
    0.0, 0.0, 0.0, 0.0, 0.0, 1.0e7
 
307
  };
 
308
  size_t index = 0;
 
309
  while (detector != detectors.end()) {
 
310
    const GlobalDS::DetectorPoint id = detector->first;
 
311
    EXPECT_EQ(ids[index], id);
 
312
    const MAUS::CovarianceMatrix uncertainties
 
313
      = detector->second.uncertainties();
 
314
    for (size_t row = 0; row < 6; ++row) {
 
315
      for (size_t column = 0; column < 6; ++column) {
 
316
        EXPECT_EQ(uncertainty_data[row*6+column],
 
317
                  uncertainties(row+1, column+1));
 
318
      }
 
319
    }
 
320
    ++index;
 
321
    ++detector;
 
322
  }
 
323
 
 
324
  // *** Test wrong size uncertainty matrices (GetJsonCovariance failure) ***
 
325
  std::stringstream bad_rows_string;
 
326
  bad_rows_string << "{"
 
327
    << "\"global_recon_detector_attributes\":["
 
328
      << "{"
 
329
        << "\"id\":2,"
 
330
        << "\"uncertainties\":["
 
331
          << "[4.0, 0.0, 0.0, 0.0, 0.0, 0.0],"
 
332
          << "[0.0, 1.44e7, 0.0, 0.0, 0.0, 0.0],"
 
333
          << "[0.0, 0.0, 300.0, 0.0, 0.0, 0.0],"
 
334
          << "[0.0, 0.0, 0.0, 0.0, 300.0, 0.0],"
 
335
          << "[0.0, 0.0, 0.0, 0.0, 0.0, 1.0e7]"
 
336
        << "]"
 
337
      << "}]}";
 
338
  success = false;
 
339
  try {
 
340
    const Json::Value configuration
 
341
      = JsonWrapper::StringToJson(bad_rows_string.str());
 
342
 
 
343
    DataStructureHelper::GetInstance().GetDetectorAttributes(
 
344
        configuration, detectors);
 
345
  } catch (MAUS::Exception& exception) {
 
346
    success = true;
 
347
  } catch (std::exception& exc) {
 
348
    success = true;
 
349
  }
 
350
  ASSERT_TRUE(success);
 
351
 
 
352
  std::stringstream bad_columns_string;
 
353
  bad_columns_string << "{"
 
354
    << "\"global_recon_detector_attributes\":["
 
355
      << "{"
 
356
        << "\"id\":2,"
 
357
        << "\"uncertainties\":["
 
358
          << "[4.0, 0.0, 0.0, 0.0, 0.0, 0.0],"
 
359
          << "[0.0, 1.44e7, 0.0, 0.0, 0.0, 0.0],"
 
360
          << "[0.0, 0.0, 300.0, 0.0, 0.0, 0.0],"
 
361
          << "[0.0, 0.0, 0.0, 1.0e7, 0.0],"
 
362
          << "[0.0, 0.0, 0.0, 0.0, 300.0, 0.0],"
 
363
          << "[0.0, 0.0, 0.0, 0.0, 0.0, 1.0e7]"
 
364
        << "]"
 
365
      << "}]}";
 
366
  success = false;
 
367
  try {
 
368
    const Json::Value configuration
 
369
      = JsonWrapper::StringToJson(bad_columns_string.str());
 
370
 
 
371
    DataStructureHelper::GetInstance().GetDetectorAttributes(
 
372
        configuration, detectors);
 
373
  } catch (MAUS::Exception& exception) {
 
374
    success = true;
 
375
  } catch (std::exception& exc) {
 
376
    success = true;
 
377
  }
 
378
  ASSERT_TRUE(success);
 
379
}
 
380
 
 
381
TEST_F(DataStructureHelperTest, TrackPoint2PhaseSpaceVector) {
 
382
  // t, E, x, Px, y, Py
 
383
  const double coordinates[6] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0};
 
384
  GlobalDS::TrackPoint track_point;
 
385
  TLorentzVector position(coordinates[2], coordinates[4], 0., coordinates[0]);
 
386
  track_point.set_position(position);
 
387
  TLorentzVector momentum(coordinates[3], coordinates[5], 0., coordinates[1]);
 
388
  track_point.set_momentum(momentum);
 
389
  MAUS::PhaseSpaceVector phase_space_vector = DataStructureHelper::GetInstance()
 
390
    .TrackPoint2PhaseSpaceVector(track_point);
 
391
  for (int index = 0; index < 6; ++index) {
 
392
    std::cout << "Coordinate " << index << "..." << std::endl;
 
393
    // Order is t, E, x, Px, y, Py
 
394
    EXPECT_EQ(coordinates[index], phase_space_vector[index]);
 
395
  }
 
396
}
 
397
 
 
398
TEST_F(DataStructureHelperTest, PhaseSpaceVector2TrackPoint) {
 
399
  // *** Test on-mass-shell conversion ***
 
400
  const double positions[4] = {1.0, 3.0, 5.0, 7.0};       // x, y, z, t
 
401
  const double momenta[4] = {2.0, 4.0, 792.979388, 800.0}; // Px, Py, Pz, E
 
402
  MAUS::PhaseSpaceVector phase_space_vector(positions[3], momenta[3],
 
403
                                            positions[0], momenta[0],
 
404
                                            positions[1], momenta[1]);
 
405
  GlobalDS::TrackPoint track_point = DataStructureHelper::GetInstance()
 
406
    .PhaseSpaceVector2TrackPoint(phase_space_vector, positions[2],
 
407
                                 GlobalDS::kMuMinus);
 
408
  TLorentzVector position = track_point.get_position();
 
409
  TLorentzVector momentum = track_point.get_momentum();
 
410
  for (int index = 0; index < 4; ++index) {
 
411
    // Order is x, y, z, t
 
412
    EXPECT_NEAR(positions[index], position[index], 1.e-6);
 
413
    EXPECT_NEAR(momenta[index], momentum[index], 1.e-6);
 
414
  }
 
415
 
 
416
  // *** Test off-mass-shell conversion (Pz set to 0.) ***
 
417
  const double bad_momenta[4] = {2.0, 4.0, 0.0, 8.0}; // Px, Py, Pz, E
 
418
  phase_space_vector = MAUS::PhaseSpaceVector(positions[3], bad_momenta[3],
 
419
                                              positions[0], bad_momenta[0],
 
420
                                              positions[1], bad_momenta[1]);
 
421
  track_point = DataStructureHelper::GetInstance().PhaseSpaceVector2TrackPoint(
 
422
    phase_space_vector, positions[2], GlobalDS::kMuMinus);
 
423
  position = track_point.get_position();
 
424
  momentum = track_point.get_momentum();
 
425
  for (int index = 0; index < 4; ++index) {
 
426
    // Order is x, y, z, t
 
427
    EXPECT_NEAR(positions[index], position[index], 1.e-6);
 
428
    EXPECT_NEAR(bad_momenta[index], momentum[index], 1.e-6);
 
429
  }
 
430
}