From aa67dab1eebe3cdba573741857051005ba72cc3b Mon Sep 17 00:00:00 2001 From: mcarrier Date: Mon, 8 May 2017 16:56:25 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@2406 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 2d857904595833667d469db97c746bbd8696eac4 --- src/Nerve_GIC/doc/Intro_graph_induced_complex.h | 130 ++++++++ src/Nerve_GIC/example/GIC.cpp | 77 +++++ src/Nerve_GIC/example/MapperDelta.cpp | 62 ++++ src/Nerve_GIC/example/Nerve.cpp | 63 ++++ src/Nerve_GIC/example/simple_GIC.cpp | 40 ++- src/Nerve_GIC/include/gudhi/GIC.h | 384 ++++++++++++++++-------- 6 files changed, 600 insertions(+), 156 deletions(-) create mode 100644 src/Nerve_GIC/doc/Intro_graph_induced_complex.h create mode 100644 src/Nerve_GIC/example/GIC.cpp create mode 100644 src/Nerve_GIC/example/MapperDelta.cpp create mode 100644 src/Nerve_GIC/example/Nerve.cpp diff --git a/src/Nerve_GIC/doc/Intro_graph_induced_complex.h b/src/Nerve_GIC/doc/Intro_graph_induced_complex.h new file mode 100644 index 00000000..0b51e345 --- /dev/null +++ b/src/Nerve_GIC/doc/Intro_graph_induced_complex.h @@ -0,0 +1,130 @@ +/* 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): Clément Maria, Pawel Dlotko, Vincent Rouvreau + * + * Copyright (C) 2016 INRIA + * + * 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 DOC_GRAPH_INDUCED_COMPLEX_INTRO_GRAPH_INDUCED_COMPLEX_H_ +#define DOC_GRAPH_INDUCED_COMPLEX_INTRO_GRAPH_INDUCED_COMPLEX_H_ + +namespace Gudhi { + +namespace graph_induced_complex { + +/** \defgroup graph_induced_complex Graph induced complex + * + * \author Mathieu Carrière + * + * @{ + * + * \section complexes Graph induced complexes (GIC) and Nerves + * + * GIC and Nerves are simplicial complexes built on top of a point cloud P. + * + * \subsection nervedefinition Nerve definition + * + * Assume you are given a cover C of your point cloud P, that is a set of subsets of P + * whose union is P itself. Then, the Nerve of this cover + * is the simplicial complex that has one k-simplex per k-fold intersection of cover elements. + * See also Wikipedia . + * + * \subsection nerveexample Example + * + * This example builds the Nerve of a point cloud sampled on a 3D human shape. + * The cover C comes from the preimages of intervals covering the height function. + * All intervals have the resolution (either the length or the number of the intervals) + * and gain (overlap percentage). + * + * \include + * + * When launching: + * + * \code $> + * \endcode + * + * the program output is: + * + * \include + * + * \section gicdefinition GIC definition + * + * Again, assume you are given a cover C of your point cloud P. Moreover, assume + * you are also given a graph G built on top of P. Then, for any clique in G + * whose nodes all belong to different elements of C, the GIC includes a corresponding + * simplex, whose dimension is the number of nodes in the clique minus one. + * + * \subsection gicexample Example + * + * This example builds the GIC of a point cloud sampled on a 3D human shape. + * The cover C comes from the preimages of intervals covering the height function, + * and the graph G comes from a Rips complex built with a threshold parameter. + * Note that if the gain is too big, the number of cliques increases a lot, + * which make the computation time much larger. + * + * \include + * + * When launching: + * + * \code $> + * \endcode + * + * the program output is: + * + * \include + * + * \subsection mapperdeltadefinition Mapper Delta + * + * If one restricts to the cliques in G whose nodes all belong to preimages of consecutive + * intervals (assuming the cover of the height function is minimal, i.e. no more than + * two intervals can intersect at a time), the GIC is of dimension one, i.e. a graph. + * We call this graph the Mapper Delta, since it is related to the usual Mapper (see + * this article ). + * + * \subsection mapperdeltaexample Example + * + * Mapper Delta comes with optimal selection for the Rips threshold, + * the resolution and the gain of the function cover. In this example, + * we compute the Mapper Delta of a point cloud sampled on a 3D human shape, + * where the graph G comes from a Rips complex with optimal threshold, + * and the cover C comes from the preimages of intervals covering the height function, + * with optimal resolution and gain. Note that optimal threshold, resolution and gain + * also exist for the Nerve of this cover. + * + * \include + * + * When launching: + * + * \code $> + * \endcode + * + * the program output is: + * + * \include + * + * + * \copyright GNU General Public License v3. + * \verbatim Contact: gudhi-users@lists.gforge.inria.fr \endverbatim + */ +/** @} */ // end defgroup graph_induced_complex + +} // namespace graph_induced_complex + +} // namespace Gudhi + +#endif // DOC_GRAPH_INDUCED_COMPLEX_INTRO_GRAPH_INDUCED_COMPLEX_H_ diff --git a/src/Nerve_GIC/example/GIC.cpp b/src/Nerve_GIC/example/GIC.cpp new file mode 100644 index 00000000..5161a46b --- /dev/null +++ b/src/Nerve_GIC/example/GIC.cpp @@ -0,0 +1,77 @@ +#include + +void usage(int nbArgs, char * const progName) { + std::cerr << "Error: Number of arguments (" << nbArgs << ") is not correct\n"; + std::cerr << "Usage: " << progName << " filename.off threshold coordinate resolution gain\n"; + std::cerr << " i.e.: " << progName << " ../../data/points/test.off 1.5 1 10 0.3 \n"; + exit(-1); // ----- >> +} + +int main(int argc, char **argv) { + if ((argc != 6) && (argc != 7)) usage(argc, (argv[0] - 1)); + + std::string off_file_name(argv[1]); + double threshold = atof(argv[2]); + //std::string func_file_name = argv[3]; + int coord = atoi(argv[3]); + double resolution = atof(argv[4]); + double gain = atof(argv[5]); + bool verb = 0; if(argc == 7) verb = 1; + + // Type definitions + using Graph_t = boost::adjacency_list < boost::vecS, boost::vecS, boost::undirectedS,\ + boost::property < vertex_filtration_t, Filtration_value >,\ + boost::property < edge_filtration_t, Filtration_value > >; + + // ---------------------------------------------------------------------------- + // Init of a graph induced complex from an OFF file + // ---------------------------------------------------------------------------- + + Gudhi::graph_induced_complex::Graph_induced_complex GIC; + GIC.set_verbose(verb); + + GIC.set_graph_from_automatic_rips(off_file_name); + //GIC.set_graph_from_rips(threshold, off_file_name); + //GIC.set_graph_from_OFF(off_file_name); + + GIC.set_function_from_coordinate(coord, off_file_name); + //GIC.set_function_from_file(func_file_name); + + GIC.set_color_from_coordinate(coord, off_file_name); + //GIC.set_color_from_file(func_file_name); + + GIC.set_automatic_resolution_for_GICMAP(); + GIC.set_gain(); + GIC.set_cover_from_function(1); + + //GIC.find_GIC_simplices(); + //GIC.find_Nerve_simplices(); + GIC.find_GICMAP_simplices_with_functional_minimal_cover(); + + GIC.plot_with_KeplerMapper(); + + Simplex_tree stree; GIC.create_complex(stree); + + std::streambuf* streambufffer = std::cout.rdbuf(); + std::ostream output_stream(streambufffer); + + // ---------------------------------------------------------------------------- + // Display information about the graph induced complex + // ---------------------------------------------------------------------------- + + if(verb){ + output_stream << "Graph induced complex is of dimension " << stree.dimension() << + " - " << stree.num_simplices() << " simplices - " << + stree.num_vertices() << " vertices." << std::endl; + + output_stream << "Iterator on graph induced complex simplices" << std::endl; + for (auto f_simplex : stree.filtration_simplex_range()) { + for (auto vertex : stree.simplex_vertex_range(f_simplex)) { + output_stream << vertex << " "; + } + output_stream << std::endl; + } + } + + return 0; +} diff --git a/src/Nerve_GIC/example/MapperDelta.cpp b/src/Nerve_GIC/example/MapperDelta.cpp new file mode 100644 index 00000000..1f8f1582 --- /dev/null +++ b/src/Nerve_GIC/example/MapperDelta.cpp @@ -0,0 +1,62 @@ +#include + +void usage(int nbArgs, char * const progName) { + std::cerr << "Error: Number of arguments (" << nbArgs << ") is not correct\n"; + std::cerr << "Usage: " << progName << " filename.off coordinate \n"; + std::cerr << " i.e.: " << progName << " ../../data/points/human.off 2 --v \n"; + exit(-1); // ----- >> +} + +int main(int argc, char **argv) { + if ((argc != 6) && (argc != 7)) usage(argc, (argv[0] - 1)); + + std::string off_file_name(argv[1]); + int coord = atoi(argv[2]); + bool verb = 0; if(argc == 4) verb = 1; + + // Type definitions + using Graph_t = boost::adjacency_list < boost::vecS, boost::vecS, boost::undirectedS,\ + boost::property < vertex_filtration_t, Filtration_value >,\ + boost::property < edge_filtration_t, Filtration_value > >; + + // --------------------------------------- + // Init of a Mapper Delta from an OFF file + // --------------------------------------- + + Gudhi::graph_induced_complex::Graph_induced_complex GIC; + GIC.set_verbose(verb); + + GIC.set_graph_from_automatic_rips(off_file_name); + GIC.set_function_from_coordinate(coord, off_file_name); + GIC.set_color_from_coordinate(coord, off_file_name); + GIC.set_automatic_resolution_for_GICMAP(); + GIC.set_gain(); + GIC.set_cover_from_function(1); + GIC.find_GICMAP_simplices_with_functional_minimal_cover(); + GIC.plot_with_KeplerMapper(); + + Simplex_tree stree; GIC.create_complex(stree); + + std::streambuf* streambufffer = std::cout.rdbuf(); + std::ostream output_stream(streambufffer); + + // ------------------------------------------ + // Display information about the Mapper Delta + // ------------------------------------------ + + if(verb){ + output_stream << "Mapper Delta is of dimension " << stree.dimension() << + " - " << stree.num_simplices() << " simplices - " << + stree.num_vertices() << " vertices." << std::endl; + + output_stream << "Iterator on Mapper Delta simplices" << std::endl; + for (auto f_simplex : stree.filtration_simplex_range()) { + for (auto vertex : stree.simplex_vertex_range(f_simplex)) { + output_stream << vertex << " "; + } + output_stream << std::endl; + } + } + + return 0; +} diff --git a/src/Nerve_GIC/example/Nerve.cpp b/src/Nerve_GIC/example/Nerve.cpp new file mode 100644 index 00000000..adcc715d --- /dev/null +++ b/src/Nerve_GIC/example/Nerve.cpp @@ -0,0 +1,63 @@ +#include + +void usage(int nbArgs, char * const progName) { + std::cerr << "Error: Number of arguments (" << nbArgs << ") is not correct\n"; + std::cerr << "Usage: " << progName << " filename.off coordinate resolution gain --v \n"; + std::cerr << " i.e.: " << progName << " ../../data/points/human.off 2 1 0.3 --v \n"; + exit(-1); // ----- >> +} + +int main(int argc, char **argv) { + if ((argc != 6) && (argc != 7)) usage(argc, (argv[0] - 1)); + + std::string off_file_name(argv[1]); + int coord = atoi(argv[2]); + double resolution = atof(argv[3]); + double gain = atof(argv[4]); + bool verb = 0; if(argc == 6) verb = 1; + + // Type definitions + using Graph_t = boost::adjacency_list < boost::vecS, boost::vecS, boost::undirectedS,\ + boost::property < vertex_filtration_t, Filtration_value >,\ + boost::property < edge_filtration_t, Filtration_value > >; + + // -------------------------------- + // Init of a Nerve from an OFF file + // -------------------------------- + + Gudhi::graph_induced_complex::Graph_induced_complex GIC; + GIC.set_verbose(verb); + + GIC.set_function_from_coordinate(coord, off_file_name); + GIC.set_color_from_coordinate(coord, off_file_name); + GIC.set_resolution_double(resolution); + GIC.set_gain(gain); + GIC.set_cover_from_function(1); + GIC.find_Nerve_simplices(); + GIC.plot_with_KeplerMapper(); + + Simplex_tree stree; GIC.create_complex(stree); + + std::streambuf* streambufffer = std::cout.rdbuf(); + std::ostream output_stream(streambufffer); + + // ---------------------------------------------------------------------------- + // Display information about the graph induced complex + // ---------------------------------------------------------------------------- + + if(verb){ + output_stream << "Nerve is of dimension " << stree.dimension() << + " - " << stree.num_simplices() << " simplices - " << + stree.num_vertices() << " vertices." << std::endl; + + output_stream << "Iterator on Nerve simplices" << std::endl; + for (auto f_simplex : stree.filtration_simplex_range()) { + for (auto vertex : stree.simplex_vertex_range(f_simplex)) { + output_stream << vertex << " "; + } + output_stream << std::endl; + } + } + + return 0; +} diff --git a/src/Nerve_GIC/example/simple_GIC.cpp b/src/Nerve_GIC/example/simple_GIC.cpp index 8b3aecb8..5161a46b 100644 --- a/src/Nerve_GIC/example/simple_GIC.cpp +++ b/src/Nerve_GIC/example/simple_GIC.cpp @@ -2,7 +2,7 @@ void usage(int nbArgs, char * const progName) { std::cerr << "Error: Number of arguments (" << nbArgs << ") is not correct\n"; - std::cerr << "Usage: " << progName << " filename.off threshold coordinate resolution gain [ouput_file.txt]\n"; + std::cerr << "Usage: " << progName << " filename.off threshold coordinate resolution gain\n"; std::cerr << " i.e.: " << progName << " ../../data/points/test.off 1.5 1 10 0.3 \n"; exit(-1); // ----- >> } @@ -16,6 +16,7 @@ int main(int argc, char **argv) { int coord = atoi(argv[3]); double resolution = atof(argv[4]); double gain = atof(argv[5]); + bool verb = 0; if(argc == 7) verb = 1; // Type definitions using Graph_t = boost::adjacency_list < boost::vecS, boost::vecS, boost::undirectedS,\ @@ -27,8 +28,9 @@ int main(int argc, char **argv) { // ---------------------------------------------------------------------------- Gudhi::graph_induced_complex::Graph_induced_complex GIC; + GIC.set_verbose(verb); - GIC.set_graph_from_automatic_rips(100, off_file_name); + GIC.set_graph_from_automatic_rips(off_file_name); //GIC.set_graph_from_rips(threshold, off_file_name); //GIC.set_graph_from_OFF(off_file_name); @@ -38,9 +40,9 @@ int main(int argc, char **argv) { GIC.set_color_from_coordinate(coord, off_file_name); //GIC.set_color_from_file(func_file_name); - resolution = GIC.set_automatic_resolution_for_GICMAP(); - - GIC.set_cover_from_function(resolution,gain,1); + GIC.set_automatic_resolution_for_GICMAP(); + GIC.set_gain(); + GIC.set_cover_from_function(1); //GIC.find_GIC_simplices(); //GIC.find_Nerve_simplices(); @@ -50,34 +52,26 @@ int main(int argc, char **argv) { Simplex_tree stree; GIC.create_complex(stree); - std::streambuf* streambufffer; - std::ofstream ouput_file_stream; - - if (argc == 7) { - ouput_file_stream.open(std::string(argv[4])); - streambufffer = ouput_file_stream.rdbuf(); - } else { - streambufffer = std::cout.rdbuf(); - } - + std::streambuf* streambufffer = std::cout.rdbuf(); std::ostream output_stream(streambufffer); // ---------------------------------------------------------------------------- // Display information about the graph induced complex // ---------------------------------------------------------------------------- - output_stream << "Graph induced complex is of dimension " << stree.dimension() << + + if(verb){ + output_stream << "Graph induced complex is of dimension " << stree.dimension() << " - " << stree.num_simplices() << " simplices - " << stree.num_vertices() << " vertices." << std::endl; - output_stream << "Iterator on graph induced complex simplices" << std::endl; - for (auto f_simplex : stree.filtration_simplex_range()) { - for (auto vertex : stree.simplex_vertex_range(f_simplex)) { - output_stream << vertex << " "; + output_stream << "Iterator on graph induced complex simplices" << std::endl; + for (auto f_simplex : stree.filtration_simplex_range()) { + for (auto vertex : stree.simplex_vertex_range(f_simplex)) { + output_stream << vertex << " "; + } + output_stream << std::endl; } - output_stream << std::endl; } - ouput_file_stream.close(); - return 0; } diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 5ee029b7..abd833e0 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -70,33 +70,36 @@ namespace graph_induced_complex { * \ingroup graph_induced_complex * * \details - * + * The data structure is a simplicial complex, representing a + * Graph Induced simplicial Complex (GIC) or a Nerve, + * and whose simplices are computed with a cover C of a point + * cloud P, which often comes from the preimages of intervals + * covering the image of a function f defined on P. + * These intervals are parameterized by their resolution + * (either their length or their number) + * and their gain (percentage of overlap). + * To compute a GIC, one also needs a graph G built on top of P, + * whose cliques with vertices belonging to different elements of C + * correspond to the simplices of the GIC. * */ class Graph_induced_complex { private: - typedef int Cover_t; - - private: + bool verbose; + typedef int Cover_t; std::vector > simplices; - - private: std::map > cover; - - private: - int maximal_dim; int data_dimension; - - private: + int maximal_dim; + int data_dimension; std::map cover_fct; - - private: - std::map > cover_color; - - private: + std::map > cover_color; Simplex_tree<> st; std::map > adjacency_matrix; + int resolution_int; + double resolution_double; + double gain; // Simplex comparator private: @@ -144,40 +147,58 @@ class Graph_induced_complex { } } + public: + void set_verbose(bool verb = 0){verbose = verb;} + // ******************************************************************************************************************* // Graphs. // ******************************************************************************************************************* public: // Set graph from file. - void set_graph_from_file(const std::string& graph_file_name){ - int neighb; int vid; std::ifstream input(graph_file_name); std::string line; std::vector edge(2); - while(std::getline(input,line)){ - std::stringstream stream(line); stream >> vid; edge[0] = vid; - while(stream >> neighb){edge[1] = neighb; st.insert_simplex_and_subfaces(edge);} - } - } + /** \brief Creates the graph G from a file containing the edges. + * + * @param[in] graph_file_name name of the input graph file. + * + */ + void set_graph_from_file(const std::string& graph_file_name){ + int neighb; int vid; std::ifstream input(graph_file_name); std::string line; std::vector edge(2); + while(std::getline(input,line)){ + std::stringstream stream(line); stream >> vid; edge[0] = vid; + while(stream >> neighb){edge[1] = neighb; st.insert_simplex_and_subfaces(edge);} + } + } public: // Set graph from OFF file. - void set_graph_from_OFF(const std::string& off_file_name){ - int numpts, numedges, numfaces, i; std::vector edge(2); double x; int num; std::vector simplex; - std::ifstream input(off_file_name); std::string line; getline(input, line); - input >> numpts; input >> numfaces; input >> numedges; - i = 0; while(i < numpts){input >> x; input >> x; input >> x; i++;} - i = 0; while(i < numfaces){ - simplex.clear(); input >> num; - for(int j = 0; j < num; j++){int k; input >> k; simplex.push_back(k);} - for(int j = 0; j < num; j++){ - for(int k = j+1; k < num; k++){ - edge[0] = simplex[j]; edge[1] = simplex[k]; - st.insert_simplex_and_subfaces(edge); - } - } - i++; - } - - } + /** \brief Creates the graph G from the triangulation given by an .OFF file. + * + * @param[in] off_file_name name of the input .OFF file. + * + */ + void set_graph_from_OFF(const std::string& off_file_name){ + int numpts, numedges, numfaces, i; std::vector edge(2); double x; int num; std::vector simplex; + std::ifstream input(off_file_name); std::string line; getline(input, line); + input >> numpts; input >> numfaces; input >> numedges; + i = 0; while(i < numpts){input >> x; input >> x; input >> x; i++;} + i = 0; while(i < numfaces){ + simplex.clear(); input >> num; + for(int j = 0; j < num; j++){int k; input >> k; simplex.push_back(k);} + for(int j = 0; j < num; j++){ + for(int k = j+1; k < num; k++){ + edge[0] = simplex[j]; edge[1] = simplex[k]; + st.insert_simplex_and_subfaces(edge); + } + } + i++; + } + } public: // Set graph from Rips complex. + /** \brief Creates the graph G from a Rips complex. + * + * @param[in] threshold threshold value for the Rips complex. + * @param[in] off_file_name name of the input .OFF file. + * + */ void set_graph_from_rips(const double& threshold, const std::string& off_file_name){ Points_off_reader off_reader(off_file_name); Rips_complex rips_complex_from_points(off_reader.get_point_cloud(), threshold, Euclidean_distance()); @@ -185,15 +206,21 @@ class Graph_induced_complex { } public: // Automatic tuning of Rips complex. - void set_graph_from_automatic_rips(const int& N, const std::string& off_file_name){ + /** \brief Creates the graph G from a Rips complex whose threshold value is automatically tuned with subsampling. + * + * @param[in] off_file_name name of the input .OFF file. + * @param[in] N number of subsampling iteration (default value 100). + * + */ + void set_graph_from_automatic_rips(const std::string& off_file_name, int N = 100){ Points_off_reader off_reader(off_file_name); int n = off_reader.get_point_cloud().size(); int m = floor(n/pow(log(n)/log(CONSTANT),1+ETA)); m = std::min(m,n-1); std::vector samples(m); double delta = 0; int dim = off_reader.get_point_cloud()[0].size(); data_dimension = dim; - std::cout << n << " points in R^" << dim << std::endl; - std::cout << "Subsampling " << m << " points" << std::endl; + if(verbose) std::cout << n << " points in R^" << dim << std::endl; + if(verbose) std::cout << "Subsampling " << m << " points" << std::endl; std::vector > dist(n); std::vector dumb(n); for(int i = 0; i < n; i++) dist[i] = dumb; @@ -203,7 +230,7 @@ class Graph_induced_complex { std::ifstream input(distances, std::ios::out | std::ios::binary); if(input.good()){ - std::cout << "Reading distances..." << std::endl; + if(verbose) std::cout << "Reading distances..." << std::endl; for(int i = 0; i < n; i++){ for (int j = 0; j < n; j++){ input.read((char*) &d,8); dist[i][j] = d; @@ -212,17 +239,17 @@ class Graph_induced_complex { input.close(); } else{ - std::cout << "Computing distances..." << std::endl; + if(verbose) std::cout << "Computing distances..." << std::endl; input.close(); std::ofstream output(distances, std::ios::out | std::ios::binary); for(int i = 0; i < n; i++){ if( (int) floor( 100*(i*1.0+1)/(n*1.0) ) %10 == 0 ) - std::cout << "\r" << floor( 100*(i*1.0+1)/(n*1.0) ) << "%" << std::flush; + if(verbose) std::cout << "\r" << floor( 100*(i*1.0+1)/(n*1.0) ) << "%" << std::flush; for (int j = 0; j < n; j++){ double dis = 0; for(int k = 0; k < dim; k++) dis += pow(off_reader.get_point_cloud()[i][k]-off_reader.get_point_cloud()[j][k],2); dis = std::sqrt(dis); dist[i][j] = dis; output.write((char*) &dis,8); } } - output.close(); std::cout << std::endl; + output.close(); if(verbose) std::cout << std::endl; } for(int i = 0; i < n; i++){ @@ -245,7 +272,7 @@ class Graph_induced_complex { } - std::cout << "delta = " << delta << std::endl; + if(verbose) std::cout << "delta = " << delta << std::endl; Rips_complex rips_complex_from_points(off_reader.get_point_cloud(), delta, Euclidean_distance()); rips_complex_from_points.create_complex(st, 1); @@ -257,6 +284,11 @@ class Graph_induced_complex { // ******************************************************************************************************************* public: // Set function from file. + /** \brief Creates the function f from a file containing the function values. + * + * @param[in] func_file_name name of the input function file. + * + */ void set_function_from_file(const std::string& func_file_name){ int vertex_id = 0; std::ifstream input(func_file_name); std::string line; double f; while(std::getline(input,line)){ @@ -266,6 +298,12 @@ class Graph_induced_complex { } public: // Set function from kth coordinate + /** \brief Creates the function f from the k-th coordinate of the point cloud P. + * + * @param[in] k coordinate to use (start at 0). + * @param[in] off_file_name name of the input .OFF file. + * + */ void set_function_from_coordinate(const int& k, const std::string& off_file_name){ Points_off_reader off_reader(off_file_name); int numpts = off_reader.get_point_cloud().size(); @@ -273,6 +311,11 @@ class Graph_induced_complex { } public: // Set function from vector. + /** \brief Creates the function f from a vector stored in memory. + * + * @param[in] function input vector of values. + * + */ void set_function_from_vector(const std::vector& function){ for(int i = 0; i < function.size(); i++) func.insert(std::pair(i, function[i])); } @@ -282,6 +325,12 @@ class Graph_induced_complex { // ******************************************************************************************************************* public: // Set cover from file. + /** \brief Creates the cover C from a file containing the cover elements of each point (the order has to be the same + * as in the input file!). + * + * @param[in] cover_file_name name of the input cover file. + * + */ void set_cover_from_file(const std::string& cover_file_name){ int vertex_id = 0; Cover_t cov; std::vector cov_elts, cov_number; std::ifstream input(cover_file_name); std::string line; @@ -296,7 +345,9 @@ class Graph_induced_complex { } public: // Automatic tuning of resolution for Mapper Delta. - double set_automatic_resolution_for_GICMAP(){ + /** \brief Computes the optimal length of intervals for a Mapper Delta. + */ + void set_automatic_resolution_for_GICMAP(){ double reso = 0; for (auto simplex : st.complex_simplex_range()) { if(st.dimension(simplex) == 1){ @@ -305,12 +356,34 @@ class Graph_induced_complex { reso = std::max(reso, std::abs(func[vertices[0]] - func[vertices[1]])); } } - std::cout << "resolution = " << reso << std::endl; - return reso; + if(verbose) std::cout << "resolution = " << reso << std::endl; + resolution_double = reso; } + public: + /** \brief Sets a length of intervals from a value stored in memory. + * + * @param[in] reso length of intervals. + * + */ + void set_resolution_double(double reso){resolution_double = reso;} + /** \brief Sets a number of intervals from a value stored in memory. + * + * @param[in] reso number of intervals. + * + */ + void set_resolution_int(int reso){resolution_int = reso;} + /** \brief Sets a gain from a value stored in memory (default value 0.3). + * + * @param[in] g gain. + * + */ + void set_gain(double g = 0.3){gain = g;} + public: // Automatic tuning of resolution for Mapper Point. - double set_automatic_resolution_for_MAP(const double& gain){ + /** \brief Computes the optimal length of intervals for a standard Mapper. + */ + void set_automatic_resolution_for_MAP(const double& gain){ double reso = 0; for (auto simplex : st.complex_simplex_range()) { if(st.dimension(simplex) == 1){ @@ -319,40 +392,49 @@ class Graph_induced_complex { reso = std::max(reso, (std::abs(func[vertices[0]] - func[vertices[1]]))/gain); } } - return reso; + if(verbose) std::cout << "resolution = " << reso << std::endl; + resolution_double = reso; } public: // Set cover with preimages of function. - void set_cover_from_function(const double& resolution, const double& gain, const bool& token){ + /** \brief Creates a cover C from the preimages of the function f. + * + * @param[in] token boolean specifying whether we use the length or the number of intervals for the cover of im(f). + * + */ + void set_cover_from_function(const bool& token){ // Read function values and compute min and max std::map::iterator it; double maxf, minf; minf = std::numeric_limits::max(); maxf = std::numeric_limits::min(); for(it = func.begin(); it != func.end(); it++){minf = std::min(minf, it->second); maxf = std::max(maxf, it->second);} - int num_pts = func.size(); std::cout << "Min function value = " << minf << " and Max function value = " << maxf << std::endl; + int num_pts = func.size(); if(verbose) std::cout << "Min function value = " << minf << " and Max function value = " << maxf << std::endl; // Compute cover of im(f) std::vector > intervals; int res; - if(!token){ - double incr = (maxf-minf)/resolution; double x = minf; double alpha = (incr*gain)/(2-2*gain); + + if(!token){ // Case we use an integer for the number of intervals. + double incr = (maxf-minf)/resolution_int; double x = minf; double alpha = (incr*gain)/(2-2*gain); double y = minf + incr + alpha; std::pair interm(x,y); intervals.push_back(interm); - for(int i = 1; i < resolution-1; i++){ + for(int i = 1; i < resolution_int-1; i++){ x = minf + i*incr - alpha; y = minf + (i+1)*incr + alpha; std::pair inter(x,y); intervals.push_back(inter); } - x = minf + (resolution-1)*incr - alpha; y = maxf; + x = minf + (resolution_int-1)*incr - alpha; y = maxf; std::pair interM(x,y); intervals.push_back(interM); res = intervals.size(); - for(int i = 0; i < res; i++) std::cout << "Interval " << i << " = [" << intervals[i].first << ", " << intervals[i].second << "]" << std::endl; + for(int i = 0; i < res; i++) + if(verbose) std::cout << "Interval " << i << " = [" << intervals[i].first << ", " << intervals[i].second << "]" << std::endl; } - else{ - double x = minf; double y = x + resolution; - while(y <= maxf && maxf - (y-gain*resolution) >= resolution){ + + else{ // Case we use a double for the length of the intervals. + double x = minf; double y = x + resolution_double; + while(y <= maxf && maxf - (y-gain*resolution_double) >= resolution_double){ std::pair inter(x,y); intervals.push_back(inter); - x = y - gain*resolution; - y = x + resolution; + x = y - gain*resolution_double; + y = x + resolution_double; } std::pair interM(x,maxf); intervals.push_back(interM); res = intervals.size(); - for(int i = 0; i < res; i++) std::cout << "Interval " << i << " = [" << intervals[i].first << ", " << intervals[i].second << "]" << std::endl; + for(int i = 0; i < res; i++) if(verbose) std::cout << "Interval " << i << " = [" << intervals[i].first << ", " << intervals[i].second << "]" << std::endl; } // Sort points according to function values @@ -420,14 +502,14 @@ class Graph_induced_complex { } // Compute the connected components with DFS - std::map visit; std::cout << "Preimage of interval " << i << std::endl; + std::map visit; if(verbose) std::cout << "Preimage of interval " << i << std::endl; for(std::map >::iterator it = prop.begin(); it != prop.end(); it++) visit.insert(std::pair(it->first, false)); if (!(prop.empty())){ for(std::map >::iterator it = prop.begin(); it != prop.end(); it++){ if ( !(visit[it->first]) ){ std::vector cc; cc.clear(); - dfs(prop,it->first,cc,visit); int cci = cc.size(); std::cout << "one CC with " << cci << " points, "; + dfs(prop,it->first,cc,visit); int cci = cc.size(); if(verbose) std::cout << "one CC with " << cci << " points, "; double average_col = 0; for(int j = 0; j < cci; j++){cover[cc[j]].push_back(id); average_col += func_color[cc[j]]/cci;} cover_fct[id] = i; cover_color[id] = std::pair(cci,average_col); @@ -435,7 +517,7 @@ class Graph_induced_complex { } } } - std::cout << std::endl; + if(verbose) std::cout << std::endl; } maximal_dim = id; @@ -447,81 +529,107 @@ class Graph_induced_complex { // ******************************************************************************************************************* public: // Set color from file. - void set_color_from_file(const std::string& color_file_name){ - int vertex_id = 0; std::ifstream input(color_file_name); std::string line; double f; - while(std::getline(input,line)){ - std::stringstream stream(line); stream >> f; - func_color.insert(std::pair(vertex_id, f)); vertex_id++; - } - } + /** \brief Computes the function used to color the nodes of the simplicial complex from a file containing the function values. + * + * @param[in] color_file_name name of the input color file. + * + */ + void set_color_from_file(const std::string& color_file_name){ + int vertex_id = 0; std::ifstream input(color_file_name); std::string line; double f; + while(std::getline(input,line)){ + std::stringstream stream(line); stream >> f; + func_color.insert(std::pair(vertex_id, f)); vertex_id++; + } + } public: // Set color from kth coordinate - void set_color_from_coordinate(const int& k, const std::string& off_file_name){ - Points_off_reader off_reader(off_file_name); - int numpts = off_reader.get_point_cloud().size(); - for(int i = 0; i < numpts; i++) func_color.insert(std::pair(i,off_reader.get_point_cloud()[i][k])); - } + /** \brief Computes the function used to color the nodes of the simplicial complex from the k-th coordinate. + * + * @param[in] k coordinate to use (start at 0). + * @param[in] off_file_name name of the input .OFF file. + * + */ + void set_color_from_coordinate(int k = 0, const std::string& off_file_name){ + Points_off_reader off_reader(off_file_name); + int numpts = off_reader.get_point_cloud().size(); + for(int i = 0; i < numpts; i++) func_color.insert(std::pair(i,off_reader.get_point_cloud()[i][k])); + } public: // Set color from vector. - void set_color_from_vector(const std::vector& color){ - for(int i = 0; i < color.size(); i++) func_color.insert(std::pair(i, color[i])); - } + /** \brief Computes the function used to color the nodes of the simplicial complex from a vector stored in memory. + * + * @param[in] color input vector of values. + * + */ + void set_color_from_vector(const std::vector& color){ + for(int i = 0; i < color.size(); i++) func_color.insert(std::pair(i, color[i])); + } public: // Create a .dot file that can be compiled with neato to produce a .pdf file - void plot_with_neato(){ - char mapp[11] = "mapper.dot"; std::ofstream graphic(mapp); graphic << "graph Mapper {" << std::endl; - double maxv, minv; maxv = std::numeric_limits::min(); minv = std::numeric_limits::max(); - for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++){ - maxv = std::max(maxv, iit->second.second); minv = std::min(minv, iit->second.second); - } - int k = 0; std::vector nodes; nodes.clear(); - for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++){ - if(iit->second.first > MASK){ - nodes.push_back(iit->first); - graphic << iit->first << "[shape=circle fontcolor=black color=black label=\"" \ - << iit->first << ":" << iit->second.first << "\" style=filled fillcolor=\"" \ - << (1-(maxv-iit->second.second)/(maxv-minv))*0.6 << ", 1, 1\"]" << std::endl; - k++; - } - } - int ke = 0; int num_simplices = simplices.size(); - for (int i = 0; i < num_simplices; i++) - if (simplices[i].size() == 2) - if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK){ - graphic << " " << simplices[i][0] << " -- " << simplices[i][1] << " [weight=15];" << std::endl; ke++;} - graphic << "}"; graphic.close(); - } + /** \brief Creates a .dot file for neato once the simplicial complex is computed to get a .pdf output. + */ + void plot_with_neato(){ + char mapp[11] = "mapper.dot"; std::ofstream graphic(mapp); graphic << "graph Mapper {" << std::endl; + double maxv, minv; maxv = std::numeric_limits::min(); minv = std::numeric_limits::max(); + for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++){ + maxv = std::max(maxv, iit->second.second); minv = std::min(minv, iit->second.second); + } + int k = 0; std::vector nodes; nodes.clear(); + for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++){ + if(iit->second.first > MASK){ + nodes.push_back(iit->first); + graphic << iit->first << "[shape=circle fontcolor=black color=black label=\"" \ + << iit->first << ":" << iit->second.first << "\" style=filled fillcolor=\"" \ + << (1-(maxv-iit->second.second)/(maxv-minv))*0.6 << ", 1, 1\"]" << std::endl; + k++; + } + } + int ke = 0; int num_simplices = simplices.size(); + for (int i = 0; i < num_simplices; i++) + if (simplices[i].size() == 2) + if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK){ + graphic << " " << simplices[i][0] << " -- " << simplices[i][1] << " [weight=15];" << std::endl; ke++;} + graphic << "}"; graphic.close(); + } public: // Create a .txt file that can be compiled with KeplerMapper to produce a .html file - void plot_with_KeplerMapper(){ - - int num_simplices = simplices.size(); int num_edges = 0; - char mapp[11] = "mapper.txt"; std::ofstream graphic(mapp); - for (int i = 0; i < num_simplices; i++) - if (simplices[i].size() == 2) - if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK) - num_edges++; - - graphic << "Cloud" << std::endl; - graphic << "Function" << std::endl; - graphic << 0 << " " << 0 << std::endl; - graphic << cover_color.size() << " " << num_edges << std::endl; - - for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++) - graphic << iit->first << " " << iit->second.second << " " << iit->second.first << std::endl; - - for (int i = 0; i < num_simplices; i++) - if (simplices[i].size() == 2) - if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK) - graphic << simplices[i][0] << " " << simplices[i][1] << std::endl; - graphic.close(); - } + /** \brief Creates a .html file for KeplerMapper once the simplicial complex is computed to get a nice visualization in browser. + */ + void plot_with_KeplerMapper(){ + + int num_simplices = simplices.size(); int num_edges = 0; + char mapp[11] = "mapper.txt"; std::ofstream graphic(mapp); + for (int i = 0; i < num_simplices; i++) + if (simplices[i].size() == 2) + if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK) + num_edges++; + + graphic << "Cloud" << std::endl; + graphic << "Function" << std::endl; + graphic << 0 << " " << 0 << std::endl; + graphic << cover_color.size() << " " << num_edges << std::endl; + + for (std::map >::iterator iit = cover_color.begin(); iit != cover_color.end(); iit++) + graphic << iit->first << " " << iit->second.second << " " << iit->second.first << std::endl; + + for (int i = 0; i < num_simplices; i++) + if (simplices[i].size() == 2) + if (cover_color[simplices[i][0]].first > MASK && cover_color[simplices[i][1]].first > MASK) + graphic << simplices[i][0] << " " << simplices[i][1] << std::endl; + graphic.close(); + } // ******************************************************************************************************************* // ******************************************************************************************************************* public: + /** \brief Creates the simplicial complex. + * + * \tparam SimplicialComplexForGIC must meet `SimplicialComplexForRips` concept. + * @param[in] complex SimplicialComplexForGIC to be created. + * + */ template void create_complex(SimplicialComplexForGIC & complex) { size_t sz = simplices.size(); int dimension = 0; @@ -533,6 +641,11 @@ class Graph_induced_complex { } public: + /** \brief Finds the maximal clique formed by different elements of the cover in a set of points. + * + * @param[in] cover_elts vector of points represented by vectors of cover elements (the ones to which they belong). + * + */ void find_all_simplices(const std::vector > & cover_elts){ int num_nodes = cover_elts.size(); std::vector simplex; @@ -541,12 +654,13 @@ class Graph_induced_complex { simplex.push_back(cover_elts[i][j]); std::sort(simplex.begin(),simplex.end()); std::vector::iterator it = std::unique(simplex.begin(),simplex.end()); simplex.resize(std::distance(simplex.begin(),it)); - simplices.push_back(simplex); //for(int i = 0; i < simplex.size(); i++) std::cout << simplex[i] << " "; std::cout << std::endl; + simplices.push_back(simplex); } public: + /** \brief Computes the simplices in the Nerve of the cover C. + */ void find_Nerve_simplices(){ - simplices.clear(); for(std::map >::iterator it = cover.begin(); it!=cover.end(); it++) simplices.push_back(it->second); std::vector >::iterator it; std::sort(simplices.begin(),simplices.end()); it = std::unique(simplices.begin(),simplices.end()); @@ -554,6 +668,8 @@ class Graph_induced_complex { } public: + /** \brief Computes the simplices in the GIC of the graph G and the cover C. + */ void find_GIC_simplices() { // Find IDs of edges to remove @@ -595,6 +711,8 @@ class Graph_induced_complex { } public: + /** \brief Computes the simplices in the Mapper Delta. + */ void find_GICMAP_simplices_with_functional_minimal_cover(){ int v1, v2; -- cgit v1.2.3