Speed up transformation by calculating the transposed inverse matrix once

This commit is contained in:
Giles Bathgate 2021-04-06 19:11:04 +01:00
parent d979e070be
commit 6ef65b8d59
9 changed files with 63 additions and 21 deletions

View File

@ -145,11 +145,16 @@ public:
transform(const Plane_3& p) const
{ return this->Ptr()->transform(p); }
Plane_3
transform(const Plane_3& p, bool is_even, const Aff_transformation_3& transposed_inverse) const
{ return this->Ptr()->transform(p, is_even, transposed_inverse); }
Plane_3
operator()(const Plane_3& p) const
{ return transform(p); } // FIXME : not compiled by the test-suite !
Aff_transformation_3 inverse() const { return this->Ptr()->inverse(); }
Aff_transformation_3 transpose() const { return this->Ptr()->transpose(); }
bool is_even() const { return this->Ptr()->is_even(); }
bool is_odd() const { return ! (this->Ptr()->is_even()); }
@ -183,8 +188,6 @@ public:
return !(*this == t);
}
protected:
Aff_transformation_3 transpose() const { return this->Ptr()->transpose(); }
};

View File

@ -40,6 +40,7 @@ public:
virtual Vector_3 transform(const Vector_3 &v) const = 0;
virtual Direction_3 transform(const Direction_3 &d) const = 0;
virtual Plane_3 transform(const Plane_3& p) const = 0;
virtual Plane_3 transform(const Plane_3& p, bool is_even, const Aff_transformation_3& transposed_inverse) const = 0;
virtual Aff_transformation_3 operator*(
const Aff_transformation_rep_baseC3 &t) const = 0;
@ -133,14 +134,19 @@ public:
t31 * v.x() + t32 * v.y() + t33 * v.z());
}
virtual Plane_3 transform(const Plane_3& p) const
virtual Plane_3 transform(const Plane_3& p, bool is_even, const Aff_transformation_3& transposed_inverse) const
{
if (is_even())
if (is_even)
return Plane_3(transform(p.point()),
transpose().inverse().transform(p.orthogonal_direction()));
transposed_inverse.transform(p.orthogonal_direction()));
else
return Plane_3(transform(p.point()),
- transpose().inverse().transform(p.orthogonal_direction()));
- transposed_inverse.transform(p.orthogonal_direction()));
}
virtual Plane_3 transform(const Plane_3& p) const
{
return transform(p, is_even(), transpose().inverse());
}

View File

@ -60,6 +60,11 @@ public:
return d;
}
virtual Plane_3 transform(const Plane_3 &p, bool, const Aff_transformation_3&) const
{
return transform(p);
}
virtual Plane_3 transform(const Plane_3 &p) const
{
// direction ( which is (p.a(), p.b(), p.c())) does not change

View File

@ -57,6 +57,11 @@ public:
return d;
}
virtual Plane_3 transform(const Plane_3 &p, bool, const Aff_transformation_3&) const
{
return transform(p);
}
virtual Plane_3 transform(const Plane_3 &p) const
{
// direction ( which is (p.a(), p.b(), p.c())) does not change

View File

@ -102,6 +102,11 @@ public:
return t.transform(*this);
}
Plane_3 transform(const Aff_transformation_3 &t, bool is_even, const Aff_transformation_3 &transposed_inverse) const
{
return t.transform(*this, is_even, transposed_inverse);
}
Plane_3 opposite() const
{
return R().construct_opposite_plane_3_object()(*this);

View File

@ -136,15 +136,19 @@ public:
Halffacet_const_iterator facets_begin() { return facet_list.begin(); }
Halffacet_const_iterator facets_end() { return facet_list.end(); }
void transform(const Aff_transformation_3& t) {
void transform(const Aff_transformation_3& t, bool is_even, const Aff_transformation_3& transposed_inverse) {
if(left_node != nullptr) {
CGAL_assertion(right_node != nullptr);
left_node->transform(t);
right_node->transform(t);
splitting_plane = splitting_plane.transform(t);
left_node->transform(t, is_even, transposed_inverse);
right_node->transform(t, is_even, transposed_inverse);
splitting_plane = splitting_plane.transform(t, is_even, transposed_inverse);
}
}
void transform(const Aff_transformation_3& t) {
transform(t, t.is_even(), t.transposed_inverse());
}
void add_facet(Halffacet_handle f, int depth) {
if(left_node == nullptr) {
facet_list.push_back(f);
@ -455,17 +459,20 @@ public:
V.post_visit(current);
}
void transform(const Aff_transformation_3& t) {
void transform(const Aff_transformation_3& t, bool is_even, const Aff_transformation_3& transposed_inverse) {
if(root == nullptr){
return;
}
root->transform(t);
root->transform(t, is_even, transposed_inverse);
BBox_updater bbup;
visit_k3tree(root, bbup);
bounding_box = bbup.box();
}
void transform(const Aff_transformation_3& t) {
transform(t, t.is_even(), t.transpose().inverse());
}
#ifdef CODE_DOES_NOT_WORK_WITH_BOTH_KERNELS_AT_THE_SAME_TIME
template <typename T>

View File

@ -111,6 +111,8 @@ public:
virtual void transform(const Aff_transformation_3& t) = 0;
virtual void transform(const Aff_transformation_3& t, bool is_even, const Aff_transformation_3& transposed_inverse) = 0;
virtual void add_facet(Halffacet_handle) {}
virtual void add_edge(Halfedge_handle) {}
@ -217,6 +219,10 @@ public:
candidate_provider->transform(t);
}
virtual void transform(const Aff_transformation_3& t, bool is_even, const Aff_transformation_3& transposed_inverse) {
candidate_provider->transform(t, is_even, transposed_inverse);
}
virtual ~SNC_point_locator_by_spatial_subdivision() noexcept(!CGAL_ASSERTIONS_ENABLED)
{
CGAL_destructor_warning(initialized ||

View File

@ -1680,7 +1680,8 @@ protected:
aff.hm(1,0), aff.hm(1,1), aff.hm(1,2),
aff.hm(2,0), aff.hm(2,1), aff.hm(2,2),
aff.hm(3,3));
const Aff_transformation_3& transposed_inverse = aff.transpose().inverse();
bool is_even = aff.is_even();
SNC_constructor cstr(snc());
std::list<Vertex_handle> vertex_list;
@ -1712,7 +1713,7 @@ protected:
vi->point() = vi->point().transform( aff);
if(! translate){
SM_decorator sdeco(&*vi);
sdeco.transform( linear);
sdeco.transform( linear, is_even, transposed_inverse);
}
}
}
@ -1722,7 +1723,7 @@ protected:
CGAL_forall_facets(fi, snc()) {
if(!is_standard(fi) || is_bounded(fi)) continue;
Plane_3 pt = fi->plane();
pt = pt.transform(aff);
pt = pt.transform(aff, is_even, transposed_inverse);
std::list<Point_3> points(Infi_box::find_points_of_box_with_plane(cstr,pt));
std::list<Vertex_handle> newVertices;
newVertices = Infi_box::create_vertices_on_infibox(cstr,
@ -1832,7 +1833,7 @@ protected:
Halffacet_iterator fi;
CGAL_forall_halffacets(fi,snc()) {
if(is_standard(fi) || ninety) {
fi->plane() = fi->plane().transform( aff);
fi->plane() = fi->plane().transform( aff, is_even, transposed_inverse);
}
}
@ -1853,7 +1854,7 @@ protected:
pl()->initialize(&snc());
delete old_pl;
}
else pl()->transform(aff);
else pl()->transform(aff, is_even, transposed_inverse);
}
SNC_constructor C(snc());

View File

@ -757,7 +757,7 @@ The macros are then |CGAL_forall_svertices_of(v,V)|,
|CGAL_forall_shalfedges_of(e,V)|, |CGAL_forall_sedges_of(e,V)|,
|CGAL_forall_sfaces_of(f,V)|, |CGAL_forall_sface_cycles_of(fc,F)|.}*/
void transform( const Aff_transformation_3& linear) {
void transform( const Aff_transformation_3& linear, bool is_even, const Aff_transformation_3& transposed_inverse) {
// CGAL_NEF_TRACEN("transform sphere map of vertex" << center_vertex()->point());
// The affine transformation is linear, i.e., no translation part.
CGAL_precondition( linear.hm(0,3) == 0 &&
@ -768,15 +768,19 @@ void transform( const Aff_transformation_3& linear) {
for (SVertex_iterator i = svertices_begin(); i != svertices_end(); ++i)
i->point() = normalized(Sphere_point( i->point().transform( linear)));
for (SHalfedge_iterator i = shalfedges_begin(); i !=shalfedges_end(); ++i)
i->circle() = Sphere_circle( i->circle().transform( linear));
i->circle() = Sphere_circle( i->circle().transform( linear, is_even, transposed_inverse));
if ( has_shalfloop()) {
shalfloop()->circle() = Sphere_circle(shalfloop()->circle()
.transform( linear));
.transform( linear, is_even, transposed_inverse));
shalfloop()->twin()->circle()
= Sphere_circle(shalfloop()->twin()->circle().transform( linear));
= Sphere_circle(shalfloop()->twin()->circle().transform( linear, is_even, transposed_inverse));
}
}
void transform( const Aff_transformation_3& linear) {
return transform(linear, linear.is_even(), linear.transpose().inverse());
}
void extract_complement() {
SVertex_handle sv;