ObjectTest.h 40 KB
Newer Older
Russell Taylor's avatar
Russell Taylor committed
1
2
3
4
5
6
7
8
9
10
#ifndef MANTID_TESTOBJECT__
#define MANTID_TESTOBJECT__

#include <cxxtest/TestSuite.h>
#include <cmath>
#include <ostream>
#include <vector>
#include <algorithm>
#include <ctime>

11
12
#include <boost/shared_ptr.hpp>

13
14
15
16
17
18
19
#include "MantidGeometry/Objects/Object.h"
#include "MantidGeometry/Surfaces/Cylinder.h"
#include "MantidGeometry/Surfaces/Sphere.h"
#include "MantidGeometry/Surfaces/Plane.h"
#include "MantidGeometry/Math/Algebra.h"
#include "MantidGeometry/Surfaces/SurfaceFactory.h"
#include "MantidGeometry/Objects/Track.h"
Nick Draper's avatar
re #843    
Nick Draper committed
20
#include "MantidGeometry/Rendering/GluGeometryHandler.h"
21
22
#include "MantidGeometry/Objects/ShapeFactory.h"

23
24
#include "MantidKernel/Material.h"

25
#include "MantidTestHelpers/ComponentCreationHelper.h"
26

Russell Taylor's avatar
Russell Taylor committed
27
28
using namespace Mantid;
using namespace Geometry;
29
using Mantid::Kernel::V3D;
Russell Taylor's avatar
Russell Taylor committed
30

31
class ObjectTest : public CxxTest::TestSuite {
Russell Taylor's avatar
Russell Taylor committed
32
33

public:
34
  void testDefaultObjectHasEmptyMaterial() {
35
36
    Object obj;

37
38
    TSM_ASSERT_DELTA("Expected a zero number density", 0.0,
                     obj.material().numberDensity(), 1e-12);
39
40
  }

41
  void testObjectSetMaterialReplacesExisting() {
42
43
44
45
46
    using Mantid::Kernel::Material;
    Object obj;

    TSM_ASSERT_DELTA("Expected a zero number density", 0.0,
                     obj.material().numberDensity(), 1e-12);
47
48
    obj.setMaterial(
        Material("arm", PhysicalConstants::getNeutronAtom(13), 45.0));
49
50
51
52
    TSM_ASSERT_DELTA("Expected a number density of 45", 45.0,
                     obj.material().numberDensity(), 1e-12);
  }

53
54
55
  void testCopyConstructorGivesObjectWithSameAttributes() {
    Object_sptr original =
        ComponentCreationHelper::createSphere(1.0, V3D(), "sphere");
56
57
58
59
60
    int objType(-1);
    double radius(-1.0), height(-1.0);
    std::vector<V3D> pts;
    original->GetObjectGeom(objType, pts, radius, height);
    TS_ASSERT_EQUALS(2, objType);
61
62
    TS_ASSERT(boost::dynamic_pointer_cast<GluGeometryHandler>(
        original->getGeometryHandler()));
63
64
65
66
67
68
69

    Object copy(*original);
    // The copy should be a primitive object with a GluGeometryHandler
    objType = -1;
    copy.GetObjectGeom(objType, pts, radius, height);

    TS_ASSERT_EQUALS(2, objType);
70
71
    TS_ASSERT(boost::dynamic_pointer_cast<GluGeometryHandler>(
        copy.getGeometryHandler()));
72
73
74
75
76
77
    TS_ASSERT_EQUALS(copy.getName(), original->getName());
    // Check the string representation is the same
    TS_ASSERT_EQUALS(copy.str(), original->str());
    TS_ASSERT_EQUALS(copy.getSurfaceIndex(), original->getSurfaceIndex());
  }

78
79
80
  void testAssignmentOperatorGivesObjectWithSameAttributes() {
    Object_sptr original =
        ComponentCreationHelper::createSphere(1.0, V3D(), "sphere");
81
82
83
84
85
    int objType(-1);
    double radius(-1.0), height(-1.0);
    std::vector<V3D> pts;
    original->GetObjectGeom(objType, pts, radius, height);
    TS_ASSERT_EQUALS(2, objType);
86
87
    TS_ASSERT(boost::dynamic_pointer_cast<GluGeometryHandler>(
        original->getGeometryHandler()));
88

89
    Object lhs;      // initialize
90
91
92
93
94
95
    lhs = *original; // assign
    // The copy should be a primitive object with a GluGeometryHandler
    objType = -1;
    lhs.GetObjectGeom(objType, pts, radius, height);

    TS_ASSERT_EQUALS(2, objType);
96
97
    TS_ASSERT(boost::dynamic_pointer_cast<GluGeometryHandler>(
        lhs.getGeometryHandler()));
98
  }
Russell Taylor's avatar
Russell Taylor committed
99

100
  void testCreateUnitCube() {
101
    Object_sptr geom_obj = createUnitCube();
102

103
    TS_ASSERT_EQUALS(geom_obj->str(), "68 -6 5 -4 3 -2 1");
104
105

    double xmin(0.0), xmax(0.0), ymin(0.0), ymax(0.0), zmin(0.0), zmax(0.0);
106
    geom_obj->getBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin);
Russell Taylor's avatar
Russell Taylor committed
107
108
  }

109
  void testIsOnSideCappedCylinder() {
110
    Object_sptr geom_obj = createCappedCylinder();
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
    // inside
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 0)), false); // origin
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 2.9, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -2.9, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -2.9)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 2.9)), false);
    // on the side
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 3)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.2, 0, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.2, 0, 0)), true);

    // on the edges
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.2, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.2, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.2, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.2, 0, 3)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.2, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.2, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.2, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.2, 0, 3)), true);
    // out side
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 3.1, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -3.1, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -3.1)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 3.1)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(1.3, 0, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(-3.3, 0, 0)), false);
Russell Taylor's avatar
Russell Taylor committed
141
142
  }

143
  void testIsValidCappedCylinder() {
144
    Object_sptr geom_obj = createCappedCylinder();
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
    // inside
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 0)), true); // origin
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 2.9, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -2.9, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -2.9)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 2.9)), true);
    // on the side
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 3)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.2, 0, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.2, 0, 0)), true);

    // on the edges
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.2, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.2, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.2, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.2, 0, 3)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.2, 3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.2, -3, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.2, 0, -3)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.2, 0, 3)), true);
    // out side
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 3.1, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -3.1, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -3.1)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 3.1)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(1.3, 0, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(-3.3, 0, 0)), false);
Russell Taylor's avatar
Russell Taylor committed
175
176
  }

177
  void testIsOnSideSphere() {
178
    Object_sptr geom_obj = createSphere();
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
    // inside
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 0)), false); // origin
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 4.0, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -4.0, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -4.0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 4.0)), false);
    // on the side
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 4.1, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -4.1, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -4.1)), true);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 4.1)), true);

    // out side
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 4.2, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, -4.2, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, -4.2)), false);
    TS_ASSERT_EQUALS(geom_obj->isOnSide(V3D(0, 0, 4.2)), false);
Russell Taylor's avatar
Russell Taylor committed
196
197
  }

198
  void testIsValidSphere() {
199
    Object_sptr geom_obj = createSphere();
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
    // inside
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 0)), true); // origin
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 4.0, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -4.0, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -4.0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 4.0)), true);
    // on the side
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 4.1, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -4.1, 0)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -4.1)), true);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 4.1)), true);

    // out side
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 4.2, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, -4.2, 0)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, -4.2)), false);
    TS_ASSERT_EQUALS(geom_obj->isValid(V3D(0, 0, 4.2)), false);
Russell Taylor's avatar
Russell Taylor committed
217
218
  }

219
  void testCalcValidTypeSphere() {
220
    Object_sptr geom_obj = createSphere();
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
    // entry on the normal
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-4.1, 0, 0), V3D(1, 0, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-4.1, 0, 0), V3D(-1, 0, 0)),
                     -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(4.1, 0, 0), V3D(1, 0, 0)), -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(4.1, 0, 0), V3D(-1, 0, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, -4.1, 0), V3D(0, 1, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, -4.1, 0), V3D(0, -1, 0)),
                     -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, 4.1, 0), V3D(0, 1, 0)), -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, 4.1, 0), V3D(0, -1, 0)), 1);

    // a glancing blow
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-4.1, 0, 0), V3D(0, 1, 0)), 0);
    // not quite on the normal
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-4.1, 0, 0), V3D(0.5, 0.5, 0)),
                     1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(4.1, 0, 0), V3D(0.5, 0.5, 0)),
                     -1);
Russell Taylor's avatar
Russell Taylor committed
240
241
  }

242
  void testGetBoundingBoxForSphere() {
243
    Object_sptr geom_obj = createSphere();
244
245
    const double tolerance(1e-10);

246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
    double xmax, ymax, zmax, xmin, ymin, zmin;
    xmax = ymax = zmax = 20;
    xmin = ymin = zmin = -20;
    geom_obj->getBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin);
    TS_ASSERT_DELTA(xmax, 4.1, tolerance);
    TS_ASSERT_DELTA(ymax, 4.1, tolerance);
    TS_ASSERT_DELTA(zmax, 4.1, tolerance);
    TS_ASSERT_DELTA(xmin, -4.1, tolerance);
    TS_ASSERT_DELTA(ymin, -4.1, tolerance);
    TS_ASSERT_DELTA(zmin, -4.1, tolerance);

    const BoundingBox &bbox = geom_obj->getBoundingBox();

    TS_ASSERT_DELTA(bbox.xMax(), 4.1, tolerance);
    TS_ASSERT_DELTA(bbox.yMax(), 4.1, tolerance);
    TS_ASSERT_DELTA(bbox.zMax(), 4.1, tolerance);
    TS_ASSERT_DELTA(bbox.xMin(), -4.1, tolerance);
    TS_ASSERT_DELTA(bbox.yMin(), -4.1, tolerance);
    TS_ASSERT_DELTA(bbox.zMin(), -4.1, tolerance);
265
266
  }

267
  void testCalcValidTypeCappedCylinder() {
268
    Object_sptr geom_obj = createCappedCylinder();
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
    // entry on the normal
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-3.2, 0, 0), V3D(1, 0, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-3.2, 0, 0), V3D(-1, 0, 0)),
                     -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(1.2, 0, 0), V3D(1, 0, 0)), -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(1.2, 0, 0), V3D(-1, 0, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, -3, 0), V3D(0, 1, 0)), 1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, -3, 0), V3D(0, -1, 0)), -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, 3, 0), V3D(0, 1, 0)), -1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(0, 3, 0), V3D(0, -1, 0)), 1);

    // a glancing blow
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-3.2, 0, 0), V3D(0, 1, 0)), 0);
    // not quite on the normal
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-3.2, 0, 0), V3D(0.5, 0.5, 0)),
                     1);
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(1.2, 0, 0), V3D(0.5, 0.5, 0)),
                     -1);
Russell Taylor's avatar
Russell Taylor committed
287
288
  }

289
  void testInterceptSurfaceSphereZ() {
290
    std::vector<Link> expectedResults;
291
    std::string S41 = "s 1 1 1 4"; // Sphere at (1,1,1) radius 4
Russell Taylor's avatar
Russell Taylor committed
292
293

    // First create some surfaces
294
295
    std::map<int, std::unique_ptr<Surface>> SphSurMap;
    SphSurMap[41] = std::unique_ptr<Sphere>(new Sphere());
Russell Taylor's avatar
Russell Taylor committed
296
297
298
    SphSurMap[41]->setSurface(S41);
    SphSurMap[41]->setName(41);

299
    // A sphere
300
    std::string ObjSphere = "-41";
Russell Taylor's avatar
Russell Taylor committed
301

302
    Object_sptr geom_obj = Object_sptr(new Object);
303
    geom_obj->setObject(41, ObjSphere);
304
    geom_obj->populate(SphSurMap);
Russell Taylor's avatar
Russell Taylor committed
305

306
    Track track(V3D(-1, 1.5, 1), V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
307

308
    // format = startPoint, endPoint, total distance so far
Russell Taylor's avatar
Russell Taylor committed
309
    // forward only intercepts means that start point should be track origin
310
311
312
    expectedResults.push_back(Link(V3D(-1, 1.5, 1),
                                   V3D(sqrt(16 - 0.25) + 1, 1.5, 1.0),
                                   sqrt(15.75) + 2, *geom_obj));
Russell Taylor's avatar
Russell Taylor committed
313

314
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
315
316
  }

317
  void testInterceptSurfaceSphereY() {
318
    std::vector<Link> expectedResults;
319
    Object_sptr geom_obj = createSphere();
320
    Track track(V3D(0, -10, 0), V3D(0, 1, 0));
Russell Taylor's avatar
Russell Taylor committed
321

322
323
324
    // format = startPoint, endPoint, total distance so far
    expectedResults.push_back(
        Link(V3D(0, -4.1, 0), V3D(0, 4.1, 0), 14.1, *geom_obj));
Russell Taylor's avatar
Russell Taylor committed
325

326
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
327
328
  }

329
  void testInterceptSurfaceSphereX() {
330
    std::vector<Link> expectedResults;
331
    Object_sptr geom_obj = createSphere();
332
    Track track(V3D(-10, 0, 0), V3D(1, 0, 0));
333

334
335
336
337
    // format = startPoint, endPoint, total distance so far
    expectedResults.push_back(
        Link(V3D(-4.1, 0, 0), V3D(4.1, 0, 0), 14.1, *geom_obj));
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
338
339
  }

340
  void testInterceptSurfaceCappedCylinderY() {
341
    std::vector<Link> expectedResults;
342
    Object_sptr geom_obj = createCappedCylinder();
343
344
    // format = startPoint, endPoint, total distance so far
    expectedResults.push_back(Link(V3D(0, -3, 0), V3D(0, 3, 0), 13, *geom_obj));
Russell Taylor's avatar
Russell Taylor committed
345

346
347
    Track track(V3D(0, -10, 0), V3D(0, 1, 0));
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
348
349
  }

350
  void testInterceptSurfaceCappedCylinderX() {
351
    std::vector<Link> expectedResults;
352
    Object_sptr geom_obj = createCappedCylinder();
353
    Track track(V3D(-10, 0, 0), V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
354

355
356
357
358
    // format = startPoint, endPoint, total distance so far
    expectedResults.push_back(
        Link(V3D(-3.2, 0, 0), V3D(1.2, 0, 0), 11.2, *geom_obj));
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
359
360
  }

361
362
363
  void testInterceptSurfaceCappedCylinderMiss() {
    std::vector<Link>
        expectedResults; // left empty as there are no expected results
364
    Object_sptr geom_obj = createCappedCylinder();
365
    Track track(V3D(-10, 0, 0), V3D(1, 1, 0));
Russell Taylor's avatar
Russell Taylor committed
366

367
    checkTrackIntercept(geom_obj, track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
368
369
  }

370
371
  void checkTrackIntercept(Track &track,
                           const std::vector<Link> &expectedResults) {
Russell Taylor's avatar
Russell Taylor committed
372
    int index = 0;
373
374
375
376
377
378
379
380
381
    for (Track::LType::const_iterator it = track.begin(); it != track.end();
         ++it) {
      TS_ASSERT_DELTA(it->distFromStart, expectedResults[index].distFromStart,
                      1e-6);
      TS_ASSERT_DELTA(it->distInsideObject,
                      expectedResults[index].distInsideObject, 1e-6);
      TS_ASSERT_EQUALS(it->componentID, expectedResults[index].componentID);
      TS_ASSERT_EQUALS(it->entryPoint, expectedResults[index].entryPoint);
      TS_ASSERT_EQUALS(it->exitPoint, expectedResults[index].exitPoint);
Russell Taylor's avatar
Russell Taylor committed
382
383
      ++index;
    }
384
    TS_ASSERT_EQUALS(index, static_cast<int>(expectedResults.size()));
Russell Taylor's avatar
Russell Taylor committed
385
386
  }

387
388
  void checkTrackIntercept(Object_sptr obj, Track &track,
                           const std::vector<Link> &expectedResults) {
389
    int unitCount = obj->interceptSurface(track);
390
391
    TS_ASSERT_EQUALS(unitCount, expectedResults.size());
    checkTrackIntercept(track, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
392
393
394
  }

  void xtestTrackTwoIsolatedCubes()
395
396
397
  /**
  Test a track going through an object
  */
Russell Taylor's avatar
Russell Taylor committed
398
  {
399
400
    std::string ObjA = "60001 -60002 60003 -60004 60005 -60006";
    std::string ObjB = "80001 -80002 60003 -60004 60005 -60006";
401

402
    createSurfaces(ObjA);
403
404
    Object object1 = Object();
    object1.setObject(3, ObjA);
Russell Taylor's avatar
Russell Taylor committed
405
406
    object1.populate(SMap);

407
    createSurfaces(ObjB);
408
409
    Object object2 = Object();
    object2.setObject(4, ObjB);
Russell Taylor's avatar
Russell Taylor committed
410
411
    object2.populate(SMap);

412
    Track TL(Kernel::V3D(-5, 0, 0), Kernel::V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
413
414

    // CARE: This CANNOT be called twice
415
416
    TS_ASSERT(object1.interceptSurface(TL) != 0);
    TS_ASSERT(object2.interceptSurface(TL) != 0);
Russell Taylor's avatar
Russell Taylor committed
417

418
    std::vector<Link> expectedResults;
419
420
421
422
    expectedResults.push_back(Link(V3D(-1, 0, 0), V3D(1, 0, 0), 6, object1));
    expectedResults.push_back(
        Link(V3D(4.5, 0, 0), V3D(6.5, 0, 0), 11.5, object2));
    checkTrackIntercept(TL, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
423
424
425
  }

  void testTrackTwoTouchingCubes()
426
427
428
  /**
  Test a track going through an object
  */
Russell Taylor's avatar
Russell Taylor committed
429
  {
430
431
    std::string ObjA = "60001 -60002 60003 -60004 60005 -60006";
    std::string ObjB = "60002 -80002 60003 -60004 60005 -60006";
Russell Taylor's avatar
Russell Taylor committed
432

433
    createSurfaces(ObjA);
434
435
    Object object1 = Object();
    object1.setObject(3, ObjA);
Russell Taylor's avatar
Russell Taylor committed
436
437
    object1.populate(SMap);

438
    createSurfaces(ObjB);
439
440
    Object object2 = Object();
    object2.setObject(4, ObjB);
Russell Taylor's avatar
Russell Taylor committed
441
442
    object2.populate(SMap);

443
    Track TL(Kernel::V3D(-5, 0, 0), Kernel::V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
444
445

    // CARE: This CANNOT be called twice
446
447
    TS_ASSERT(object1.interceptSurface(TL) != 0);
    TS_ASSERT(object2.interceptSurface(TL) != 0);
Russell Taylor's avatar
Russell Taylor committed
448

449
    std::vector<Link> expectedResults;
450
451
452
    expectedResults.push_back(Link(V3D(-1, 0, 0), V3D(1, 0, 0), 6, object1));
    expectedResults.push_back(
        Link(V3D(1, 0, 0), V3D(6.5, 0, 0), 11.5, object2));
Russell Taylor's avatar
Russell Taylor committed
453

454
    checkTrackIntercept(TL, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
455
456
457
  }

  void testTrackCubeWithInternalSphere()
458
459
460
  /**
  Test a track going through an object
  */
Russell Taylor's avatar
Russell Taylor committed
461
  {
462
463
    std::string ObjA = "60001 -60002 60003 -60004 60005 -60006 71";
    std::string ObjB = "-71";
Russell Taylor's avatar
Russell Taylor committed
464

465
    createSurfaces(ObjA);
466
467
    Object object1 = Object();
    object1.setObject(3, ObjA);
Russell Taylor's avatar
Russell Taylor committed
468
469
    object1.populate(SMap);

470
    createSurfaces(ObjB);
471
472
    Object object2 = Object();
    object2.setObject(4, ObjB);
Russell Taylor's avatar
Russell Taylor committed
473
474
    object2.populate(SMap);

475
    Track TL(Kernel::V3D(-5, 0, 0), Kernel::V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
476
477

    // CARE: This CANNOT be called twice
478
479
    TS_ASSERT(object1.interceptSurface(TL) != 0);
    TS_ASSERT(object2.interceptSurface(TL) != 0);
Russell Taylor's avatar
Russell Taylor committed
480

481
    std::vector<Link> expectedResults;
482
483
484
485
486
487
    expectedResults.push_back(
        Link(V3D(-1, 0, 0), V3D(-0.8, 0, 0), 4.2, object1));
    expectedResults.push_back(
        Link(V3D(-0.8, 0, 0), V3D(0.8, 0, 0), 5.8, object1));
    expectedResults.push_back(Link(V3D(0.8, 0, 0), V3D(1, 0, 0), 6, object2));
    checkTrackIntercept(TL, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
488
489
490
  }

  void testTrack_CubePlusInternalEdgeTouchSpheres()
491
492
493
  /**
  Test a track going through an object
  */
Russell Taylor's avatar
Russell Taylor committed
494
  {
495
496
    std::string ObjA = "60001 -60002 60003 -60004 60005 -60006 72 73";
    std::string ObjB = "(-72 : -73)";
497

498
    createSurfaces(ObjA);
499
500
    Object object1 = Object();
    object1.setObject(3, ObjA);
Russell Taylor's avatar
Russell Taylor committed
501
502
    object1.populate(SMap);

503
    createSurfaces(ObjB);
504
505
    Object object2 = Object();
    object2.setObject(4, ObjB);
Russell Taylor's avatar
Russell Taylor committed
506
507
    object2.populate(SMap);

508
    Track TL(Kernel::V3D(-5, 0, 0), Kernel::V3D(1, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
509
510

    // CARE: This CANNOT be called twice
511
512
    TS_ASSERT(object1.interceptSurface(TL) != 0);
    TS_ASSERT(object2.interceptSurface(TL) != 0);
Russell Taylor's avatar
Russell Taylor committed
513

514
    std::vector<Link> expectedResults;
515
516
517
518
519
520
    expectedResults.push_back(
        Link(V3D(-1, 0, 0), V3D(-0.4, 0, 0), 4.6, object1));
    expectedResults.push_back(
        Link(V3D(-0.4, 0, 0), V3D(0.2, 0, 0), 5.2, object1));
    expectedResults.push_back(Link(V3D(0.2, 0, 0), V3D(1, 0, 0), 6, object2));
    checkTrackIntercept(TL, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
521
522
523
  }

  void testTrack_CubePlusInternalEdgeTouchSpheresMiss()
524
525
526
  /**
  Test a track missing an object
  */
Russell Taylor's avatar
Russell Taylor committed
527
  {
528
529
    std::string ObjA = "60001 -60002 60003 -60004 60005 -60006 72 73";
    std::string ObjB = "(-72 : -73)";
530

531
    createSurfaces(ObjA);
532
533
    Object object1 = Object();
    object1.setObject(3, ObjA);
Russell Taylor's avatar
Russell Taylor committed
534
535
    object1.populate(SMap);

536
    createSurfaces(ObjB);
537
538
    Object object2 = Object();
    object2.setObject(4, ObjB);
Russell Taylor's avatar
Russell Taylor committed
539
540
    object2.populate(SMap);

541
    Track TL(Kernel::V3D(-5, 0, 0), Kernel::V3D(0, 1, 0));
Russell Taylor's avatar
Russell Taylor committed
542
543

    // CARE: This CANNOT be called twice
544
545
    TS_ASSERT_EQUALS(object1.interceptSurface(TL), 0);
    TS_ASSERT_EQUALS(object2.interceptSurface(TL), 0);
Russell Taylor's avatar
Russell Taylor committed
546

547
548
    std::vector<Link> expectedResults; // left empty as this should miss
    checkTrackIntercept(TL, expectedResults);
Russell Taylor's avatar
Russell Taylor committed
549
550
551
  }

  void testFindPointInCube()
552
553
554
  /**
  Test find point in cube
  */
Russell Taylor's avatar
Russell Taylor committed
555
  {
556
    Object_sptr geom_obj = createUnitCube();
Russell Taylor's avatar
Russell Taylor committed
557
    // initial guess in object
558
    Kernel::V3D pt;
559
560
    TS_ASSERT_EQUALS(geom_obj->getPointInObject(pt), 1);
    TS_ASSERT_EQUALS(pt, V3D(0, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
561
562
    // initial guess not in object, but on x-axis
    std::vector<std::string> planes;
563
564
565
566
567
568
569
570
571
    planes.push_back("px 10");
    planes.push_back("px 11");
    planes.push_back("py -0.5");
    planes.push_back("py 0.5");
    planes.push_back("pz -0.5");
    planes.push_back("pz 0.5");
    Object_sptr B = createCuboid(planes);
    TS_ASSERT_EQUALS(B->getPointInObject(pt), 1);
    TS_ASSERT_EQUALS(pt, V3D(10, 0, 0));
Russell Taylor's avatar
Russell Taylor committed
572
573
    // on y axis
    planes.clear();
574
575
576
577
578
579
580
581
582
    planes.push_back("px -0.5");
    planes.push_back("px 0.5");
    planes.push_back("py -22");
    planes.push_back("py -21");
    planes.push_back("pz -0.5");
    planes.push_back("pz 0.5");
    Object_sptr C = createCuboid(planes);
    TS_ASSERT_EQUALS(C->getPointInObject(pt), 1);
    TS_ASSERT_EQUALS(pt, V3D(0, -21, 0));
Russell Taylor's avatar
Russell Taylor committed
583
584
    // not on principle axis, now works using getBoundingBox
    planes.clear();
585
586
587
588
589
590
591
592
593
594
595
    planes.push_back("px 0.5");
    planes.push_back("px 1.5");
    planes.push_back("py -22");
    planes.push_back("py -21");
    planes.push_back("pz -0.5");
    planes.push_back("pz 0.5");
    Object_sptr D = createCuboid(planes);
    TS_ASSERT_EQUALS(D->getPointInObject(pt), 1);
    TS_ASSERT_DELTA(pt.X(), 1.0, 1e-6);
    TS_ASSERT_DELTA(pt.Y(), -21.5, 1e-6);
    TS_ASSERT_DELTA(pt.Z(), 0.0, 1e-6);
Russell Taylor's avatar
Russell Taylor committed
596
    planes.clear();
597
598
599
600
601
602
603
604
    // Test non axis aligned (AA) case - getPointInObject works because the
    // object is on a principle axis
    // However, if not on a principle axis then the getBoundingBox fails to find
    // correct minima (maxima are OK)
    // This is related to use of the complement for -ve surfaces and might be
    // avoided by only using +ve surfaces
    // for defining non-AA objects. However, BoundingBox is poor for non-AA and
    // needs improvement if these are
Russell Taylor's avatar
Russell Taylor committed
605
    // common
606
607
608
609
610
611
612
613
614
615
616
    planes.push_back("p 1 0 0 -0.5");
    planes.push_back("p 1 0 0 0.5");
    planes.push_back("p 0 .70710678118 .70710678118 -1.1");
    planes.push_back("p 0 .70710678118 .70710678118 -0.1");
    planes.push_back("p 0 -.70710678118 .70710678118 -0.5");
    planes.push_back("p 0 -.70710678118 .70710678118 0.5");
    Object_sptr E = createCuboid(planes);
    TS_ASSERT_EQUALS(E->getPointInObject(pt), 1);
    TS_ASSERT_DELTA(pt.X(), 0.0, 1e-6);
    TS_ASSERT_DELTA(pt.Y(), -0.1414213562373, 1e-6);
    TS_ASSERT_DELTA(pt.Z(), 0.0, 1e-6);
Russell Taylor's avatar
Russell Taylor committed
617
    planes.clear();
618
619
    // This test fails to find a point in object, as object not on a principle
    // axis
Russell Taylor's avatar
Russell Taylor committed
620
    // and getBoundingBox does not give a useful result in this case.
621
622
    // Object is unit cube located at +-0.5 in x but centred on z=y=-1.606.. and
    // rotated 45deg
Russell Taylor's avatar
Russell Taylor committed
623
    // to these two axes
624
625
626
627
628
629
630
631
632
633
634
635
636
637
    planes.push_back("p 1 0 0 -0.5");
    planes.push_back("p 1 0 0 0.5");
    planes.push_back("p 0  .70710678118 .70710678118 -2");
    planes.push_back("p 0  .70710678118 .70710678118 -1");
    planes.push_back("p 0 -.70710678118 .70710678118 -0.5");
    planes.push_back("p 0 -.70710678118 .70710678118 0.5");
    Object_sptr F = createCuboid(planes);
    TS_ASSERT_EQUALS(F->getPointInObject(pt), 0);
    // Test use of defineBoundingBox to explictly set the bounding box, when the
    // automatic method fails
    F->defineBoundingBox(0.5, -1 / (2.0 * sqrt(2.0)), -1.0 / (2.0 * sqrt(2.0)),
                         -0.5, -sqrt(2.0) - 1.0 / (2.0 * sqrt(2.0)),
                         -sqrt(2.0) - 1.0 / (2.0 * sqrt(2.0)));
    TS_ASSERT_EQUALS(F->getPointInObject(pt), 1);
638
    Object_sptr S = createSphere();
639
640
    TS_ASSERT_EQUALS(S->getPointInObject(pt), 1);
    TS_ASSERT_EQUALS(pt, V3D(0.0, 0.0, 0));
Russell Taylor's avatar
Russell Taylor committed
641
642
643
  }

  void testSolidAngleSphere()
644
645
646
  /**
  Test solid angle calculation for a sphere
  */
Russell Taylor's avatar
Russell Taylor committed
647
  {
648
    Object_sptr geom_obj = createSphere();
649
    double satol = 2e-2; // tolerance for solid angle
Russell Taylor's avatar
Russell Taylor committed
650
651
652
653
654

    // Solid angle at distance 8.1 from centre of sphere radius 4.1 x/y/z
    // Expected solid angle calculated values from sa=2pi(1-cos(arcsin(R/r))
    // where R is sphere radius and r is distance of observer from sphere centre
    // Intercept for track in reverse direction now worked round
655
656
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(8.1, 0, 0)), 0.864364,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
657
    // internal point (should be 4pi)
658
659
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(0, 0, 0)), 4 * M_PI,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
660
    // surface point
661
662
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(4.1, 0, 0)), 2 * M_PI,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
663
664
665
  }

  void testSolidAngleCappedCylinder()
666
667
668
  /**
  Test solid angle calculation for a capped cylinder
  */
Russell Taylor's avatar
Russell Taylor committed
669
  {
670
    Object_sptr geom_obj = createSmallCappedCylinder();
671
    // Want to test triangulation so setup a geometry handler
672
673
674
675
    boost::shared_ptr<GluGeometryHandler> h =
        boost::shared_ptr<GluGeometryHandler>(
            new GluGeometryHandler(geom_obj.get()));
    h->setCylinder(V3D(-0.0015, 0.0, 0.0), V3D(1., 0.0, 0.0), 0.005, 0.003);
676
    geom_obj->setGeometryHandler(h);
677

678
    double satol(1e-8); // tolerance for solid angle
Russell Taylor's avatar
Russell Taylor committed
679

680
681
682
683
684
685
    // solid angle at point -0.5 from capped cyl -1.0 -0.997 in x, rad 0.005 -
    // approx WISH cylinder
    // We intentionally exclude the cylinder end caps so they this should
    // produce 0
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-0.5, 0.0, 0.0)), 0.0,
                    satol);
686
    // Other end
687
688
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-1.497, 0.0, 0.0)), 0.0,
                    satol);
689

690
    // Side values
691
692
693
694
695
696
697
698
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, 0.1)), 0.00301186,
                    satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, -0.1)), 0.00301186,
                    satol);
    // Sweep in the axis of the cylinder angle to see if the solid angle
    // decreases (as we are excluding the end caps)
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0.1, 0.0, 0.1)),
                    0.00100267, satol);
699

Russell Taylor's avatar
Russell Taylor committed
700
    // internal point (should be 4pi)
701
702
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-0.999, 0.0, 0.0)),
                    4 * M_PI, satol);
703
704

    // surface points
705
706
707
708
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-1.0, 0.0, 0.0)), 2 * M_PI,
                    satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-0.997, 0.0, 0.0)),
                    2 * M_PI, satol);
Russell Taylor's avatar
Russell Taylor committed
709
710
711
  }

  void testSolidAngleCubeTriangles()
712
713
714
715
  /**
  Test solid angle calculation for a cube using triangles
  - test for using Open Cascade surface triangulation for all solid angles.
  */
Russell Taylor's avatar
Russell Taylor committed
716
  {
717
    Object_sptr geom_obj = createUnitCube();
718
    double satol = 1e-3; // tolerance for solid angle
Russell Taylor's avatar
Russell Taylor committed
719
720
721
722
723

    // solid angle at distance 0.5 should be 4pi/6 by symmetry
    //
    // tests for Triangulated cube
    //
724
725
726
727
728
729
730
731
732
733
734
735
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(1.0, 0, 0)),
                    M_PI * 2.0 / 3.0, satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-1.0, 0, 0)),
                    M_PI * 2.0 / 3.0, satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 1.0, 0)),
                    M_PI * 2.0 / 3.0, satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, -1.0, 0)),
                    M_PI * 2.0 / 3.0, satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, 1.0)),
                    M_PI * 2.0 / 3.0, satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, -1.0)),
                    M_PI * 2.0 / 3.0, satol);
Russell Taylor's avatar
Russell Taylor committed
736
737
  }

738
  /** Add a scale factor */
739
  void testSolidAngleCubeTriangles_WithScaleFactor() {
740
    Object_sptr geom_obj = createUnitCube();
741
    double satol = 1e-3; // tolerance for solid angle
742
    // solid angle at distance 0.5 should be 4pi/6 by symmetry
743
744
745
746
    double expected = M_PI * 2.0 / 3.0;
    V3D scaleFactor(2.0, 2.0, 2.0);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(2.0, 0, 0), scaleFactor),
                    expected, satol);
747
748
  }

749
  void testGetBoundingBoxForCylinder()
750
751
752
  /**
  Test bounding box for a object capped cylinder
  */
Russell Taylor's avatar
Russell Taylor committed
753
  {
754
    Object_sptr geom_obj = createCappedCylinder();
755
756
757
758
759
760
761
762
763
764
    double xmax, ymax, zmax, xmin, ymin, zmin;
    xmax = ymax = zmax = 100;
    xmin = ymin = zmin = -100;
    geom_obj->getBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin);
    TS_ASSERT_DELTA(xmax, 1.2, 0.0001);
    TS_ASSERT_DELTA(ymax, 3.0, 0.0001);
    TS_ASSERT_DELTA(zmax, 3.0, 0.0001);
    TS_ASSERT_DELTA(xmin, -3.2, 0.0001);
    TS_ASSERT_DELTA(ymin, -3.0, 0.0001);
    TS_ASSERT_DELTA(zmin, -3.0, 0.0001);
Russell Taylor's avatar
Russell Taylor committed
765
766
  }

767
  void testdefineBoundingBox()
768
769
770
  /**
  Test use of defineBoundingBox
  */
Russell Taylor's avatar
Russell Taylor committed
771
  {
772
    Object_sptr geom_obj = createCappedCylinder();
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
    double xmax, ymax, zmax, xmin, ymin, zmin;
    xmax = 1.2;
    ymax = 3.0;
    zmax = 3.0;
    xmin = -3.2;
    ymin = -3.0;
    zmin = -3.0;

    TS_ASSERT_THROWS_NOTHING(
        geom_obj->defineBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin));

    const BoundingBox &boundBox = geom_obj->getBoundingBox();

    TS_ASSERT_EQUALS(boundBox.xMax(), 1.2);
    TS_ASSERT_EQUALS(boundBox.yMax(), 3.0);
    TS_ASSERT_EQUALS(boundBox.zMax(), 3.0);
    TS_ASSERT_EQUALS(boundBox.xMin(), -3.2);
    TS_ASSERT_EQUALS(boundBox.yMin(), -3.0);
    TS_ASSERT_EQUALS(boundBox.zMin(), -3.0);

    // Inconsistent bounding box
    xmax = 1.2;
    xmin = 3.0;
    TS_ASSERT_THROWS(
        geom_obj->defineBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin),
        std::invalid_argument);
Russell Taylor's avatar
Russell Taylor committed
799
800
  }
  void testSurfaceTriangulation()
801
802
803
  /**
  Test triangle solid angle calc
  */
Russell Taylor's avatar
Russell Taylor committed
804
  {
805
    Object_sptr geom_obj = createCappedCylinder();
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
    double xmax, ymax, zmax, xmin, ymin, zmin;
    xmax = 20;
    ymax = 20.0;
    zmax = 20.0;
    xmin = -20.0;
    ymin = -20.0;
    zmin = -20.0;
    geom_obj->getBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin);
    double saTri, saRay;
    V3D observer(4.2, 0, 0);

    double satol = 1e-3; // typical result tolerance

    //    if(timeTest)
    //    {
    //      // block to test time of solid angle methods
    //      // change false to true to include
    //      int iter=4000;
    //      int starttime=clock();
    //      for (int i=0;i<iter;i++)
    //        saTri=geom_obj->triangleSolidAngle(observer);
    //      int endtime=clock();
    //      std::cout << std::endl << "Cyl tri time=" <<
    //      (endtime-starttime)/(static_cast<double>(CLOCKS_PER_SEC*iter)) <<
    //      std::endl;
    //      iter=50;
    //      starttime=clock();
    //      for (int i=0;i<iter;i++)
    //        saRay=geom_obj->rayTraceSolidAngle(observer);
    //      endtime=clock();
    //      std::cout << "Cyl ray time=" <<
    //      (endtime-starttime)/(static_cast<double>(CLOCKS_PER_SEC*iter)) <<
    //      std::endl;
    //    }

    saTri = geom_obj->triangleSolidAngle(observer);
    saRay = geom_obj->rayTraceSolidAngle(observer);
    TS_ASSERT_DELTA(saTri, 1.840302, 0.001);
    TS_ASSERT_DELTA(saRay, 1.840302, 0.01);

    observer = V3D(-7.2, 0, 0);
    saTri = geom_obj->triangleSolidAngle(observer);
    saRay = geom_obj->rayTraceSolidAngle(observer);

    TS_ASSERT_DELTA(saTri, 1.25663708, 0.001);
    TS_ASSERT_DELTA(saRay, 1.25663708, 0.001);
Russell Taylor's avatar
Russell Taylor committed
852
853

    // No analytic value for side on SA, using hi-res value
854
855
856
857
858
859
860
861
862
863
864
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, 7)), 0.7531,
                    0.753 * satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 7, 0)), 0.7531,
                    0.753 * satol);

    saTri = geom_obj->triangleSolidAngle(V3D(20, 0, 0));
    TS_ASSERT_DELTA(saTri, 0.07850147, satol * 0.0785);
    saTri = geom_obj->triangleSolidAngle(V3D(200, 0, 0));
    TS_ASSERT_DELTA(saTri, 0.000715295, satol * 0.000715);
    saTri = geom_obj->triangleSolidAngle(V3D(2000, 0, 0));
    TS_ASSERT_DELTA(saTri, 7.08131e-6, satol * 7.08e-6);
Russell Taylor's avatar
Russell Taylor committed
865
866
  }
  void testSolidAngleSphereTri()
867
868
869
  /**
  Test solid angle calculation for a sphere from triangulation
  */
Russell Taylor's avatar
Russell Taylor committed
870
  {
871
    Object_sptr geom_obj = createSphere();
872
    double satol = 1e-3; // tolerance for solid angle
Russell Taylor's avatar
Russell Taylor committed
873
874
875
876
877

    // Solid angle at distance 8.1 from centre of sphere radius 4.1 x/y/z
    // Expected solid angle calculated values from sa=2pi(1-cos(arcsin(R/r))
    // where R is sphere radius and r is distance of observer from sphere centre
    // Intercept for track in reverse direction now worked round
878
879
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(8.1, 0, 0)), 0.864364,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
880
    // internal point (should be 4pi)
881
882
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0, 0, 0)), 4 * M_PI,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
883
    // surface point
884
885
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(4.1, 0, 0)), 2 * M_PI,
                    satol);
Russell Taylor's avatar
Russell Taylor committed
886
887
888
889
  }

private:
  /// Surface type
890
  typedef std::map<int, std::unique_ptr<Surface>> STYPE;
Russell Taylor's avatar
Russell Taylor committed
891
892

  /// set timeTest true to get time comparisons of soild angle methods
893
  const static bool timeTest = false;
894

895
  STYPE SMap; ///< Surface Map
Russell Taylor's avatar
Russell Taylor committed
896

897
898
899
900
  Object_sptr createCappedCylinder() {
    std::string C31 = "cx 3.0"; // cylinder x-axis radius 3
    std::string C32 = "px 1.2";
    std::string C33 = "px -3.2";
Russell Taylor's avatar
Russell Taylor committed
901
902

    // First create some surfaces
903
904
905
906
    std::map<int, std::unique_ptr<Surface>> CylSurMap;
    CylSurMap[31] = std::unique_ptr<Surface>(new Cylinder());
    CylSurMap[32] = std::unique_ptr<Surface>(new Plane());
    CylSurMap[33] = std::unique_ptr<Surface>(new Plane());
Russell Taylor's avatar
Russell Taylor committed
907
908

    CylSurMap[31]->setSurface(C31);
909
910
911
912
913
914
    CylSurMap[32]->setSurface(C32);
    CylSurMap[33]->setSurface(C33);
    CylSurMap[31]->setName(31);
    CylSurMap[32]->setName(32);
    CylSurMap[33]->setName(33);

915
    // Capped cylinder (id 21)
916
    // using surface ids: 31 (cylinder) 32 (plane (top) ) and 33 (plane (base))
917
    std::string ObjCapCylinder = "-31 -32 33";
918

919
    Object_sptr retVal = Object_sptr(new Object);
920
    retVal->setObject(21, ObjCapCylinder);
921
922
    retVal->populate(CylSurMap);

923
    TS_ASSERT(retVal.get());
924
925
926

    return retVal;
  }
927

928
929
  // This creates a cylinder to test the solid angle that is more realistic in
  // size
930
  // for a detector cylinder
931
932
933
934
935
  Object_sptr createSmallCappedCylinder() {
    std::string C31 =
        "cx 0.005"; // cylinder x-axis radius 0.005 and height 0.003
    std::string C32 = "px -0.997";
    std::string C33 = "px -1.0";
936
937

    // First create some surfaces
938
939
940
941
    std::map<int, std::unique_ptr<Surface>> CylSurMap;
    CylSurMap[31] = std::unique_ptr<Cylinder>(new Cylinder());
    CylSurMap[32] = std::unique_ptr<Plane>(new Plane());
    CylSurMap[33] = std::unique_ptr<Plane>(new Plane());
942
943

    CylSurMap[31]->setSurface(C31);
Russell Taylor's avatar
Russell Taylor committed
944
945
946
947
948
949
    CylSurMap[32]->setSurface(C32);
    CylSurMap[33]->setSurface(C33);
    CylSurMap[31]->setName(31);
    CylSurMap[32]->setName(32);
    CylSurMap[33]->setName(33);

950
    // Capped cylinder (id 21)
Russell Taylor's avatar
Russell Taylor committed
951
    // using surface ids: 31 (cylinder) 32 (plane (top) ) and 33 (plane (base))
952
    std::string ObjCapCylinder = "-31 -32 33";
Russell Taylor's avatar
Russell Taylor committed
953

954
    Object_sptr retVal = Object_sptr(new Object);
955
    retVal->setObject(21, ObjCapCylinder);
956
    retVal->populate(CylSurMap);
Russell Taylor's avatar
Russell Taylor committed
957
958
959
960

    return retVal;
  }

961
962
  Object_sptr createSphere() {
    std::string S41 = "so 4.1"; // Sphere at origin radius 4.1
Russell Taylor's avatar
Russell Taylor committed
963
964

    // First create some surfaces
965
966
    std::map<int, std::unique_ptr<Surface>> SphSurMap;
    SphSurMap[41] = std::unique_ptr<Sphere>(new Sphere());
Russell Taylor's avatar
Russell Taylor committed
967
968
969
    SphSurMap[41]->setSurface(S41);
    SphSurMap[41]->setName(41);

970
    // A sphere
971
    std::string ObjSphere = "-41";
Russell Taylor's avatar
Russell Taylor committed
972

973
    Object_sptr retVal = Object_sptr(new Object);
974
    retVal->setObject(41, ObjSphere);
975
    retVal->populate(SphSurMap);
Russell Taylor's avatar
Russell Taylor committed
976
977
978
979
980

    return retVal;
  }

  void clearSurfMap()
981
982
983
984
  /**
  Clears the surface map for a new test
  or destruction.
  */
Russell Taylor's avatar
Russell Taylor committed
985
  {
986
    SMap.clear();
Russell Taylor's avatar
Russell Taylor committed
987
988
989
    return;
  }

990
991
992
993
994
  void createSurfaces(const std::string &desired)
  /**
  Creates a list of surfaces for used in the objects
  and populates the MObj layers.
  */
Russell Taylor's avatar
Russell Taylor committed
995
996
997
998
999
  {
    clearSurfMap();

    // PLANE SURFACES:

1000
    typedef std::pair<int, std::string> SCompT;