From 56b5f54a72e23a1feb7d12ad8af30cf7baf25060 Mon Sep 17 00:00:00 2001 From: mcarrier Date: Thu, 4 May 2017 15:40:34 +0000 Subject: git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/Nerve_GIC@2397 636b058d-ea47-450e-bf9e-a15bfbe3eedb Former-commit-id: a58d8b1539aeb53957e843dcfa3deb87cc59ebbd --- src/Nerve_GIC/include/gudhi/GIC.h | 107 ++++++++++++++++++-------------------- 1 file changed, 51 insertions(+), 56 deletions(-) (limited to 'src/Nerve_GIC/include/gudhi/GIC.h') diff --git a/src/Nerve_GIC/include/gudhi/GIC.h b/src/Nerve_GIC/include/gudhi/GIC.h index 93a7e5dc..d73a6d86 100644 --- a/src/Nerve_GIC/include/gudhi/GIC.h +++ b/src/Nerve_GIC/include/gudhi/GIC.h @@ -82,6 +82,9 @@ class Graph_induced_complex { private: int maximal_dim; + private: + std::map cover_fct; + private: Simplex_tree<> st; std::map > adjacency_matrix; @@ -202,8 +205,8 @@ class Graph_induced_complex { void set_cover_from_function(const double& resolution, const double& gain, const bool& token){ // Read function values and compute min and max - std::map::iterator it; /*std::vector range;*/ double maxf, minf; minf = std::numeric_limits::max(); maxf = std::numeric_limits::min(); - for(it = func.begin(); it != func.end(); it++){/*range.push_back(it->second);*/ minf = std::min(minf, it->second); maxf = std::max(maxf, it->second);} + 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 << minf << " " << maxf << std::endl; // Compute cover of im(f) @@ -218,6 +221,7 @@ class Graph_induced_complex { } x = minf + (resolution-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 << intervals[i].first << " " << intervals[i].second << std::endl; } else{ double x = minf; double y = x + resolution; @@ -244,27 +248,31 @@ class Graph_induced_complex { } } - int id = 0; int pos = 0; - for(int i = 0; i < res; i++){ + int id = 0; int pos = 0; double min_prop_int; double max_prop_int; + + for(int i = 0; i < res; i++){ // Find points in the preimage std::map > prop; prop.clear(); - std::pair inter1 = intervals[i]; + std::pair inter1 = intervals[i]; min_prop_int = inter1.first; int tmp = pos; if(i != res-1){ + if(i != 0){ - std::pair inter3 = intervals[i-1]; + std::pair inter3 = intervals[i-1]; min_prop_int = inter3.second; while(func[points[tmp]] < inter3.second && tmp != num_pts){ prop.insert(std::pair >(points[tmp],adjacency_matrix[points[tmp]])); tmp++; } } - std::pair inter2 = intervals[i+1]; + + std::pair inter2 = intervals[i+1]; max_prop_int = inter2.first; while(func[points[tmp]] < inter2.first && tmp != num_pts){ prop.insert(std::pair >(points[tmp],adjacency_matrix[points[tmp]])); tmp++; } + pos = tmp; while(func[points[tmp]] < inter1.second && tmp != num_pts){ prop.insert(std::pair >(points[tmp],adjacency_matrix[points[tmp]])); @@ -272,12 +280,16 @@ class Graph_induced_complex { } } + else{ + std::pair inter3 = intervals[i-1]; + min_prop_int = inter3.second; max_prop_int = inter1.second; while(func[points[tmp]] < inter3.second && tmp != num_pts){ prop.insert(std::pair >(points[tmp],adjacency_matrix[points[tmp]])); tmp++; } + while(tmp != num_pts){ prop.insert(std::pair >(points[tmp],adjacency_matrix[points[tmp]])); tmp++; @@ -286,20 +298,21 @@ class Graph_induced_complex { } // Compute the connected components with DFS - std::map visit; + std::map visit; //std::cout << 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 << cci << std::endl; - for(int i = 0; i < cci; i++) cover[cc[i]].push_back(id); + dfs(prop,it->first,cc,visit); int cci = cc.size(); //std::cout << cci << " "; + for(int j = 0; j < cci; j++) cover[cc[j]].push_back(id); + cover_fct[id] = i; id++; } } } - + //std::cout << std::endl; } maximal_dim = id; @@ -323,45 +336,22 @@ class Graph_induced_complex { } public: - void find_all_simplices(const std::vector > & cover_elts\ - //int token, std::vector & simplex_tmp - ){ + void find_all_simplices(const std::vector > & cover_elts){ int num_nodes = cover_elts.size(); - - // Old method. - - /*if(token == num_nodes-1){ - int num_clus = cover_elts[token].size(); - for(int i = 0; i < num_clus; i++){ - std::vector simplex = simplex_tmp; simplex.push_back(cover_elts[token][i]); - std::vector::iterator it; - std::sort(simplex.begin(),simplex.end()); it = std::unique(simplex.begin(),simplex.end()); - simplex.resize(std::distance(simplex.begin(),it)); - simplices.push_back(simplex); - } - } - else{ - int num_clus = cover_elts[token].size(); - for(int i = 0; i < num_clus; i++){ - std::vector simplex = simplex_tmp; simplex.push_back(cover_elts[token][i]); - find_all_simplices(cover_elts, token+1, simplex); - } - }*/ - std::vector simplex; for(int i = 0; i < num_nodes; i++) for(int j = 0; j < cover_elts[i].size(); j++) 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); //for(int i = 0; i < simplex.size(); i++) std::cout << simplex[i] << " "; std::cout << std::endl; } public: void find_Nerve_simplices(){ simplices.clear(); for(std::map >::iterator it = cover.begin(); it!=cover.end(); it++){ - simplices.push_back(it->second); + simplices.push_back(it->second); //std::cout << it->second[0] << std::endl; } std::vector >::iterator it; std::sort(simplices.begin(),simplices.end()); it = std::unique(simplices.begin(),simplices.end()); @@ -399,9 +389,9 @@ class Graph_induced_complex { simplices.clear(); for (auto simplex : st.complex_simplex_range()) { if(!st.has_children(simplex)){ - std::vector > cover_elts; //std::vector sim; + std::vector > cover_elts; for (auto vertex : st.simplex_vertex_range(simplex)) cover_elts.push_back(cover[vertex]); - find_all_simplices(cover_elts/*,0,sim*/); + find_all_simplices(cover_elts); } } std::vector >::iterator it; @@ -411,26 +401,31 @@ class Graph_induced_complex { public: void find_GIC_simplices_with_functional_minimal_cover(){ - for(std::map >::iterator it = cover.begin(); it != cover.end(); it++){ + + int v1, v2; + + // Loop on all points. + for(std::map >::iterator it = cover.begin(); it != cover.end(); it++){ + int vid = it->first; std::vector neighbors = adjacency_matrix[vid]; int num_neighb = neighbors.size(); - std::vector vids(1); vids[0] = cover[vid][0]; simplices.push_back(vids); + //std::cout << vid << " " << num_neighb << std::endl; + + // Find cover of current point (vid). + if(cover[vid].size() == 2) v1 = std::min(cover[vid][0],cover[vid][1]); else v1 = cover[vid][0]; + + // Loop on neighbors. for(int i = 0; i < num_neighb; i++){ - int neighb = neighbors[i]; int v1, v2; - if(func[vid] > func[neighb]){ - if(cover[vid].size() == 2) v1 = std::max(cover[vid][0],cover[vid][1]); else v1 = cover[vid][0]; - if(cover[neighb].size() == 2) v2 = std::min(cover[neighb][0],cover[neighb][1]); else v2 = cover[neighb][0]; - std::vector edge(2); edge[0] = v1; edge[1] = v2; - //std::cout << v1 << " " << v2 << std::endl; - if(v1 != v2) simplices.push_back(edge); + + int neighb = neighbors[i]; + + // Find cover of neighbor (neighb). + if(cover[neighb].size() == 2) v2 = std::max(cover[neighb][0],cover[neighb][1]); else v2 = cover[neighb][0]; + + // If neighbor is in next interval, add edge. + if(cover_fct[v2] == cover_fct[v1] + 1){ + std::vector edge(2); edge[0] = v1; edge[1] = v2; //std::cout << v1 << " " << v2 << std::endl; + simplices.push_back(edge); } - /*else{ - if( func[neighb]-func[vid] <= resolution*(2-gain) ){ - if(cover[vid].size() == 2) v1 = std::min(cover[vid][0],cover[vid][1]); else v1 = cover[vid][0]; - if(cover[neighb].size() == 2) v2 = std::max(cover[neighb][0],cover[neighb][1]); else v2 = cover[neighb][0]; - std::vector edge(2); edge[0] = v1; edge[1] = v2; - if(v2 > v1) simplices.push_back(edge); - } - }*/ } } std::vector >::iterator it; -- cgit v1.2.3