Belle II Software development
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
13using namespace std;
14using namespace Belle2;
15using namespace ECL;
16
17namespace {
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::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
std
STL namespace.
Belle2::ECL::shape_t
shape
Definition: shapes.h:29
Belle2::ECL::zr_t
simple struct with z and r coordinates
Definition: BelleLathe.h:25