GeoLib::SurfaceGrid Class Reference

Detailed Description

Definition at line 30 of file SurfaceGrid.h.

#include <SurfaceGrid.h>

Inheritance diagram for GeoLib::SurfaceGrid:

Collaboration diagram for GeoLib::SurfaceGrid:

Public Member Functions
	SurfaceGrid (GeoLib::Surface const *const sfc)
bool	isPointInSurface (MathLib::Point3d const &pnt, double eps=std::numeric_limits< double >::epsilon()) const
Public Member Functions inherited from GeoLib::AABB
template<typename PNT_TYPE >
	AABB (std::vector< PNT_TYPE * > const &pnts, std::vector< std::size_t > const &ids)
template<typename InputIterator >
	AABB (InputIterator first, InputIterator last)
template<typename PNT_TYPE >
bool	update (PNT_TYPE const &p)
template<typename T >
bool	containsPoint (T const &pnt, double eps) const
template<typename T >
bool	containsPointXY (T const &pnt) const
MinMaxPoints	getMinMaxPoints () const
Eigen::Vector3d const &	getMinPoint () const
Eigen::Vector3d const &	getMaxPoint () const
bool	containsAABB (AABB const &other_aabb) const

Private Member Functions
void	sortTrianglesInGridCells (GeoLib::Surface const *const sfc)
bool	sortTriangleInGridCells (GeoLib::Triangle const *const triangle)
std::optional< std::array< std::size_t, 3 > >	getGridCellCoordinates (MathLib::Point3d const &p) const

Private Attributes
std::array< double, 3 >	_step_sizes {}
std::array< double, 3 >	_inverse_step_sizes {}
std::array< std::size_t, 3 >	_n_steps
std::vector< std::vector< GeoLib::Triangle const * > >	_triangles_in_grid_box

Constructor & Destructor Documentation

SurfaceGrid()

GeoLib::SurfaceGrid::SurfaceGrid ( GeoLib::Surface const *const sfc )

explicit

Definition at line 25 of file SurfaceGrid.cpp.

    : AABB(sfc->getAABB()), _n_steps({{1, 1, 1}})
{
    auto [min_point, max_point] = getMinMaxPoints();
    // enlarge the bounding, such that the points with maximal coordinates
    // fits into the grid
    for (std::size_t k(0); k < 3; ++k)
    {
        max_point[k] += std::abs(max_point[k]) * 1e-6;
        if (std::abs(max_point[k]) < std::numeric_limits<double>::epsilon())
        {
            max_point[k] = (max_point[k] - min_point[k]) * (1.0 + 1e-6);
        }
    }
 
    Eigen::Vector3d delta = max_point - min_point;
 
    if (delta[0] < std::numeric_limits<double>::epsilon())
    {
        const double max_delta(std::max(delta[1], delta[2]));
        min_point[0] -= max_delta * 0.5e-3;
        max_point[0] += max_delta * 0.5e-3;
        delta[0] = max_point[0] - min_point[0];
    }
 
    if (delta[1] < std::numeric_limits<double>::epsilon())
    {
        const double max_delta(std::max(delta[0], delta[2]));
        min_point[1] -= max_delta * 0.5e-3;
        max_point[1] += max_delta * 0.5e-3;
        delta[1] = max_point[1] - min_point[1];
    }
 
    if (delta[2] < std::numeric_limits<double>::epsilon())
    {
        const double max_delta(std::max(delta[0], delta[1]));
        min_point[2] -= max_delta * 0.5e-3;
        max_point[2] += max_delta * 0.5e-3;
        delta[2] = max_point[2] - min_point[2];
    }
 
    update(min_point);
    update(max_point);
 
    const std::size_t n_tris(sfc->getNumberOfTriangles());
    const std::size_t n_tris_per_cell(5);
 
    Eigen::Matrix<bool, 3, 1> dim =
        delta.array() >= std::numeric_limits<double>::epsilon();
 
    // *** condition: n_tris / n_cells < n_tris_per_cell
    //                where n_cells = _n_steps[0] * _n_steps[1] * _n_steps[2]
    // *** with _n_steps[0] =
    // ceil(pow(n_tris*delta[0]*delta[0]/(n_tris_per_cell*delta[1]*delta[2]),
    // 1/3.)));
    //          _n_steps[1] = _n_steps[0] * delta[1]/delta[0],
    //          _n_steps[2] = _n_steps[0] * delta[2]/delta[0]
    auto sc_ceil = [](double v)
    { return static_cast<std::size_t>(std::ceil(v)); };
    switch (dim.count())
    {
        case 3:  // 3d case
            _n_steps[0] =
                sc_ceil(std::cbrt(n_tris * delta[0] * delta[0] /
                                  (n_tris_per_cell * delta[1] * delta[2])));
            _n_steps[1] =
                sc_ceil(_n_steps[0] * std::min(delta[1] / delta[0], 100.0));
            _n_steps[2] =
                sc_ceil(_n_steps[0] * std::min(delta[2] / delta[0], 100.0));
            break;
        case 2:  // 2d cases
            if (dim[0] && dim[2])
            {  // 2d case: xz plane, y = const
                _n_steps[0] = sc_ceil(std::sqrt(n_tris * delta[0] /
                                                (n_tris_per_cell * delta[2])));
                _n_steps[2] = sc_ceil(_n_steps[0] * delta[2] / delta[0]);
            }
            else if (dim[0] && dim[1])
            {  // 2d case: xy plane, z = const
                _n_steps[0] = sc_ceil(std::sqrt(n_tris * delta[0] /
                                                (n_tris_per_cell * delta[1])));
                _n_steps[1] = sc_ceil(_n_steps[0] * delta[1] / delta[0]);
            }
            else if (dim[1] && dim[2])
            {  // 2d case: yz plane, x = const
                _n_steps[1] = sc_ceil(std::sqrt(n_tris * delta[1] /
                                                (n_tris_per_cell * delta[2])));
                _n_steps[2] =
                    sc_ceil(n_tris * delta[2] / (n_tris_per_cell * delta[1]));
            }
            break;
        case 1:  // 1d cases
            for (std::size_t k(0); k < 3; ++k)
            {
                if (dim[k])
                {
                    _n_steps[k] =
                        sc_ceil(static_cast<double>(n_tris) / n_tris_per_cell);
                }
            }
    }
 
    // some frequently used expressions to fill the grid vectors
    for (std::size_t k(0); k < 3; k++)
    {
        _step_sizes[k] = delta[k] / _n_steps[k];
        if (delta[k] > std::numeric_limits<double>::epsilon())
        {
            _inverse_step_sizes[k] = 1.0 / _step_sizes[k];
        }
        else
        {
            _inverse_step_sizes[k] = 0;
        }
    }
 
    _triangles_in_grid_box.resize(_n_steps[0] * _n_steps[1] * _n_steps[2]);
    sortTrianglesInGridCells(sfc);
}

Member Function Documentation

getGridCellCoordinates()

std::optional< std::array< std::size_t, 3 > > GeoLib::SurfaceGrid::getGridCellCoordinates ( MathLib::Point3d const & p ) const

private

Definition at line 219 of file SurfaceGrid.cpp.

{
    auto const& min_point{getMinPoint()};
    std::array<std::size_t, 3> coords{
        {static_cast<std::size_t>((p[0] - min_point[0]) *
                                  _inverse_step_sizes[0]),
         static_cast<std::size_t>((p[1] - min_point[1]) *
                                  _inverse_step_sizes[1]),
         static_cast<std::size_t>((p[2] - min_point[2]) *
                                  _inverse_step_sizes[2])}};
 
    if (coords[0] >= _n_steps[0] || coords[1] >= _n_steps[1] ||
        coords[2] >= _n_steps[2])
    {
        DBUG(
            "Computed indices ({:d},{:d},{:d}), max grid cell indices "
            "({:d},{:d},{:d})",
            coords[0], coords[1], coords[2], _n_steps[0], _n_steps[1],
            _n_steps[2]);
        return {};
    }
    return coords;
}

References _inverse_step_sizes, _n_steps, DBUG(), and GeoLib::AABB::getMinPoint().

Referenced by isPointInSurface(), and sortTriangleInGridCells().

isPointInSurface()

bool GeoLib::SurfaceGrid::isPointInSurface	(	MathLib::Point3d const &	pnt,
		double	eps = std::numeric_limits<double>::epsilon() ) const

Definition at line 244 of file SurfaceGrid.cpp.

{
    auto const grid_cell_coordinates = getGridCellCoordinates(pnt);
    if (!grid_cell_coordinates)
    {
        return false;
    }
 
    std::size_t const grid_cell_index =
        (*grid_cell_coordinates)[0] +
        (*grid_cell_coordinates)[1] * _n_steps[0] +
        (*grid_cell_coordinates)[2] * _n_steps[0] * _n_steps[1];
 
    std::vector<Triangle const*> const& triangles(
        _triangles_in_grid_box[grid_cell_index]);
    return std::any_of(triangles.begin(), triangles.end(),
                       [eps, pnt](auto const* triangle)
                       { return triangle->containsPoint(pnt, eps); });
}

References _n_steps, _triangles_in_grid_box, and getGridCellCoordinates().

sortTriangleInGridCells()

bool GeoLib::SurfaceGrid::sortTriangleInGridCells ( GeoLib::Triangle const *const triangle )

private

Definition at line 165 of file SurfaceGrid.cpp.

{
    // compute grid coordinates for each triangle point
    std::optional<std::array<std::size_t, 3> const> c_p0(
        getGridCellCoordinates(*(triangle->getPoint(0))));
    if (!c_p0)
    {
        return false;
    }
    std::optional<std::array<std::size_t, 3> const> c_p1(
        getGridCellCoordinates(*(triangle->getPoint(1))));
    if (!c_p1)
    {
        return false;
    }
    std::optional<std::array<std::size_t, 3> const> c_p2(
        getGridCellCoordinates(*(triangle->getPoint(2))));
    if (!c_p2)
    {
        return false;
    }
 
    // determine interval in grid (grid cells the triangle will be inserted)
    std::size_t const i_min(
        std::min(std::min((*c_p0)[0], (*c_p1)[0]), (*c_p2)[0]));
    std::size_t const i_max(
        std::max(std::max((*c_p0)[0], (*c_p1)[0]), (*c_p2)[0]));
    std::size_t const j_min(
        std::min(std::min((*c_p0)[1], (*c_p1)[1]), (*c_p2)[1]));
    std::size_t const j_max(
        std::max(std::max((*c_p0)[1], (*c_p1)[1]), (*c_p2)[1]));
    std::size_t const k_min(
        std::min(std::min((*c_p0)[2], (*c_p1)[2]), (*c_p2)[2]));
    std::size_t const k_max(
        std::max(std::max((*c_p0)[2], (*c_p1)[2]), (*c_p2)[2]));
 
    const std::size_t n_plane(_n_steps[0] * _n_steps[1]);
 
    // insert the triangle into the grid cells
    for (std::size_t i(i_min); i <= i_max; i++)
    {
        for (std::size_t j(j_min); j <= j_max; j++)
        {
            for (std::size_t k(k_min); k <= k_max; k++)
            {
                _triangles_in_grid_box[i + j * _n_steps[0] + k * n_plane]
                    .push_back(triangle);
            }
        }
    }
 
    return true;
}

References _n_steps, _triangles_in_grid_box, getGridCellCoordinates(), and GeoLib::Triangle::getPoint().

Referenced by sortTrianglesInGridCells().

sortTrianglesInGridCells()

void GeoLib::SurfaceGrid::sortTrianglesInGridCells ( GeoLib::Surface const *const sfc )

private

Definition at line 145 of file SurfaceGrid.cpp.

{
    for (std::size_t l(0); l < sfc->getNumberOfTriangles(); l++)
    {
        if (!sortTriangleInGridCells((*sfc)[l]))
        {
            Point const& p0(*((*sfc)[l]->getPoint(0)));
            Point const& p1(*((*sfc)[l]->getPoint(1)));
            Point const& p2(*((*sfc)[l]->getPoint(2)));
            auto const [min, max] = getMinMaxPoints();
            OGS_FATAL(
                "Sorting triangle {:d} [({:f},{:f},{:f}), ({:f},{:f},{:f}), "
                "({:f},{:f},{:f}) into grid. Bounding box is [{:f}, {:f}] x "
                "[{:f}, {:f}] x [{:f}, {:f}].",
                l, p0[0], p0[1], p0[2], p1[0], p1[1], p1[2], p2[0], p2[1],
                p2[2], min[0], max[0], min[1], max[1], min[2], max[2]);
        }
    }
}

References GeoLib::AABB::getMinMaxPoints(), GeoLib::Surface::getNumberOfTriangles(), OGS_FATAL, and sortTriangleInGridCells().

Member Data Documentation

_inverse_step_sizes

std::array<double, 3> GeoLib::SurfaceGrid::_inverse_step_sizes {}

private

Definition at line 44 of file SurfaceGrid.h.

{};

Referenced by getGridCellCoordinates().

_n_steps

std::array<std::size_t, 3> GeoLib::SurfaceGrid::_n_steps

private

Definition at line 45 of file SurfaceGrid.h.

Referenced by getGridCellCoordinates(), isPointInSurface(), and sortTriangleInGridCells().

_step_sizes

std::array<double, 3> GeoLib::SurfaceGrid::_step_sizes {}

private

Definition at line 43 of file SurfaceGrid.h.

{};

_triangles_in_grid_box

std::vector<std::vector<GeoLib::Triangle const*> > GeoLib::SurfaceGrid::_triangles_in_grid_box

private

Definition at line 46 of file SurfaceGrid.h.

Referenced by isPointInSurface(), and sortTriangleInGridCells().

---

The documentation for this class was generated from the following files:

• GeoLib/SurfaceGrid.h
• GeoLib/SurfaceGrid.cpp