From 991a1b06c258f21845dab412174e25e05ebb362c Mon Sep 17 00:00:00 2001 From: mcarrier Date: Thu, 4 May 2017 16:50:43 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@2399 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: 021dad464dcacbb4c24690e05770339a0f7f0d4f --- src/Nerve_GIC/include/gudhi/GIC.h | 98 +++++++++++++++++++++++++++++++++++---- 1 file changed, 90 insertions(+), 8 deletions(-) diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index d73a6d86..62b5a237 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -44,6 +44,10 @@ #include #include // for std::max #include // for std::uint32_t +#include + +#define CONSTANT 10 +#define ETA 0.001 using Simplex_tree = Gudhi::Simplex_tree<>; using Filtration_value = Simplex_tree::Filtration_value; @@ -116,6 +120,25 @@ class Graph_induced_complex { dfs(G,G[p][i],cc,visit); } + // Find random number in [0,1]. + double GetUniform(){ + static std::default_random_engine re; + static std::uniform_real_distribution Dist(0,1); + return Dist(re); + } + + // Subsample points. + void SampleWithoutReplacement(int populationSize, int sampleSize, std::vector & samples){ + int& n = sampleSize; int& N = populationSize; + int t = 0; int m = 0; double u; + while (m < n){ + u = GetUniform(); + if ( (N - t)*u >= n - m ) + t++; + else{samples[m] = t; t++; m++;} + } + } + // ******************************************************************************************************************* // Graphs. // ******************************************************************************************************************* @@ -129,7 +152,7 @@ class Graph_induced_complex { } } - public: + 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); @@ -156,6 +179,42 @@ class Graph_induced_complex { rips_complex_from_points.create_complex(st, 1); } + public: // Automatic tuning of Rips complex. + void set_graph_from_automatic_rips(const int& N, const std::string& off_file_name){ + + 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(); + + std::vector > dist; std::vector dumb(n); + for(int i = 0; i < n; i++) dist.push_back(dumb); + + for(int i = 0; i < n; i++){ + for(int j = i; 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); + dist[i][j] = std::sqrt(dis); dist[j][i] = std::sqrt(dis); + } + } + + #pragma omp parallel for + for (int i = 0; i < N; i++){ + + SampleWithoutReplacement(n,m,samples); + double hausdorff_dist = 0; + for (int j = 0; j < n; j++){ + double mj = dist[j][samples[0]]; for (int k = 1; k < m; k++) mj = std::min(mj, dist[j][samples[k]]); + hausdorff_dist = std::max(hausdorff_dist, mj); + } + delta += hausdorff_dist/N; + + } + + Rips_complex rips_complex_from_points(off_reader.get_point_cloud(), delta, Euclidean_distance()); + rips_complex_from_points.create_complex(st, 1); + + } + // ******************************************************************************************************************* // Functions. @@ -173,9 +232,8 @@ class Graph_induced_complex { public: // Set function from kth coordinate void set_function_from_coordinate(const int& k, const std::string& off_file_name){ Points_off_reader off_reader(off_file_name); - //std::vector cloud = off_reader.get_point_cloud(); - int numpts = off_reader.get_point_cloud().size(); //cloud.size(); - for(int i = 0; i < numpts; i++) func.insert(std::pair(i,off_reader.get_point_cloud()[i][k])); //std::cout << cloud[i][k] << std::endl; + int numpts = off_reader.get_point_cloud().size(); + for(int i = 0; i < numpts; i++) func.insert(std::pair(i,off_reader.get_point_cloud()[i][k])); } public: // Set function from vector. @@ -201,6 +259,32 @@ class Graph_induced_complex { cov_number.resize(std::distance(cov_number.begin(),it)); maximal_dim = cov_number.size(); } + public: // Automatic tuning of resolution for Mapper Delta. + double set_automatic_resolution_for_GICMAP(){ + double reso = 0; + for (auto simplex : st.complex_simplex_range()) { + if(st.dimension(simplex) == 1){ + std::vector vertices; + for(auto vertex : st.simplex_vertex_range(simplex)) vertices.push_back(vertex); + reso = std::max(reso, std::abs(func[vertices[0]] - func[vertices[1]])); + } + } + return reso; + } + + public: // Automatic tuning of resolution for Mapper Point. + double set_automatic_resolution_for_MAP(const double& gain){ + double reso = 0; + for (auto simplex : st.complex_simplex_range()) { + if(st.dimension(simplex) == 1){ + std::vector vertices; + for(auto vertex : st.simplex_vertex_range(simplex)) vertices.push_back(vertex); + reso = std::max(reso, (std::abs(func[vertices[0]] - func[vertices[1]]))/gain); + } + } + return reso; + } + public: // Set cover with preimages of function. void set_cover_from_function(const double& resolution, const double& gain, const bool& token){ @@ -350,9 +434,7 @@ class Graph_induced_complex { public: void find_Nerve_simplices(){ simplices.clear(); - for(std::map >::iterator it = cover.begin(); it!=cover.end(); it++){ - simplices.push_back(it->second); //std::cout << it->second[0] << std::endl; - } + 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()); simplices.resize(std::distance(simplices.begin(),it)); @@ -400,7 +482,7 @@ class Graph_induced_complex { } public: - void find_GIC_simplices_with_functional_minimal_cover(){ + void find_GICMAP_simplices_with_functional_minimal_cover(){ int v1, v2; -- cgit v1.2.3