OPAL (Object Oriented Parallel Accelerator Library) 2022.1
OPAL
Ellipse.cpp
Go to the documentation of this file.
4#include "Physics/Physics.h"
5#include "Physics/Units.h"
6
7#include <boost/regex.hpp>
8
9namespace mslang {
10 void Ellipse::print(int indentwidth) {
11 std::string indent(indentwidth, ' ');
12 std::string indent2(indentwidth + 8, ' ');
13 Vector_t origin = trafo_m.getOrigin();
14 double angle = trafo_m.getAngle() * Units::rad2deg;
15 std::cout << indent << "ellipse, \n"
16 << indent2 << "w: " << width_m << ", \n"
17 << indent2 << "h: " << height_m << ", \n"
18 << indent2 << "origin: " << origin[0] << ", " << origin[1] << ",\n"
19 << indent2 << "angle: " << angle << "\n"
20 << indent2 << std::setw(14) << trafo_m(0, 0) << std::setw(14) << trafo_m(0, 1) << std::setw(14) << trafo_m(0, 2) << "\n"
21 << indent2 << std::setw(14) << trafo_m(1, 0) << std::setw(14) << trafo_m(1, 1) << std::setw(14) << trafo_m(1, 2) << "\n"
22 << indent2 << std::setw(14) << trafo_m(2, 0) << std::setw(14) << trafo_m(2, 1) << std::setw(14) << trafo_m(2, 2)
23 << std::endl;
24 }
25
26 void Ellipse::writeGnuplot(std::ofstream &out) const {
27 const unsigned int N = 101;
28 const double dp = Physics::two_pi / (N - 1);
29 const unsigned int colwidth = out.precision() + 8;
30
31 double phi = 0;
32 for (unsigned int i = 0; i < N; ++ i, phi += dp) {
33 Vector_t pt(0.0);
34 pt[0] = std::copysign(sqrt(std::pow(height_m * width_m * 0.25, 2) /
35 (std::pow(height_m * 0.5, 2) +
36 std::pow(width_m * 0.5 * tan(phi), 2))),
37 cos(phi));
38 pt[1] = pt[0] * tan(phi);
39 pt = trafo_m.transformFrom(pt);
40
41 out << std::setw(colwidth) << pt[0]
42 << std::setw(colwidth) << pt[1]
43 << std::endl;
44 }
45 out << std::endl;
46
47 for (auto item: divisor_m) {
48 item->writeGnuplot(out);
49 }
50
51 // bb_m.writeGnuplot(out);
52 }
53
54 void Ellipse::apply(std::vector<std::shared_ptr<Base> > &bfuncs) {
55 bfuncs.emplace_back(this->clone());
56 }
57
58 std::shared_ptr<Base> Ellipse::clone() const{
59 std::shared_ptr<Ellipse> elps(new Ellipse);
60 elps->width_m = width_m;
61 elps->height_m = height_m;
62 elps->trafo_m = trafo_m;
63 elps->bb_m = bb_m;
64
65 for (auto item: divisor_m) {
66 elps->divisor_m.emplace_back(item->clone());
67 }
68
69 return std::static_pointer_cast<Base>(elps);
70 }
71
73 Vector_t llc(0.0), urc(0.0);
74 const Vector_t e_x(1.0, 0.0, 0.0), e_y(0.0, 1.0, 0.0);
75 const Vector_t center = trafo_m.transformFrom(Vector_t(0.0));
76 const Vector_t e_xp = trafo_m.transformFrom(e_x) - center;
77 const Vector_t e_yp = trafo_m.transformFrom(e_y) - center;
78 const double &M11 = e_xp[0];
79 const double &M12 = e_yp[0];
80 const double &M21 = e_xp[1];
81 const double &M22 = e_yp[1];
82
83 double t = atan2(height_m * M12, width_m * M11);
84 double halfwidth = 0.5 * (M11 * width_m * cos(t) +
85 M12 * height_m * sin(t));
86 llc[0] = center[0] - std::abs(halfwidth);
87 urc[0] = center[0] + std::abs(halfwidth);
88
89 t = atan2(height_m * M22, width_m * M21);
90
91 double halfheight = 0.5 * (M21 * width_m * cos(t) +
92 M22 * height_m * sin(t));
93
94 llc[1] = center[1] - std::abs(halfheight);
95 urc[1] = center[1] + std::abs(halfheight);
96
97 bb_m = BoundingBox2D(llc, urc);
98
99 for (auto item: divisor_m) {
100 item->computeBoundingBox();
101 }
102 }
103
104 bool Ellipse::isInside(const Vector_t &R) const {
105 if (!bb_m.isInside(R)) return false;
106
108 if (4 * (std::pow(X[0] / width_m, 2) + std::pow(X[1] / height_m, 2)) <= 1) {
109
110 for (auto item: divisor_m) {
111 if (item->isInside(R)) return false;
112 }
113
114 return true;
115 }
116
117 return false;
118 }
119
121
122 ArgumentExtractor arguments(std::string(it, end));
123 double width, height;
124 try {
125 width = parseMathExpression(arguments.get(0));
126 height = parseMathExpression(arguments.get(1));
127 } catch (std::runtime_error &e) {
128 std::cout << e.what() << std::endl;
129 return false;
130 }
131
132 Ellipse *elps = static_cast<Ellipse*>(fun);
133 elps->width_m = width;
134 elps->height_m = height;
135
136 if (elps->width_m < 0.0) {
137 std::cout << "Ellipse: a negative width provided '"
138 << arguments.get(0) << " = " << elps->width_m << "'"
139 << std::endl;
140 return false;
141 }
142 if (elps->height_m < 0.0) {
143 std::cout << "Ellipse: a negative height provided '"
144 << arguments.get(1) << " = " << elps->height_m << "'"
145 << std::endl;
146 return false;
147 }
148
149 it += (arguments.getLengthConsumed() + 1);
150
151 return true;
152 }
153}
PartBunchBase< T, Dim >::ConstIterator end(PartBunchBase< T, Dim > const &bunch)
Tps< T > cos(const Tps< T > &x)
Cosine.
Definition: TpsMath.h:129
Tps< T > pow(const Tps< T > &x, int y)
Integer power.
Definition: TpsMath.h:76
Tps< T > tan(const Tps< T > &x)
Tangent.
Definition: TpsMath.h:147
Tps< T > sin(const Tps< T > &x)
Sine.
Definition: TpsMath.h:111
Tps< T > sqrt(const Tps< T > &x)
Square root.
Definition: TpsMath.h:91
#define X(arg)
Definition: fftpack.cpp:112
PETE_TBTree< FnCopysign, PETE_Scalar< Vektor< T1, Dim > >, typename T2::PETE_Expr_t > copysign(const Vektor< T1, Dim > &l, const PETE_Expr< T2 > &r)
PETE_TUTree< FnAbs, typename T::PETE_Expr_t > abs(const PETE_Expr< T > &l)
PETE_TBTree< FnArcTan2, PETE_Scalar< Vektor< T1, Dim > >, typename T2::PETE_Expr_t > atan2(const Vektor< T1, Dim > &l, const PETE_Expr< T2 > &r)
Inform & endl(Inform &inf)
Definition: Inform.cpp:42
constexpr double two_pi
The value of.
Definition: Physics.h:33
constexpr double e
The value of.
Definition: Physics.h:39
constexpr double rad2deg
Definition: Units.h:146
double parseMathExpression(const std::string &str)
Definition: matheval.cpp:4
std::string::iterator iterator
Definition: MSLang.h:16
std::vector< std::shared_ptr< Base > > divisor_m
Definition: MSLang.h:43
AffineTransformation trafo_m
Definition: MSLang.h:41
BoundingBox2D bb_m
Definition: MSLang.h:42
Vector_t transformTo(const Vector_t &v) const
Vector_t transformFrom(const Vector_t &v) const
std::string get(unsigned int i) const
unsigned int getLengthConsumed() const
bool isInside(const Vector_t &X) const
virtual void computeBoundingBox()
Definition: Ellipse.cpp:72
virtual void print(int indentwidth)
Definition: Ellipse.cpp:10
virtual bool isInside(const Vector_t &R) const
Definition: Ellipse.cpp:104
double width_m
Definition: Ellipse.h:8
virtual std::shared_ptr< Base > clone() const
Definition: Ellipse.cpp:58
virtual void writeGnuplot(std::ofstream &out) const
Definition: Ellipse.cpp:26
virtual void apply(std::vector< std::shared_ptr< Base > > &bfuncs)
Definition: Ellipse.cpp:54
double height_m
Definition: Ellipse.h:9
static bool parse_detail(iterator &it, const iterator &end, Function *fun)
Definition: Ellipse.cpp:120
Vektor< double, 3 > Vector_t
Definition: Vektor.h:6