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

#include <cxxtest/TestSuite.h>
#include <cmath>
#include <ostream>
#include <vector>
#include <algorithm>
#include <ctime>
10
#include <iomanip>
Russell Taylor's avatar
Russell Taylor committed
11

12
13
#include <boost/shared_ptr.hpp>

Russell Taylor's avatar
Russell Taylor committed
14
#include "MantidGeometry/V3D.h" 
Nick Draper's avatar
re #843    
Nick Draper committed
15
16
17
18
19
20
21
22
#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" 
#include "MantidGeometry/Rendering/GluGeometryHandler.h"
23
#include "MantidGeometry/Objects/BoundingBox.h"
24

Russell Taylor's avatar
Russell Taylor committed
25
26
27
using namespace Mantid;
using namespace Geometry;

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

public:


Nick Draper's avatar
re #843    
Nick Draper committed
34
  void testCreateUnitCube()
Russell Taylor's avatar
Russell Taylor committed
35
  {
36
    Object_sptr geom_obj = createUnitCube();
37

38
    TS_ASSERT_EQUALS(geom_obj->str(),"68 -6 5 -4 3 -2 1");
39
40

    double xmin(0.0), xmax(0.0), ymin(0.0), ymax(0.0), zmin(0.0), zmax(0.0);
41
    geom_obj->getBoundingBox(xmax, ymax, zmax, xmin, ymin, zmin);
42

Russell Taylor's avatar
Russell Taylor committed
43
44
  }

Nick Draper's avatar
re #843    
Nick Draper committed
45
  void testIsOnSideCappedCylinder()
Russell Taylor's avatar
Russell Taylor committed
46
  {
47
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
48
    //inside
49
50
51
52
53
    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);
Russell Taylor's avatar
Russell Taylor committed
54
    //on the side
55
56
57
58
59
60
    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);
Russell Taylor's avatar
Russell Taylor committed
61
62

    //on the edges
63
64
65
66
67
68
69
70
    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);
Russell Taylor's avatar
Russell Taylor committed
71
    //out side
72
73
74
75
76
77
    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
78
79
80
81
  }

  void testIsValidCappedCylinder()
  {
82
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
83
    //inside
84
85
86
87
88
    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);
Russell Taylor's avatar
Russell Taylor committed
89
    //on the side
90
91
92
93
94
95
    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);
Russell Taylor's avatar
Russell Taylor committed
96
97

    //on the edges
98
99
100
101
102
103
104
105
    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);
Russell Taylor's avatar
Russell Taylor committed
106
    //out side
107
108
109
110
111
112
    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
113
114
115
116
  }

  void testIsOnSideSphere()
  {
117
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
118
    //inside
119
120
121
122
123
    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);
Russell Taylor's avatar
Russell Taylor committed
124
    //on the side
125
126
127
128
    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);
Russell Taylor's avatar
Russell Taylor committed
129
130

    //out side
131
132
133
134
    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
135
136
137
138
  }

  void testIsValidSphere()
  {
139
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
140
    //inside
141
142
143
144
145
    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);
Russell Taylor's avatar
Russell Taylor committed
146
    //on the side
147
148
149
150
    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);
Russell Taylor's avatar
Russell Taylor committed
151
152

    //out side
153
154
155
156
    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
157
158
159
160
  }

  void testCalcValidTypeSphere()
  {
161
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
162
    //entry on the normal
163
164
165
166
167
168
169
170
    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);
Russell Taylor's avatar
Russell Taylor committed
171
172

    //a glancing blow
173
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-4.1,0,0),V3D(0,1,0)),0);
Russell Taylor's avatar
Russell Taylor committed
174
    //not quite on the normal
175
176
    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
177
178
  }

179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
  void testGetBoundingBoxForSphere()
  {
    Object_sptr geom_obj = createSphere();    
    const double tolerance(1e-10);

    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);

195
196
197
198
199
200
201
202
    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);
203
204
  }

Russell Taylor's avatar
Russell Taylor committed
205
206
  void testCalcValidTypeCappedCylinder()
  {
207
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
208
    //entry on the normal
209
210
211
212
213
214
215
216
    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);
Russell Taylor's avatar
Russell Taylor committed
217
218

    //a glancing blow
219
    TS_ASSERT_EQUALS(geom_obj->calcValidType(V3D(-3.2,0,0),V3D(0,1,0)),0);
Russell Taylor's avatar
Russell Taylor committed
220
    //not quite on the normal
221
222
    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
223
224
225
226
  }

  void testInterceptSurfaceSphereZ()
  {
227
    std::vector<Link> expectedResults;
Russell Taylor's avatar
Russell Taylor committed
228
229
230
231
232
233
234
235
236
237
238
    std::string S41="s 1 1 1 4";         // Sphere at (1,1,1) radius 4

    // First create some surfaces
    std::map<int,Surface*> SphSurMap;
    SphSurMap[41]=new Sphere();
    SphSurMap[41]->setSurface(S41);
    SphSurMap[41]->setName(41);

    // A sphere 
    std::string ObjSphere="-41" ;

239
240
241
    Object_sptr geom_obj = Object_sptr(new Object); 
    geom_obj->setObject(41,ObjSphere);
    geom_obj->populate(SphSurMap);
Russell Taylor's avatar
Russell Taylor committed
242
243
244
245


    Track track(V3D(-1,1.5,1),V3D(1,0,0));

246
    // format = startPoint, endPoint, total distance so far
Russell Taylor's avatar
Russell Taylor committed
247
    // forward only intercepts means that start point should be track origin
248
249
    expectedResults.push_back(Link(V3D(-1,1.5,1),
      V3D(sqrt(16-0.25)+1,1.5,1.0),sqrt(15.75)+2));
Russell Taylor's avatar
Russell Taylor committed
250

251
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
252
253
254
255
  }

  void testInterceptSurfaceSphereY()
  {
256
    std::vector<Link> expectedResults;
257
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
258
259
    Track track(V3D(0,-10,0),V3D(0,1,0));

260
261
    //format = startPoint, endPoint, total distance so far
    expectedResults.push_back(Link(V3D(0,-4.1,0),V3D(0,4.1,0),14.1));
Russell Taylor's avatar
Russell Taylor committed
262

263
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
264
265
266
267
  }

  void testInterceptSurfaceSphereX()
  {
268
    std::vector<Link> expectedResults;
269
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
270
    Track track(V3D(-10,0,0),V3D(1,0,0));
271

272
273
    //format = startPoint, endPoint, total distance so far
    expectedResults.push_back(Link(V3D(-4.1,0,0),V3D(4.1,0,0),14.1));
274
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
275
276
277
278
  }

  void testInterceptSurfaceCappedCylinderY()
  {
279
    std::vector<Link> expectedResults;
280
    Object_sptr geom_obj = createCappedCylinder();
281
282
    //format = startPoint, endPoint, total distance so far
    expectedResults.push_back(Link(V3D(0,-3,0),V3D(0,3,0),13));
Russell Taylor's avatar
Russell Taylor committed
283
284

    Track track(V3D(0,-10,0),V3D(0,1,0));
285
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
286
287
288
289
  }

  void testInterceptSurfaceCappedCylinderX()
  {
290
    std::vector<Link> expectedResults;
291
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
292
293
    Track track(V3D(-10,0,0),V3D(1,0,0));

294
295
    //format = startPoint, endPoint, total distance so far
    expectedResults.push_back(Link(V3D(-3.2,0,0),V3D(1.2,0,0),11.2));
296
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
297
298
299
300
  }

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

305
    checkTrackIntercept(geom_obj,track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
306
307
  }

308
  void checkTrackIntercept(Track& track, const std::vector<Link>& expectedResults)
Russell Taylor's avatar
Russell Taylor committed
309
310
311
312
  {
    int index = 0;
    for (Track::LType::const_iterator it = track.begin(); it!=track.end();++it)
    {
313
314
      TS_ASSERT_DELTA(it->distFromStart,expectedResults[index].distFromStart,1e-6);
      TS_ASSERT_DELTA(it->distInsideObject,expectedResults[index].distInsideObject,1e-6);
315
      TS_ASSERT_EQUALS(it->componentID,expectedResults[index].componentID);
316
317
      TS_ASSERT_EQUALS(it->entryPoint,expectedResults[index].entryPoint);
      TS_ASSERT_EQUALS(it->exitPoint,expectedResults[index].exitPoint);
Russell Taylor's avatar
Russell Taylor committed
318
319
320
321
322
      ++index;
    }
    TS_ASSERT_EQUALS(index,static_cast<int>(expectedResults.size()));
  }

323
  void checkTrackIntercept(Object_sptr obj, Track& track, const std::vector<Link>& expectedResults)
Russell Taylor's avatar
Russell Taylor committed
324
  {
325
    int unitCount = obj->interceptSurface(track);
326
327
    TS_ASSERT_EQUALS(unitCount,expectedResults.size());
    checkTrackIntercept(track,expectedResults);
Russell Taylor's avatar
Russell Taylor committed
328
329
330
  }

  void xtestTrackTwoIsolatedCubes()
331
    /**
Russell Taylor's avatar
Russell Taylor committed
332
333
334
335
336
    Test a track going through an object
    */
  {
    std::string ObjA="60001 -60002 60003 -60004 60005 -60006";
    std::string ObjB="80001 -80002 60003 -60004 60005 -60006";
337

338
    createSurfaces(ObjA);
Russell Taylor's avatar
Russell Taylor committed
339
340
341
342
    Object object1=Object();
    object1.setObject(3,ObjA);
    object1.populate(SMap);

343
    createSurfaces(ObjB);
Russell Taylor's avatar
Russell Taylor committed
344
345
346
347
348
349
350
351
    Object object2=Object();
    object2.setObject(4,ObjB);
    object2.populate(SMap);

    Track TL(Geometry::V3D(-5,0,0),
      Geometry::V3D(1,0,0));

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

355
356
357
    std::vector<Link> expectedResults;
    expectedResults.push_back(Link(V3D(-1,0,0),V3D(1,0,0),6));
    expectedResults.push_back(Link(V3D(4.5,0,0),V3D(6.5,0,0),11.5));
Russell Taylor's avatar
Russell Taylor committed
358
359
360
361
362
    checkTrackIntercept(TL,expectedResults);

  }

  void testTrackTwoTouchingCubes()
363
    /**
Russell Taylor's avatar
Russell Taylor committed
364
365
366
367
368
369
    Test a track going through an object
    */
  {
    std::string ObjA="60001 -60002 60003 -60004 60005 -60006";
    std::string ObjB="60002 -80002 60003 -60004 60005 -60006";

370
    createSurfaces(ObjA);
Russell Taylor's avatar
Russell Taylor committed
371
372
373
374
    Object object1=Object();
    object1.setObject(3,ObjA);
    object1.populate(SMap);

375
    createSurfaces(ObjB);
Russell Taylor's avatar
Russell Taylor committed
376
377
378
379
    Object object2=Object();
    object2.setObject(4,ObjB);
    object2.populate(SMap);

380
    Track TL(Geometry::V3D(-5,0,0), Geometry::V3D(1,0,0));
Russell Taylor's avatar
Russell Taylor committed
381
382

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

386
387
388
    std::vector<Link> expectedResults;
    expectedResults.push_back(Link(V3D(-1,0,0),V3D(1,0,0),6));
    expectedResults.push_back(Link(V3D(1,0,0),V3D(6.5,0,0),11.5));
Russell Taylor's avatar
Russell Taylor committed
389
390
391
392
393
394

    checkTrackIntercept(TL,expectedResults);

  }

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

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

407
    createSurfaces(ObjB);
Russell Taylor's avatar
Russell Taylor committed
408
409
410
411
412
413
414
415
416
417
418
    Object object2=Object();
    object2.setObject(4,ObjB);
    object2.populate(SMap);

    Track TL(Geometry::V3D(-5,0,0),
      Geometry::V3D(1,0,0));

    // CARE: This CANNOT be called twice
    TS_ASSERT(object1.interceptSurface(TL)!=0);
    TS_ASSERT(object2.interceptSurface(TL)!=0);

419
420
421
422
    std::vector<Link> expectedResults;
    expectedResults.push_back(Link(V3D(-1,0,0),V3D(-0.8,0,0),4.2));
    expectedResults.push_back(Link(V3D(-0.8,0,0),V3D(0.8,0,0),5.8));
    expectedResults.push_back(Link(V3D(0.8,0,0),V3D(1,0,0),6));
Russell Taylor's avatar
Russell Taylor committed
423
424
425
426
    checkTrackIntercept(TL,expectedResults);
  }

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

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

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

444
    Track TL(Geometry::V3D(-5,0,0), Geometry::V3D(1,0,0));
Russell Taylor's avatar
Russell Taylor committed
445
446
447
448
449
450


    // CARE: This CANNOT be called twice
    TS_ASSERT(object1.interceptSurface(TL)!=0);
    TS_ASSERT(object2.interceptSurface(TL)!=0);

451
452
453
454
    std::vector<Link> expectedResults;
    expectedResults.push_back(Link(V3D(-1,0,0),V3D(-0.4,0,0),4.6));
    expectedResults.push_back(Link(V3D(-0.4,0,0),V3D(0.2,0,0),5.2));
    expectedResults.push_back(Link(V3D(0.2,0,0),V3D(1,0,0),6));
Russell Taylor's avatar
Russell Taylor committed
455
456
457
458
    checkTrackIntercept(TL,expectedResults);
  }

  void testTrack_CubePlusInternalEdgeTouchSpheresMiss()
459
    /**
Russell Taylor's avatar
Russell Taylor committed
460
461
462
463
464
    Test a track missing an object
    */
  {
    std::string ObjA="60001 -60002 60003 -60004 60005 -60006 72 73";
    std::string ObjB="(-72 : -73)";
465

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

471
    createSurfaces(ObjB);
Russell Taylor's avatar
Russell Taylor committed
472
473
474
475
476
477
478
479
480
481
482
483
    Object object2=Object();
    object2.setObject(4,ObjB);
    object2.populate(SMap);

    Track TL(Geometry::V3D(-5,0,0),
      Geometry::V3D(0,1,0));


    // CARE: This CANNOT be called twice
    TS_ASSERT_EQUALS(object1.interceptSurface(TL),0);
    TS_ASSERT_EQUALS(object2.interceptSurface(TL),0);

484
    std::vector<Link> expectedResults; //left empty as this should miss
Russell Taylor's avatar
Russell Taylor committed
485
486
487
488
    checkTrackIntercept(TL,expectedResults);
  }

  void testFindPointInCube()
489
    /**
Russell Taylor's avatar
Russell Taylor committed
490
491
492
    Test find point in cube
    */
  {
493
    Object_sptr geom_obj = createUnitCube();
Russell Taylor's avatar
Russell Taylor committed
494
495
    // initial guess in object
    Geometry::V3D pt;
496
    TS_ASSERT_EQUALS(geom_obj->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
497
498
499
500
501
502
    TS_ASSERT_EQUALS(pt,V3D(0,0,0));
    // initial guess not in object, but on x-axis
    std::vector<std::string> planes;
    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");
503
504
    Object_sptr B =createCuboid(planes);
    TS_ASSERT_EQUALS(B->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
505
506
507
508
509
510
    TS_ASSERT_EQUALS(pt,V3D(10,0,0));
    // on y axis
    planes.clear();
    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");
511
512
    Object_sptr C =createCuboid(planes);
    TS_ASSERT_EQUALS(C->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
513
514
515
516
517
518
    TS_ASSERT_EQUALS(pt,V3D(0,-21,0));
    // not on principle axis, now works using getBoundingBox
    planes.clear();
    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");
519
520
    Object_sptr D =createCuboid(planes);
    TS_ASSERT_EQUALS(D->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
521
522
523
524
525
526
527
528
529
530
531
532
    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);
    planes.clear();
    // 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
    // common
    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");
533
534
    Object_sptr E =createCuboid(planes);
    TS_ASSERT_EQUALS(E->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
535
536
537
538
539
540
541
542
543
544
545
    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);
    planes.clear();
    // This test fails to find a point in object, as object not on a principle axis
    // and getBoundingBox does not give a useful result in this case.
    // Object is unit cube located at +-0.5 in x but centred on z=y=-1.606.. and rotated 45deg
    // to these two axes
    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");
546
547
    Object_sptr F =createCuboid(planes);
    TS_ASSERT_EQUALS(F->getPointInObject(pt),0);
Russell Taylor's avatar
Russell Taylor committed
548
    // Test use of defineBoundingBox to explictly set the bounding box, when the automatic method fails
549
    F->defineBoundingBox(0.5,-1/(2.0*sqrt(2.0)),-1.0/(2.0*sqrt(2.0)),
Russell Taylor's avatar
Russell Taylor committed
550
      -0.5,-sqrt(2.0)-1.0/(2.0*sqrt(2.0)),-sqrt(2.0)-1.0/(2.0*sqrt(2.0)));
551
552
553
    TS_ASSERT_EQUALS(F->getPointInObject(pt),1);
    Object_sptr S = createSphere();
    TS_ASSERT_EQUALS(S->getPointInObject(pt),1);
Russell Taylor's avatar
Russell Taylor committed
554
555
556
557
558
    TS_ASSERT_EQUALS(pt,V3D(0.0,0.0,0));
  }


  void testSolidAngleSphere()
559
    /**
Russell Taylor's avatar
Russell Taylor committed
560
561
562
    Test solid angle calculation for a sphere
    */
  {
563
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
564
565
566
567
568
569
    double satol=2e-2; // tolerance for solid angle

    // 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
570
571
572
573
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(8.1,0,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(0,8.1,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(0,0,8.1)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(0,0,-8.1)),0.864364,satol);
Russell Taylor's avatar
Russell Taylor committed
574
    // internal point (should be 4pi)
575
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(0,0,0)),4*M_PI,satol);
Russell Taylor's avatar
Russell Taylor committed
576
    // surface point
577
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(4.1,0,0)),2*M_PI,satol);
Russell Taylor's avatar
Russell Taylor committed
578
    // distant points
579
580
581
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(20,0,0)),0.133442,satol);
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(200,0,0)),0.0013204,satol);
    TS_ASSERT_DELTA(geom_obj->rayTraceSolidAngle(V3D(2000,0,0)),1.32025e-5,satol);
Russell Taylor's avatar
Russell Taylor committed
582
583
584
    //
    // test solidAngle interface, which will be main method to solid angle
    //
585
586
587
588
    TS_ASSERT_DELTA(geom_obj->solidAngle(V3D(8.1,0,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->solidAngle(V3D(0,8.1,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->solidAngle(V3D(0,0,8.1)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->solidAngle(V3D(0,0,-8.1)),0.864364,satol);
Russell Taylor's avatar
Russell Taylor committed
589
590
591
    //
  }

592

Russell Taylor's avatar
Russell Taylor committed
593
  void testSolidAngleCappedCylinder()
594
    /**
Russell Taylor's avatar
Russell Taylor committed
595
596
597
    Test solid angle calculation for a capped cylinder
    */
  {
598
    Object_sptr geom_obj = createSmallCappedCylinder();
599
    // Want to test triangulation so setup a geometry handler
600
    boost::shared_ptr<GluGeometryHandler> h = boost::shared_ptr<GluGeometryHandler>(new GluGeometryHandler(geom_obj.get()));
601
    h->setCylinder(V3D(-0.0015,0.0,0.0), V3D(1., 0.0, 0.0), 0.005, 0.003);
602
    geom_obj->setGeometryHandler(h);
603

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

606
    // solid angle at point -0.5 from capped cyl -1.0 -0.997 in x, rad 0.005 - approx WISH cylinder
607
608
    // 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);
609
    // Other end
610
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-1.497, 0.0, 0.0)), 0.0, satol);
611

612
613
614
615
616
    // Side values
    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);
617

Russell Taylor's avatar
Russell Taylor committed
618
    // internal point (should be 4pi)
619
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(-0.999, 0.0, 0.0)),4*M_PI,satol);
620
621

    // surface points
622
623
    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);
624

Russell Taylor's avatar
Russell Taylor committed
625
626
  }

627

Russell Taylor's avatar
Russell Taylor committed
628
  void testSolidAngleCubeTriangles()
629
    /**
Russell Taylor's avatar
Russell Taylor committed
630
631
632
633
    Test solid angle calculation for a cube using triangles
    - test for using Open Cascade surface triangulation for all solid angles.
    */
  {
634
    Object_sptr geom_obj = createUnitCube();
Russell Taylor's avatar
Russell Taylor committed
635
636
637
638
639
640
    double satol=1e-3; // tolerance for solid angle

    // solid angle at distance 0.5 should be 4pi/6 by symmetry
    //
    // tests for Triangulated cube
    //
641
642
643
644
645
646
    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
647

648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
//    if(timeTest)
//    {
//      // block to test time of solid angle methods
//      // change false to true to include
//      double saRay,saTri;
//      V3D observer(1.0,0,0);
//      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 << "Cube 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 << "Cube ray time=" << (endtime-starttime)/(static_cast<double>(CLOCKS_PER_SEC*iter)) << std::endl;
//    }
Russell Taylor's avatar
Russell Taylor committed
667
668
669

  }

670
671
672
673
674
675
676
677
678
679
680
681
682
683


  /** Add a scale factor */
  void testSolidAngleCubeTriangles_WithScaleFactor()
  {
    Object_sptr geom_obj = createUnitCube();
    double satol=1e-3; // tolerance for solid angle
    // solid angle at distance 0.5 should be 4pi/6 by symmetry
    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);
  }


684
  void testGetBoundingBoxForCylinder()
685
    /**
Russell Taylor's avatar
Russell Taylor committed
686
687
688
    Test bounding box for a object capped cylinder
    */
  {
689
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
690
691
692
    double xmax,ymax,zmax,xmin,ymin,zmin;
    xmax=ymax=zmax=100;
    xmin=ymin=zmin=-100;
693
    geom_obj->getBoundingBox(xmax,ymax,zmax,xmin,ymin,zmin);
Russell Taylor's avatar
Russell Taylor committed
694
695
696
697
698
699
700
701
    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);
  }

702
  void testdefineBoundingBox()
703
    /**
Russell Taylor's avatar
Russell Taylor committed
704
705
706
    Test use of defineBoundingBox
    */
  {
707
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
708
709
710
    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;
711
712

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

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

716
717
718
719
720
721
    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);
722
723

    //Inconsistent bounding box
Russell Taylor's avatar
Russell Taylor committed
724
    xmax=1.2;xmin=3.0;
725
    TS_ASSERT_THROWS(geom_obj->defineBoundingBox(xmax,ymax,zmax,xmin,ymin,zmin),std::invalid_argument);
Russell Taylor's avatar
Russell Taylor committed
726
727
728

  }
  void testSurfaceTriangulation()
729
    /**
Russell Taylor's avatar
Russell Taylor committed
730
731
732
    Test triangle solid angle calc
    */
  {
733
    Object_sptr geom_obj = createCappedCylinder();
Russell Taylor's avatar
Russell Taylor committed
734
735
736
    double xmax,ymax,zmax,xmin,ymin,zmin;
    xmax=20;ymax=20.0;zmax=20.0;
    xmin=-20.0;ymin=-20.0;zmin=-20.0;
737
    geom_obj->getBoundingBox(xmax,ymax,zmax,xmin,ymin,zmin);
Russell Taylor's avatar
Russell Taylor committed
738
739
    double saTri,saRay;
    V3D observer(4.2,0,0);
740

Russell Taylor's avatar
Russell Taylor committed
741
742
    double satol=1e-3; // typical result tolerance

743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
//    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;
//    }
Russell Taylor's avatar
Russell Taylor committed
760

761
762
    saTri=geom_obj->triangleSolidAngle(observer);
    saRay=geom_obj->rayTraceSolidAngle(observer);
Russell Taylor's avatar
Russell Taylor committed
763
764
    TS_ASSERT_DELTA(saTri,1.840302,0.001);
    TS_ASSERT_DELTA(saRay,1.840302,0.01);
765

Russell Taylor's avatar
Russell Taylor committed
766
    observer=V3D(-7.2,0,0);
767
768
    saTri=geom_obj->triangleSolidAngle(observer);
    saRay=geom_obj->rayTraceSolidAngle(observer);
769

Russell Taylor's avatar
Russell Taylor committed
770
771
772
773
    TS_ASSERT_DELTA(saTri,1.25663708,0.001);
    TS_ASSERT_DELTA(saRay,1.25663708,0.001);

    // No analytic value for side on SA, using hi-res value
774
775
    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);
Russell Taylor's avatar
Russell Taylor committed
776

777
    saTri=geom_obj->triangleSolidAngle(V3D(20,0,0));
Russell Taylor's avatar
Russell Taylor committed
778
    TS_ASSERT_DELTA(saTri,0.07850147,satol*0.0785);
779
    saTri=geom_obj->triangleSolidAngle(V3D(200,0,0));
Russell Taylor's avatar
Russell Taylor committed
780
    TS_ASSERT_DELTA(saTri,0.000715295,satol*0.000715);
781
    saTri=geom_obj->triangleSolidAngle(V3D(2000,0,0));
Russell Taylor's avatar
Russell Taylor committed
782
    TS_ASSERT_DELTA(saTri,7.08131e-6,satol*7.08e-6);
783

Russell Taylor's avatar
Russell Taylor committed
784
785
  }
  void testSolidAngleSphereTri()
786
    /**
Russell Taylor's avatar
Russell Taylor committed
787
788
789
    Test solid angle calculation for a sphere from triangulation
    */
  {
790
    Object_sptr geom_obj = createSphere();
Russell Taylor's avatar
Russell Taylor committed
791
792
793
794
795
796
    double satol=1e-3; // tolerance for solid angle

    // 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
797
798
799
800
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(8.1,0,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0,8.1,0)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0,0,8.1)),0.864364,satol);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0,0,-8.1)),0.864364,satol);
Russell Taylor's avatar
Russell Taylor committed
801
    // internal point (should be 4pi)
802
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(0,0,0)),4*M_PI,satol);
Russell Taylor's avatar
Russell Taylor committed
803
    // surface point
804
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(4.1,0,0)),2*M_PI,satol);
Russell Taylor's avatar
Russell Taylor committed
805
    // distant points
806
807
808
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(20,0,0)),0.133442,satol*0.133);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(200,0,0)),0.0013204,satol*0.00132);
    TS_ASSERT_DELTA(geom_obj->triangleSolidAngle(V3D(2000,0,0)),1.32025e-5,satol*1.32e-5);
Russell Taylor's avatar
Russell Taylor committed
809

810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
//    if(timeTest)
//    {
//      // block to test time of solid angle methods
//      // change false to true to include
//      double saTri,saRay;
//      int iter=400;
//      V3D observer(8.1,0,0);
//      int starttime=clock();
//      for (int i=0;i<iter;i++)
//        saTri=geom_obj->triangleSolidAngle(observer);
//      int endtime=clock();
//      std::cout << std::endl << "Sphere tri time =" << (endtime-starttime)/(static_cast<double>(CLOCKS_PER_SEC*iter)) << std::endl;
//      iter=40;
//      starttime=clock();
//      for (int i=0;i<iter;i++)
//        saRay=geom_obj->rayTraceSolidAngle(observer);
//      endtime=clock();
//      std::cout << "Sphere ray time =" << (endtime-starttime)/(static_cast<double>(CLOCKS_PER_SEC*iter)) << std::endl;
//    }
Russell Taylor's avatar
Russell Taylor committed
829
830
831
832
833
834
835
836
837
838

  }

private:

  /// Surface type
  typedef std::map<int,Surface*> STYPE ; 

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

Russell Taylor's avatar
Russell Taylor committed
840
841
  STYPE SMap;   ///< Surface Map

842
  Object_sptr createCappedCylinder()
Russell Taylor's avatar
Russell Taylor committed
843
844
845
846
847
848
849
850
851
852
853
854
  {
    std::string C31="cx 3.0";         // cylinder x-axis radius 3
    std::string C32="px 1.2";
    std::string C33="px -3.2";

    // First create some surfaces
    std::map<int,Surface*> CylSurMap;
    CylSurMap[31]=new Cylinder();
    CylSurMap[32]=new Plane();
    CylSurMap[33]=new Plane();

    CylSurMap[31]->setSurface(C31);
855
856
857
858
859
860
861
862
863
864
    CylSurMap[32]->setSurface(C32);
    CylSurMap[33]->setSurface(C33);
    CylSurMap[31]->setName(31);
    CylSurMap[32]->setName(32);
    CylSurMap[33]->setName(33);

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

865
    Object_sptr retVal = Object_sptr(new Object); 
866
867
868
    retVal->setObject(21,ObjCapCylinder);
    retVal->populate(CylSurMap);

869
    TS_ASSERT(retVal.get());
870
871
872

    return retVal;
  }
873

874
875
  // This creates a cylinder to test the solid angle that is more realistic in size
  // for a detector cylinder
876
  Object_sptr createSmallCappedCylinder()
877
878
879
880
881
882
883
884
885
886
887
888
  {
    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";

    // First create some surfaces
    std::map<int,Surface*> CylSurMap;
    CylSurMap[31]=new Cylinder();
    CylSurMap[32]=new Plane();
    CylSurMap[33]=new Plane();

    CylSurMap[31]->setSurface(C31);
Russell Taylor's avatar
Russell Taylor committed
889
890
891
892
893
894
895
896
897
898
    CylSurMap[32]->setSurface(C32);
    CylSurMap[33]->setSurface(C33);
    CylSurMap[31]->setName(31);
    CylSurMap[32]->setName(32);
    CylSurMap[33]->setName(33);

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

899
900
901
    Object_sptr retVal = Object_sptr(new Object); 
    retVal->setObject(21,ObjCapCylinder);
    retVal->populate(CylSurMap);
Russell Taylor's avatar
Russell Taylor committed
902
903
904
905

    return retVal;
  }

906
  Object_sptr createSphere()
Russell Taylor's avatar
Russell Taylor committed
907
908
909
910
911
912
913
914
915
916
917
918
  {
    std::string S41="so 4.1";         // Sphere at origin radius 4.1

    // First create some surfaces
    std::map<int,Surface*> SphSurMap;
    SphSurMap[41]=new Sphere();
    SphSurMap[41]->setSurface(S41);
    SphSurMap[41]->setName(41);

    // A sphere 
    std::string ObjSphere="-41" ;

919
920
921
    Object_sptr retVal = Object_sptr(new Object); 
    retVal->setObject(41,ObjSphere);
    retVal->populate(SphSurMap);
Russell Taylor's avatar
Russell Taylor committed
922
923
924
925
926

    return retVal;
  }

  void clearSurfMap()
927
    /**
Russell Taylor's avatar
Russell Taylor committed
928
929
930
931
    Clears the surface map for a new test
    or destruction.
    */
  {
932
    SMap.clear();
Russell Taylor's avatar
Russell Taylor committed
933
934
935
    return;
  }

936
  void createSurfaces(const std::string& desired)
937
    /**
Russell Taylor's avatar
Russell Taylor committed
938
939
940
941
942
943
944
945
946
947
    Creates a list of surfaces for used in the objects
    and populates the MObj layers.
    */
  {
    clearSurfMap();

    // PLANE SURFACES:

    typedef std::pair<int,std::string> SCompT;
    std::vector<SCompT> SurfLine;
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
    if (desired.find("60001") != std::string::npos)
      SurfLine.push_back(SCompT(60001,"px -1"));
    if (desired.find("60002") != std::string::npos)
      SurfLine.push_back(SCompT(60002,"px 1"));
    if (desired.find("60003") != std::string::npos)
      SurfLine.push_back(SCompT(60003,"py -2"));
    if (desired.find("60004") != std::string::npos)
      SurfLine.push_back(SCompT(60004,"py 2"));
    if (desired.find("60005") != std::string::npos)
      SurfLine.push_back(SCompT(60005,"pz -3"));
    if (desired.find("60006") != std::string::npos)
      SurfLine.push_back(SCompT(60006,"pz 3"));

    if (desired.find("80001") != std::string::npos)
      SurfLine.push_back(SCompT(80001,"px 4.5"));
    if (desired.find("80002") != std::string::npos)
      SurfLine.push_back(SCompT(80002,"px 6.5"));

    if (desired.find("71") != std::string::npos)
      SurfLine.push_back(SCompT(71,"so 0.8"));
    if (desired.find("72") != std::string::npos)
      SurfLine.push_back(SCompT(72,"s -0.7 0 0 0.3"));
    if (desired.find("73") != std::string::npos)
      SurfLine.push_back(SCompT(73,"s 0.6 0 0 0.4"));
Russell Taylor's avatar
Russell Taylor committed
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992

    std::vector<SCompT>::const_iterator vc;

    // Note that the testObject now manages the "new Plane"
    Geometry::Surface* A;
    for(vc=SurfLine.begin();vc!=SurfLine.end();vc++)
    {  
      A=Geometry::SurfaceFactory::Instance()->processLine(vc->second);
      if (!A)
      {
        std::cerr<<"Failed to process line "<<vc->second<<std::endl;
        exit(1);
      }
      A->setName(vc->first);
      SMap.insert(STYPE::value_type(vc->first,A));
    }

    return;
  }


993
  Object_sptr createUnitCube()
Russell Taylor's avatar
Russell Taylor committed
994
995
996
997
998
999
1000
  {
    std::string C1="px -0.5";         // cube +/-0.5
    std::string C2="px 0.5";
    std::string C3="py -0.5";
    std::string C4="py 0.5";
    std::string C5="pz -0.5";
    std::string C6="pz 0.5";