From 6bc84d6e645cc4aada745fe779e985a2f62718e7 Mon Sep 17 00:00:00 2001 From: fgodi Date: Mon, 29 Jun 2015 18:16:07 +0000 Subject: Ajout de Grid_cell.h git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/bottleneckDistance@660 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: e1b4ccdcb482c992c444adc91d67531bfa3e90dc --- src/Bottleneck/include/gudhi/Graph_matching.h | 4 +- src/Bottleneck/include/gudhi/Grid_cell.h | 188 +++++++++++++++++++++ .../include/gudhi/Layered_neighbors_finder.h | 2 +- src/Bottleneck/include/gudhi/Neighbors_finder.h | 2 +- .../include/gudhi/Persistence_diagrams_graph.h | 59 ++++--- .../include/gudhi/Planar_neighbors_finder.h | 4 +- 6 files changed, 231 insertions(+), 28 deletions(-) create mode 100644 src/Bottleneck/include/gudhi/Grid_cell.h (limited to 'src') diff --git a/src/Bottleneck/include/gudhi/Graph_matching.h b/src/Bottleneck/include/gudhi/Graph_matching.h index 17412e6c..2bcc6a61 100644 --- a/src/Bottleneck/include/gudhi/Graph_matching.h +++ b/src/Bottleneck/include/gudhi/Graph_matching.h @@ -74,13 +74,13 @@ private: void update(std::deque & path); }; -Graph_matching::Graph_matching() +inline Graph_matching::Graph_matching() : r(0.), v_to_u(G::size(), null_point_index()), unmatched_in_u() { for (int u_point_index = 0; u_point_index < G::size(); ++u_point_index) unmatched_in_u.emplace_back(u_point_index); } -Graph_matching& Graph_matching::operator=(const Graph_matching& m) { +inline Graph_matching& Graph_matching::operator=(const Graph_matching& m) { r = m.r; v_to_u = m.v_to_u; unmatched_in_u = m.unmatched_in_u; diff --git a/src/Bottleneck/include/gudhi/Grid_cell.h b/src/Bottleneck/include/gudhi/Grid_cell.h new file mode 100644 index 00000000..eee938cb --- /dev/null +++ b/src/Bottleneck/include/gudhi/Grid_cell.h @@ -0,0 +1,188 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Francois Godi + * + * Copyright (C) 2015 INRIA Sophia-Antipolis (France) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#ifndef SRC_BOTTLENECK_INCLUDE_GUDHI_GRID_CELL_H_ +#define SRC_BOTTLENECK_INCLUDE_GUDHI_GRID_CELL_H_ + +#include +#include +#include + +#include "Persistence_diagrams_graph.h" + +namespace Gudhi { + +namespace bottleneck { + +/** \internal \brief TODO + * + * \ingroup bottleneck_distance + */ +class Grid_cell { +public: + Grid_cell(double r); + void add(int v_point_index); + void remove(int v_point_index); + bool contains(int v_point_index) const; + int pull_center(); + int pull_xi(int u_point_index); + int pull_xd(int u_point_index); + int pull_yi(int u_point_index); + int pull_yd(int u_point_index); + int pull_xi_yi(int u_point_index); + int pull_xi_yd(int u_point_index); + int pull_xd_yi(int u_point_index); + int pull_xd_yd(int u_point_index); + +private: + double r; + std::set xi_order; + std::set yi_order; + struct Hidden_points_tree{int point; std::list hidden;}; + typedef std::map, G::Compare_x> Corner_tree; + Corner_tree xi_yi_order; + Corner_tree xi_yd_order; + Corner_tree xd_yi_order; + Corner_tree xd_yd_order; + void remove_aux(Corner_tree t, int v_point_index); + void build_xi_yi(); + void build_xi_yd(); + void build_xd_yi(); + void build_xd_yd(); +}; + + +inline Grid_cell::Grid_cell(double r) + : r(r), xi_order(G::Compare_x(r)), yi_order(G::Compare_y(r)), xi_yi_order(G::Compare_x(r)), + xi_yd_order(G::Compare_x(r)), xd_yi_order(G::Compare_x(r)), xd_yd_order(G::Compare_x(r)) {} + +inline void Grid_cell::add(int v_point_index){ + xi_order.emplace(v_point_index); +} + +inline bool Grid_cell::contains(int v_point_index) const{ + return (xi_order.count(v_point_index) > 0); +} + +inline void Grid_cell::remove_aux(Corner_tree t, int v_point_index){ + if(t.empty()) + return; + std::list hidden_points = t.at(v_point_index); + t.erase(v_point_index); + for(auto it = hidden_points.begin(); it != hidden_points.end(); ++it) + t.emplace(it->point,it->hidden); + +} + +inline void Grid_cell::remove(int v_point_index){ + xi_order.erase(v_point_index); + yi_order.erase(v_point_index); + remove_aux(xi_yi_order,v_point_index); + remove_aux(xi_yd_order,v_point_index); + remove_aux(xd_yi_order,v_point_index); + remove_aux(xd_yd_order,v_point_index); +} + +//factorization needed \/ \/ \/ + +inline int Grid_cell::pull_center(){ + if(xi_order.empty()) + return null_point_index(); + int v_point_index = *xi_order.begin(); + remove(v_point_index); + return v_point_index; +} + +inline int Grid_cell::pull_xi(int u_point_index){ + if(xi_order.empty()) + return null_point_index(); + int v_point_index = *xi_order.begin(); //! + if(G::distance(u_point_index,v_point_index)<=r){ + remove(v_point_index); + return v_point_index; + } + return null_point_index(); +} + +inline int Grid_cell::pull_xd(int u_point_index){ + if(xi_order.empty()) + return null_point_index(); + int v_point_index = *xi_order.rbegin(); //! + if(G::distance(u_point_index,v_point_index)<=r){ + remove(v_point_index); + return v_point_index; + } + return null_point_index(); +} + +inline int Grid_cell::pull_yi(int u_point_index){ + if(xi_order.empty()) + return null_point_index(); + if(yi_order.empty()) + for(auto it = xi_order.begin(); it!= xi_order.end(); ++it) + yi_order.emplace(*it); + int v_point_index = *yi_order.begin(); //! + if(G::distance(u_point_index,v_point_index)<=r){ + remove(v_point_index); + return v_point_index; + } + return null_point_index(); +} + +inline int Grid_cell::pull_yd(int u_point_index){ + if(xi_order.empty()) + return null_point_index(); + if(yi_order.empty()) + for(auto it = xi_order.begin(); it!= xi_order.end(); ++it) + yi_order.emplace(*it); + int v_point_index = *yi_order.rbegin(); //! + if(G::distance(u_point_index,v_point_index)<=r){ + remove(v_point_index); + return v_point_index; + } + return null_point_index(); +} + +inline int Grid_cell::pull_xi_yi(int u_point_index){ + if(xi_order.empty()) + return null_point_index(); + if(xi_yi_order.empty()) + build_xi_yi(); + auto it = xi_yi_order.upper_bound(u_point_index+G::size()); + if(it==xi_yi_order.cbegin()) //! + return null_point_index(); + it--; //! + int v_point_index = it->first; + for(auto it2=it->second.begin();it2!=it->second.end();it2++) + xi_yi_order.emplace(it2->point,it2->hidden); + if(G::distance(u_point_index,v_point_index)<=r){ + remove(v_point_index); + return v_point_index; + } + return null_point_index(); +} + +} // namespace bottleneck + +} // namespace Gudhi + +#endif // SRC_BOTTLENECK_INCLUDE_GUDHI_GRID_CELL_H_ diff --git a/src/Bottleneck/include/gudhi/Layered_neighbors_finder.h b/src/Bottleneck/include/gudhi/Layered_neighbors_finder.h index e6b7f30d..58805b86 100644 --- a/src/Bottleneck/include/gudhi/Layered_neighbors_finder.h +++ b/src/Bottleneck/include/gudhi/Layered_neighbors_finder.h @@ -54,7 +54,7 @@ private: std::vector neighbors_finder; }; -Layered_neighbors_finder::Layered_neighbors_finder(double r) : +inline Layered_neighbors_finder::Layered_neighbors_finder(double r) : r(r), neighbors_finder() { } inline void Layered_neighbors_finder::add(int v_point_index, int vlayer) { diff --git a/src/Bottleneck/include/gudhi/Neighbors_finder.h b/src/Bottleneck/include/gudhi/Neighbors_finder.h index ac3f8600..be81877a 100644 --- a/src/Bottleneck/include/gudhi/Neighbors_finder.h +++ b/src/Bottleneck/include/gudhi/Neighbors_finder.h @@ -58,7 +58,7 @@ private: bool contains(int v_point_index); }; -Neighbors_finder::Neighbors_finder(double r) : +inline Neighbors_finder::Neighbors_finder(double r) : r(r), planar_neighbors_f(r), projections_f() { } inline void Neighbors_finder::add(int v_point_index) { diff --git a/src/Bottleneck/include/gudhi/Persistence_diagrams_graph.h b/src/Bottleneck/include/gudhi/Persistence_diagrams_graph.h index 66e43145..aed328e2 100644 --- a/src/Bottleneck/include/gudhi/Persistence_diagrams_graph.h +++ b/src/Bottleneck/include/gudhi/Persistence_diagrams_graph.h @@ -28,6 +28,7 @@ #include #include #include +#include namespace Gudhi { @@ -60,10 +61,10 @@ public: static int size(); /** \internal \brief Returns the O(n^2) sorted distances between the points. */ static std::unique_ptr< std::vector > sorted_distances(); - /** \internal \brief Compare points regarding x coordinate. Use v_point_index for V points and u_point_index + G::size() for U points. */ - struct Compare_x{bool operator() (const int point_index_1, const int point_index_2) const;}; - /** \internal \brief Compare points regarding y coordinate. Use v_point_index for V points and u_point_index + G::size() for U points. */ - struct Compare_y{bool operator() (const int point_index_1, const int point_index_2) const;}; + /** \internal \brief Compare points regarding x%r coordinate. Use v_point_index for V points and u_point_index + G::size() for U points. */ + struct Compare_x{double r; Compare_x(double r); bool operator()(const int point_index_1, const int point_index_2) const;}; + /** \internal \brief Compare points regarding y%r coordinate. Use v_point_index for V points and u_point_index + G::size() for U points. */ + struct Compare_y{double r; Compare_y(double r);bool operator()(const int point_index_1, const int point_index_2) const;}; private: /** \internal \typedef \brief Internal_point is the internal points representation, indexes used outside. */ @@ -74,20 +75,20 @@ private: static Internal_point get_v_point(int v_point_index); }; -/** \internal \typedef \brief Alias used outside. */ +/** \internal \typedef \brief Shorter alias */ typedef Persistence_diagrams_graph G; // static initialization, seems to work but strange -std::vector Persistence_diagrams_graph::u = [] {return std::vector();}(); -std::vector Persistence_diagrams_graph::v = [] {return std::vector();}(); +std::vector G::u = [] {return std::vector();}(); +std::vector G::v = [] {return std::vector();}(); inline int null_point_index() { return -1; } template -void Persistence_diagrams_graph::initialize(const Persistence_diagram1 &diag1, - const Persistence_diagram2 &diag2, double e){ +inline void G::initialize(const Persistence_diagram1 &diag1, + const Persistence_diagram2 &diag2, double e){ u.clear(); v.clear(); for (auto it = diag1.cbegin(); it != diag1.cend(); ++it) @@ -100,25 +101,25 @@ void Persistence_diagrams_graph::initialize(const Persistence_diagram1 &diag1, swap(u, v); } -inline bool Persistence_diagrams_graph::on_the_u_diagonal(int u_point_index) { +inline bool G::on_the_u_diagonal(int u_point_index) { return u_point_index >= static_cast (u.size()); } -inline bool Persistence_diagrams_graph::on_the_v_diagonal(int v_point_index) { +inline bool G::on_the_v_diagonal(int v_point_index) { return v_point_index >= static_cast (v.size()); } -inline int Persistence_diagrams_graph::corresponding_point_in_u(int v_point_index) { +inline int G::corresponding_point_in_u(int v_point_index) { return on_the_v_diagonal(v_point_index) ? v_point_index - static_cast (v.size()) : v_point_index + static_cast (u.size()); } -inline int Persistence_diagrams_graph::corresponding_point_in_v(int u_point_index) { +inline int G::corresponding_point_in_v(int u_point_index) { return on_the_u_diagonal(u_point_index) ? u_point_index - static_cast (u.size()) : u_point_index + static_cast (v.size()); } -inline double Persistence_diagrams_graph::distance(int u_point_index, int v_point_index) { +inline double G::distance(int u_point_index, int v_point_index) { if (on_the_u_diagonal(u_point_index) && on_the_v_diagonal(v_point_index)) return 0; Internal_point p_u = get_u_point(u_point_index); @@ -126,11 +127,11 @@ inline double Persistence_diagrams_graph::distance(int u_point_index, int v_poin return std::max(std::fabs(p_u.first - p_v.first), std::fabs(p_u.second - p_v.second)); } -inline int Persistence_diagrams_graph::size() { +inline int G::size() { return static_cast (u.size() + v.size()); } -inline std::unique_ptr< std::vector > Persistence_diagrams_graph::sorted_distances() { +inline std::unique_ptr< std::vector > G::sorted_distances() { // could be optimized std::set sorted_distances; for (int u_point_index = 0; u_point_index < size(); ++u_point_index) @@ -140,7 +141,7 @@ inline std::unique_ptr< std::vector > Persistence_diagrams_graph::sorted return sd_up; } -inline Persistence_diagrams_graph::Internal_point Persistence_diagrams_graph::get_u_point(int u_point_index) { +inline G::Internal_point G::get_u_point(int u_point_index) { if (!on_the_u_diagonal(u_point_index)) return u.at(u_point_index); Internal_point projector = v.at(corresponding_point_in_v(u_point_index)); @@ -148,7 +149,7 @@ inline Persistence_diagrams_graph::Internal_point Persistence_diagrams_graph::ge return Internal_point(x, x); } -inline Persistence_diagrams_graph::Internal_point Persistence_diagrams_graph::get_v_point(int v_point_index) { +inline G::Internal_point G::get_v_point(int v_point_index) { if (!on_the_v_diagonal(v_point_index)) return v.at(v_point_index); Internal_point projector = u.at(corresponding_point_in_u(v_point_index)); @@ -156,16 +157,30 @@ inline Persistence_diagrams_graph::Internal_point Persistence_diagrams_graph::ge return Internal_point(x, x); } -inline bool Persistence_diagrams_graph::Compare_x::operator() (const int point_index_1, const int point_index_2) const{ +G::Compare_x::Compare_x(double r) + : r(r){ } + +G::Compare_y::Compare_y(double r) + : r(r){ } + +inline bool G::Compare_x::operator()(const int point_index_1, const int point_index_2) const{ G::Internal_point p1 = point_index_1 < G::size() ? G::get_v_point(point_index_1) : G::get_u_point(point_index_1 - G::size()); G::Internal_point p2 = point_index_2 < G::size() ? G::get_v_point(point_index_2) : G::get_u_point(point_index_2 - G::size()); - return p1.first < p2.first; + double x1 = fmod(p1.first,r); + double x2 = fmod(p2.first,r); + if(x1 == x2) + return point_index_1 > point_index_2; + return x1 < x2; } -inline bool Persistence_diagrams_graph::Compare_y::operator() (const int point_index_1, const int point_index_2) const{ +inline bool G::Compare_y::operator()(const int point_index_1, const int point_index_2) const{ G::Internal_point p1 = point_index_1 < G::size() ? G::get_v_point(point_index_1) : G::get_u_point(point_index_1 - G::size()); G::Internal_point p2 = point_index_2 < G::size() ? G::get_v_point(point_index_2) : G::get_u_point(point_index_2 - G::size()); - return p1.second < p2.second; + double y1 = fmod(p1.second,r); + double y2 = fmod(p2.second,r); + if(y1 == y2) + return point_index_1 > point_index_2; + return y1 < y2; } } // namespace bottleneck diff --git a/src/Bottleneck/include/gudhi/Planar_neighbors_finder.h b/src/Bottleneck/include/gudhi/Planar_neighbors_finder.h index 4d820c7f..e403735c 100644 --- a/src/Bottleneck/include/gudhi/Planar_neighbors_finder.h +++ b/src/Bottleneck/include/gudhi/Planar_neighbors_finder.h @@ -85,7 +85,7 @@ private: typedef Naive_pnf Planar_neighbors_finder; -Abstract_planar_neighbors_finder::Abstract_planar_neighbors_finder(double r) : +inline Abstract_planar_neighbors_finder::Abstract_planar_neighbors_finder(double r) : r(r) { } inline Abstract_planar_neighbors_finder::~Abstract_planar_neighbors_finder() {} @@ -100,7 +100,7 @@ inline std::unique_ptr< std::list > Abstract_planar_neighbors_finder::pull_ return all_pull; } -Naive_pnf::Naive_pnf(double r) : +inline Naive_pnf::Naive_pnf(double r) : Abstract_planar_neighbors_finder(r), candidates() { } inline void Naive_pnf::add(int v_point_index) { -- cgit v1.2.3