cgal/Kernel_d/test/Kernel_d/afftrafo-test.cpp

354 lines
12 KiB
C++

#include <CGAL/Homogeneous_d.h>
#include <CGAL/Cartesian_d.h>
#include <CGAL/simple_transformations_d.h>
#include <CGAL/test_macros.h>
#include <CGAL/Exact_integer.h>
#include <CGAL/Exact_rational.h>
typedef CGAL::Exact_integer RT_;
typedef CGAL::Exact_rational FT_;
int main()
{
CGAL_TEST_START;
{
typedef CGAL::Homogeneous_d<RT_> Kernel;
typedef Kernel::RT RT;
typedef Kernel::FT FT;
typedef Kernel::LA LA;
typedef CGAL::Point_d<Kernel> Point;
typedef CGAL::Vector_d<Kernel> Vector;
typedef CGAL::Direction_d<Kernel> Direction;
typedef CGAL::Hyperplane_d<Kernel> Hyperplane;
typedef CGAL::Sphere_d<Kernel> Sphere;
typedef CGAL::Segment_d<Kernel> Segment;
typedef CGAL::Ray_d<Kernel> Ray;
typedef CGAL::Line_d<Kernel> Line;
typedef LA::Vector IVector;
typedef LA::Matrix IMatrix;
typedef CGAL::Aff_transformation_d<Kernel> Afftrafo;
{
Vector e1 = Vector(3,Vector::Base_vector(),0);
int IV1[]={1,0,1,1};
int IV2[]={0,0,1};
IVector iv1(IV1,IV1+4), iv2(IV2,IV2+3);
Point p0(3), p1(p0 + e1), p2(3,iv1.begin(),iv1.end()),
p3(3,iv2.begin(),iv2.end(),1), p4(3);
std::vector<Point> V =
make_vector(Point(-1,-1,-1,1),Point(1,-1,-1,1),
Point(-1,1,-1,1),Point(-1,-1,1,1));
Sphere S(3,V.begin(),V.end());
Line l1(p0,p1), l2(p0,Point(1,1,0,1));
Ray r(p0,p1); Segment s(p0,p1);
CGAL_TEST(CGAL::reflect(p1,p0)==Point(-1,0,0,1)){}
CGAL_TEST(CGAL::reflect(p2,l1)==Point(1,0,-1,1)){}
CGAL_TEST(CGAL::reflect(S,p0)==S.opposite()){}
CGAL_TEST(CGAL::reflect(S,l1)==S){}
CGAL_TEST(CGAL::reflect(l1,p0)==l1.opposite()){}
CGAL_TEST(CGAL::reflect(l1,l2)==Line(p0,Point(0,1,0,1))){}
CGAL_TEST(CGAL::reflect(r,p0)==r.opposite()){}
CGAL_TEST(CGAL::reflect(r,l2)==Ray(p0,Point(0,1,0,1))){}
CGAL_TEST(CGAL::reflect(s,p0)==Segment(p0,Point(-1,0,0,1))){}
CGAL_TEST(CGAL::reflect(s,l2)==Segment(p0,Point(0,1,0,1))){}
}
{
int IV[]={2,3,1};
IVector iv(IV,IV+3);
IMatrix M1(3),M2(3);
M1(0,0) = -1; M1(1,1) = -1; M1(2,2) = 1; // rot pi
M2(0,0) = 0; M2(0,1) = -1; M2(0,2) = 1;
M2(1,0) = 1; M2(1,1) = 0; M2(1,2) = 1; M2(2,2) = 1;
Vector rv(5,0);
Direction rd(1,1);
Afftrafo at0(2); // whole zero in 2d
Afftrafo at1(M1); // rot(pi) in 2d
Afftrafo at2(CGAL::SCALING,iv.begin(),iv.end()); // scale by (2,3) in 2d
Afftrafo at3(CGAL::TRANSLATION,rv); // translate by (5,0) in 2d
Afftrafo at4(at1); // rot(pi) in 2d
Afftrafo at5(M2); // rot(-pi/2) + trans(1,1) in 2d
Afftrafo at6(2,CGAL::IDENTITY); // id trafo in 2d
Afftrafo at9(2,CGAL::SCALING,3,1); // scale(3/1) in 2d
Afftrafo at10(2,CGAL::ROTATION,RT(1),RT(0),RT(1)); // rot(pi/2) in 2d
Afftrafo at11(2,CGAL::ROTATION,rd,1,1000); // d2_rot approx
CGAL_TEST(at1==at4){}
CGAL_TEST(at1!=at5){}
CGAL_TEST(at1.dimension()==2){}
CGAL_TEST(at1.matrix() == M1){}
CGAL_IO_TEST(at4,at0,CGAL::IO::ASCII);
Afftrafo at12 = at2 * at2.inverse();
Point p(1,1);
CGAL_TEST(p.transform(at12)==p){}
CGAL_TEST(p.transform(at1)==Point(-1,-1)){}
CGAL_TEST(p.transform(at2)==Point(2,3)){}
CGAL_TEST(p.transform(at3)==Point(6,1)){}
CGAL_TEST(p.transform(at4)==Point(-1,-1)){}
CGAL_TEST(p.transform(at5)==Point(0,2)){}
CGAL_TEST(p.transform(at6)==p){}
CGAL_TEST(p.transform(at9)==Point(3,3)){}
CGAL_TEST(p.transform(at10)==Point(-1,1)){}
Vector v_approx =
Point(1,0).transform(at11) - CGAL::ORIGIN;
/* we have v_approx = (x,y) building an angle alpha with
the x-axis then sin alpha = y/length(v_approx) and we
want to know that sin alpha <= num/den of at11.
we started with a direction rd(1,1) with angle pi/4
and thus sin pi/4 = 1/ sqrt(2)
Thus we check that the squared difference of the sine
values is bounded by sqr(num/den):
sqr(y/sqrt((sqr(x) + sqr(y))) - 1/sqrt(2)) < sqr(num/den)
sqr(a - b) < sqr(c)
<=> sqr(a) + sqr(b) - sqr(c) < 2ab
<=> sqr( sqr(a) + sqr(b) - sqr(c) ) < 4 sqr(a) sqr(b)
where sqr(b) = 1/2
<=> sqr( sqr(a) + sqr(b) - sqr(c) ) < 2 sqr(a)
in at11 we had num = 1 and den = 1000
*/
FT sqr_sin_alpha = CGAL_NTS square(v_approx.y())/
v_approx.squared_length();
// sqr(a)
FT sqr_sin_piforth(1,2);
// sqr(b)
FT sqr_bound = Kernel::make_FT(RT(1),CGAL_NTS square(RT(1000)));
// sqr(c)
FT sqr_sum = sqr_sin_alpha + sqr_sin_piforth - sqr_bound;
CGAL_TEST( (sqr_sum*sqr_sum) < (FT(2)*sqr_sin_alpha) ){}
// and analogously:
FT sqr_cos_alpha = CGAL_NTS square(v_approx.x())/
v_approx.squared_length();
FT sqr_cos_piforth(1,2);
sqr_sum = sqr_cos_alpha + sqr_cos_piforth - sqr_bound;
CGAL_TEST((sqr_sum*sqr_sum) < FT(2)*sqr_cos_alpha){}
Vector v(2,2);
CGAL_TEST(v.transform(at9)==3*v){}
CGAL_TEST(v.transform(at3)==v){}
CGAL_TEST(v.transform(at1)==-v){}
Direction dir(5,1);
CGAL_TEST(dir.transform(at9)==dir){}
CGAL_TEST(dir.transform(at3)==dir){}
CGAL_TEST(dir.transform(at1)==-dir){}
// 3 dim hyperplane:
Afftrafo at20 = Afftrafo(CGAL::TRANSLATION,Vector(5,0,0,1));
Afftrafo at21 = Afftrafo(3,CGAL::ROTATION,RT(1),RT(0),RT(1),0,2);
// pi/2-rotation in x-z-plane in 3-space
Afftrafo at22 = at21*at20; // composition of both
Point p1(1,0,0,1), p2(0,1,0,1), p3(0,0,1,1), o(3);
CGAL_TEST(p1.transform(at22)==p1.transform(at20).transform(at21)){}
std::vector<Point> U = make_vector(p1,p2,p3);
Hyperplane hyp(U.begin(),U.end(),o,CGAL::ON_NEGATIVE_SIDE);
Point pp1 = p1.transform(at22),
pp2 = p2.transform(at22),
pp3 = p3.transform(at22),
o_trafo = o.transform(at22);
std::vector<Point> V = make_vector(pp1,pp2,pp3);
Hyperplane hyp_trafo(V.begin(),V.end(),o_trafo,
CGAL::ON_NEGATIVE_SIDE);
CGAL_TEST(hyp.transform(at22)==hyp_trafo){}
Segment s(p1,p2);
CGAL_TEST(s.transform(at22)==Segment(pp1,pp2)){}
Ray r(p1,p2);
CGAL_TEST(r.transform(at22)==Ray(pp1,pp2)){}
Line l(p1,p2);
CGAL_TEST(l.transform(at22)==Line(pp1,pp2)){}
std::vector<Point> W = make_vector(p1,p2,p3,o);
Sphere S(3,W.begin(),W.end());
W = make_vector(pp1,pp2,pp3,o_trafo);
CGAL_TEST(S.transform(at22)==Sphere(3,W.begin(),W.end())){}
}
}
{
typedef CGAL::Cartesian_d<FT_> Kernel;
typedef Kernel::RT RT;
typedef Kernel::FT FT;
typedef Kernel::LA LA;
typedef CGAL::Point_d<Kernel> Point;
typedef CGAL::Vector_d<Kernel> Vector;
typedef CGAL::Direction_d<Kernel> Direction;
typedef CGAL::Hyperplane_d<Kernel> Hyperplane;
typedef CGAL::Sphere_d<Kernel> Sphere;
typedef CGAL::Segment_d<Kernel> Segment;
typedef CGAL::Ray_d<Kernel> Ray;
typedef CGAL::Line_d<Kernel> Line;
typedef LA::Vector IVector;
typedef LA::Matrix IMatrix;
typedef CGAL::Aff_transformation_d<Kernel> Afftrafo;
{
Vector e1 = Vector(3,Vector::Base_vector(),0);
int IV1[]={1,0,1,1};
int IV2[]={0,0,1};
IVector iv1(IV1,IV1+4), iv2(IV2,IV2+3);
Point p0(3), p1(p0 + e1), p2(3,iv1.begin(),iv1.end()),
p3(3,iv2.begin(),iv2.end(),1), p4(3);
std::vector<Point> V =
make_vector(Point(-1,-1,-1,1),Point(1,-1,-1,1),
Point(-1,1,-1,1),Point(-1,-1,1,1));
Sphere S(3,V.begin(),V.end());
Line l1(p0,p1), l2(p0,Point(1,1,0,1));
Ray r(p0,p1); Segment s(p0,p1);
CGAL_TEST(CGAL::reflect(p1,p0)==Point(-1,0,0,1)){}
CGAL_TEST(CGAL::reflect(p2,l1)==Point(1,0,-1,1)){}
CGAL_TEST(CGAL::reflect(S,p0)==S.opposite()){}
CGAL_TEST(CGAL::reflect(S,l1)==S){}
CGAL_TEST(CGAL::reflect(l1,p0)==l1.opposite()){}
CGAL_TEST(CGAL::reflect(l1,l2)==Line(p0,Point(0,1,0,1))){}
CGAL_TEST(CGAL::reflect(r,p0)==r.opposite()){}
CGAL_TEST(CGAL::reflect(r,l2)==Ray(p0,Point(0,1,0,1))){}
CGAL_TEST(CGAL::reflect(s,p0)==Segment(p0,Point(-1,0,0,1))){}
CGAL_TEST(CGAL::reflect(s,l2)==Segment(p0,Point(0,1,0,1))){}
}
{
int IV[]={2,3,1};
IVector iv(IV,IV+3);
IMatrix M1(3),M2(3);
M1(0,0) = -1; M1(1,1) = -1; M1(2,2) = 1; // rot pi
M2(0,0) = 0; M2(0,1) = -1; M2(0,2) = 1;
M2(1,0) = 1; M2(1,1) = 0; M2(1,2) = 1; M2(2,2) = 1;
Vector rv(5,0);
Direction rd(1,1);
Afftrafo at0(2); // whole zero in 2d
Afftrafo at1(M1); // rot(pi) in 2d
Afftrafo at2(CGAL::SCALING,iv.begin(),iv.end()); // scale by (2,3) in 2d
Afftrafo at3(CGAL::TRANSLATION,rv); // translate by (5,0) in 2d
Afftrafo at4(at1); // rot(pi) in 2d
Afftrafo at5(M2); // rot(-pi/2) + trans(1,1) in 2d
Afftrafo at6(2,CGAL::IDENTITY); // id trafo in 2d
Afftrafo at9(2,CGAL::SCALING,3,1); // scale(3/1) in 2d
Afftrafo at10(2,CGAL::ROTATION,RT(1),RT(0),RT(1)); // rot(pi/2) in 2d
Afftrafo at11(2,CGAL::ROTATION,rd,1,1000); // d2_rot approx
CGAL_TEST(at1==at4){}
CGAL_TEST(at1!=at5){}
CGAL_TEST(at1.dimension()==2){}
CGAL_TEST(at1.matrix() == M1){}
CGAL_IO_TEST(at4,at0,CGAL::IO::ASCII);
Afftrafo at12 = at2 * at2.inverse();
Point p(1,1);
CGAL_TEST(p.transform(at12)==p){}
CGAL_TEST(p.transform(at1)==Point(-1,-1)){}
CGAL_TEST(p.transform(at2)==Point(2,3)){}
CGAL_TEST(p.transform(at3)==Point(6,1)){}
CGAL_TEST(p.transform(at4)==Point(-1,-1)){}
CGAL_TEST(p.transform(at5)==Point(0,2)){}
CGAL_TEST(p.transform(at6)==p){}
CGAL_TEST(p.transform(at9)==Point(3,3)){}
CGAL_TEST(p.transform(at10)==Point(-1,1)){}
Vector v_approx =
Point(1,0).transform(at11) - CGAL::ORIGIN;
/* we have v_approx = (x,y) building an angle alpha with
the x-axis then sin alpha = y/length(v_approx) and we
want to know that sin alpha <= num/den of at11.
we started with a direction rd(1,1) with angle pi/4
and thus sin pi/4 = 1/ sqrt(2)
Thus we check that the squared difference of the sine
values is bounded by sqr(num/den):
sqr(y/sqrt((sqr(x) + sqr(y))) - 1/sqrt(2)) < sqr(num/den)
sqr(a - b) < sqr(c)
<=> sqr(a) + sqr(b) - sqr(c) < 2ab
<=> sqr( sqr(a) + sqr(b) - sqr(c) ) < 4 sqr(a) sqr(b)
where sqr(b) = 1/2
<=> sqr( sqr(a) + sqr(b) - sqr(c) ) < 2 sqr(a)
in at11 we had num = 1 and den = 1000
*/
FT sqr_sin_alpha = CGAL_NTS square(v_approx.y())/
v_approx.squared_length();
// sqr(a)
FT sqr_sin_piforth(1,2);
// sqr(b)
FT sqr_bound = Kernel::make_FT(RT(1),CGAL_NTS square(RT(1000)));
// sqr(c)
FT sqr_sum = sqr_sin_alpha + sqr_sin_piforth - sqr_bound;
CGAL_TEST( (sqr_sum*sqr_sum) < (FT(2)*sqr_sin_alpha) ){}
// and analogously:
FT sqr_cos_alpha = CGAL_NTS square(v_approx.x())/
v_approx.squared_length();
FT sqr_cos_piforth(1,2);
sqr_sum = sqr_cos_alpha + sqr_cos_piforth - sqr_bound;
CGAL_TEST((sqr_sum*sqr_sum) < FT(2)*sqr_cos_alpha){}
Vector v(2,2);
CGAL_TEST(v.transform(at9)==3*v){}
CGAL_TEST(v.transform(at3)==v){}
CGAL_TEST(v.transform(at1)==-v){}
Direction dir(5,1);
CGAL_TEST(dir.transform(at9)==dir){}
CGAL_TEST(dir.transform(at3)==dir){}
CGAL_TEST(dir.transform(at1)==-dir){}
// 3 dim hyperplane:
Afftrafo at20 = Afftrafo(CGAL::TRANSLATION,Vector(5,0,0,1));
Afftrafo at21 = Afftrafo(3,CGAL::ROTATION,RT(1),RT(0),RT(1),0,2);
// pi/2-rotation in x-z-plane in 3-space
Afftrafo at22 = at21*at20; // composition of both
Point p1(1,0,0,1), p2(0,1,0,1), p3(0,0,1,1), o(3);
CGAL_TEST(p1.transform(at22)==p1.transform(at20).transform(at21)){}
std::vector<Point> U = make_vector(p1,p2,p3);
Hyperplane hyp(U.begin(),U.end(),o,CGAL::ON_NEGATIVE_SIDE);
Point pp1 = p1.transform(at22),
pp2 = p2.transform(at22),
pp3 = p3.transform(at22),
o_trafo = o.transform(at22);
std::vector<Point> V = make_vector(pp1,pp2,pp3);
Hyperplane hyp_trafo(V.begin(),V.end(),o_trafo,
CGAL::ON_NEGATIVE_SIDE);
CGAL_TEST(hyp.transform(at22)==hyp_trafo){}
Segment s(p1,p2);
CGAL_TEST(s.transform(at22)==Segment(pp1,pp2)){}
Ray r(p1,p2);
CGAL_TEST(r.transform(at22)==Ray(pp1,pp2)){}
Line l(p1,p2);
CGAL_TEST(l.transform(at22)==Line(pp1,pp2)){}
std::vector<Point> W = make_vector(p1,p2,p3,o);
Sphere S(3,W.begin(),W.end());
W = make_vector(pp1,pp2,pp3,o_trafo);
CGAL_TEST(S.transform(at22)==Sphere(3,W.begin(),W.end())){}
}
}
CGAL_TEST_END;
}