2 files changed, 6 insertions, 7 deletions

diff --git a/src/Persistence_representations/doc/Persistence_representations_doc.h b/src/Persistence_representations/doc/Persistence_representations_doc.h
index a8f8ba8b..091612ff 100644
--- a/src/Persistence_representations/doc/Persistence_representations_doc.h
+++ b/src/Persistence_representations/doc/Persistence_representations_doc.h
@@ -230,7 +230,7 @@ namespace Persistence_representations {
  In addition to the previous method, we also provide two more methods to perform exact calculations, in the sense that we use functions
  instead of matrices to define the kernel between the points of the diagrams.
  Indeed, in both of these exact methods, the kernel is no longer provided as a square matrix, or a filter (see parameters above), but rather as
- a function assigning a real value to a 2D point in the plane.
+ a function assigning a real value to a pair of points in the plane.
 
  In the first of these exact methods, we aim at obtaining a finite-dimensional representation of the diagram, so we still use a grid of pixels.
  On the other hand, in the second exact method, we represent diagrams implicitly as functions (i.e. infinite-dimensional representations). This way, we no longer require grids,
diff --git a/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h b/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h
index 2e23f69e..8bdcef65 100644
--- a/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h
+++ b/src/Persistence_representations/include/gudhi/Sliced_Wasserstein.h
@@ -148,19 +148,19 @@ class Sliced_Wasserstein {
     if(this->approx == -1){
 
       // Add projections onto diagonal.
-      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();
+      int n1, n2; n1 = diagram1.size(); n2 = diagram2.size(); double max_ordinate = 0; double max_abscissa = 0;
       for (int i = 0; i < n2; i++){
-        min_ordinate = std::min(min_ordinate, diagram2[i].second); min_abscissa = std::min(min_abscissa, diagram2[i].first);
+        max_ordinate = std::max(max_ordinate, std::abs(diagram2[i].second)); max_abscissa = std::max(max_abscissa, std::abs(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++){
-        min_ordinate = std::min(min_ordinate, diagram1[i].second); min_abscissa = std::min(min_abscissa, diagram1[i].first);
+        max_ordinate = std::max(max_ordinate, std::abs(diagram1[i].second)); max_abscissa = std::max(max_abscissa, std::abs(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.
-      double thresh_y = pow(10,log10(min_ordinate)-5); double thresh_x = pow(10,log10(min_abscissa)-5);
+      double thresh_y = max_ordinate * 0.00001; double thresh_x = max_abscissa * 0.00001;
       srand(time(NULL));
       for (int i = 0; i < num_pts_dgm; i++){
         diagram1[i].first += thresh_x*(1.0-2.0*rand()/RAND_MAX); diagram1[i].second += thresh_y*(1.0-2.0*rand()/RAND_MAX);
@@ -189,8 +189,7 @@ class Sliced_Wasserstein {
 
       // Find the inverses of the orders.
       std::vector<int> order1(num_pts_dgm); std::vector<int> order2(num_pts_dgm);
-      for(int i = 0; i < num_pts_dgm; i++)  for (int j = 0; j < num_pts_dgm; j++)  if(orderp1[j] == i){  order1[i] = j; break;  }
-      for(int i = 0; i < num_pts_dgm; i++)  for (int j = 0; j < num_pts_dgm; j++)  if(orderp2[j] == i){  order2[i] = j; break;  }
+      for(int i = 0; i < num_pts_dgm; i++){ order1[orderp1[i]] = i; order2[orderp2[i]] = i; }
 
       // Record all inversions of points in the orders as theta varies along the positive half-disk.
       std::vector<std::vector<std::pair<int,double> > > anglePerm1(num_pts_dgm);