Exact closest neighbor not working correctly

classic Classic list List threaded Threaded
8 messages Options
Reply | Threaded
Open this post in threaded view
|

Exact closest neighbor not working correctly

Juan Jose Casafranca
I am trying to find the closest neighbor in a point cloud to a given point. I have followed the documentation example in which the points are stored in a user vector and the Tree only contains indices. Since the results were not what I was expecting, I am also computing the closest neighbor by brute force. Ill explain my code:

I have a vector of points (coarsePoints) and I have a vector of indices (xC_vertices). The indices are the points in coarsePoints that can be marked as closest neighbor. 

This is how I compute the closest point using a brute force algorithm
```
for (int i = 0; i < d_vertices.size(); ++i) {
        Real minDistance = std::numeric_limits<Real>::infinity();
        int vertexId = -1;
        const auto &dPoint = highResPoints[d_vertices[i]];
        const CPF::Point_3 pD{dPoint[0], dPoint[1], dPoint[2]};
        for (int j = 0; j < xC_vertices.size(); ++j) {
            const auto &xPoint = coarsePoints[xC_vertices[j]];
            const CPF::Point_3 pC{xPoint[0], xPoint[1], xPoint[2]};

            if (CGAL::squared_distance(pC, pD) < minDistance) {
                vertexId = j;
                minDistance = CGAL::squared_distance(pC, pD);
            }
        }
```

When using CGAL, I have the following (based on the documentation example):

```
class My_point_property_map
{
    const sofa::helper::vector<CPF::SofaTypes::Point> &points;
    const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices;
public:
    typedef CPF::Point_3 value_type;
    typedef value_type &reference;
    typedef std::size_t key_type;
    typedef boost::lvalue_property_map_tag category;
    My_point_property_map(const sofa::helper::vector<CPF::SofaTypes::Point> &pts,
                          const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices_)
        : points(pts)
        , xC_vertices(xC_vertices_)
    {
    }
    value_type operator[](key_type k) const
    {
        const CPF::SofaTypes::Point &sPoint = points[xC_vertices[k]];
        return CPF::Point_3{sPoint[0], sPoint[1], sPoint[2]};
    }
    friend value_type get(const My_point_property_map &ppmap, key_type i) { return ppmap[i]; }
};

using SearchTraitsBase = CGAL::Search_traits_3<CPF::K>;
using SearchTraits = CGAL::Search_traits_adapter<std::size_t, My_point_property_map, SearchTraitsBase>;
using NeighborSearch = CGAL::Orthogonal_k_neighbor_search<SearchTraits>;
using Tree = NeighborSearch::Tree;
using Splitter = typename Tree::Splitter;
using SearchDistance = typename NeighborSearch::Distance;

int searchNeighbor() {
    My_point_property_map ppmap(coarsePoints, xC_vertices);
    Tree searchTree(boost::counting_iterator<int>(0),
                    boost::counting_iterator<int>(xC_vertices.size()),
                    Splitter(),
                    SearchTraits(ppmap));
    SearchDistance trDistance(ppmap);

for (int i = 0; i < d_vertices.size(); ++i) {
        NeighborSearch search(searchTree, CPF::Point_3{dPoint[0], dPoint[1], dPoint[2]}, 1, 0, true, trDistance);
        for (const auto &searchResult : search) {
            if (searchResult.first != vertexId) {
                spdlog::get("cpf")->info("Brute force id = {}", vertexId);
                spdlog::get("cpf")->info("CGAL id = {}", searchResult.first);
                spdlog::get("cpf")->info("Brute force distance = {}", minDistance);
                spdlog::get("cpf")->info("CGAL distance = {}", searchResult.second);
            }
        }
} }
```

This return a los of vertices that are actually different:

[2020-04-17 18:10:31.988] [cpf] [info] Brute force id = 386

[2020-04-17 18:10:31.988] [cpf] [info] CGAL id = 84

[2020-04-17 18:10:31.988] [cpf] [info] Brute force distance = 0.0

[2020-04-17 18:10:31.988] [cpf] [info] CGAL distance = 0.041339367896971506


I am not sure if I am using incorrectly the NeighborSearch or what is happening. Any ideas what might be wrong?




Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Sebastien Loriot (GeometryFactory)
Could you please provide a minimal example that we could compile, run
and that would be showing the error?

At first sight, I don't see anything suspicious.

Thanks,

Sebastien.

On 4/17/20 6:25 PM, Juan Jose Casafranca wrote:

> I am trying to find the closest neighbor in a point cloud to a given
> point. I have followed the documentation example in which the points are
> stored in a user vector and the Tree only contains indices. Since the
> results were not what I was expecting, I am also computing the closest
> neighbor by brute force. Ill explain my code:
>
> I have a vector of points (coarsePoints) and I have a vector of indices
> (xC_vertices). The indices are the points in coarsePoints that can be
> marked as closest neighbor.
>
> This is how I compute the closest point using a brute force algorithm
> ```
>
> for (int i = 0; i < d_vertices.size(); ++i) {
>
> Real minDistance = std::numeric_limits<Real>::infinity();
>
> int vertexId = -1;
>
> const auto &dPoint = highResPoints[d_vertices[i]];
>
> const CPF::Point_3 pD{dPoint[0], dPoint[1], dPoint[2]};
>
> for (int j = 0; j < xC_vertices.size(); ++j) {
>
> const auto &xPoint = coarsePoints[xC_vertices[j]];
>
> const CPF::Point_3 pC{xPoint[0], xPoint[1], xPoint[2]};
>
>
> if (CGAL::squared_distance(pC, pD) < minDistance) {
>
> vertexId = j;
>
> minDistance = CGAL::squared_distance(pC, pD);
>
> }
>
> }
>
> ```
>
>
> When using CGAL, I have the following (based on the documentation example):
>
>
> ```
>
> class My_point_property_map
>
> {
>
> const sofa::helper::vector<CPF::SofaTypes::Point> &points;
>
> const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices;
>
> public:
>
> typedef CPF::Point_3 value_type;
>
> typedef value_type &reference;
>
> typedef std::size_t key_type;
>
> typedef boost::lvalue_property_map_tag category;
>
> My_point_property_map(const sofa::helper::vector<CPF::SofaTypes::Point>
> &pts,
>
> const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices_)
>
> : points(pts)
>
> , xC_vertices(xC_vertices_)
>
> {
>
> }
>
> value_type operator[](key_type k) const
>
> {
>
> const CPF::SofaTypes::Point &sPoint = points[xC_vertices[k]];
>
> return CPF::Point_3{sPoint[0], sPoint[1], sPoint[2]};
>
> }
>
> friend value_type get(const My_point_property_map &ppmap, key_type i) {
> return ppmap[i]; }
>
> };
>
>
> using SearchTraitsBase = CGAL::Search_traits_3<CPF::K>;
>
> using SearchTraits = CGAL::Search_traits_adapter<std::size_t,
> My_point_property_map, SearchTraitsBase>;
>
> using NeighborSearch = CGAL::Orthogonal_k_neighbor_search<SearchTraits>;
>
> using Tree = NeighborSearch::Tree;
>
> using Splitter = typename Tree::Splitter;
>
> using SearchDistance = typename NeighborSearch::Distance;
>
>
> int searchNeighbor() {
>
> My_point_property_map ppmap(coarsePoints, xC_vertices);
>
> Tree searchTree(boost::counting_iterator<int>(0),
>
> boost::counting_iterator<int>(xC_vertices.size()),
>
> Splitter(),
>
> SearchTraits(ppmap));
>
> SearchDistance trDistance(ppmap);
>
>
> for (int i = 0; i < d_vertices.size(); ++i) {
>
> NeighborSearch search(searchTree, CPF::Point_3{dPoint[0], dPoint[1],
> dPoint[2]}, 1, 0, true, trDistance);
>
> for (const auto &searchResult : search) {
>
> if (searchResult.first != vertexId) {
>
> spdlog::get("cpf")->info("Brute force id = {}", vertexId);
>
> spdlog::get("cpf")->info("CGAL id = {}", searchResult.first);
>
> spdlog::get("cpf")->info("Brute force distance = {}", minDistance);
>
> spdlog::get("cpf")->info("CGAL distance = {}", searchResult.second);
>
> }
>
> }
>
> } }
>
> ```
>
>
> This return a los of vertices that are actually different:
>
>
> [2020-04-17 18:10:31.988] [cpf] [info] Brute force id = 386
>
> [2020-04-17 18:10:31.988] [cpf] [info] CGAL id = 84
>
> [2020-04-17 18:10:31.988] [cpf] [info] Brute force distance = 0.0
>
> [2020-04-17 18:10:31.988] [cpf] [info] CGAL distance = 0.041339367896971506
>
>
> I am not sure if I am using incorrectly the NeighborSearch or what is
> happening. Any ideas what might be wrong?
>
>
>
>

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss


Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Juan Jose Casafranca
Sure, find attached a minimal example and the data I am using.
Thanks


El vie., 17 abr. 2020 a las 18:41, Sebastien Loriot (GeometryFactory) (<[hidden email]>) escribió:
Could you please provide a minimal example that we could compile, run
and that would be showing the error?

At first sight, I don't see anything suspicious.

Thanks,

Sebastien.

On 4/17/20 6:25 PM, Juan Jose Casafranca wrote:
> I am trying to find the closest neighbor in a point cloud to a given
> point. I have followed the documentation example in which the points are
> stored in a user vector and the Tree only contains indices. Since the
> results were not what I was expecting, I am also computing the closest
> neighbor by brute force. Ill explain my code:
>
> I have a vector of points (coarsePoints) and I have a vector of indices
> (xC_vertices). The indices are the points in coarsePoints that can be
> marked as closest neighbor.
>
> This is how I compute the closest point using a brute force algorithm
> ```
>
> for (int i = 0; i < d_vertices.size(); ++i) {
>
> Real minDistance = std::numeric_limits<Real>::infinity();
>
> int vertexId = -1;
>
> const auto &dPoint = highResPoints[d_vertices[i]];
>
> const CPF::Point_3 pD{dPoint[0], dPoint[1], dPoint[2]};
>
> for (int j = 0; j < xC_vertices.size(); ++j) {
>
> const auto &xPoint = coarsePoints[xC_vertices[j]];
>
> const CPF::Point_3 pC{xPoint[0], xPoint[1], xPoint[2]};
>
>
> if (CGAL::squared_distance(pC, pD) < minDistance) {
>
> vertexId = j;
>
> minDistance = CGAL::squared_distance(pC, pD);
>
> }
>
> }
>
> ```
>
>
> When using CGAL, I have the following (based on the documentation example):
>
>
> ```
>
> class My_point_property_map
>
> {
>
> const sofa::helper::vector<CPF::SofaTypes::Point> &points;
>
> const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices;
>
> public:
>
> typedef CPF::Point_3 value_type;
>
> typedef value_type &reference;
>
> typedef std::size_t key_type;
>
> typedef boost::lvalue_property_map_tag category;
>
> My_point_property_map(const sofa::helper::vector<CPF::SofaTypes::Point>
> &pts,
>
> const std::vector<CPF::SofaTypes::Tetra::value_type> &xC_vertices_)
>
> : points(pts)
>
> , xC_vertices(xC_vertices_)
>
> {
>
> }
>
> value_type operator[](key_type k) const
>
> {
>
> const CPF::SofaTypes::Point &sPoint = points[xC_vertices[k]];
>
> return CPF::Point_3{sPoint[0], sPoint[1], sPoint[2]};
>
> }
>
> friend value_type get(const My_point_property_map &ppmap, key_type i) {
> return ppmap[i]; }
>
> };
>
>
> using SearchTraitsBase = CGAL::Search_traits_3<CPF::K>;
>
> using SearchTraits = CGAL::Search_traits_adapter<std::size_t,
> My_point_property_map, SearchTraitsBase>;
>
> using NeighborSearch = CGAL::Orthogonal_k_neighbor_search<SearchTraits>;
>
> using Tree = NeighborSearch::Tree;
>
> using Splitter = typename Tree::Splitter;
>
> using SearchDistance = typename NeighborSearch::Distance;
>
>
> int searchNeighbor() {
>
> My_point_property_map ppmap(coarsePoints, xC_vertices);
>
> Tree searchTree(boost::counting_iterator<int>(0),
>
> boost::counting_iterator<int>(xC_vertices.size()),
>
> Splitter(),
>
> SearchTraits(ppmap));
>
> SearchDistance trDistance(ppmap);
>
>
> for (int i = 0; i < d_vertices.size(); ++i) {
>
> NeighborSearch search(searchTree, CPF::Point_3{dPoint[0], dPoint[1],
> dPoint[2]}, 1, 0, true, trDistance);
>
> for (const auto &searchResult : search) {
>
> if (searchResult.first != vertexId) {
>
> spdlog::get("cpf")->info("Brute force id = {}", vertexId);
>
> spdlog::get("cpf")->info("CGAL id = {}", searchResult.first);
>
> spdlog::get("cpf")->info("Brute force distance = {}", minDistance);
>
> spdlog::get("cpf")->info("CGAL distance = {}", searchResult.second);
>
> }
>
> }
>
> } }
>
> ```
>
>
> This return a los of vertices that are actually different:
>
>
> [2020-04-17 18:10:31.988] [cpf] [info] Brute force id = 386
>
> [2020-04-17 18:10:31.988] [cpf] [info] CGAL id = 84
>
> [2020-04-17 18:10:31.988] [cpf] [info] Brute force distance = 0.0
>
> [2020-04-17 18:10:31.988] [cpf] [info] CGAL distance = 0.041339367896971506
>
>
> I am not sure if I am using incorrectly the NeighborSearch or what is
> happening. Any ideas what might be wrong?
>
>
>
>

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss



minimalExample.cpp (5K) Download Attachment
test.bin (173K) Download Attachment
Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Marc Glisse
On Fri, 17 Apr 2020, Juan Jose Casafranca wrote:

> Sure, find attached a minimal example and the data I am using.

It is far from minimal, there are several unnecessary #includes at least.

Your property map looks suspicious, with operator[] that returns a value
but reference that is a typedef for a reference type, and a category that
implies a reference.

Did you try -fsanitize=address?

--
Marc Glisse

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss


Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Juan Jose Casafranca
So if I understand correctly, the problem here is that the property map must model an LValuePropertyMap and that concept requires th operator[] to return a reference, right? I have tried storing a vector of K::Point_3 points and works fine. Is there anyway I can build the K::Point_3 when needed, instead, similarly to how I am doing it right now?

Thanks

El vie., 17 abr. 2020 a las 21:03, Marc Glisse (<[hidden email]>) escribió:
On Fri, 17 Apr 2020, Juan Jose Casafranca wrote:

> Sure, find attached a minimal example and the data I am using.

It is far from minimal, there are several unnecessary #includes at least.

Your property map looks suspicious, with operator[] that returns a value
but reference that is a typedef for a reference type, and a category that
implies a reference.

Did you try -fsanitize=address?

--
Marc Glisse

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss


Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Lokesh Kumar
Hi
Are You using any query algorithm for find the closest neighbor point.
Thanks

On Sat, Apr 18, 2020 at 1:08 AM Juan Jose Casafranca <[hidden email]> wrote:
So if I understand correctly, the problem here is that the property map must model an LValuePropertyMap and that concept requires th operator[] to return a reference, right? I have tried storing a vector of K::Point_3 points and works fine. Is there anyway I can build the K::Point_3 when needed, instead, similarly to how I am doing it right now?

Thanks

El vie., 17 abr. 2020 a las 21:03, Marc Glisse (<[hidden email]>) escribió:
On Fri, 17 Apr 2020, Juan Jose Casafranca wrote:

> Sure, find attached a minimal example and the data I am using.

It is far from minimal, there are several unnecessary #includes at least.

Your property map looks suspicious, with operator[] that returns a value
but reference that is a typedef for a reference type, and a category that
implies a reference.

Did you try -fsanitize=address?

--
Marc Glisse

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss




--
With Warm regards
Lokesh kumar
Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Marc Glisse
In reply to this post by Juan Jose Casafranca
On Fri, 17 Apr 2020, Juan Jose Casafranca wrote:

> So if I understand correctly, the problem here is that the property map
> must model an LValuePropertyMap and that concept requires th operator[] to
> return a reference, right? I have tried storing a vector of K::Point_3
> points and works fine. Is there anyway I can build the K::Point_3 when
> needed, instead, similarly to how I am doing it right now?

I am not sure, but LValuePropertyMap might be ok with some proxy type
(like vector<bool> uses), or not. Also, the doc says LValuePropertyMap is
required, but the code actually contains a fallback for more general
property maps, so you could try that.

It may be easier to work with a single point type. If you cannot use
CGAL's points everywhere, you could instead provide traits that define
suitable functions for MyPoint.

--
Marc Glisse

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss


Reply | Threaded
Open this post in threaded view
|

Re: Exact closest neighbor not working correctly

Juan Jose Casafranca
I have finally decide to use an LValuePropertyMap and create a std::vector<K::Point_3>, at least for now. Maybe on the future I will create the appropriate traits for the other Point type I am using. Thanks, it would had been really complex to debug this without your help.

El sáb., 18 abr. 2020 a las 9:05, Marc Glisse (<[hidden email]>) escribió:
On Fri, 17 Apr 2020, Juan Jose Casafranca wrote:

> So if I understand correctly, the problem here is that the property map
> must model an LValuePropertyMap and that concept requires th operator[] to
> return a reference, right? I have tried storing a vector of K::Point_3
> points and works fine. Is there anyway I can build the K::Point_3 when
> needed, instead, similarly to how I am doing it right now?

I am not sure, but LValuePropertyMap might be ok with some proxy type
(like vector<bool> uses), or not. Also, the doc says LValuePropertyMap is
required, but the code actually contains a fallback for more general
property maps, so you could try that.

It may be easier to work with a single point type. If you cannot use
CGAL's points everywhere, you could instead provide traits that define
suitable functions for MyPoint.

--
Marc Glisse

--
You are currently subscribed to cgal-discuss.
To unsubscribe or access the archives, go to
https://sympa.inria.fr/sympa/info/cgal-discuss