1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
|
// ************************************************************************************************
//
// BornAgain: simulate and fit reflection and scattering
//
//! @file Sample/Shape/IShape3D.cpp
//! @brief Implements default methods of interface IShape3D.
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
// ************************************************************************************************
#include "Sample/Shape/IShape3D.h"
#include <cmath>
#include <numbers>
using std::numbers::pi;
// Value of 24 ensures that real points stick out of the convex hull at most
// 1% of the radius
const size_t IShape3D::N_Circle = 24;
std::vector<R3> RectangleVertices(double length, double width, double z)
{
std::vector<R3> result = {{length / 2.0, width / 2.0, z},
{-length / 2.0, width / 2.0, z},
{-length / 2.0, -width / 2.0, z},
{length / 2.0, -width / 2.0, z}};
return result;
}
std::vector<R3> EllipseVerticesZ(double r_x, double r_y, double z)
{
static constexpr double delta_angle = (2 * pi) / IShape3D::N_Circle;
std::vector<R3> result(IShape3D::N_Circle);
for (size_t i = 0; i < IShape3D::N_Circle; ++i) {
double angle = i * delta_angle;
double x = r_x * std::cos(angle);
double y = r_y * std::sin(angle);
result[i] = R3(x, y, z);
}
return result;
}
std::vector<R3> EllipseVerticesXtrunc(double x, double r_y, double r_z, double z_b, double z_t)
{
static constexpr double delta_angle = (2 * pi) / IShape3D::N_Circle;
// for full ellipse
// std::vector<R3> result(IShape3D::N_Circle);
// for (size_t i = 0; i < IShape3D::N_Circle; ++i) {
// double angle = i * delta_angle;
// double y = r_y * std::cos(angle);
// double z = r_z * std::sin(angle) + r_z; // bottom is at 0, top at 2*r_z
// result[i] = R3(x, y, z);
// }
// for truncated ellipse
std::vector<R3> result;
result.reserve(IShape3D::N_Circle + 4);
for (size_t i = 0; i < IShape3D::N_Circle; ++i) {
double angle = i * delta_angle;
double y = r_y * std::cos(angle);
double z = r_z * std::sin(angle) - z_b; // bottom is at 0, top at z_t-z_b
if (0 <= z && z <= (z_t - z_b))
result.emplace_back(x, y, z);
}
double alpha_t = atan(z_t / sqrt(r_z * r_z - z_t * z_t) * r_z / r_y);
double alpha_b = atan(z_b / sqrt(r_z * r_z - z_b * z_b) * r_z / r_y);
double y_t = r_y * std::cos(alpha_t);
double y_b = r_y * std::cos(alpha_b);
result.emplace_back(x, y_t, z_t - z_b);
result.emplace_back(x, -y_t, z_t - z_b);
result.emplace_back(x, y_b, 0);
result.emplace_back(x, -y_b, 0);
return result;
}
|
