Commit 885ed95b authored by Julius Nehring-Wirxel's avatar Julius Nehring-Wirxel
Browse files

Merge origin/develop into jn/develop

parents 7d442a02 e208eb52
......@@ -19,10 +19,14 @@
#include <typed-geometry/functions/objects/intersection.hh>
#include <typed-geometry/functions/objects/normal.hh>
#include <typed-geometry/functions/objects/perimeter.hh>
#include <typed-geometry/functions/objects/plane.hh>
#include <typed-geometry/functions/objects/project.hh>
#include <typed-geometry/functions/objects/rasterize.hh>
#include <typed-geometry/functions/objects/segmentize.hh>
#include <typed-geometry/functions/objects/size.hh>
#include <typed-geometry/functions/objects/tangent.hh>
#include <typed-geometry/functions/objects/triangle.hh>
#include <typed-geometry/functions/objects/triangulate.hh>
#include <typed-geometry/functions/objects/triangulation.hh>
#include <typed-geometry/functions/objects/vertices.hh>
#include <typed-geometry/functions/objects/volume.hh>
......@@ -5,9 +5,11 @@
#include <typed-geometry/types/objects/capsule.hh>
#include <typed-geometry/types/objects/cylinder.hh>
#include <typed-geometry/types/objects/ellipse.hh>
#include <typed-geometry/types/objects/halfspace.hh>
#include <typed-geometry/types/objects/hemisphere.hh>
#include <typed-geometry/types/objects/inf_cone.hh>
#include <typed-geometry/types/objects/inf_cylinder.hh>
#include <typed-geometry/types/objects/plane.hh>
#include <typed-geometry/types/objects/pyramid.hh>
#include <typed-geometry/types/objects/sphere.hh>
......@@ -75,6 +77,11 @@ template <class ScalarT, class TraitsT>
{
return {v.center, v.radius, v.normal};
}
template <int D, class ScalarT>
[[nodiscard]] constexpr plane<D, ScalarT> boundary_of(halfspace<D, ScalarT> const& v)
{
return {v.normal, v.dis};
}
// === no caps versions ===
......@@ -134,6 +141,7 @@ template <class ScalarT, class TraitsT>
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr pos<1, ScalarT> caps_of(hemisphere<1, ScalarT, TraitsT> const& v)
{
static_assert(!std::is_same_v<TraitsT, boundary_no_caps_tag> && "1D hemisphere_boundary_no_caps does not have any caps");
return v.center;
}
......
......@@ -2,13 +2,9 @@
#include <typed-geometry/detail/utility.hh>
#include <typed-geometry/types/objects/line.hh>
#include <typed-geometry/types/objects/plane.hh>
#include <typed-geometry/types/objects/segment.hh>
#include <typed-geometry/types/pos.hh>
#include <typed-geometry/types/quadric.hh>
#include <typed-geometry/functions/basic/mix.hh>
#include "coordinates.hh"
#include "project.hh"
......@@ -45,7 +41,7 @@ template <class ScalarT>
auto d0d1 = dot(l0.dir, l1.dir);
auto b0 = dot(l1.pos - l0.pos, l0.dir);
auto b1 = dot(l1.pos - l0.pos, l1.dir);
auto [t0, t1] = inverse(tg::mat<2, 2, ScalarT>::from_cols({1, d0d1}, {-d0d1, -1})) * tg::vec2(b0, b1);
auto [t0, t1] = inverse(mat<2, 2, ScalarT>::from_cols({1, d0d1}, {-d0d1, -1})) * vec2(b0, b1);
return {t0, t1};
}
......@@ -56,6 +52,17 @@ template <class ScalarT>
return {l0[t0], l1[t1]};
}
template <class ScalarT>
[[nodiscard]] constexpr pair<ScalarT, ScalarT> closest_points_parameters(segment<3, ScalarT> const& s, line<3, ScalarT> const& l)
{
auto ls = inf_of(s);
auto len = length(s);
auto [ts, tl] = closest_points_parameters(ls, l);
auto tClamped = clamp(ts, ScalarT(0), len);
return {tClamped / len, coordinates(l, ls[tClamped])};
}
// =========== Other Implementations ===========
......
......@@ -46,6 +46,20 @@ template <int D, class ScalarT>
return b == o;
}
// default implementation for contains(objA, objB) that works when objA is solid and vertices_of(objB) is defined
template <class A, class B>
[[nodiscard]] constexpr auto contains(A const& a, B const& b, dont_deduce<typename B::scalar_t> eps = static_cast<typename B::scalar_t>(0))
-> enable_if<std::is_same_v<typename object_traits<A>::tag_t, default_object_tag>, decltype(vertices_of(b), false)>
{
for (auto const& vertex : vertices_of(b))
if (!contains(a, vertex, eps))
return false;
return true;
}
// object specific implementations for contains(obj, pos)
template <class ScalarT>
[[nodiscard]] constexpr bool contains(aabb<1, ScalarT> const& b, ScalarT const& o, dont_deduce<ScalarT> eps = ScalarT(0))
{
......@@ -120,8 +134,7 @@ template <int D, class ScalarT>
if (ri > bi + eps)
return false; // False if outside of the aabb in any dimension
if (!onSomeBoundary && (ri >= bi - eps))
onSomeBoundary = true;
onSomeBoundary = onSomeBoundary || (ri >= bi - eps);
}
return onSomeBoundary; // True, if at on the boundary in at least one dimension
}
......
......@@ -91,45 +91,72 @@ template <class ScalarT>
}
template <class ScalarT>
[[nodiscard]] constexpr fractional_result<ScalarT> distance(segment<2, ScalarT> const& s0, segment<2, ScalarT> const& s1)
[[nodiscard]] constexpr fractional_result<ScalarT> distance_sqr(segment<2, ScalarT> const& s0, segment<2, ScalarT> const& s1)
{
auto l0 = line<2, ScalarT>::from_points(s0.pos0, s0.pos1);
auto l1 = line<2, ScalarT>::from_points(s1.pos0, s1.pos1);
auto sl0 = distance(s0.pos0, s0.pos1);
auto sl1 = distance(s1.pos0, s1.pos1);
auto l0 = inf_of(s0);
auto l1 = inf_of(s1);
auto len0 = length(s0);
auto len1 = length(s1);
auto [t0, t1] = intersection_parameters(l0, l1);
if (ScalarT(0) <= t0 && t0 <= ScalarT(sl0) && //
ScalarT(0) <= t1 && t1 <= ScalarT(sl1))
if (ScalarT(0) <= t0 && t0 <= len0 && //
ScalarT(0) <= t1 && t1 <= len1)
return ScalarT(0); // intersects
auto p0 = t0 < ScalarT(0) ? s0.pos0 : s0.pos1;
auto p1 = t1 < ScalarT(0) ? s1.pos0 : s1.pos1;
auto p0 = t0 * ScalarT(2) < len0 ? s0.pos0 : s0.pos1;
auto p1 = t1 * ScalarT(2) < len1 ? s1.pos0 : s1.pos1;
return min(distance(p0, s1), distance(p1, s0));
return min(distance_sqr(p0, s1), distance_sqr(p1, s0));
}
template <class ScalarT>
[[nodiscard]] constexpr fractional_result<ScalarT> distance(segment<3, ScalarT> const& s0, segment<3, ScalarT> const& s1)
[[nodiscard]] constexpr fractional_result<ScalarT> distance_sqr(segment<3, ScalarT> const& s0, segment<3, ScalarT> const& s1)
{
auto l0 = line<3, ScalarT>::from_points(s0.pos0, s0.pos1);
auto l1 = line<3, ScalarT>::from_points(s1.pos0, s1.pos1);
auto sl0 = distance(s0.pos0, s0.pos1);
auto sl1 = distance(s1.pos0, s1.pos1);
auto l0 = inf_of(s0);
auto l1 = inf_of(s1);
auto len0 = length(s0);
auto len1 = length(s1);
auto [t0, t1] = closest_points_parameters(l0, l1);
if (ScalarT(0) <= t0 && t0 <= ScalarT(sl0) && //
ScalarT(0) <= t1 && t1 <= ScalarT(sl1))
return distance(l0[t0], l1[t1]); // closest points is inside segments
if (ScalarT(0) <= t0 && t0 <= len0 && //
ScalarT(0) <= t1 && t1 <= len1)
return distance_sqr(l0[t0], l1[t1]); // closest points is inside segments
auto p0 = t0 * ScalarT(2) < len0 ? s0.pos0 : s0.pos1;
auto p1 = t1 * ScalarT(2) < len1 ? s1.pos0 : s1.pos1;
return min(distance_sqr(p0, s1), distance_sqr(p1, s0));
}
template <class ScalarT>
[[nodiscard]] constexpr fractional_result<ScalarT> distance_sqr(segment<2, ScalarT> const& s, line<2, ScalarT> const& l)
{
auto ls = inf_of(s);
auto len = length(s);
auto [ts, tl] = intersection_parameters(ls, l);
if (ScalarT(0) <= ts && ts <= len)
return ScalarT(0); // intersects
auto p0 = t0 < ScalarT(0) ? s0.pos0 : s0.pos1;
auto p1 = t1 < ScalarT(0) ? s1.pos0 : s1.pos1;
auto p = ts * ScalarT(2) < len ? s.pos0 : s.pos1;
return distance_sqr(p, l);
}
template <class ScalarT>
[[nodiscard]] constexpr fractional_result<ScalarT> distance_sqr(segment<3, ScalarT> const& s, line<3, ScalarT> const& l)
{
auto ls = inf_of(s);
auto len = length(s);
return min(distance(p0, s1), distance(p1, s0));
auto [ts, tl] = closest_points_parameters(ls, l);
auto tClamped = clamp(ts, ScalarT(0), len);
return distance_sqr(ls[tClamped], l);
}
template <int D, class ScalarT>
[[nodiscard]] constexpr fractional_result<ScalarT> distance_sqr(line<D, ScalarT> const& l, segment<D, ScalarT> const& s)
{
return distance_sqr(s, l);
}
// TODO: use GJK or something?
......
......@@ -17,7 +17,7 @@ template <int D, class ScalarT>
}
template <int D, class ScalarT>
[[nodiscard]] constexpr array<segment<2, ScalarT>, 4> edges_of(quad<D, ScalarT> const& q)
[[nodiscard]] constexpr array<segment<D, ScalarT>, 4> edges_of(quad<D, ScalarT> const& q)
{
return {{{q.pos00, q.pos10}, {q.pos10, q.pos11}, {q.pos11, q.pos01}, {q.pos01, q.pos00}}};
}
......@@ -60,14 +60,14 @@ template <class ScalarT, class TraitsT>
}
template <class ScalarT, int DomainD, class TraitsT>
[[nodiscard]] constexpr array<segment<2, ScalarT>, 4> edges_of(box<2, ScalarT, DomainD, TraitsT> const& b)
[[nodiscard]] constexpr array<segment<DomainD, ScalarT>, 4> edges_of(box<2, ScalarT, DomainD, TraitsT> const& b)
{
const auto vs = vertices_of(b);
return {{{vs[0], vs[1]}, {vs[1], vs[2]}, {vs[2], vs[3]}, {vs[3], vs[0]}}};
}
template <class ScalarT, int DomainD, class TraitsT>
[[nodiscard]] constexpr array<segment<3, ScalarT>, 12> edges_of(box<3, ScalarT, DomainD, TraitsT> const& b)
[[nodiscard]] constexpr array<segment<DomainD, ScalarT>, 12> edges_of(box<3, ScalarT, DomainD, TraitsT> const& b)
{
const auto vs = vertices_of(b);
return {{
......@@ -79,12 +79,10 @@ template <class ScalarT, int DomainD, class TraitsT>
}};
}
template <class BaseT, class TraitsT>
template <class BaseT, class TraitsT, typename = std::enable_if_t<!std::is_same_v<BaseT, sphere<2, typename BaseT::scalar_t, 3>>>>
[[nodiscard]] constexpr auto edges_of(pyramid<BaseT, TraitsT> const& py)
{
using ScalarT = typename BaseT::scalar_t;
static_assert(!std::is_same_v<BaseT, sphere<2, ScalarT, 3>>, "not possible for cones");
const auto apex = apex_of(py);
const auto edgesBase = edges_of(py.base);
auto res = array<segment<3, ScalarT>, edgesBase.size() * 2>();
......
......@@ -62,12 +62,11 @@ template <class ScalarT, class TraitsT>
return {{faceX0, faceX1, faceY0, faceY1, faceZ0, faceZ1}};
}
template <class BaseT>
template <class BaseT, typename = std::enable_if_t<!std::is_same_v<BaseT, sphere<2, typename BaseT::scalar_t, 3>>>>
[[nodiscard]] constexpr auto faces_of(pyramid<BaseT, boundary_no_caps_tag> const& py)
{
using ScalarT = typename BaseT::scalar_t;
static_assert(is_floating_point<ScalarT>, "cannot be guaranteed for integers");
static_assert(!std::is_same_v<BaseT, sphere<2, ScalarT, 3>>, "not possible for cones");
const auto apex = apex_of(py);
const auto verts = vertices_of(py.base);
......@@ -77,7 +76,7 @@ template <class BaseT>
return triangles;
}
template <class BaseT, class TraitsT>
template <class BaseT, class TraitsT, typename = std::enable_if_t<!std::is_same_v<BaseT, sphere<2, typename BaseT::scalar_t, 3>>>>
[[nodiscard]] constexpr auto faces_of(pyramid<BaseT, TraitsT> const& py)
{
auto mantle = faces_of(boundary_no_caps_of(py));
......
#pragma once
#include <typed-geometry/functions/objects/normal.hh>
#include <typed-geometry/types/objects/box.hh>
#include <typed-geometry/types/objects/ellipse.hh>
#include <typed-geometry/types/objects/plane.hh>
#include <typed-geometry/types/objects/sphere.hh>
#include <typed-geometry/types/objects/triangle.hh>
// tg::plane_of(obj) returns the supporting plane of an object
// (when unambiguously defined)
// TODO: maybe also for "base planes" as in cone / hemisphere / pyramid?
namespace tg
{
template <int D, class ScalarT>
[[nodiscard]] constexpr plane<D, ScalarT> plane_of(plane<D, ScalarT> const& p)
{
return p;
}
template <int D, class ScalarT>
[[nodiscard]] constexpr plane<D, ScalarT> plane_of(halfspace<D, ScalarT> const& p)
{
return plane<3, ScalarT>(p.normal, p.dis);
}
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr plane<3, ScalarT> plane_of(box<2, ScalarT, 3, TraitsT> const& b)
{
return plane<3, ScalarT>(normal_of(b), b.center);
}
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr plane<3, ScalarT> plane_of(ellipse<2, ScalarT, 3, TraitsT> const& e)
{
return plane<3, ScalarT>(normal_of(e), e.center);
}
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr plane<3, ScalarT> plane_of(sphere<2, ScalarT, 3, TraitsT> const& s)
{
return plane<3, ScalarT>(s.normal, s.center);
}
template <class ScalarT>
[[nodiscard]] constexpr plane<3, ScalarT> plane_of(triangle<3, ScalarT> const& t)
{
return plane<3, ScalarT>(normal_of(t), t.pos0);
}
}
#pragma once
#include <typed-geometry/feature/basic.hh>
#include <typed-geometry/functions/objects/edges.hh>
namespace tg
{
/// calls on_segment for each edge of the object
/// on_segment: (tg::segment) -> void
template <class Obj, class OnSegment, std::enable_if_t<has_edges_of<Obj>, int> = 0>
constexpr void segmentize(Obj const& obj, OnSegment&& on_segment)
{
for (auto&& s : edges_of(obj))
on_segment(s);
}
/// calls on_segment for a discretized version of the sphere
template <class ScalarT, class TraitsT, class OnSegment>
void segmentize(sphere<2, ScalarT, 2, TraitsT> const& s, int segs, OnSegment&& on_segment)
{
TG_ASSERT(segs >= 3);
using dir_t = dir<2, ScalarT>;
auto const dir_of = [&](int i) -> dir_t {
if (i == 0 || i == segs)
return dir_t(1, 0);
auto [sa, ca] = sin_cos(tau<ScalarT> * i / ScalarT(segs));
return dir_t(ca, sa);
};
auto const pos_of = [&](int i) {
auto d = dir_of(i);
return s.center + d * s.radius;
};
auto prev = pos_of(0);
for (auto i = 1; i <= segs; ++i)
{
auto curr = pos_of(i);
on_segment(tg::segment<2, ScalarT>(prev, curr));
prev = curr;
}
}
template <class ScalarT, class TraitsT, class OnSegment>
void segmentize(sphere<2, ScalarT, 3, TraitsT> const& s, int segs, OnSegment&& on_segment)
{
TG_ASSERT(segs >= 3);
auto const dx = any_normal(s.normal);
auto const dy = cross(s.normal, dx);
using dir_t = dir<2, ScalarT>;
auto const dir_of = [&](int i) -> dir_t {
if (i == 0 || i == segs)
return dir_t(1, 0);
auto [sa, ca] = sin_cos(tau<ScalarT> * i / ScalarT(segs));
return dir_t(ca, sa);
};
auto const pos_of = [&](int i) {
auto d = dir_of(i) * s.radius;
return s.center + dx * d.x + dy * d.y;
};
auto prev = pos_of(0);
for (auto i = 1; i <= segs; ++i)
{
auto curr = pos_of(i);
on_segment(tg::segment<3, ScalarT>(prev, curr));
prev = curr;
}
}
}
#pragma once
#include <typed-geometry/feature/basic.hh>
#include <typed-geometry/functions/objects/triangulation.hh>
namespace tg
{
/// calls on_triangle for each triangle of the objects triangulation
/// on_triangle: (tg::triangle) -> void
template <class Obj, class OnTriangle, std::enable_if_t<has_triangulation_of<Obj>, int> = 0>
constexpr void triangulate(Obj const& obj, OnTriangle&& on_triangle)
{
for (auto&& t : triangulation_of(obj))
on_triangle(t);
}
/// computes a uv triangulation of the given sphere and calls on_triangle with each triangle
/// on_triangle: (tg::triangle) -> void
/// NOTE: currently recomputes a lot of sin/cos and thus is not the fastest
/// NOTE: normal_of(t) points outwards
/// NOTE: segs_u must be >= 3, segs_v must be >= 2
template <class ScalarT, class TraitsT, class OnTriangle>
void triangulate_uv(sphere<3, ScalarT, 3, TraitsT> const& s, int segs_u, int segs_v, OnTriangle&& on_triangle)
{
TG_ASSERT(segs_u >= 3);
TG_ASSERT(segs_v >= 2);
// TODO: some caching of sin/cos
using dir_t = dir<3, ScalarT>;
auto const dir_of = [&](int u, int v) -> dir_t {
if (v == 0)
return dir_t(0, 1, 0);
else if (v == segs_v)
return dir_t(0, -1, 0);
auto [su, cu] = sin_cos(tau<ScalarT> * (u == segs_u ? 0 : u) / ScalarT(segs_u));
auto [sv, cv] = sin_cos(pi<ScalarT> * v / ScalarT(segs_v));
return dir_t(sv * su, cv, sv * cu);
};
auto const pos_of = [&](int u, int v) {
auto d = dir_of(u, v);
return s.center + d * s.radius;
};
// cap u
{
auto p0 = pos_of(0, 0);
for (auto i = 0; i < segs_u; ++i)
{
auto p1 = pos_of(i, 1);
auto p2 = pos_of(i + 1, 1);
on_triangle(tg::triangle<3, ScalarT>(p0, p1, p2));
}
}
// inner grid
for (auto j = 1; j < segs_v - 1; ++j)
{
for (auto i = 0; i < segs_u; ++i)
{
auto p00 = pos_of(i + 0, j + 0);
auto p01 = pos_of(i + 0, j + 1);
auto p10 = pos_of(i + 1, j + 0);
auto p11 = pos_of(i + 1, j + 1);
on_triangle(tg::triangle<3, ScalarT>(p00, p01, p11));
on_triangle(tg::triangle<3, ScalarT>(p00, p11, p10));
}
}
// cap v
{
auto p0 = pos_of(0, segs_v);
for (auto i = 0; i < segs_u; ++i)
{
auto p1 = pos_of(i, segs_v - 1);
auto p2 = pos_of(i + 1, segs_v - 1);
on_triangle(tg::triangle<3, ScalarT>(p0, p2, p1));
}
}
}
}
#pragma once
#include <typed-geometry/functions/objects/faces.hh>
// triangulation_of(obj) returns a collection of triangles approximating the surface of obj
// sometimes, a second parameter can be provided for fine control
namespace tg
{
template <int D, class ScalarT>
[[nodiscard]] constexpr array<triangle<D, ScalarT>, 1> triangulation_of(triangle<D, ScalarT> const& t)
{
return {{t}};
}
template <int D, class ScalarT>
[[nodiscard]] constexpr array<triangle<D, ScalarT>, 2> triangulation_of(quad<D, ScalarT> const& q)
{
// TODO: orientation guarantees?
return {{
{q.pos00, q.pos10, q.pos11}, //
{q.pos00, q.pos11, q.pos01},
}};
}
template <int D, class ScalarT, class TraitsT>
[[nodiscard]] constexpr array<triangle<D, ScalarT>, 2> triangulation_of(box<2, ScalarT, D, TraitsT> const& b)
{
auto const p00 = b.center - b.half_extents[0] - b.half_extents[1];
auto const p10 = b.center + b.half_extents[0] - b.half_extents[1];
auto const p01 = b.center - b.half_extents[0] + b.half_extents[1];
auto const p11 = b.center + b.half_extents[0] + b.half_extents[1];
// TODO: orientation guarantees?
return {{
{p00, p10, p11}, //
{p00, p11, p01},
}};
}
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr array<triangle<3, ScalarT>, 12> triangulation_of(box<3, ScalarT, 3, TraitsT> const& b)
{
int idx = 0;
array<triangle<3, ScalarT>, 12> r;
for (auto const& f : faces_of(b))
for (auto const& t : triangulation_of(f))
r[idx++] = t;
return r;
}
// normal_of(t) points outwards
template <class ScalarT, class TraitsT>
[[nodiscard]] constexpr array<triangle<3, ScalarT>, 12> triangulation_of(aabb<3, ScalarT, TraitsT> const& b)
{
int idx = 0;
array<triangle<3, ScalarT>, 12> r;
auto dir = b.max - b.min;
auto add_face = [&](pos<3, ScalarT> const& p, vec<3, ScalarT> const& d0, vec<3, ScalarT> const& d1) {
auto const p00 = p;
auto const p10 = p + d0;
auto const p01 = p + d1;
auto const p11 = p + d0 + d1;