diff options
author | mcarrier <mcarrier@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2018-08-23 21:11:47 +0000 |
---|---|---|
committer | mcarrier <mcarrier@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2018-08-23 21:11:47 +0000 |
commit | 6d00273077db54d609262e79702cbd5a94491105 (patch) | |
tree | b5e8a915c8ececfdc8ebb5ef4ffebfcbd9d973ad /src | |
parent | 5f5a7a21e9db73eaf9dc2604cb0de3066f7a4fb6 (diff) |
git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/kernels@3827 636b058d-ea47-450e-bf9e-a15bfbe3eedb
Former-commit-id: a6a297a7f14703e55954706328072acdd447484c
Diffstat (limited to 'src')
3 files changed, 12 insertions, 34 deletions
diff --git a/src/Persistence_representations/include/gudhi/Persistence_heat_maps.h b/src/Persistence_representations/include/gudhi/Persistence_heat_maps.h index 12188526..43f10b8c 100644 --- a/src/Persistence_representations/include/gudhi/Persistence_heat_maps.h +++ b/src/Persistence_representations/include/gudhi/Persistence_heat_maps.h @@ -1002,8 +1002,6 @@ double Persistence_heat_maps<Scalling_of_kernels>::compute_scalar_product(const } else{ - GUDHI_CHECK(this->approx != second.approx || this->f != second.f, std::invalid_argument("Error: different values for representations")); - int num_pts1 = this->d.size(); int num_pts2 = second.d.size(); double kernel_val = 0; for(int i = 0; i < num_pts1; i++) for(int j = 0; j < num_pts2; j++) diff --git a/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h b/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h index a3c0dc2f..6f67f7bc 100644 --- a/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h +++ b/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h @@ -97,23 +97,7 @@ class Sliced_Wasserstein { // Compute the angle formed by two points of a PD double compute_angle(const Persistence_diagram & diag, int i, int j) const { - std::pair<double,double> vect; double x1,y1, x2,y2; - x1 = diag[i].first; y1 = diag[i].second; - x2 = diag[j].first; y2 = diag[j].second; - if (y1 - y2 > 0){ - vect.first = y1 - y2; - vect.second = x2 - x1;} - else{ - if(y1 - y2 < 0){ - vect.first = y2 - y1; - vect.second = x1 - x2; - } - else{ - vect.first = 0; - vect.second = abs(x1 - x2);} - } - double norm = std::sqrt(vect.first*vect.first + vect.second*vect.second); - return asin(vect.second/norm); + if(diag[i].second == diag[j].second) return pi/2; else return atan((diag[j].first-diag[i].first)/(diag[i].second-diag[j].second)); } // Compute the integral of |cos()| between alpha and beta, valid only if alpha is in [-pi,pi] and beta-alpha is in [0,pi] @@ -145,10 +129,7 @@ class Sliced_Wasserstein { double compute_int(double theta1, double theta2, int p, int q, const Persistence_diagram & diag1, const Persistence_diagram & diag2) const { double norm = std::sqrt( (diag1[p].first-diag2[q].first)*(diag1[p].first-diag2[q].first) + (diag1[p].second-diag2[q].second)*(diag1[p].second-diag2[q].second) ); double angle1; - if (diag1[p].first > diag2[q].first) - angle1 = theta1 - asin( (diag1[p].second-diag2[q].second)/norm ); - else - angle1 = theta1 - asin( (diag2[q].second-diag1[p].second)/norm ); + if (diag1[p].first == diag2[q].first) angle1 = theta1 - pi/2; else angle1 = theta1 - atan((diag1[p].second-diag2[q].second)/(diag1[p].first-diag2[q].first)); double angle2 = angle1 + theta2 - theta1; double integral = compute_int_cos(angle1,angle2); return norm*integral; @@ -164,24 +145,23 @@ class Sliced_Wasserstein { if(this->approx == -1){ // Add projections onto diagonal. - int n1, n2; n1 = diagram1.size(); n2 = diagram2.size(); double max_ordinate = std::numeric_limits<double>::lowest(); + int n1, n2; n1 = diagram1.size(); n2 = diagram2.size(); double min_ordinate = std::numeric_limits<double>::max(); double min_abscissa = std::numeric_limits<double>::max(); for (int i = 0; i < n2; i++){ - max_ordinate = std::max(max_ordinate, diagram2[i].second); + min_ordinate = std::min(min_ordinate, diagram2[i].second); min_abscissa = std::min(min_abscissa, diagram2[i].first); diagram1.emplace_back( (diagram2[i].first+diagram2[i].second)/2, (diagram2[i].first+diagram2[i].second)/2 ); } for (int i = 0; i < n1; i++){ - max_ordinate = std::max(max_ordinate, diagram1[i].second); + min_ordinate = std::min(min_ordinate, diagram1[i].second); min_abscissa = std::min(min_abscissa, diagram1[i].first); diagram2.emplace_back( (diagram1[i].first+diagram1[i].second)/2, (diagram1[i].first+diagram1[i].second)/2 ); } int num_pts_dgm = diagram1.size(); // Slightly perturb the points so that the PDs are in generic positions. - int mag = 0; while(max_ordinate > 10){mag++; max_ordinate/=10;} - double thresh = pow(10,-5+mag); + double thresh_y = pow(10,log10(min_ordinate)-5); double thresh_x = pow(10,log10(min_abscissa)-5); srand(time(NULL)); for (int i = 0; i < num_pts_dgm; i++){ - diagram1[i].first += thresh*(1.0-2.0*rand()/RAND_MAX); diagram1[i].second += thresh*(1.0-2.0*rand()/RAND_MAX); - diagram2[i].first += thresh*(1.0-2.0*rand()/RAND_MAX); diagram2[i].second += thresh*(1.0-2.0*rand()/RAND_MAX); + diagram1[i].first += thresh_x*(1.0-2.0*rand()/RAND_MAX); diagram1[i].second += thresh_y*(1.0-2.0*rand()/RAND_MAX); + diagram2[i].first += thresh_x*(1.0-2.0*rand()/RAND_MAX); diagram2[i].second += thresh_y*(1.0-2.0*rand()/RAND_MAX); } // Compute all angles in both PDs. @@ -201,8 +181,8 @@ class Sliced_Wasserstein { // Initialize orders of the points of both PDs (given by ordinates when theta = -pi/2). std::vector<int> orderp1, orderp2; for (int i = 0; i < num_pts_dgm; i++){ orderp1.push_back(i); orderp2.push_back(i); } - std::sort( orderp1.begin(), orderp1.end(), [=](int i, int j){ if(diagram1[i].second != diagram1[j].second) return (diagram1[i].second < diagram1[j].second); else return (diagram1[i].first > diagram1[j].first); } ); - std::sort( orderp2.begin(), orderp2.end(), [=](int i, int j){ if(diagram2[i].second != diagram2[j].second) return (diagram2[i].second < diagram2[j].second); else return (diagram2[i].first > diagram2[j].first); } ); + std::sort( orderp1.begin(), orderp1.end(), [&](int i, int j){ if(diagram1[i].second != diagram1[j].second) return (diagram1[i].second < diagram1[j].second); else return (diagram1[i].first > diagram1[j].first); } ); + std::sort( orderp2.begin(), orderp2.end(), [&](int i, int j){ if(diagram2[i].second != diagram2[j].second) return (diagram2[i].second < diagram2[j].second); else return (diagram2[i].first > diagram2[j].first); } ); // Find the inverses of the orders. std::vector<int> order1(num_pts_dgm); std::vector<int> order2(num_pts_dgm); @@ -274,6 +254,7 @@ class Sliced_Wasserstein { public: /** \brief Sliced Wasserstein kernel constructor. + * \implements Topological_data_with_distances, Real_valued_topological_data, Topological_data_with_scalar_product * \ingroup Sliced_Wasserstein * * @param[in] _diagram persistence diagram. @@ -282,7 +263,6 @@ class Sliced_Wasserstein { * points on all directions are stored in memory to reduce computation time. * */ - // This class implements the following concepts: Topological_data_with_distances, Real_valued_topological_data, Topological_data_with_scalar_product Sliced_Wasserstein(const Persistence_diagram & _diagram, double _sigma = 1.0, int _approx = 10):diagram(_diagram), approx(_approx), sigma(_sigma) {build_rep();} /** \brief Evaluation of the kernel on a pair of diagrams. diff --git a/src/cmake/modules/GUDHI_modules.cmake b/src/cmake/modules/GUDHI_modules.cmake index 276fb2cc..f95d0c34 100644 --- a/src/cmake/modules/GUDHI_modules.cmake +++ b/src/cmake/modules/GUDHI_modules.cmake @@ -17,7 +17,7 @@ function(add_gudhi_module file_path) endfunction(add_gudhi_module) option(WITH_GUDHI_BENCHMARK "Activate/desactivate benchmark compilation" OFF) -option(WITH_GUDHI_EXAMPLE "Activate/desactivate examples compilation and installation" ON) +option(WITH_GUDHI_EXAMPLE "Activate/desactivate examples compilation and installation" OFF) option(WITH_GUDHI_PYTHON "Activate/desactivate python module compilation and installation" ON) option(WITH_GUDHI_TEST "Activate/desactivate examples compilation and installation" ON) option(WITH_GUDHI_UTILITIES "Activate/desactivate utilities compilation and installation" ON) |