Belle II Software  release-08-01-10
eclGeometry.cc
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 #include <ecl/geometry/BelleCrystal.h>
9 #include <ecl/geometry/BelleLathe.h>
10 #include <gtest/gtest.h>
11 #include "ecl/geometry/shapes.h"
12 
13 using namespace std;
14 using namespace Belle2;
15 using namespace ECL;
16 
17 namespace {
18 
20  class BelleCrystalTest : public ::testing::Test {};
21 
23  class BelleLatheTest : public ::testing::Test {
24  protected:
25  double absolute_tolerance = 0.00001;
26  };
27 
28  TEST_F(BelleCrystalTest, BoundingBoxFourSides)
29  {
30  // Create box size 2 2 2
31  G4ThreeVector vertices[] = {G4ThreeVector(-1, -1, -1), G4ThreeVector(1, -1, -1), G4ThreeVector(1, 1, -1), G4ThreeVector(-1, 1, -1),
32  G4ThreeVector(-1, -1, 1), G4ThreeVector(1, -1, 1), G4ThreeVector(1, 1, 1), G4ThreeVector(-1, 1, 1)
33  };
34  BelleCrystal* crystal = new BelleCrystal("solid4", 4, vertices);
35 
36  // Get bounding box
37  G4ThreeVector pMin; G4ThreeVector pMax;
38  crystal->BoundingLimits(pMin, pMax);
39 
40  // Manually calculated bounding box
41  G4ThreeVector pMin_real(-1, -1, -1); G4ThreeVector pMax_real(1, 1, 1);
42 
43  EXPECT_DOUBLE_EQ(pMin_real.x(), pMin.x());
44  EXPECT_DOUBLE_EQ(pMin_real.y(), pMin.y());
45  EXPECT_DOUBLE_EQ(pMin_real.z(), pMin.z());
46 
47  EXPECT_DOUBLE_EQ(pMax_real.x(), pMax.x());
48  EXPECT_DOUBLE_EQ(pMax_real.y(), pMax.y());
49  EXPECT_DOUBLE_EQ(pMax_real.z(), pMax.z());
50  }
51 
52  TEST_F(BelleCrystalTest, BoundingBoxFiveSides)
53  {
54  // Create pentagon side
55  double z = 1;
56  double c1 = 0.25 * (sqrt(5) - 1), c2 = 0.25 * (sqrt(5) + 1);
57  double s1 = 0.25 * sqrt(10 + 2 * sqrt(5)), s2 = 0.25 * sqrt(10 - 2 * sqrt(5));
58  G4ThreeVector vertices[] = {G4ThreeVector(0, 1, -z), G4ThreeVector(-s1, c1, -z), G4ThreeVector(-s2, -c2, -z), G4ThreeVector(s2, -c2, -z), G4ThreeVector(s1, c1, -z),
59  G4ThreeVector(0, 1, z), G4ThreeVector(-s1, c1, z), G4ThreeVector(-s2, -c2, z), G4ThreeVector(s2, -c2, z), G4ThreeVector(s1, c1, z)
60  };
61  BelleCrystal* crystal = new BelleCrystal("solid5", 5, vertices);
62 
63  // Get bounding box
64  G4ThreeVector pMin; G4ThreeVector pMax;
65  crystal->BoundingLimits(pMin, pMax);
66 
67  // Manually calculated bounding box
68  G4ThreeVector pMin_real(-s1, -c2, -z); G4ThreeVector pMax_real(s1, 1, z);
69 
70  EXPECT_DOUBLE_EQ(pMin_real.x(), pMin.x());
71  EXPECT_DOUBLE_EQ(pMin_real.y(), pMin.y());
72  EXPECT_DOUBLE_EQ(pMin_real.z(), pMin.z());
73 
74  EXPECT_DOUBLE_EQ(pMax_real.x(), pMax.x());
75  EXPECT_DOUBLE_EQ(pMax_real.y(), pMax.y());
76  EXPECT_DOUBLE_EQ(pMax_real.z(), pMax.z());
77  }
78 
79  TEST_F(BelleCrystalTest, PentagonSurface)
80  {
81  // Create pentagon side
82  double z = 1;
83  double c1 = 0.25 * (sqrt(5) - 1), c2 = 0.25 * (sqrt(5) + 1);
84  double s1 = 0.25 * sqrt(10 + 2 * sqrt(5)), s2 = 0.25 * sqrt(10 - 2 * sqrt(5));
85  G4ThreeVector p0m(0, 1, -z), p1m(-s1, c1, -z), p2m(-s2, -c2, -z), p3m(s2, -c2, -z), p4m(s1, c1, -z);
86  G4ThreeVector p0p(0, 1, z), p1p(-s1, c1, z), p2p(-s2, -c2, z), p3p(s2, -c2, z), p4p(s1, c1, z);
87  G4ThreeVector vertices[] = {p0m, p1m, p2m, p3m, p4m, p0p, p1p, p2p, p3p, p4p};
88  BelleCrystal* c = new BelleCrystal("solid5", 5, vertices);
89 
90  for (int i = 0; i < 1000 * 100; i++) {
91  G4ThreeVector p = c->GetPointOnSurface();
92  EXPECT_EQ(c->Inside(p), kSurface);
93  }
94  delete c;
95  }
96 
97  TEST_F(BelleCrystalTest, BarrelSurfaces)
98  {
99  vector<shape_t*> cryst = load_shapes("/ecl/data/crystal_shape_barrel.dat");
100  double wrapthickness = 0.17;
101  for (auto it = cryst.begin(); it != cryst.end(); it++) {
102  shape_t* s = *it;
103  std::string prefix("sv_"); prefix += "barrel"; prefix += "_wrap";
104  G4Translate3D tw;
105  G4VSolid* c = s->get_solid(prefix, wrapthickness, tw);
106  for (int i = 0; i < 1000 * 100; i++) {
107  G4ThreeVector p = c->GetPointOnSurface();
108  EInside I = c->Inside(p);
109  EXPECT_EQ(I, kSurface);
110  }
111  delete c;
112  }
113  for (auto it = cryst.begin(); it != cryst.end(); it++) delete *it;
114  }
115 
116  TEST_F(BelleCrystalTest, ForwardSurfaces)
117  {
118  vector<shape_t*> cryst = load_shapes("/ecl/data/crystal_shape_forward.dat");
119  double wrapthickness = 0.17;
120  for (auto it = cryst.begin(); it != cryst.end(); it++) {
121  shape_t* s = *it;
122  std::string prefix("sv_"); prefix += "barrel"; prefix += "_wrap";
123  G4Translate3D tw;
124  G4VSolid* c = s->get_solid(prefix, wrapthickness, tw);
125  for (int i = 0; i < 1000 * 100; i++) {
126  G4ThreeVector p = c->GetPointOnSurface();
127  EInside I = c->Inside(p);
128  EXPECT_EQ(I, kSurface);
129  }
130  delete c;
131  }
132  for (auto it = cryst.begin(); it != cryst.end(); it++) delete *it;
133  }
134 
135  TEST_F(BelleCrystalTest, BackwardSurfaces)
136  {
137  vector<shape_t*> cryst = load_shapes("/ecl/data/crystal_shape_backward.dat");
138  double wrapthickness = 0.17;
139  for (auto it = cryst.begin(); it != cryst.end(); it++) {
140  shape_t* s = *it;
141  std::string prefix("sv_"); prefix += "barrel"; prefix += "_wrap";
142  G4Translate3D tw;
143  G4VSolid* c = s->get_solid(prefix, wrapthickness, tw);
144  for (int i = 0; i < 1000 * 100; i++) {
145  G4ThreeVector p = c->GetPointOnSurface();
146  EInside I = c->Inside(p);
147  EXPECT_EQ(I, kSurface);
148  }
149  delete c;
150  }
151  for (auto it = cryst.begin(); it != cryst.end(); it++) delete *it;
152  }
153 
154  // Bounding limits for half circle starting at pi/2
155  TEST_F(BelleLatheTest, BoundingBoxSect0)
156  {
157  // Create the Belle Lathe
158  zr_t bint[] = {{ -1, 0}, {1, 0}, {1, 1}, { -1, 1}}; //z, r
159  std::vector<zr_t> contourb(bint, bint + sizeof(bint) / sizeof(zr_t));
160  BelleLathe* sect = new BelleLathe("sect", M_PI / 2, M_PI, contourb);
161 
162  // Get bounding box
163  G4ThreeVector pMin; G4ThreeVector pMax;
164  sect->BoundingLimits(pMin, pMax);
165 
166  // Manually calculated bounding box
167  G4ThreeVector pMin_real(-1, -1, -1);
168  G4ThreeVector pMax_real(0, 1, 1);
169 
170  EXPECT_NEAR(pMin_real.x(), pMin.x(), absolute_tolerance);
171  EXPECT_NEAR(pMin_real.y(), pMin.y(), absolute_tolerance);
172  EXPECT_NEAR(pMin_real.z(), pMin.z(), absolute_tolerance);
173 
174  EXPECT_NEAR(pMax_real.x(), pMax.x(), absolute_tolerance);
175  EXPECT_NEAR(pMax_real.y(), pMax.y(), absolute_tolerance);
176  EXPECT_NEAR(pMax_real.z(), pMax.z(), absolute_tolerance);
177  }
178 
179  // Bouding limits for half circle starting at 3pi/2
180  TEST_F(BelleLatheTest, BoundingBoxSect1)
181  {
182  // Create the Belle Lathe
183  zr_t bint[] = {{ -1, 0}, {1, 0}, {1, 1}, { -1, 1}}; //z, r
184  std::vector<zr_t> contourb(bint, bint + sizeof(bint) / sizeof(zr_t));
185  BelleLathe* sect = new BelleLathe("sect", 3 * M_PI / 2, M_PI, contourb);
186 
187  // Get bounding box
188  G4ThreeVector pMin; G4ThreeVector pMax;
189  sect->BoundingLimits(pMin, pMax);
190 
191  // Manually calculated bounding box
192  G4ThreeVector pMin_real(0, -1, -1);
193  G4ThreeVector pMax_real(1, 1, 1);
194 
195  EXPECT_NEAR(pMin_real.x(), pMin.x(), absolute_tolerance);
196  EXPECT_NEAR(pMin_real.y(), pMin.y(), absolute_tolerance);
197  EXPECT_NEAR(pMin_real.z(), pMin.z(), absolute_tolerance);
198 
199  EXPECT_NEAR(pMax_real.x(), pMax.x(), absolute_tolerance);
200  EXPECT_NEAR(pMax_real.y(), pMax.y(), absolute_tolerance);
201  EXPECT_NEAR(pMax_real.z(), pMax.z(), absolute_tolerance);
202  }
203 
204  // Bouding limits for quarter starting at pi
205  TEST_F(BelleLatheTest, BoundingBoxSect2)
206  {
207  // Create the Belle Lathe
208  zr_t bint[] = {{ -1, 0}, {1, 0}, {1, 1}, { -1, 1}}; //z, r
209  std::vector<zr_t> contourb(bint, bint + sizeof(bint) / sizeof(zr_t));
210  BelleLathe* sect = new BelleLathe("sect", M_PI, M_PI / 2, contourb);
211 
212  // Get bounding box
213  G4ThreeVector pMin; G4ThreeVector pMax;
214  sect->BoundingLimits(pMin, pMax);
215 
216  // Manually calculated bounding box
217  G4ThreeVector pMin_real(-1, -1, -1);
218  G4ThreeVector pMax_real(0, 0, 1);
219 
220  EXPECT_NEAR(pMin_real.x(), pMin.x(), absolute_tolerance);
221  EXPECT_NEAR(pMin_real.y(), pMin.y(), absolute_tolerance);
222  EXPECT_NEAR(pMin_real.z(), pMin.z(), absolute_tolerance);
223 
224  EXPECT_NEAR(pMax_real.x(), pMax.x(), absolute_tolerance);
225  EXPECT_NEAR(pMax_real.y(), pMax.y(), absolute_tolerance);
226  EXPECT_NEAR(pMax_real.z(), pMax.z(), absolute_tolerance);
227  }
228 
229  // Bouding limits for three quarter starting at 0
230  TEST_F(BelleLatheTest, BoundingBoxSect3)
231  {
232  // Create the Belle Lathe
233  zr_t bint[] = {{ -1, 0}, {1, 0}, {1, 1}, { -1, 1}}; //z, r
234  std::vector<zr_t> contourb(bint, bint + sizeof(bint) / sizeof(zr_t));
235  BelleLathe* sect = new BelleLathe("sect", 0, 3 * M_PI / 2, contourb);
236 
237  // Get bounding box
238  G4ThreeVector pMin; G4ThreeVector pMax;
239  sect->BoundingLimits(pMin, pMax);
240 
241  // Manually calculated bounding box
242  G4ThreeVector pMin_real(-1, -1, -1);
243  G4ThreeVector pMax_real(1, 1, 1);
244 
245  EXPECT_NEAR(pMin_real.x(), pMin.x(), absolute_tolerance);
246  EXPECT_NEAR(pMin_real.y(), pMin.y(), absolute_tolerance);
247  EXPECT_NEAR(pMin_real.z(), pMin.z(), absolute_tolerance);
248 
249  EXPECT_NEAR(pMax_real.x(), pMax.x(), absolute_tolerance);
250  EXPECT_NEAR(pMax_real.y(), pMax.y(), absolute_tolerance);
251  EXPECT_NEAR(pMax_real.z(), pMax.z(), absolute_tolerance);
252  }
253 
254  // Bouding limits for quarter starting at pi/4 w/ hole
255  TEST_F(BelleLatheTest, BoundingBoxSect4)
256  {
257  // Create the Belle Lathe
258  zr_t bint[] = {{ -1, 0.5}, {1, 0.5}, {1, 1}, { -1, 1}}; //z, r, rmin is 0.5
259  std::vector<zr_t> contourb(bint, bint + sizeof(bint) / sizeof(zr_t));
260  BelleLathe* sect = new BelleLathe("sect", M_PI / 4, M_PI / 2, contourb);
261 
262  // Get bounding box
263  G4ThreeVector pMin; G4ThreeVector pMax;
264  sect->BoundingLimits(pMin, pMax);
265 
266  // Manually calculated bounding box
267  G4ThreeVector pMin_real(cos(3 * M_PI / 4), 0.5 * sin(M_PI / 4), -1);
268  G4ThreeVector pMax_real(cos(M_PI / 4), 1, 1);
269 
270  EXPECT_NEAR(pMin_real.x(), pMin.x(), absolute_tolerance);
271  EXPECT_NEAR(pMin_real.y(), pMin.y(), absolute_tolerance);
272  EXPECT_NEAR(pMin_real.z(), pMin.z(), absolute_tolerance);
273 
274  EXPECT_NEAR(pMax_real.x(), pMax.x(), absolute_tolerance);
275  EXPECT_NEAR(pMax_real.y(), pMax.y(), absolute_tolerance);
276  EXPECT_NEAR(pMax_real.z(), pMax.z(), absolute_tolerance);
277  }
278 }
Belle2::ECL::BelleCrystal
a Belle crystal in Geant4
Definition: BelleCrystal.h:37
Belle2::ECL::BelleCrystal::BoundingLimits
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Two vectors define an axis-parallel bounding box for the shape.
Definition: BelleCrystal.cc:733
Belle2::ECL::BelleLathe
BelleLathe class.
Definition: BelleLathe.h:67
Belle2::ECL::BelleLathe::BoundingLimits
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Two vectors define an axis-parallel bounding box for the shape.
Definition: BelleLathe.cc:1855
Belle2::TEST_F
TEST_F(GlobalLabelTest, LargeNumberOfTimeDependentParameters)
Test large number of time-dep params for registration and retrieval.
Definition: globalLabel.cc:72
Belle2::sqrt
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:17
Belle2::ECL::shape_t
shape
Definition: shapes.h:29
Belle2::ECL::zr_t
simple struct with z and r coordinates
Definition: BelleLathe.h:25