From f911ed3882c03b049a18f2a638758cb5ef994bcc Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 25 May 2022 11:59:37 +0200 Subject: Add tests for silhouettes --- src/python/test/test_representations.py | 29 +++++++++++++++++++++++++++++ 1 file changed, 29 insertions(+) (limited to 'src') diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index d219ce7a..8be3a1db 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -168,3 +168,32 @@ def test_kernel_empty_diagrams(): # PersistenceFisherKernel(bandwidth_fisher=1., bandwidth=1.)(empty_diag, empty_diag) # PersistenceFisherKernel(bandwidth_fisher=1., bandwidth=1., kernel_approx=RBFSampler(gamma=1./2, n_components=100000).fit(np.ones([1,2])))(empty_diag, empty_diag) + +def test_silhouette_permutation_invariance(): + dgm = _n_diags(1)[0] + dgm_permuted = dgm[np.random.permutation(dgm.shape[0]).astype(int)] + random_resolution = random.randint(50, 100) * 10 + slt = Silhouette(resolution=random_resolution, weight=pow(2)) + + assert np.all(np.isclose(slt(dgm), slt(dgm_permuted))) + + +def test_silhouette_multiplication_invariance(): + dgm = _n_diags(1)[0] + n_repetitions = np.random.randint(2, high=10) + dgm_augmented = np.repeat(dgm, repeats=n_repetitions, axis=0) + + random_resolution = random.randint(50, 100) * 10 + slt = Silhouette(resolution=random_resolution, weight=pow(2)) + assert np.all(np.isclose(slt(dgm), slt(dgm_augmented))) + + +def test_silhouette_numeric(): + dgm = np.array([[2., 3.], [5., 6.]]) + slt = Silhouette(resolution=9, weight=pow(1)) + #slt.fit([dgm]) + # x_values = array([2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.]) + + expected_silhouette = np.array([0., 0.5, 0., 0., 0., 0., 0., 0.5, 0.])/np.sqrt(2) + output_silhouette = slt(dgm) + assert np.all(np.isclose(output_silhouette, expected_silhouette)) -- cgit v1.2.3 From 3aa89676d1dc2cafcc692480bbf424a97dbbd501 Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 25 May 2022 14:20:12 +0200 Subject: Vectorize Silhouette implementation --- src/python/gudhi/representations/vector_methods.py | 48 ++++++---------------- 1 file changed, 13 insertions(+), 35 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index f8078d03..62b35389 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -85,7 +85,7 @@ class PersistenceImage(BaseEstimator, TransformerMixin): Xfit.append(image.flatten()[np.newaxis,:]) - Xfit = np.concatenate(Xfit,0) + Xfit = np.concatenate(Xfit, 0) return Xfit @@ -235,6 +235,7 @@ class Silhouette(BaseEstimator, TransformerMixin): sample_range ([double, double]): minimum and maximum for the weighted average domain, of the form [x_min, x_max] (default [numpy.nan, numpy.nan]). It is the interval on which samples will be drawn evenly. If one of the values is numpy.nan, it can be computed from the persistence diagrams with the fit() method. """ self.weight, self.resolution, self.sample_range = weight, resolution, sample_range + self.im_range = None def fit(self, X, y=None): """ @@ -245,6 +246,7 @@ class Silhouette(BaseEstimator, TransformerMixin): y (n x 1 array): persistence diagram labels (unused). """ self.sample_range = _automatic_sample_range(np.array(self.sample_range), X, y) + self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.resolution) return self def transform(self, X): @@ -257,44 +259,20 @@ class Silhouette(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of diagrams) x (**resolution**): output persistence silhouettes. """ - num_diag, Xfit = len(X), [] - x_values = np.linspace(self.sample_range[0], self.sample_range[1], self.resolution) - step_x = x_values[1] - x_values[0] + Xfit = [] + x_values = self.im_range - for i in range(num_diag): - - diagram, num_pts_in_diag = X[i], X[i].shape[0] - - sh, weights = np.zeros(self.resolution), np.zeros(num_pts_in_diag) - for j in range(num_pts_in_diag): - weights[j] = self.weight(diagram[j,:]) + for i, diag in enumerate(X): + midpoints, heights = (diag[:, 0] + diag[:, 1])/2., (diag[:, 1] - diag[:, 0])/2. + weights = np.array([self.weight(point) for point in diag]) total_weight = np.sum(weights) - for j in range(num_pts_in_diag): - - [px,py] = diagram[j,:2] - weight = weights[j] / total_weight - min_idx = np.clip(np.ceil((px - self.sample_range[0]) / step_x).astype(int), 0, self.resolution) - mid_idx = np.clip(np.ceil((0.5*(py+px) - self.sample_range[0]) / step_x).astype(int), 0, self.resolution) - max_idx = np.clip(np.ceil((py - self.sample_range[0]) / step_x).astype(int), 0, self.resolution) - - if min_idx < self.resolution and max_idx > 0: - - silhouette_value = self.sample_range[0] + min_idx * step_x - px - for k in range(min_idx, mid_idx): - sh[k] += weight * silhouette_value - silhouette_value += step_x - - silhouette_value = py - self.sample_range[0] - mid_idx * step_x - for k in range(mid_idx, max_idx): - sh[k] += weight * silhouette_value - silhouette_value -= step_x - - Xfit.append(np.reshape(np.sqrt(2) * sh, [1,-1])) - - Xfit = np.concatenate(Xfit, 0) + tent_functions = heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]) + tent_functions[tent_functions < 0.] = 0. + silhouette = np.sum(weights[None, :]/total_weight * tent_functions, axis=1) + Xfit.append(silhouette * np.sqrt(2)) - return Xfit + return np.stack(Xfit, axis=0) def __call__(self, diag): """ -- cgit v1.2.3 From 1a76ecc3e7459e3461e1f182004362dcb663addd Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 25 May 2022 14:34:11 +0200 Subject: Compactify --- src/python/gudhi/representations/vector_methods.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 62b35389..e6289a37 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -267,8 +267,7 @@ class Silhouette(BaseEstimator, TransformerMixin): weights = np.array([self.weight(point) for point in diag]) total_weight = np.sum(weights) - tent_functions = heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]) - tent_functions[tent_functions < 0.] = 0. + tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) silhouette = np.sum(weights[None, :]/total_weight * tent_functions, axis=1) Xfit.append(silhouette * np.sqrt(2)) -- cgit v1.2.3 From e8d0cbc3311765900e098b472608dc40b84d07d8 Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 25 May 2022 15:14:15 +0200 Subject: Optimize using keops --- src/python/gudhi/representations/vector_methods.py | 27 +++++++++++++++------- 1 file changed, 19 insertions(+), 8 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index e6289a37..55dc2c5b 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -10,6 +10,7 @@ # - 2021/11 Vincent Rouvreau: factorize _automatic_sample_range import numpy as np +from pykeops.numpy import Genred from sklearn.base import BaseEstimator, TransformerMixin from sklearn.exceptions import NotFittedError from sklearn.preprocessing import MinMaxScaler, MaxAbsScaler @@ -259,19 +260,29 @@ class Silhouette(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of diagrams) x (**resolution**): output persistence silhouettes. """ - Xfit = [] - x_values = self.im_range + silhouette_formula = "normalized_weights * ReLU(heights - Abs(x_values - midpoints))" + variables = [ + "normalized_weights = Vi(1)", + "heights = Vi(1)", + "midpoints = Vi(1)", + "x_values = Vj(1)", + ] + silhouette = Genred(silhouette_formula, variables, reduction_op="Sum", axis=0) + + silhouettes_list = [] + x_values = self.im_range for i, diag in enumerate(X): - midpoints, heights = (diag[:, 0] + diag[:, 1])/2., (diag[:, 1] - diag[:, 0])/2. + midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. weights = np.array([self.weight(point) for point in diag]) - total_weight = np.sum(weights) + weights /= np.sum(weights) - tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) - silhouette = np.sum(weights[None, :]/total_weight * tent_functions, axis=1) - Xfit.append(silhouette * np.sqrt(2)) + silhouettes_list.append( + np.sqrt(2) * silhouette(weights[:, None], heights[:, None], + midpoints[:, None], x_values[:, None])[:, 0] + ) - return np.stack(Xfit, axis=0) + return np.stack(silhouettes_list, axis=0) def __call__(self, diag): """ -- cgit v1.2.3 From 912156b36da1dce1f73f8d2a63cc18e67c173d54 Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 25 May 2022 16:43:15 +0200 Subject: Move the initialisation of the Genred method to the constructor of Silhouette --- src/python/gudhi/representations/vector_methods.py | 21 ++++++++++----------- 1 file changed, 10 insertions(+), 11 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 55dc2c5b..c250c98c 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -238,6 +238,15 @@ class Silhouette(BaseEstimator, TransformerMixin): self.weight, self.resolution, self.sample_range = weight, resolution, sample_range self.im_range = None + silhouette_formula = "normalized_weights * ReLU(heights - Abs(x_values - midpoints))" + variables = [ + "normalized_weights = Vi(1)", + "heights = Vi(1)", + "midpoints = Vi(1)", + "x_values = Vj(1)", + ] + self.silhouette = Genred(silhouette_formula, variables, reduction_op="Sum", axis=0) + def fit(self, X, y=None): """ Fit the Silhouette class on a list of persistence diagrams: if any of the values in **sample_range** is numpy.nan, replace it with the corresponding value computed on the given list of persistence diagrams. @@ -260,16 +269,6 @@ class Silhouette(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of diagrams) x (**resolution**): output persistence silhouettes. """ - - silhouette_formula = "normalized_weights * ReLU(heights - Abs(x_values - midpoints))" - variables = [ - "normalized_weights = Vi(1)", - "heights = Vi(1)", - "midpoints = Vi(1)", - "x_values = Vj(1)", - ] - silhouette = Genred(silhouette_formula, variables, reduction_op="Sum", axis=0) - silhouettes_list = [] x_values = self.im_range for i, diag in enumerate(X): @@ -278,7 +277,7 @@ class Silhouette(BaseEstimator, TransformerMixin): weights /= np.sum(weights) silhouettes_list.append( - np.sqrt(2) * silhouette(weights[:, None], heights[:, None], + np.sqrt(2) * self.silhouette(weights[:, None], heights[:, None], midpoints[:, None], x_values[:, None])[:, 0] ) -- cgit v1.2.3 From 500de3b265c227eadd81fb860fade4f6c416cfbc Mon Sep 17 00:00:00 2001 From: wreise Date: Thu, 2 Jun 2022 10:49:07 +0200 Subject: Add tests for Landscapes --- src/python/test/test_representations.py | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) (limited to 'src') diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index 8be3a1db..0324b026 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -197,3 +197,29 @@ def test_silhouette_numeric(): expected_silhouette = np.array([0., 0.5, 0., 0., 0., 0., 0., 0.5, 0.])/np.sqrt(2) output_silhouette = slt(dgm) assert np.all(np.isclose(output_silhouette, expected_silhouette)) + + +def test_landscape_small_persistence_invariance(): + dgm = np.array([[2., 6.], [2., 5.], [3., 7.]]) + small_persistence_pts = np.random.rand(10, 2) + small_persistence_pts[:, 1] += small_persistence_pts[:, 0] + small_persistence_pts += np.min(dgm) + dgm_augmented = np.concatenate([dgm, small_persistence_pts], axis=0) + + lds = Landscape(num_landscapes=2, resolution=5) + lds_dgm, lds_dgm_augmented = lds(dgm), lds(dgm_augmented) + + assert np.all(np.isclose(lds_dgm, lds_dgm_augmented)) + + +def test_landscape_numeric(): + dgm = np.array([[2., 6.], [3., 5.]]) + lds_ref = np.array([ + 0., 0.5, 1., 1.5, 2., 1.5, 1., 0.5, 0., # tent of [2, 6] + 0., 0., 0., 0.5, 1., 0.5, 0., 0., 0. + ]) + lds_ref *= np.sqrt(2) + lds = Landscape(num_landscapes=2, resolution=9, sample_range=[2., 6.]) + lds_dgm = lds(dgm) + print(lds.sample_range, lds.new_resolution, lds_dgm.shape) + assert np.all(np.isclose(lds_dgm, lds_ref)) -- cgit v1.2.3 From 775cd158c869f75a4925ac6f7a94a1847f7d7788 Mon Sep 17 00:00:00 2001 From: wreise Date: Thu, 2 Jun 2022 15:38:20 +0200 Subject: Remove print from test --- src/python/test/test_representations.py | 1 - 1 file changed, 1 deletion(-) (limited to 'src') diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index 0324b026..5f29740f 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -221,5 +221,4 @@ def test_landscape_numeric(): lds_ref *= np.sqrt(2) lds = Landscape(num_landscapes=2, resolution=9, sample_range=[2., 6.]) lds_dgm = lds(dgm) - print(lds.sample_range, lds.new_resolution, lds_dgm.shape) assert np.all(np.isclose(lds_dgm, lds_ref)) -- cgit v1.2.3 From a70ad064ac3e0aee5c9b8084d35b7ce329c6bddc Mon Sep 17 00:00:00 2001 From: wreise Date: Thu, 2 Jun 2022 15:38:38 +0200 Subject: Start landscape optimization --- src/python/gudhi/representations/vector_methods.py | 64 ++++++++-------------- 1 file changed, 22 insertions(+), 42 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index c250c98c..acc62943 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -141,6 +141,16 @@ class Landscape(BaseEstimator, TransformerMixin): self.nan_in_range = np.isnan(np.array(self.sample_range)) self.new_resolution = self.resolution + self.nan_in_range.sum() + landscape_formula = "(-1)*ReLU(heights - Abs(x_values - midpoints))" + variables = [ + "heights = Vi(1)", + "midpoints = Vi(1)", + "x_values = Vj(1)", + ] + from pykeops.numpy import Genred + self.landscape = Genred(landscape_formula, variables, reduction_op="KMin", + axis=0, opt_arg=self.num_landscapes) + def fit(self, X, y=None): """ Fit the Landscape class on a list of persistence diagrams: if any of the values in **sample_range** is numpy.nan, replace it with the corresponding value computed on the given list of persistence diagrams. @@ -162,53 +172,23 @@ class Landscape(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of diagrams) x (number of samples = **num_landscapes** x **resolution**): output persistence landscapes. """ - num_diag, Xfit = len(X), [] - x_values = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) - step_x = x_values[1] - x_values[0] - - for i in range(num_diag): - - diagram, num_pts_in_diag = X[i], X[i].shape[0] - ls = np.zeros([self.num_landscapes, self.new_resolution]) - - events = [] - for j in range(self.new_resolution): - events.append([]) - - for j in range(num_pts_in_diag): - [px,py] = diagram[j,:2] - min_idx = np.clip(np.ceil((px - self.sample_range[0]) / step_x).astype(int), 0, self.new_resolution) - mid_idx = np.clip(np.ceil((0.5*(py+px) - self.sample_range[0]) / step_x).astype(int), 0, self.new_resolution) - max_idx = np.clip(np.ceil((py - self.sample_range[0]) / step_x).astype(int), 0, self.new_resolution) - - if min_idx < self.new_resolution and max_idx > 0: - - landscape_value = self.sample_range[0] + min_idx * step_x - px - for k in range(min_idx, mid_idx): - events[k].append(landscape_value) - landscape_value += step_x - - landscape_value = py - self.sample_range[0] - mid_idx * step_x - for k in range(mid_idx, max_idx): - events[k].append(landscape_value) - landscape_value -= step_x - - for j in range(self.new_resolution): - events[j].sort(reverse=True) - for k in range( min(self.num_landscapes, len(events[j])) ): - ls[k,j] = events[j][k] + Xfit = [] + x_values = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) + for i, diag in enumerate(X): + midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. + tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) + tent_functions.partition(diag.shape[0] - self.num_landscapes, axis=1) + landscapes = np.sort(tent_functions[-self.num_landscapes:, :])[::-1].T if self.nan_in_range[0]: - ls = ls[:,1:] + landscapes = landscapes[:,1:] if self.nan_in_range[1]: - ls = ls[:,:-1] - ls = np.sqrt(2)*np.reshape(ls,[1,-1]) - Xfit.append(ls) - - Xfit = np.concatenate(Xfit,0) + landscapes = landscapes[:,:-1] + landscapes = np.sqrt(2)*np.reshape(landscapes, [1, -1]) + Xfit.append(landscapes) - return Xfit + return np.stack(Xfit, axis=0) def __call__(self, diag): """ -- cgit v1.2.3 From 42b18e60e418f4078cd6406dcc202b696798c844 Mon Sep 17 00:00:00 2001 From: wreise Date: Sun, 5 Jun 2022 16:59:02 +0200 Subject: Import pykeops locally in Silhouette and Lanscapes --- src/python/gudhi/representations/vector_methods.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index acc62943..b0843120 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -10,7 +10,6 @@ # - 2021/11 Vincent Rouvreau: factorize _automatic_sample_range import numpy as np -from pykeops.numpy import Genred from sklearn.base import BaseEstimator, TransformerMixin from sklearn.exceptions import NotFittedError from sklearn.preprocessing import MinMaxScaler, MaxAbsScaler @@ -218,6 +217,7 @@ class Silhouette(BaseEstimator, TransformerMixin): self.weight, self.resolution, self.sample_range = weight, resolution, sample_range self.im_range = None + from pykeops.numpy import Genred silhouette_formula = "normalized_weights * ReLU(heights - Abs(x_values - midpoints))" variables = [ "normalized_weights = Vi(1)", -- cgit v1.2.3 From 60e57f9c86a7aae67c2931200066aba059ec2721 Mon Sep 17 00:00:00 2001 From: wreise Date: Fri, 5 Aug 2022 22:19:30 +0200 Subject: Test the numpy version --- src/python/gudhi/representations/vector_methods.py | 42 ++++++---------------- 1 file changed, 11 insertions(+), 31 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index b0843120..7f311b3b 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -140,16 +140,6 @@ class Landscape(BaseEstimator, TransformerMixin): self.nan_in_range = np.isnan(np.array(self.sample_range)) self.new_resolution = self.resolution + self.nan_in_range.sum() - landscape_formula = "(-1)*ReLU(heights - Abs(x_values - midpoints))" - variables = [ - "heights = Vi(1)", - "midpoints = Vi(1)", - "x_values = Vj(1)", - ] - from pykeops.numpy import Genred - self.landscape = Genred(landscape_formula, variables, reduction_op="KMin", - axis=0, opt_arg=self.num_landscapes) - def fit(self, X, y=None): """ Fit the Landscape class on a list of persistence diagrams: if any of the values in **sample_range** is numpy.nan, replace it with the corresponding value computed on the given list of persistence diagrams. @@ -178,13 +168,13 @@ class Landscape(BaseEstimator, TransformerMixin): midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) tent_functions.partition(diag.shape[0] - self.num_landscapes, axis=1) - landscapes = np.sort(tent_functions[-self.num_landscapes:, :])[::-1].T + landscapes = np.sort(tent_functions, axis=1)[:, -self.num_landscapes:][:, ::-1].T if self.nan_in_range[0]: - landscapes = landscapes[:,1:] + landscapes = landscapes[:, 1:] if self.nan_in_range[1]: - landscapes = landscapes[:,:-1] - landscapes = np.sqrt(2)*np.reshape(landscapes, [1, -1]) + landscapes = landscapes[:, :-1] + landscapes = np.sqrt(2) * np.ravel(landscapes) Xfit.append(landscapes) return np.stack(Xfit, axis=0) @@ -217,16 +207,6 @@ class Silhouette(BaseEstimator, TransformerMixin): self.weight, self.resolution, self.sample_range = weight, resolution, sample_range self.im_range = None - from pykeops.numpy import Genred - silhouette_formula = "normalized_weights * ReLU(heights - Abs(x_values - midpoints))" - variables = [ - "normalized_weights = Vi(1)", - "heights = Vi(1)", - "midpoints = Vi(1)", - "x_values = Vj(1)", - ] - self.silhouette = Genred(silhouette_formula, variables, reduction_op="Sum", axis=0) - def fit(self, X, y=None): """ Fit the Silhouette class on a list of persistence diagrams: if any of the values in **sample_range** is numpy.nan, replace it with the corresponding value computed on the given list of persistence diagrams. @@ -249,19 +229,19 @@ class Silhouette(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of diagrams) x (**resolution**): output persistence silhouettes. """ - silhouettes_list = [] + Xfit = [] x_values = self.im_range + for i, diag in enumerate(X): midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. weights = np.array([self.weight(point) for point in diag]) - weights /= np.sum(weights) + total_weight = np.sum(weights) - silhouettes_list.append( - np.sqrt(2) * self.silhouette(weights[:, None], heights[:, None], - midpoints[:, None], x_values[:, None])[:, 0] - ) + tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) + silhouette = np.sum(weights[None, :] / total_weight * tent_functions, axis=1) + Xfit.append(silhouette * np.sqrt(2)) - return np.stack(silhouettes_list, axis=0) + return np.stack(Xfit, axis=0) def __call__(self, diag): """ -- cgit v1.2.3 From 96f1dc0310d146e2bc368624683f46f1ca9a5116 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Wed, 14 Sep 2022 10:29:22 +0200 Subject: c++17 is the new standard and MSVC 17 also --- src/cmake/modules/GUDHI_compilation_flags.cmake | 2 +- src/common/doc/installation.h | 4 ++-- src/python/doc/installation.rst | 2 +- 3 files changed, 4 insertions(+), 4 deletions(-) (limited to 'src') diff --git a/src/cmake/modules/GUDHI_compilation_flags.cmake b/src/cmake/modules/GUDHI_compilation_flags.cmake index 567fbc40..e2d3d872 100644 --- a/src/cmake/modules/GUDHI_compilation_flags.cmake +++ b/src/cmake/modules/GUDHI_compilation_flags.cmake @@ -11,7 +11,7 @@ macro(add_cxx_compiler_flag _flag) endif() endmacro() -set (CMAKE_CXX_STANDARD 14) +set (CMAKE_CXX_STANDARD 17) enable_testing() diff --git a/src/common/doc/installation.h b/src/common/doc/installation.h index 28526498..a05ac588 100644 --- a/src/common/doc/installation.h +++ b/src/common/doc/installation.h @@ -5,9 +5,9 @@ * Examples of GUDHI headers inclusion can be found in \ref utilities. * * \section compiling Compiling - * The library uses c++14 and requires Boost ≥ 1.66.0 + * The library uses c++17 and requires Boost ≥ 1.66.0 * and CMake ≥ 3.5. - * It is a multi-platform library and compiles on Linux, Mac OSX and Visual Studio 2015. + * It is a multi-platform library and compiles on Linux, Mac OSX and Visual Studio 2017. * * \subsection utilities Utilities and examples * To build the utilities, run the following commands in a terminal: diff --git a/src/python/doc/installation.rst b/src/python/doc/installation.rst index 4eefd415..b704f778 100644 --- a/src/python/doc/installation.rst +++ b/src/python/doc/installation.rst @@ -39,7 +39,7 @@ If you are instead using a git checkout, beware that the paths are a bit different, and in particular the `python/` subdirectory is actually `src/python/` there. -The library uses c++14 and requires `Boost `_ :math:`\geq` 1.66.0, +The library uses c++17 and requires `Boost `_ :math:`\geq` 1.66.0, `CMake `_ :math:`\geq` 3.5 to generate makefiles, Python :math:`\geq` 3.5, `NumPy `_ :math:`\geq` 1.15.0, `Cython `_ :math:`\geq` 0.27 and `pybind11 `_ to compile the GUDHI Python module. -- cgit v1.2.3 From eee4180428d97ffef99e86123134835682181594 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Wed, 14 Sep 2022 10:30:19 +0200 Subject: std::ofstream::streampos is deprecated with c++17 --- src/Tangential_complex/benchmark/XML_exporter.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/Tangential_complex/benchmark/XML_exporter.h b/src/Tangential_complex/benchmark/XML_exporter.h index 16b62eb6..38fe049f 100644 --- a/src/Tangential_complex/benchmark/XML_exporter.h +++ b/src/Tangential_complex/benchmark/XML_exporter.h @@ -157,7 +157,7 @@ class Streaming_XML_exporter { m_xml_fstream << " " << std::endl; // Save current pointer position - std::ofstream::streampos pos = m_xml_fstream.tellp(); + auto pos = m_xml_fstream.tellp(); // Close the XML file (temporarily) so that the XML file is always correct m_xml_fstream << "" << std::endl; // Restore the pointer position so that the next "add_element" will overwrite -- cgit v1.2.3 From a034f774490ea06c018726ed58a0cfbe9d12d97c Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Mon, 26 Sep 2022 13:40:57 +0200 Subject: Cech_persistence utility for safe/exact/fast mode --- src/Cech_complex/utilities/CMakeLists.txt | 20 +++++++- src/Cech_complex/utilities/cech_persistence.cpp | 67 ++++++++++++++++++------- src/Cech_complex/utilities/cechcomplex.md | 8 ++- 3 files changed, 74 insertions(+), 21 deletions(-) (limited to 'src') diff --git a/src/Cech_complex/utilities/CMakeLists.txt b/src/Cech_complex/utilities/CMakeLists.txt index e80a698e..64557cee 100644 --- a/src/Cech_complex/utilities/CMakeLists.txt +++ b/src/Cech_complex/utilities/CMakeLists.txt @@ -9,8 +9,24 @@ if (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 5.0.1) target_link_libraries(cech_persistence ${TBB_LIBRARIES}) endif() - add_test(NAME Cech_complex_utility_from_rips_on_tore_3D COMMAND $ - "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" "-r" "0.25" "-m" "0.5" "-d" "3" "-p" "3") + add_test(NAME Cech_complex_utility_from_rips_on_tore_3D_safe COMMAND $ + "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" "-r" "0.25" "-m" "0.5" "-d" "3" "-p" "3" "-o" "safe.pers") + add_test(NAME Cech_complex_utility_from_rips_on_tore_3D_fast COMMAND $ + "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" "-r" "0.25" "-m" "0.5" "-d" "3" "-p" "3" "-o" "fast.pers" "-f") + add_test(NAME Cech_complex_utility_from_rips_on_tore_3D_exact COMMAND $ + "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" "-r" "0.25" "-m" "0.5" "-d" "3" "-p" "3" "-o" "exact.pers" "-e") + + if (DIFF_PATH) + add_test(Cech_complex_utilities_diff_exact ${DIFF_PATH} + "exact.pers" "safe.pers") + set_tests_properties(Cech_complex_utilities_diff_exact PROPERTIES DEPENDS + "Cech_complex_utility_from_rips_on_tore_3D_safe;Cech_complex_utility_from_rips_on_tore_3D_exact") + + add_test(Cech_complex_utilities_diff_fast ${DIFF_PATH} + "fast.pers" "safe.pers") + set_tests_properties(Cech_complex_utilities_diff_fast PROPERTIES DEPENDS + "Cech_complex_utility_from_rips_on_tore_3D_safe;Cech_complex_utility_from_rips_on_tore_3D_fast") + endif() install(TARGETS cech_persistence DESTINATION bin) endif() diff --git a/src/Cech_complex/utilities/cech_persistence.cpp b/src/Cech_complex/utilities/cech_persistence.cpp index 75d10c0f..a07ba212 100644 --- a/src/Cech_complex/utilities/cech_persistence.cpp +++ b/src/Cech_complex/utilities/cech_persistence.cpp @@ -16,6 +16,7 @@ #include #include // For EXACT or SAFE version +#include // For FAST version #include #include @@ -25,41 +26,66 @@ using Simplex_tree = Gudhi::Simplex_tree; using Filtration_value = Simplex_tree::Filtration_value; -using Kernel = CGAL::Epeck_d; -using Point = typename Kernel::Point_d; -using Points_off_reader = Gudhi::Points_off_reader; -using Cech_complex = Gudhi::cech_complex::Cech_complex; using Field_Zp = Gudhi::persistent_cohomology::Field_Zp; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology; -void program_options(int argc, char* argv[], std::string& off_file_points, std::string& filediag, - Filtration_value& max_radius, int& dim_max, int& p, Filtration_value& min_persistence); +void program_options(int argc, char* argv[], std::string& off_file_points, bool& exact, bool& fast, + std::string& filediag, Filtration_value& max_radius, int& dim_max, int& p, + Filtration_value& min_persistence); + +template +Simplex_tree create_simplex_tree(const std::string &off_file_points, bool exact_version, + Filtration_value max_radius, int dim_max) { + using Point = typename Kernel::Point_d; + using Points_off_reader = Gudhi::Points_off_reader; + using Cech_complex = Gudhi::cech_complex::Cech_complex; + + Simplex_tree stree; + + Points_off_reader off_reader(off_file_points); + Cech_complex cech_complex_from_file(off_reader.get_point_cloud(), max_radius, exact_version); + cech_complex_from_file.create_complex(stree, dim_max); + + return stree; +} int main(int argc, char* argv[]) { std::string off_file_points; std::string filediag; + bool exact_version = false; + bool fast_version = false; Filtration_value max_radius; int dim_max; int p; Filtration_value min_persistence; - program_options(argc, argv, off_file_points, filediag, max_radius, dim_max, p, min_persistence); + program_options(argc, argv, off_file_points, exact_version, fast_version, filediag, max_radius, dim_max, p, + min_persistence); - Points_off_reader off_reader(off_file_points); - Cech_complex cech_complex_from_file(off_reader.get_point_cloud(), max_radius); + if ((exact_version) && (fast_version)) { + std::cerr << "You cannot set the exact and the fast version." << std::endl; + exit(-1); + } - // Construct the Cech complex in a Simplex Tree - Simplex_tree simplex_tree; + Simplex_tree stree; + if (fast_version) { + // WARNING : CGAL::Epick_d is fast but not safe (unlike CGAL::Epeck_d) + // (i.e. when the points are on a grid) + using Fast_kernel = CGAL::Epick_d; + stree = create_simplex_tree(off_file_points, exact_version, max_radius, dim_max); + } else { + using Kernel = CGAL::Epeck_d; + stree = create_simplex_tree(off_file_points, exact_version, max_radius, dim_max); + } - cech_complex_from_file.create_complex(simplex_tree, dim_max); - std::clog << "The complex contains " << simplex_tree.num_simplices() << " simplices \n"; - std::clog << " and has dimension " << simplex_tree.dimension() << " \n"; + std::clog << "The complex contains " << stree.num_simplices() << " simplices \n"; + std::clog << " and has dimension " << stree.dimension() << " \n"; // Sort the simplices in the order of the filtration - simplex_tree.initialize_filtration(); + stree.initialize_filtration(); // Compute the persistence diagram of the complex - Persistent_cohomology pcoh(simplex_tree); + Persistent_cohomology pcoh(stree); // initializes the coefficient field for homology pcoh.init_coefficients(p); @@ -77,8 +103,9 @@ int main(int argc, char* argv[]) { return 0; } -void program_options(int argc, char* argv[], std::string& off_file_points, std::string& filediag, - Filtration_value& max_radius, int& dim_max, int& p, Filtration_value& min_persistence) { +void program_options(int argc, char* argv[], std::string& off_file_points, bool& exact, bool& fast, + std::string& filediag, Filtration_value& max_radius, int& dim_max, int& p, + Filtration_value& min_persistence) { namespace po = boost::program_options; po::options_description hidden("Hidden options"); hidden.add_options()("input-file", po::value(&off_file_points), @@ -86,6 +113,10 @@ void program_options(int argc, char* argv[], std::string& off_file_points, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( + "exact,e", po::bool_switch(&exact), + "To activate exact version of Cech complex (default is false, not available if fast is set)")( + "fast,f", po::bool_switch(&fast), + "To activate fast version of Cech complex (default is false, not available if exact is set)")( "output-file,o", po::value(&filediag)->default_value(std::string()), "Name of file in which the persistence diagram is written. Default print in std::clog")( "max-radius,r", diff --git a/src/Cech_complex/utilities/cechcomplex.md b/src/Cech_complex/utilities/cechcomplex.md index 821e4dad..0e82674d 100644 --- a/src/Cech_complex/utilities/cechcomplex.md +++ b/src/Cech_complex/utilities/cechcomplex.md @@ -26,7 +26,11 @@ a prime number). **Usage** -`cech_persistence [options] ` +`cech_persistence [options] ` + +where +`` is the path to the input point cloud in +[nOFF ASCII format]({{ site.officialurl }}/doc/latest/fileformats.html#FileFormatsOFF). **Allowed options** @@ -36,6 +40,8 @@ a prime number). * `-d [ --cpx-dimension ]` (default = 1) Maximal dimension of the Čech complex we want to compute. * `-p [ --field-charac ]` (default = 11) Characteristic p of the coefficient field Z/pZ for computing homology. * `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals. +* `-e [ --exact ]` for the exact computation version. +* `-f [ --fast ]` for the fast computation version. Beware: this program may use a lot of RAM and take a lot of time if `max-edge-length` is set to a large value. -- cgit v1.2.3 From 565cc1fb0c81b627fa6bec0e9ca76571663a7834 Mon Sep 17 00:00:00 2001 From: wreise Date: Fri, 7 Oct 2022 17:59:51 +0200 Subject: Test: resolution when nan in the range --- src/python/test/test_representations.py | 8 ++++++++ 1 file changed, 8 insertions(+) (limited to 'src') diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index 2ca72d07..823e8620 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -241,3 +241,11 @@ def test_landscape_numeric(): lds = Landscape(num_landscapes=2, resolution=9, sample_range=[2., 6.]) lds_dgm = lds(dgm) assert np.all(np.isclose(lds_dgm, lds_ref)) + + +def test_landscape_nan_range(): + dgm = np.array([[2., 6.], [3., 5.]]) + lds = Landscape(num_landscapes=2, resolution=9, sample_range=[np.nan, 6.]) + lds_dgm = lds(dgm) + assert (lds.sample_range[0] == 2) & (lds.sample_range[1] == 6) + assert lds.new_resolution == 10 -- cgit v1.2.3 From d9dfffdb580ab865a829fce851779f33fa47e4f7 Mon Sep 17 00:00:00 2001 From: wreise Date: Fri, 7 Oct 2022 18:13:26 +0200 Subject: Add triming of the range; Marcs' comments --- src/python/gudhi/representations/vector_methods.py | 31 ++++++++++++++-------- src/python/test/test_representations.py | 2 +- 2 files changed, 21 insertions(+), 12 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 8c8b46db..8114885e 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -123,6 +123,15 @@ def _automatic_sample_range(sample_range, X, y): pass return sample_range + +def trim_on_edges(x, are_endpoints_nan): + if are_endpoints_nan[0]: + x = x[1:] + if are_endpoints_nan[1]: + x = x[:-1] + return x + + class Landscape(BaseEstimator, TransformerMixin): """ This is a class for computing persistence landscapes from a list of persistence diagrams. A persistence landscape is a collection of 1D piecewise-linear functions computed from the rank function associated to the persistence diagram. These piecewise-linear functions are then sampled evenly on a given range and the corresponding vectors of samples are concatenated and returned. See http://jmlr.org/papers/v16/bubenik15a.html for more details. @@ -149,6 +158,8 @@ class Landscape(BaseEstimator, TransformerMixin): y (n x 1 array): persistence diagram labels (unused). """ self.sample_range = _automatic_sample_range(np.array(self.sample_range), X, y) + self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) + self.im_range = trim_on_edges(self.im_range, self.nan_in_range) return self def transform(self, X): @@ -163,17 +174,13 @@ class Landscape(BaseEstimator, TransformerMixin): """ Xfit = [] - x_values = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) + x_values = self.im_range for i, diag in enumerate(X): midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) tent_functions.partition(diag.shape[0] - self.num_landscapes, axis=1) - landscapes = np.sort(tent_functions, axis=1)[:, -self.num_landscapes:][:, ::-1].T + landscapes = tent_functions[:, -self.num_landscapes:][:, ::-1].T - if self.nan_in_range[0]: - landscapes = landscapes[:, 1:] - if self.nan_in_range[1]: - landscapes = landscapes[:, :-1] landscapes = np.sqrt(2) * np.ravel(landscapes) Xfit.append(landscapes) @@ -189,7 +196,7 @@ class Landscape(BaseEstimator, TransformerMixin): Returns: numpy array with shape (number of samples = **num_landscapes** x **resolution**): output persistence landscape. """ - return self.fit_transform([diag])[0,:] + return self.fit_transform([diag])[0, :] class Silhouette(BaseEstimator, TransformerMixin): """ @@ -205,7 +212,8 @@ class Silhouette(BaseEstimator, TransformerMixin): sample_range ([double, double]): minimum and maximum for the weighted average domain, of the form [x_min, x_max] (default [numpy.nan, numpy.nan]). It is the interval on which samples will be drawn evenly. If one of the values is numpy.nan, it can be computed from the persistence diagrams with the fit() method. """ self.weight, self.resolution, self.sample_range = weight, resolution, sample_range - self.im_range = None + self.nan_in_range = np.isnan(np.array(self.sample_range)) + self.new_resolution = self.resolution + self.nan_in_range.sum() def fit(self, X, y=None): """ @@ -216,7 +224,8 @@ class Silhouette(BaseEstimator, TransformerMixin): y (n x 1 array): persistence diagram labels (unused). """ self.sample_range = _automatic_sample_range(np.array(self.sample_range), X, y) - self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.resolution) + self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) + self.im_range = trim_on_edges(self.im_range, self.nan_in_range) return self def transform(self, X): @@ -232,9 +241,9 @@ class Silhouette(BaseEstimator, TransformerMixin): Xfit = [] x_values = self.im_range - for i, diag in enumerate(X): + for diag in X: midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. - weights = np.array([self.weight(point) for point in diag]) + weights = np.array([self.weight(pt) for pt in diag]) total_weight = np.sum(weights) tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index 823e8620..948d7804 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -209,7 +209,7 @@ def test_silhouette_multiplication_invariance(): def test_silhouette_numeric(): dgm = np.array([[2., 3.], [5., 6.]]) - slt = Silhouette(resolution=9, weight=pow(1)) + slt = Silhouette(resolution=9, weight=pow(1), sample_range=[2., 6.]) #slt.fit([dgm]) # x_values = array([2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.]) -- cgit v1.2.3 From 059ff0c42a069c744ed121c948bc3d39b5cc7f10 Mon Sep 17 00:00:00 2001 From: wreise Date: Fri, 7 Oct 2022 18:24:32 +0200 Subject: Remove enumerate --- src/python/gudhi/representations/vector_methods.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 8114885e..5ea4ea48 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -175,7 +175,7 @@ class Landscape(BaseEstimator, TransformerMixin): Xfit = [] x_values = self.im_range - for i, diag in enumerate(X): + for diag in X: midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) tent_functions.partition(diag.shape[0] - self.num_landscapes, axis=1) -- cgit v1.2.3 From 4aac9e03c400bd43f237504cf4ff9d25f041e473 Mon Sep 17 00:00:00 2001 From: wreise Date: Wed, 12 Oct 2022 15:55:22 +0200 Subject: Clean argpartition --- src/python/gudhi/representations/vector_methods.py | 20 ++++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 5ea4ea48..3a91eccd 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -124,7 +124,7 @@ def _automatic_sample_range(sample_range, X, y): return sample_range -def trim_on_edges(x, are_endpoints_nan): +def _trim_on_edges(x, are_endpoints_nan): if are_endpoints_nan[0]: x = x[1:] if are_endpoints_nan[1]: @@ -159,7 +159,7 @@ class Landscape(BaseEstimator, TransformerMixin): """ self.sample_range = _automatic_sample_range(np.array(self.sample_range), X, y) self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) - self.im_range = trim_on_edges(self.im_range, self.nan_in_range) + self.im_range = _trim_on_edges(self.im_range, self.nan_in_range) return self def transform(self, X): @@ -178,9 +178,17 @@ class Landscape(BaseEstimator, TransformerMixin): for diag in X: midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) - tent_functions.partition(diag.shape[0] - self.num_landscapes, axis=1) - landscapes = tent_functions[:, -self.num_landscapes:][:, ::-1].T - + n_points = diag.shape[0] + # Get indices of largest elements: can't take more than n_points - 1 (the last ones are in the right position) + argpartition = np.argpartition(-tent_functions, min(self.num_landscapes, n_points-1), axis=1) + landscapes = np.take_along_axis(tent_functions, argpartition, axis=1) + landscapes = landscapes[:, :min(self.num_landscapes, n_points)].T + + # Complete the array with zeros to get the right number of landscapes + if self.num_landscapes > n_points: + landscapes = np.concatenate([ + landscapes, np.zeros((self.num_landscapes-n_points, *landscapes.shape[1:])) + ], axis=0) landscapes = np.sqrt(2) * np.ravel(landscapes) Xfit.append(landscapes) @@ -225,7 +233,7 @@ class Silhouette(BaseEstimator, TransformerMixin): """ self.sample_range = _automatic_sample_range(np.array(self.sample_range), X, y) self.im_range = np.linspace(self.sample_range[0], self.sample_range[1], self.new_resolution) - self.im_range = trim_on_edges(self.im_range, self.nan_in_range) + self.im_range = _trim_on_edges(self.im_range, self.nan_in_range) return self def transform(self, X): -- cgit v1.2.3 From e94e490554ce56040b998893977255a704d35e59 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Sat, 15 Oct 2022 01:18:56 +0200 Subject: Fix badly centered barcode --- src/python/gudhi/persistence_graphical_tools.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/persistence_graphical_tools.py b/src/python/gudhi/persistence_graphical_tools.py index 21275cdd..e438aa66 100644 --- a/src/python/gudhi/persistence_graphical_tools.py +++ b/src/python/gudhi/persistence_graphical_tools.py @@ -194,19 +194,21 @@ def plot_persistence_barcode( y=[(death - birth) if death != float("inf") else (infinity - birth) for (dim,(birth,death)) in persistence] c=[colormap[dim] for (dim,(birth,death)) in persistence] - axes.barh(list(reversed(range(len(x)))), y, height=0.8, left=x, alpha=alpha, color=c, linewidth=0) + axes.barh(range(len(x)), y, left=x, alpha=alpha, color=c, linewidth=0) if legend: - dimensions = list(set(item[0] for item in persistence)) + dimensions = set(item[0] for item in persistence) axes.legend( handles=[mpatches.Patch(color=colormap[dim], label=str(dim)) for dim in dimensions], loc="lower right", ) axes.set_title("Persistence barcode", fontsize=fontsize) + axes.set_yticks([]) + axes.invert_yaxis() # Ends plot on infinity value and starts a little bit before min_birth if len(x) != 0: - axes.axis([axis_start, infinity, 0, len(x)]) + axes.set_xlim((axis_start, infinity)) return axes except ImportError as import_error: -- cgit v1.2.3 From e3d4dd49c848eedecd70bd743ae4a94b98b7f086 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Sat, 15 Oct 2022 15:38:01 +0200 Subject: Remove unnecessary includes --- src/Alpha_complex/test/Alpha_complex_unit_test.cpp | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/test/Alpha_complex_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_unit_test.cpp index f74ad217..1ac0093d 100644 --- a/src/Alpha_complex/test/Alpha_complex_unit_test.cpp +++ b/src/Alpha_complex/test/Alpha_complex_unit_test.cpp @@ -13,18 +13,14 @@ #include #include -#include #include #include -#include // float comparison -#include +#include // std::out_of_range #include #include #include -// to construct a simplex_tree from Delaunay_triangulation -#include #include #include -- cgit v1.2.3 From b3b5d1da1ad441c8a76e8137b25e752e5654b938 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Sat, 15 Oct 2022 15:40:41 +0200 Subject: Duplicate Alpha_complex_unit_test.cpp history in Alpha_complex_dim3_unit_test.cpp history. --- .../test/Alpha_complex_dim3_unit_test.cpp | 308 +++++++++++++++++++++ 1 file changed, 308 insertions(+) create mode 100644 src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp (limited to 'src') diff --git a/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp new file mode 100644 index 00000000..1ac0093d --- /dev/null +++ b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp @@ -0,0 +1,308 @@ +/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT. + * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. + * Author(s): Vincent Rouvreau + * + * Copyright (C) 2015 Inria + * + * Modification(s): + * - YYYY/MM Author: Description of the modification + */ + +#define BOOST_TEST_DYN_LINK +#define BOOST_TEST_MODULE "alpha_complex" +#include +#include + +#include +#include + +#include // std::out_of_range +#include +#include + +#include +#include +#include + +// Use dynamic_dimension_tag for the user to be able to set dimension +typedef CGAL::Epeck_d< CGAL::Dynamic_dimension_tag > Exact_kernel_d; +// Use static dimension_tag for the user not to be able to set dimension +typedef CGAL::Epeck_d< CGAL::Dimension_tag<3> > Exact_kernel_s; +// Use dynamic_dimension_tag for the user to be able to set dimension +typedef CGAL::Epick_d< CGAL::Dynamic_dimension_tag > Inexact_kernel_d; +// Use static dimension_tag for the user not to be able to set dimension +typedef CGAL::Epick_d< CGAL::Dimension_tag<3> > Inexact_kernel_s; +// The triangulation uses the default instantiation of the TriangulationDataStructure template parameter + +typedef boost::mpl::list list_of_kernel_variants; + +BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_OFF_file, TestedKernel, list_of_kernel_variants) { + // ---------------------------------------------------------------------------- + // + // Init of an alpha-complex from a OFF file + // + // ---------------------------------------------------------------------------- + std::string off_file_name("alphacomplexdoc.off"); + double max_alpha_square_value = 60.0; + std::clog << "========== OFF FILE NAME = " << off_file_name << " - alpha²=" << + max_alpha_square_value << "==========" << std::endl; + + Gudhi::alpha_complex::Alpha_complex alpha_complex_from_file(off_file_name); + + Gudhi::Simplex_tree<> simplex_tree_60; + BOOST_CHECK(alpha_complex_from_file.create_complex(simplex_tree_60, max_alpha_square_value)); + + std::clog << "alpha_complex_from_file.num_vertices()=" << alpha_complex_from_file.num_vertices() << std::endl; + BOOST_CHECK(alpha_complex_from_file.num_vertices() == 7); + + std::clog << "simplex_tree_60.dimension()=" << simplex_tree_60.dimension() << std::endl; + BOOST_CHECK(simplex_tree_60.dimension() == 2); + + std::clog << "simplex_tree_60.num_vertices()=" << simplex_tree_60.num_vertices() << std::endl; + BOOST_CHECK(simplex_tree_60.num_vertices() == 7); + + std::clog << "simplex_tree_60.num_simplices()=" << simplex_tree_60.num_simplices() << std::endl; + BOOST_CHECK(simplex_tree_60.num_simplices() == 25); + + max_alpha_square_value = 59.0; + std::clog << "========== OFF FILE NAME = " << off_file_name << " - alpha²=" << + max_alpha_square_value << "==========" << std::endl; + + Gudhi::Simplex_tree<> simplex_tree_59; + BOOST_CHECK(alpha_complex_from_file.create_complex(simplex_tree_59, max_alpha_square_value)); + + std::clog << "alpha_complex_from_file.num_vertices()=" << alpha_complex_from_file.num_vertices() << std::endl; + BOOST_CHECK(alpha_complex_from_file.num_vertices() == 7); + + std::clog << "simplex_tree_59.dimension()=" << simplex_tree_59.dimension() << std::endl; + BOOST_CHECK(simplex_tree_59.dimension() == 2); + + std::clog << "simplex_tree_59.num_vertices()=" << simplex_tree_59.num_vertices() << std::endl; + BOOST_CHECK(simplex_tree_59.num_vertices() == 7); + + std::clog << "simplex_tree_59.num_simplices()=" << simplex_tree_59.num_simplices() << std::endl; + BOOST_CHECK(simplex_tree_59.num_simplices() == 23); +} + +// Use static dimension_tag for the user not to be able to set dimension +typedef CGAL::Epeck_d< CGAL::Dimension_tag<4> > Kernel_4; +typedef Kernel_4::Point_d Point_4; +typedef std::vector Vector_4_Points; + +bool is_point_in_list(Vector_4_Points points_list, Point_4 point) { + for (auto& point_in_list : points_list) { + if (point_in_list == point) { + return true; // point found + } + } + return false; // point not found +} + +BOOST_AUTO_TEST_CASE(Alpha_complex_from_points) { + // ---------------------------------------------------------------------------- + // Init of a list of points + // ---------------------------------------------------------------------------- + Vector_4_Points points; + std::vector coords = { 0.0, 0.0, 0.0, 1.0 }; + points.push_back(Point_4(coords.begin(), coords.end())); + coords = { 0.0, 0.0, 1.0, 0.0 }; + points.push_back(Point_4(coords.begin(), coords.end())); + coords = { 0.0, 1.0, 0.0, 0.0 }; + points.push_back(Point_4(coords.begin(), coords.end())); + coords = { 1.0, 0.0, 0.0, 0.0 }; + points.push_back(Point_4(coords.begin(), coords.end())); + + // ---------------------------------------------------------------------------- + // Init of an alpha complex from the list of points + // ---------------------------------------------------------------------------- + Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); + + std::clog << "========== Alpha_complex_from_points ==========" << std::endl; + + Gudhi::Simplex_tree<> simplex_tree; + BOOST_CHECK(alpha_complex_from_points.create_complex(simplex_tree)); + + std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; + BOOST_CHECK(alpha_complex_from_points.num_vertices() == points.size()); + + // Another way to check num_simplices + std::clog << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; + int num_simplices = 0; + for (auto f_simplex : simplex_tree.filtration_simplex_range()) { + num_simplices++; + std::clog << " ( "; + for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex)) { + std::clog << vertex << " "; + } + std::clog << ") -> " << "[" << simplex_tree.filtration(f_simplex) << "] "; + std::clog << std::endl; + } + BOOST_CHECK(num_simplices == 15); + std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; + BOOST_CHECK(simplex_tree.num_simplices() == 15); + + std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; + BOOST_CHECK(simplex_tree.dimension() == 3); + std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; + BOOST_CHECK(simplex_tree.num_vertices() == points.size()); + + for (auto f_simplex : simplex_tree.filtration_simplex_range()) { + switch (simplex_tree.dimension(f_simplex)) { + case 0: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 0.0); + break; + case 1: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 1.0/2.0); + break; + case 2: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 2.0/3.0); + break; + case 3: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 3.0/4.0); + break; + default: + BOOST_CHECK(false); // Shall not happen + break; + } + } + + Point_4 p0 = alpha_complex_from_points.get_point(0); + std::clog << "alpha_complex_from_points.get_point(0)=" << p0 << std::endl; + BOOST_CHECK(4 == p0.dimension()); + BOOST_CHECK(is_point_in_list(points, p0)); + + Point_4 p1 = alpha_complex_from_points.get_point(1); + std::clog << "alpha_complex_from_points.get_point(1)=" << p1 << std::endl; + BOOST_CHECK(4 == p1.dimension()); + BOOST_CHECK(is_point_in_list(points, p1)); + + Point_4 p2 = alpha_complex_from_points.get_point(2); + std::clog << "alpha_complex_from_points.get_point(2)=" << p2 << std::endl; + BOOST_CHECK(4 == p2.dimension()); + BOOST_CHECK(is_point_in_list(points, p2)); + + Point_4 p3 = alpha_complex_from_points.get_point(3); + std::clog << "alpha_complex_from_points.get_point(3)=" << p3 << std::endl; + BOOST_CHECK(4 == p3.dimension()); + BOOST_CHECK(is_point_in_list(points, p3)); + + // Test to the limit + BOOST_CHECK_THROW (alpha_complex_from_points.get_point(4), std::out_of_range); + BOOST_CHECK_THROW (alpha_complex_from_points.get_point(-1), std::out_of_range); + BOOST_CHECK_THROW (alpha_complex_from_points.get_point(1234), std::out_of_range); + + // Test after prune_above_filtration + bool modified = simplex_tree.prune_above_filtration(0.6); + BOOST_CHECK(modified); + + // Another way to check num_simplices + std::clog << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; + num_simplices = 0; + for (auto f_simplex : simplex_tree.filtration_simplex_range()) { + num_simplices++; + std::clog << " ( "; + for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex)) { + std::clog << vertex << " "; + } + std::clog << ") -> " << "[" << simplex_tree.filtration(f_simplex) << "] "; + std::clog << std::endl; + } + BOOST_CHECK(num_simplices == 10); + std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; + BOOST_CHECK(simplex_tree.num_simplices() == 10); + + std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; + BOOST_CHECK(simplex_tree.dimension() == 1); + std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; + BOOST_CHECK(simplex_tree.num_vertices() == 4); + + for (auto f_simplex : simplex_tree.filtration_simplex_range()) { + switch (simplex_tree.dimension(f_simplex)) { + case 0: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 0.0); + break; + case 1: + GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 1.0/2.0); + break; + default: + BOOST_CHECK(false); // Shall not happen + break; + } + } + +} + +BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_empty_points, TestedKernel, list_of_kernel_variants) { + std::clog << "========== Alpha_complex_from_empty_points ==========" << std::endl; + + // ---------------------------------------------------------------------------- + // Init of an empty list of points + // ---------------------------------------------------------------------------- + std::vector points; + + // ---------------------------------------------------------------------------- + // Init of an alpha complex from the list of points + // ---------------------------------------------------------------------------- + Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); + + std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; + BOOST_CHECK(alpha_complex_from_points.num_vertices() == points.size()); + + // Test to the limit + BOOST_CHECK_THROW (alpha_complex_from_points.get_point(0), std::out_of_range); + + Gudhi::Simplex_tree<> simplex_tree; + BOOST_CHECK(!alpha_complex_from_points.create_complex(simplex_tree)); + + std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; + BOOST_CHECK(simplex_tree.num_simplices() == 0); + + std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; + BOOST_CHECK(simplex_tree.dimension() == -1); + + std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; + BOOST_CHECK(simplex_tree.num_vertices() == points.size()); +} + +using Inexact_kernel_2 = CGAL::Epick_d< CGAL::Dimension_tag<2> >; +using Exact_kernel_2 = CGAL::Epeck_d< CGAL::Dimension_tag<2> >; +using list_of_kernel_2_variants = boost::mpl::list; + +BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_with_duplicated_points, TestedKernel, list_of_kernel_2_variants) { + std::clog << "========== Alpha_complex_with_duplicated_points ==========" << std::endl; + + using Point = typename TestedKernel::Point_d; + using Vector_of_points = std::vector; + + // ---------------------------------------------------------------------------- + // Init of a list of points + // ---------------------------------------------------------------------------- + Vector_of_points points; + points.push_back(Point(1.0, 1.0)); + points.push_back(Point(7.0, 0.0)); + points.push_back(Point(4.0, 6.0)); + points.push_back(Point(9.0, 6.0)); + points.push_back(Point(0.0, 14.0)); + points.push_back(Point(2.0, 19.0)); + points.push_back(Point(9.0, 17.0)); + // duplicated points + points.push_back(Point(1.0, 1.0)); + points.push_back(Point(7.0, 0.0)); + + // ---------------------------------------------------------------------------- + // Init of an alpha complex from the list of points + // ---------------------------------------------------------------------------- + std::clog << "Init" << std::endl; + Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); + + Gudhi::Simplex_tree<> simplex_tree; + std::clog << "create_complex" << std::endl; + BOOST_CHECK(alpha_complex_from_points.create_complex(simplex_tree)); + + std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; + BOOST_CHECK(alpha_complex_from_points.num_vertices() < points.size()); + + std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() + << std::endl; + BOOST_CHECK(simplex_tree.num_vertices() < points.size()); +} -- cgit v1.2.3 From b380380aecb601fbc601710da1a62cca202304db Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Sat, 15 Oct 2022 15:46:09 +0200 Subject: Split Alpha_complex_unit_test.cpp into 2 files It was taking too much memory to compile --- .../test/Alpha_complex_dim3_unit_test.cpp | 191 --------------------- src/Alpha_complex/test/Alpha_complex_unit_test.cpp | 91 ---------- src/Alpha_complex/test/CMakeLists.txt | 6 +- 3 files changed, 5 insertions(+), 283 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp index 1ac0093d..0085ae67 100644 --- a/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp +++ b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp @@ -22,7 +22,6 @@ #include #include -#include // Use dynamic_dimension_tag for the user to be able to set dimension typedef CGAL::Epeck_d< CGAL::Dynamic_dimension_tag > Exact_kernel_d; @@ -84,153 +83,6 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_OFF_file, TestedKernel, list_of BOOST_CHECK(simplex_tree_59.num_simplices() == 23); } -// Use static dimension_tag for the user not to be able to set dimension -typedef CGAL::Epeck_d< CGAL::Dimension_tag<4> > Kernel_4; -typedef Kernel_4::Point_d Point_4; -typedef std::vector Vector_4_Points; - -bool is_point_in_list(Vector_4_Points points_list, Point_4 point) { - for (auto& point_in_list : points_list) { - if (point_in_list == point) { - return true; // point found - } - } - return false; // point not found -} - -BOOST_AUTO_TEST_CASE(Alpha_complex_from_points) { - // ---------------------------------------------------------------------------- - // Init of a list of points - // ---------------------------------------------------------------------------- - Vector_4_Points points; - std::vector coords = { 0.0, 0.0, 0.0, 1.0 }; - points.push_back(Point_4(coords.begin(), coords.end())); - coords = { 0.0, 0.0, 1.0, 0.0 }; - points.push_back(Point_4(coords.begin(), coords.end())); - coords = { 0.0, 1.0, 0.0, 0.0 }; - points.push_back(Point_4(coords.begin(), coords.end())); - coords = { 1.0, 0.0, 0.0, 0.0 }; - points.push_back(Point_4(coords.begin(), coords.end())); - - // ---------------------------------------------------------------------------- - // Init of an alpha complex from the list of points - // ---------------------------------------------------------------------------- - Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); - - std::clog << "========== Alpha_complex_from_points ==========" << std::endl; - - Gudhi::Simplex_tree<> simplex_tree; - BOOST_CHECK(alpha_complex_from_points.create_complex(simplex_tree)); - - std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; - BOOST_CHECK(alpha_complex_from_points.num_vertices() == points.size()); - - // Another way to check num_simplices - std::clog << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; - int num_simplices = 0; - for (auto f_simplex : simplex_tree.filtration_simplex_range()) { - num_simplices++; - std::clog << " ( "; - for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex)) { - std::clog << vertex << " "; - } - std::clog << ") -> " << "[" << simplex_tree.filtration(f_simplex) << "] "; - std::clog << std::endl; - } - BOOST_CHECK(num_simplices == 15); - std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; - BOOST_CHECK(simplex_tree.num_simplices() == 15); - - std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; - BOOST_CHECK(simplex_tree.dimension() == 3); - std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; - BOOST_CHECK(simplex_tree.num_vertices() == points.size()); - - for (auto f_simplex : simplex_tree.filtration_simplex_range()) { - switch (simplex_tree.dimension(f_simplex)) { - case 0: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 0.0); - break; - case 1: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 1.0/2.0); - break; - case 2: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 2.0/3.0); - break; - case 3: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 3.0/4.0); - break; - default: - BOOST_CHECK(false); // Shall not happen - break; - } - } - - Point_4 p0 = alpha_complex_from_points.get_point(0); - std::clog << "alpha_complex_from_points.get_point(0)=" << p0 << std::endl; - BOOST_CHECK(4 == p0.dimension()); - BOOST_CHECK(is_point_in_list(points, p0)); - - Point_4 p1 = alpha_complex_from_points.get_point(1); - std::clog << "alpha_complex_from_points.get_point(1)=" << p1 << std::endl; - BOOST_CHECK(4 == p1.dimension()); - BOOST_CHECK(is_point_in_list(points, p1)); - - Point_4 p2 = alpha_complex_from_points.get_point(2); - std::clog << "alpha_complex_from_points.get_point(2)=" << p2 << std::endl; - BOOST_CHECK(4 == p2.dimension()); - BOOST_CHECK(is_point_in_list(points, p2)); - - Point_4 p3 = alpha_complex_from_points.get_point(3); - std::clog << "alpha_complex_from_points.get_point(3)=" << p3 << std::endl; - BOOST_CHECK(4 == p3.dimension()); - BOOST_CHECK(is_point_in_list(points, p3)); - - // Test to the limit - BOOST_CHECK_THROW (alpha_complex_from_points.get_point(4), std::out_of_range); - BOOST_CHECK_THROW (alpha_complex_from_points.get_point(-1), std::out_of_range); - BOOST_CHECK_THROW (alpha_complex_from_points.get_point(1234), std::out_of_range); - - // Test after prune_above_filtration - bool modified = simplex_tree.prune_above_filtration(0.6); - BOOST_CHECK(modified); - - // Another way to check num_simplices - std::clog << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; - num_simplices = 0; - for (auto f_simplex : simplex_tree.filtration_simplex_range()) { - num_simplices++; - std::clog << " ( "; - for (auto vertex : simplex_tree.simplex_vertex_range(f_simplex)) { - std::clog << vertex << " "; - } - std::clog << ") -> " << "[" << simplex_tree.filtration(f_simplex) << "] "; - std::clog << std::endl; - } - BOOST_CHECK(num_simplices == 10); - std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; - BOOST_CHECK(simplex_tree.num_simplices() == 10); - - std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; - BOOST_CHECK(simplex_tree.dimension() == 1); - std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; - BOOST_CHECK(simplex_tree.num_vertices() == 4); - - for (auto f_simplex : simplex_tree.filtration_simplex_range()) { - switch (simplex_tree.dimension(f_simplex)) { - case 0: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 0.0); - break; - case 1: - GUDHI_TEST_FLOAT_EQUALITY_CHECK(simplex_tree.filtration(f_simplex), 1.0/2.0); - break; - default: - BOOST_CHECK(false); // Shall not happen - break; - } - } - -} BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_empty_points, TestedKernel, list_of_kernel_variants) { std::clog << "========== Alpha_complex_from_empty_points ==========" << std::endl; @@ -263,46 +115,3 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_empty_points, TestedKernel, lis std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; BOOST_CHECK(simplex_tree.num_vertices() == points.size()); } - -using Inexact_kernel_2 = CGAL::Epick_d< CGAL::Dimension_tag<2> >; -using Exact_kernel_2 = CGAL::Epeck_d< CGAL::Dimension_tag<2> >; -using list_of_kernel_2_variants = boost::mpl::list; - -BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_with_duplicated_points, TestedKernel, list_of_kernel_2_variants) { - std::clog << "========== Alpha_complex_with_duplicated_points ==========" << std::endl; - - using Point = typename TestedKernel::Point_d; - using Vector_of_points = std::vector; - - // ---------------------------------------------------------------------------- - // Init of a list of points - // ---------------------------------------------------------------------------- - Vector_of_points points; - points.push_back(Point(1.0, 1.0)); - points.push_back(Point(7.0, 0.0)); - points.push_back(Point(4.0, 6.0)); - points.push_back(Point(9.0, 6.0)); - points.push_back(Point(0.0, 14.0)); - points.push_back(Point(2.0, 19.0)); - points.push_back(Point(9.0, 17.0)); - // duplicated points - points.push_back(Point(1.0, 1.0)); - points.push_back(Point(7.0, 0.0)); - - // ---------------------------------------------------------------------------- - // Init of an alpha complex from the list of points - // ---------------------------------------------------------------------------- - std::clog << "Init" << std::endl; - Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); - - Gudhi::Simplex_tree<> simplex_tree; - std::clog << "create_complex" << std::endl; - BOOST_CHECK(alpha_complex_from_points.create_complex(simplex_tree)); - - std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; - BOOST_CHECK(alpha_complex_from_points.num_vertices() < points.size()); - - std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() - << std::endl; - BOOST_CHECK(simplex_tree.num_vertices() < points.size()); -} diff --git a/src/Alpha_complex/test/Alpha_complex_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_unit_test.cpp index 1ac0093d..b474917f 100644 --- a/src/Alpha_complex/test/Alpha_complex_unit_test.cpp +++ b/src/Alpha_complex/test/Alpha_complex_unit_test.cpp @@ -24,66 +24,6 @@ #include #include -// Use dynamic_dimension_tag for the user to be able to set dimension -typedef CGAL::Epeck_d< CGAL::Dynamic_dimension_tag > Exact_kernel_d; -// Use static dimension_tag for the user not to be able to set dimension -typedef CGAL::Epeck_d< CGAL::Dimension_tag<3> > Exact_kernel_s; -// Use dynamic_dimension_tag for the user to be able to set dimension -typedef CGAL::Epick_d< CGAL::Dynamic_dimension_tag > Inexact_kernel_d; -// Use static dimension_tag for the user not to be able to set dimension -typedef CGAL::Epick_d< CGAL::Dimension_tag<3> > Inexact_kernel_s; -// The triangulation uses the default instantiation of the TriangulationDataStructure template parameter - -typedef boost::mpl::list list_of_kernel_variants; - -BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_OFF_file, TestedKernel, list_of_kernel_variants) { - // ---------------------------------------------------------------------------- - // - // Init of an alpha-complex from a OFF file - // - // ---------------------------------------------------------------------------- - std::string off_file_name("alphacomplexdoc.off"); - double max_alpha_square_value = 60.0; - std::clog << "========== OFF FILE NAME = " << off_file_name << " - alpha²=" << - max_alpha_square_value << "==========" << std::endl; - - Gudhi::alpha_complex::Alpha_complex alpha_complex_from_file(off_file_name); - - Gudhi::Simplex_tree<> simplex_tree_60; - BOOST_CHECK(alpha_complex_from_file.create_complex(simplex_tree_60, max_alpha_square_value)); - - std::clog << "alpha_complex_from_file.num_vertices()=" << alpha_complex_from_file.num_vertices() << std::endl; - BOOST_CHECK(alpha_complex_from_file.num_vertices() == 7); - - std::clog << "simplex_tree_60.dimension()=" << simplex_tree_60.dimension() << std::endl; - BOOST_CHECK(simplex_tree_60.dimension() == 2); - - std::clog << "simplex_tree_60.num_vertices()=" << simplex_tree_60.num_vertices() << std::endl; - BOOST_CHECK(simplex_tree_60.num_vertices() == 7); - - std::clog << "simplex_tree_60.num_simplices()=" << simplex_tree_60.num_simplices() << std::endl; - BOOST_CHECK(simplex_tree_60.num_simplices() == 25); - - max_alpha_square_value = 59.0; - std::clog << "========== OFF FILE NAME = " << off_file_name << " - alpha²=" << - max_alpha_square_value << "==========" << std::endl; - - Gudhi::Simplex_tree<> simplex_tree_59; - BOOST_CHECK(alpha_complex_from_file.create_complex(simplex_tree_59, max_alpha_square_value)); - - std::clog << "alpha_complex_from_file.num_vertices()=" << alpha_complex_from_file.num_vertices() << std::endl; - BOOST_CHECK(alpha_complex_from_file.num_vertices() == 7); - - std::clog << "simplex_tree_59.dimension()=" << simplex_tree_59.dimension() << std::endl; - BOOST_CHECK(simplex_tree_59.dimension() == 2); - - std::clog << "simplex_tree_59.num_vertices()=" << simplex_tree_59.num_vertices() << std::endl; - BOOST_CHECK(simplex_tree_59.num_vertices() == 7); - - std::clog << "simplex_tree_59.num_simplices()=" << simplex_tree_59.num_simplices() << std::endl; - BOOST_CHECK(simplex_tree_59.num_simplices() == 23); -} - // Use static dimension_tag for the user not to be able to set dimension typedef CGAL::Epeck_d< CGAL::Dimension_tag<4> > Kernel_4; typedef Kernel_4::Point_d Point_4; @@ -232,37 +172,6 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_from_points) { } -BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_from_empty_points, TestedKernel, list_of_kernel_variants) { - std::clog << "========== Alpha_complex_from_empty_points ==========" << std::endl; - - // ---------------------------------------------------------------------------- - // Init of an empty list of points - // ---------------------------------------------------------------------------- - std::vector points; - - // ---------------------------------------------------------------------------- - // Init of an alpha complex from the list of points - // ---------------------------------------------------------------------------- - Gudhi::alpha_complex::Alpha_complex alpha_complex_from_points(points); - - std::clog << "alpha_complex_from_points.num_vertices()=" << alpha_complex_from_points.num_vertices() << std::endl; - BOOST_CHECK(alpha_complex_from_points.num_vertices() == points.size()); - - // Test to the limit - BOOST_CHECK_THROW (alpha_complex_from_points.get_point(0), std::out_of_range); - - Gudhi::Simplex_tree<> simplex_tree; - BOOST_CHECK(!alpha_complex_from_points.create_complex(simplex_tree)); - - std::clog << "simplex_tree.num_simplices()=" << simplex_tree.num_simplices() << std::endl; - BOOST_CHECK(simplex_tree.num_simplices() == 0); - - std::clog << "simplex_tree.dimension()=" << simplex_tree.dimension() << std::endl; - BOOST_CHECK(simplex_tree.dimension() == -1); - - std::clog << "simplex_tree.num_vertices()=" << simplex_tree.num_vertices() << std::endl; - BOOST_CHECK(simplex_tree.num_vertices() == points.size()); -} using Inexact_kernel_2 = CGAL::Epick_d< CGAL::Dimension_tag<2> >; using Exact_kernel_2 = CGAL::Epeck_d< CGAL::Dimension_tag<2> >; diff --git a/src/Alpha_complex/test/CMakeLists.txt b/src/Alpha_complex/test/CMakeLists.txt index 0595ca92..dd2c235f 100644 --- a/src/Alpha_complex/test/CMakeLists.txt +++ b/src/Alpha_complex/test/CMakeLists.txt @@ -8,14 +8,18 @@ if (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 4.11.0) add_executable ( Alpha_complex_test_unit Alpha_complex_unit_test.cpp ) target_link_libraries(Alpha_complex_test_unit ${CGAL_LIBRARY}) + add_executable ( Alpha_complex_dim3_test_unit Alpha_complex_dim3_unit_test.cpp ) + target_link_libraries(Alpha_complex_dim3_test_unit ${CGAL_LIBRARY}) add_executable ( Delaunay_complex_test_unit Delaunay_complex_unit_test.cpp ) target_link_libraries(Delaunay_complex_test_unit ${CGAL_LIBRARY}) if (TBB_FOUND) target_link_libraries(Alpha_complex_test_unit ${TBB_LIBRARIES}) + target_link_libraries(Alpha_complex_dim3_test_unit ${TBB_LIBRARIES}) target_link_libraries(Delaunay_complex_test_unit ${TBB_LIBRARIES}) endif() gudhi_add_boost_test(Alpha_complex_test_unit) + gudhi_add_boost_test(Alpha_complex_dim3_test_unit) gudhi_add_boost_test(Delaunay_complex_test_unit) endif (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 4.11.0) @@ -73,4 +77,4 @@ if (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 5.1.0) endif() gudhi_add_boost_test(Zero_weighted_alpha_complex_test_unit) -endif (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 5.1.0) \ No newline at end of file +endif (NOT CGAL_WITH_EIGEN3_VERSION VERSION_LESS 5.1.0) -- cgit v1.2.3 From 80edbbaccc344abfa18ebc48a8493f747a775476 Mon Sep 17 00:00:00 2001 From: wreise Date: Sat, 15 Oct 2022 16:54:04 +0200 Subject: Simplify sorting in landscapes --- src/python/test/test_representations.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/python/test/test_representations.py b/src/python/test/test_representations.py index 948d7804..58caab21 100755 --- a/src/python/test/test_representations.py +++ b/src/python/test/test_representations.py @@ -235,10 +235,12 @@ def test_landscape_numeric(): dgm = np.array([[2., 6.], [3., 5.]]) lds_ref = np.array([ 0., 0.5, 1., 1.5, 2., 1.5, 1., 0.5, 0., # tent of [2, 6] - 0., 0., 0., 0.5, 1., 0.5, 0., 0., 0. + 0., 0., 0., 0.5, 1., 0.5, 0., 0., 0., + 0., 0., 0., 0., 0., 0., 0., 0., 0., + 0., 0., 0., 0., 0., 0., 0., 0., 0., ]) lds_ref *= np.sqrt(2) - lds = Landscape(num_landscapes=2, resolution=9, sample_range=[2., 6.]) + lds = Landscape(num_landscapes=4, resolution=9, sample_range=[2., 6.]) lds_dgm = lds(dgm) assert np.all(np.isclose(lds_dgm, lds_ref)) -- cgit v1.2.3 From 74617f0673aa13bce47833c366321a8838a7d123 Mon Sep 17 00:00:00 2001 From: wreise Date: Sat, 15 Oct 2022 17:00:25 +0200 Subject: Remove old lines --- src/python/gudhi/representations/vector_methods.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 3a91eccd..6267e077 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -179,10 +179,8 @@ class Landscape(BaseEstimator, TransformerMixin): midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) n_points = diag.shape[0] - # Get indices of largest elements: can't take more than n_points - 1 (the last ones are in the right position) - argpartition = np.argpartition(-tent_functions, min(self.num_landscapes, n_points-1), axis=1) - landscapes = np.take_along_axis(tent_functions, argpartition, axis=1) - landscapes = landscapes[:, :min(self.num_landscapes, n_points)].T + tent_functions.partition(n_points-self.num_landscapes, axis=1) + landscapes = np.sort(tent_functions[:, -self.num_landscapes:], axis=1)[:, ::-1].T # Complete the array with zeros to get the right number of landscapes if self.num_landscapes > n_points: -- cgit v1.2.3 From cd7dea8627f4b1c624e88d5ff28b32d1602f5e39 Mon Sep 17 00:00:00 2001 From: wreise Date: Sat, 15 Oct 2022 18:45:42 +0200 Subject: Treat the case when there are less points than landscape layers --- src/python/gudhi/representations/vector_methods.py | 13 +++++++------ 1 file changed, 7 insertions(+), 6 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/representations/vector_methods.py b/src/python/gudhi/representations/vector_methods.py index 6267e077..a169aee8 100644 --- a/src/python/gudhi/representations/vector_methods.py +++ b/src/python/gudhi/representations/vector_methods.py @@ -179,14 +179,15 @@ class Landscape(BaseEstimator, TransformerMixin): midpoints, heights = (diag[:, 0] + diag[:, 1]) / 2., (diag[:, 1] - diag[:, 0]) / 2. tent_functions = np.maximum(heights[None, :] - np.abs(x_values[:, None] - midpoints[None, :]), 0) n_points = diag.shape[0] - tent_functions.partition(n_points-self.num_landscapes, axis=1) - landscapes = np.sort(tent_functions[:, -self.num_landscapes:], axis=1)[:, ::-1].T - # Complete the array with zeros to get the right number of landscapes if self.num_landscapes > n_points: - landscapes = np.concatenate([ - landscapes, np.zeros((self.num_landscapes-n_points, *landscapes.shape[1:])) - ], axis=0) + tent_functions = np.concatenate( + [tent_functions, np.zeros((tent_functions.shape[0], self.num_landscapes-n_points))], + axis=1 + ) + tent_functions.partition(tent_functions.shape[1]-self.num_landscapes, axis=1) + landscapes = np.sort(tent_functions[:, -self.num_landscapes:], axis=1)[:, ::-1].T + landscapes = np.sqrt(2) * np.ravel(landscapes) Xfit.append(landscapes) -- cgit v1.2.3 From d18ef465b79fe53b103bb05d75d7f71b37792276 Mon Sep 17 00:00:00 2001 From: hschreiber Date: Tue, 18 Oct 2022 16:50:12 +0200 Subject: added possibility of negative values in persistance diagram --- .../include/gudhi/Persistence_graph.h | 57 +++++++--- .../test/bottleneck_unit_test.cpp | 120 +++++++++++++++++---- .../example/example_spatial_searching.cpp | 4 +- src/Spatial_searching/test/test_Kd_tree_search.cpp | 4 +- 4 files changed, 144 insertions(+), 41 deletions(-) (limited to 'src') diff --git a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h index 33f03b9c..9b663cc2 100644 --- a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h +++ b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h @@ -20,6 +20,7 @@ #include #include #include // for numeric_limits +#include namespace Gudhi { @@ -70,27 +71,51 @@ Persistence_graph::Persistence_graph(const Persistence_diagram1 &diag1, : u(), v(), b_alive(0.) { std::vector u_alive; std::vector v_alive; + std::vector u_nalive; + std::vector v_nalive; + int u_inf = 0; + int v_inf = 0; + double inf = std::numeric_limits::infinity(); + double neginf = -1 * inf; + for (auto it = std::begin(diag1); it != std::end(diag1); ++it) { - if (std::get<1>(*it) == std::numeric_limits::infinity()) - u_alive.push_back(std::get<0>(*it)); - else if (std::get<1>(*it) - std::get<0>(*it) > e) - u.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), u.size())); + if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ + if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) + u_inf++; + else if (std::get<0>(*it) == neginf) + u_nalive.push_back(std::get<1>(*it)); + else if (std::get<1>(*it) == inf) + u_alive.push_back(std::get<0>(*it)); + else if (std::get<1>(*it) - std::get<0>(*it) > e) + u.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), u.size())); + } } for (auto it = std::begin(diag2); it != std::end(diag2); ++it) { - if (std::get<1>(*it) == std::numeric_limits::infinity()) - v_alive.push_back(std::get<0>(*it)); - else if (std::get<1>(*it) - std::get<0>(*it) > e) - v.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), v.size())); + if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ + if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) + v_inf++; + else if (std::get<0>(*it) == neginf) + v_nalive.push_back(std::get<1>(*it)); + else if (std::get<1>(*it) == inf) + v_alive.push_back(std::get<0>(*it)); + else if (std::get<1>(*it) - std::get<0>(*it) > e) + v.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), v.size())); + } } if (u.size() < v.size()) swap(u, v); - std::sort(u_alive.begin(), u_alive.end()); - std::sort(v_alive.begin(), v_alive.end()); - if (u_alive.size() != v_alive.size()) { + + if (u_alive.size() != v_alive.size() || u_nalive.size() != v_nalive.size() || u_inf != v_inf) { b_alive = std::numeric_limits::infinity(); } else { + std::sort(u_alive.begin(), u_alive.end()); + std::sort(v_alive.begin(), v_alive.end()); + std::sort(u_nalive.begin(), u_nalive.end()); + std::sort(v_nalive.begin(), v_nalive.end()); for (auto it_u = u_alive.cbegin(), it_v = v_alive.cbegin(); it_u != u_alive.cend(); ++it_u, ++it_v) b_alive = (std::max)(b_alive, std::fabs(*it_u - *it_v)); + for (auto it_u = u_nalive.cbegin(), it_v = v_nalive.cbegin(); it_u != u_nalive.cend(); ++it_u, ++it_v) + b_alive = (std::max)(b_alive, std::fabs(*it_u - *it_v)); } } @@ -114,7 +139,7 @@ inline int Persistence_graph::corresponding_point_in_v(int u_point_index) const inline double Persistence_graph::distance(int u_point_index, int v_point_index) const { if (on_the_u_diagonal(u_point_index) && on_the_v_diagonal(v_point_index)) - return 0.; + return 0.; Internal_point p_u = get_u_point(u_point_index); Internal_point p_v = get_v_point(v_point_index); return (std::max)(std::fabs(p_u.x() - p_v.x()), std::fabs(p_u.y() - p_v.y())); @@ -132,9 +157,9 @@ inline std::vector Persistence_graph::sorted_distances() const { std::vector distances; distances.push_back(0.); // for empty diagrams for (int u_point_index = 0; u_point_index < size(); ++u_point_index) { - distances.push_back(distance(u_point_index, corresponding_point_in_v(u_point_index))); - for (int v_point_index = 0; v_point_index < size(); ++v_point_index) - distances.push_back(distance(u_point_index, v_point_index)); + distances.push_back(distance(u_point_index, corresponding_point_in_v(u_point_index))); + for (int v_point_index = 0; v_point_index < size(); ++v_point_index) + distances.push_back(distance(u_point_index, v_point_index)); } #ifdef GUDHI_USE_TBB tbb::parallel_sort(distances.begin(), distances.end()); @@ -146,7 +171,7 @@ inline std::vector Persistence_graph::sorted_distances() const { inline Internal_point Persistence_graph::get_u_point(int u_point_index) const { if (!on_the_u_diagonal(u_point_index)) - return u.at(u_point_index); + return u.at(u_point_index); Internal_point projector = v.at(corresponding_point_in_v(u_point_index)); double m = (projector.x() + projector.y()) / 2.; return Internal_point(m, m, u_point_index); diff --git a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp index 44141baa..79ee9c2c 100644 --- a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp +++ b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp @@ -31,14 +31,14 @@ std::vector< std::pair > v1, v2; BOOST_AUTO_TEST_CASE(persistence_graph) { // Random construction for (int i = 0; i < n1; i++) { - double a = unif(re); - double b = unif(re); - v1.emplace_back(std::min(a, b), std::max(a, b)); + double a = unif(re); + double b = unif(re); + v1.emplace_back(std::min(a, b), std::max(a, b)); } for (int i = 0; i < n2; i++) { - double a = unif(re); - double b = unif(re); - v2.emplace_back(std::min(a, b), std::max(a, b)); + double a = unif(re); + double b = unif(re); + v2.emplace_back(std::min(a, b), std::max(a, b)); } Persistence_graph g(v1, v2, 0.); std::vector d(g.sorted_distances()); @@ -84,14 +84,14 @@ BOOST_AUTO_TEST_CASE(neighbors_finder) { Persistence_graph g(v1, v2, 0.); Neighbors_finder nf(g, 1.); for (int v_point_index = 1; v_point_index < ((n2 + n1)*9 / 10); v_point_index += 2) - nf.add(v_point_index); + nf.add(v_point_index); // int v_point_index_1 = nf.pull_near(n2 / 2); BOOST_CHECK((v_point_index_1 == -1) || (g.distance(n2 / 2, v_point_index_1) <= 1.)); std::vector l = nf.pull_all_near(n2 / 2); bool v = true; for (auto it = l.cbegin(); it != l.cend(); ++it) - v = v && (g.distance(n2 / 2, *it) > 1.); + v = v && (g.distance(n2 / 2, *it) > 1.); BOOST_CHECK(v); int v_point_index_2 = nf.pull_near(n2 / 2); BOOST_CHECK(v_point_index_2 == -1); @@ -101,7 +101,7 @@ BOOST_AUTO_TEST_CASE(layered_neighbors_finder) { Persistence_graph g(v1, v2, 0.); Layered_neighbors_finder lnf(g, 1.); for (int v_point_index = 1; v_point_index < ((n2 + n1)*9 / 10); v_point_index += 2) - lnf.add(v_point_index, v_point_index % 7); + lnf.add(v_point_index, v_point_index % 7); // int v_point_index_1 = lnf.pull_near(n2 / 2, 6); BOOST_CHECK((v_point_index_1 == -1) || (g.distance(n2 / 2, v_point_index_1) <= 1.)); @@ -119,7 +119,7 @@ BOOST_AUTO_TEST_CASE(graph_matching) { m1.set_r(0.); int e = 0; while (m1.multi_augment()) - ++e; + ++e; BOOST_CHECK(e > 0); BOOST_CHECK(e <= 2 * sqrt(2 * (n1 + n2))); Graph_matching m2 = m1; @@ -127,7 +127,7 @@ BOOST_AUTO_TEST_CASE(graph_matching) { m2.set_r(upper_bound); e = 0; while (m2.multi_augment()) - ++e; + ++e; BOOST_CHECK(e <= 2 * sqrt(2 * (n1 + n2))); BOOST_CHECK(m2.perfect()); BOOST_CHECK(!m1.perfect()); @@ -139,16 +139,16 @@ BOOST_AUTO_TEST_CASE(global) { std::default_random_engine re; std::vector< std::pair > v1, v2; for (int i = 0; i < n1; i++) { - double a = unif1(re); - double b = unif1(re); - double x = unif2(re); - double y = unif2(re); - v1.emplace_back(std::min(a, b), std::max(a, b)); - v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); - if (i % 5 == 0) - v1.emplace_back(std::min(a, b), std::min(a, b) + x); - if (i % 3 == 0) - v2.emplace_back(std::max(a, b), std::max(a, b) + y); + double a = unif1(re); + double b = unif1(re); + double x = unif2(re); + double y = unif2(re); + v1.emplace_back(std::min(a, b), std::max(a, b)); + v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); + if (i % 5 == 0) + v1.emplace_back(std::min(a, b), std::min(a, b) + x); + if (i % 3 == 0) + v2.emplace_back(std::max(a, b), std::max(a, b) + y); } BOOST_CHECK(bottleneck_distance(v1, v2, 0.) <= upper_bound / 100.); BOOST_CHECK(bottleneck_distance(v1, v2, upper_bound / 10000.) <= upper_bound / 100. + upper_bound / 10000.); @@ -159,3 +159,81 @@ BOOST_AUTO_TEST_CASE(global) { BOOST_CHECK(bottleneck_distance(empty, empty) == 0); BOOST_CHECK(bottleneck_distance(empty, one) == 1); } + +BOOST_AUTO_TEST_CASE(neg_global) { + std::uniform_real_distribution unif1(0., upper_bound); + std::uniform_real_distribution unif2(upper_bound / 10000., upper_bound / 100.); + std::default_random_engine re; + std::vector< std::pair > v1, v2; + for (int i = 0; i < n1; i++) { + double a = std::log(unif1(re)); + double b = std::log(unif1(re)); + double x = std::log(unif2(re)); + double y = std::log(unif2(re)); + v1.emplace_back(std::min(a, b), std::max(a, b)); + v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); + if (i % 5 == 0) + v1.emplace_back(std::min(a, b), std::min(a, b) + x); + if (i % 3 == 0) + v2.emplace_back(std::max(a, b), std::max(a, b) + y); + } + BOOST_CHECK(bottleneck_distance(v1, v2, 0.) <= upper_bound / 100.); + BOOST_CHECK(bottleneck_distance(v1, v2, upper_bound / 10000.) <= upper_bound / 100. + upper_bound / 10000.); + BOOST_CHECK(std::abs(bottleneck_distance(v1, v2, 0.) - bottleneck_distance(v1, v2, upper_bound / 10000.)) <= upper_bound / 10000.); + + std::vector< std::pair > empty; + std::vector< std::pair > one = {{8, 10}}; + BOOST_CHECK(bottleneck_distance(empty, empty) == 0); + BOOST_CHECK(bottleneck_distance(empty, one) == 1); +} + +BOOST_AUTO_TEST_CASE(bottleneck_simple_test) { + std::vector< std::pair > v1, v2; + double inf = std::numeric_limits::infinity(); + double neginf = -1 * inf; + double b; + + v1.emplace_back(9.6, 14.); + v2.emplace_back(9.5, 14.1); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.09 && b < 0.11); + + v1.emplace_back(-34.974, -34.2); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.386 && b < 0.388); + + v1.emplace_back(neginf, 3.7); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK_EQUAL(b, inf); + + v2.emplace_back(neginf, 4.45); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.74 && b < 0.76); + + v1.emplace_back(-60.6, 52.1); + v2.emplace_back(-61.5, 53.); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.89 && b < 0.91); + + v1.emplace_back(3., inf); + v2.emplace_back(3.2, inf); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.89 && b < 0.91); + + v1.emplace_back(neginf, inf); + v2.emplace_back(neginf, inf); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK(b > 0.89 && b < 0.91); + + v2.emplace_back(6, inf); + + b = Gudhi::persistence_diagram::bottleneck_distance(v1, v2, 0.); + BOOST_CHECK_EQUAL(b, inf); +} diff --git a/src/Spatial_searching/example/example_spatial_searching.cpp b/src/Spatial_searching/example/example_spatial_searching.cpp index 8f9151fc..09c2dabf 100644 --- a/src/Spatial_searching/example/example_spatial_searching.cpp +++ b/src/Spatial_searching/example/example_spatial_searching.cpp @@ -25,7 +25,7 @@ int main(void) { // 10-nearest neighbor query std::clog << "10 nearest neighbors from points[20]:\n"; auto knn_range = points_ds.k_nearest_neighbors(points[20], 10, true); - for (auto const& nghb : knn_range) + for (auto const nghb : knn_range) std::clog << nghb.first << " (sq. dist. = " << nghb.second << ")\n"; // Incremental nearest neighbor query @@ -38,7 +38,7 @@ int main(void) { // 10-furthest neighbor query std::clog << "10 furthest neighbors from points[20]:\n"; auto kfn_range = points_ds.k_furthest_neighbors(points[20], 10, true); - for (auto const& nghb : kfn_range) + for (auto const nghb : kfn_range) std::clog << nghb.first << " (sq. dist. = " << nghb.second << ")\n"; // Incremental furthest neighbor query diff --git a/src/Spatial_searching/test/test_Kd_tree_search.cpp b/src/Spatial_searching/test/test_Kd_tree_search.cpp index d6c6fba3..e9acfaa7 100644 --- a/src/Spatial_searching/test/test_Kd_tree_search.cpp +++ b/src/Spatial_searching/test/test_Kd_tree_search.cpp @@ -45,7 +45,7 @@ BOOST_AUTO_TEST_CASE(test_Kd_tree_search) { std::vector knn_result; FT last_dist = -1.; - for (auto const& nghb : kns_range) { + for (auto const nghb : kns_range) { BOOST_CHECK(nghb.second > last_dist); knn_result.push_back(nghb.second); last_dist = nghb.second; @@ -76,7 +76,7 @@ BOOST_AUTO_TEST_CASE(test_Kd_tree_search) { std::vector kfn_result; last_dist = kfn_range.begin()->second; - for (auto const& nghb : kfn_range) { + for (auto const nghb : kfn_range) { BOOST_CHECK(nghb.second <= last_dist); kfn_result.push_back(nghb.second); last_dist = nghb.second; -- cgit v1.2.3 From ff59ebb1a3417701a9282783adeb254af14b856c Mon Sep 17 00:00:00 2001 From: hschreiber Date: Tue, 18 Oct 2022 17:31:34 +0200 Subject: syntax correction --- src/Bottleneck_distance/include/gudhi/Persistence_graph.h | 2 +- src/Bottleneck_distance/test/bottleneck_unit_test.cpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h index 9b663cc2..de989e1b 100644 --- a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h +++ b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h @@ -76,7 +76,7 @@ Persistence_graph::Persistence_graph(const Persistence_diagram1 &diag1, int u_inf = 0; int v_inf = 0; double inf = std::numeric_limits::infinity(); - double neginf = -1 * inf; + double neginf = -inf; for (auto it = std::begin(diag1); it != std::end(diag1); ++it) { if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ diff --git a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp index 79ee9c2c..38ed89a8 100644 --- a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp +++ b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp @@ -190,7 +190,7 @@ BOOST_AUTO_TEST_CASE(neg_global) { BOOST_AUTO_TEST_CASE(bottleneck_simple_test) { std::vector< std::pair > v1, v2; double inf = std::numeric_limits::infinity(); - double neginf = -1 * inf; + double neginf = -inf; double b; v1.emplace_back(9.6, 14.); -- cgit v1.2.3 From 0a183262e27e133bc8a77c91f5c97506497cd380 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Wed, 19 Oct 2022 08:28:02 +0200 Subject: persistence diagram default output is std::cout, not std::clog. Modify utilities help to be in agreement with this --- src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp | 2 +- src/Alpha_complex/utilities/alpha_complex_persistence.cpp | 2 +- src/Cech_complex/utilities/cech_persistence.cpp | 2 +- .../utilities/distance_matrix_edge_collapse_rips_persistence.cpp | 2 +- src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp | 2 +- src/Persistent_cohomology/example/persistence_from_file.cpp | 2 +- src/Persistent_cohomology/example/rips_multifield_persistence.cpp | 2 +- src/Persistent_cohomology/example/rips_persistence_step_by_step.cpp | 2 +- .../example/rips_persistence_via_boundary_matrix.cpp | 2 +- src/Rips_complex/utilities/rips_correlation_matrix_persistence.cpp | 2 +- src/Rips_complex/utilities/rips_distance_matrix_persistence.cpp | 2 +- src/Rips_complex/utilities/rips_persistence.cpp | 2 +- src/Rips_complex/utilities/sparse_rips_persistence.cpp | 2 +- src/Witness_complex/utilities/strong_witness_persistence.cpp | 2 +- src/Witness_complex/utilities/weak_witness_persistence.cpp | 2 +- src/Witness_complex/utilities/witnesscomplex.md | 4 ++-- 16 files changed, 17 insertions(+), 17 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp index 91899040..e65d8c6f 100644 --- a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp @@ -263,7 +263,7 @@ void program_options(int argc, char *argv[], std::string &off_file_points, bool "cuboid-file,c", po::value(&cuboid_file), "Name of file describing the periodic domain. Format is:\n min_hx min_hy min_hz\n max_hx max_hy max_hz")( "output-file,o", po::value(&output_file_diag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-alpha-square-value,r", po::value(&alpha_square_max_value) ->default_value(std::numeric_limits::infinity()), diff --git a/src/Alpha_complex/utilities/alpha_complex_persistence.cpp b/src/Alpha_complex/utilities/alpha_complex_persistence.cpp index e86b34e2..29edbd8e 100644 --- a/src/Alpha_complex/utilities/alpha_complex_persistence.cpp +++ b/src/Alpha_complex/utilities/alpha_complex_persistence.cpp @@ -163,7 +163,7 @@ void program_options(int argc, char *argv[], std::string &off_file_points, bool "weight-file,w", po::value(&weight_file)->default_value(std::string()), "Name of file containing a point weights. Format is one weight per line:\n W1\n ...\n Wn ")( "output-file,o", po::value(&output_file_diag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-alpha-square-value,r", po::value(&alpha_square_max_value) ->default_value(std::numeric_limits::infinity()), "Maximal alpha square value for the Alpha complex construction.")( diff --git a/src/Cech_complex/utilities/cech_persistence.cpp b/src/Cech_complex/utilities/cech_persistence.cpp index a07ba212..e6419f3d 100644 --- a/src/Cech_complex/utilities/cech_persistence.cpp +++ b/src/Cech_complex/utilities/cech_persistence.cpp @@ -118,7 +118,7 @@ void program_options(int argc, char* argv[], std::string& off_file_points, bool& "fast,f", po::bool_switch(&fast), "To activate fast version of Cech complex (default is false, not available if exact is set)")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-radius,r", po::value(&max_radius)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Cech complex construction.")( diff --git a/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp b/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp index 38efb9e6..70b489b5 100644 --- a/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp +++ b/src/Collapse/utilities/distance_matrix_edge_collapse_rips_persistence.cpp @@ -111,7 +111,7 @@ void program_options(int argc, char* argv[], std::string& csv_matrix_file, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::cout")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.")( diff --git a/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp b/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp index d8f42ab6..a8fd6f14 100644 --- a/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp +++ b/src/Collapse/utilities/point_cloud_edge_collapse_rips_persistence.cpp @@ -140,7 +140,7 @@ void program_options(int argc, char* argv[], std::string& off_file_points, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::cout")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.")( diff --git a/src/Persistent_cohomology/example/persistence_from_file.cpp b/src/Persistent_cohomology/example/persistence_from_file.cpp index 38c44514..7f89c001 100644 --- a/src/Persistent_cohomology/example/persistence_from_file.cpp +++ b/src/Persistent_cohomology/example/persistence_from_file.cpp @@ -93,7 +93,7 @@ void program_options(int argc, char * argv[] visible.add_options() ("help,h", "produce help message") ("output-file,o", po::value(&output_file)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog") + "Name of file in which the persistence diagram is written. Default print in standard output") ("field-charac,p", po::value(&p)->default_value(11), "Characteristic p of the coefficient field Z/pZ for computing homology.") ("min-persistence,m", po::value(&min_persistence), diff --git a/src/Persistent_cohomology/example/rips_multifield_persistence.cpp b/src/Persistent_cohomology/example/rips_multifield_persistence.cpp index ca26a5b9..84453898 100644 --- a/src/Persistent_cohomology/example/rips_multifield_persistence.cpp +++ b/src/Persistent_cohomology/example/rips_multifield_persistence.cpp @@ -96,7 +96,7 @@ void program_options(int argc, char * argv[] visible.add_options() ("help,h", "produce help message") ("output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog") + "Name of file in which the persistence diagram is written. Default print in standard output") ("max-edge-length,r", po::value(&threshold)->default_value(0), "Maximal length of an edge for the Rips complex construction.") ("cpx-dimension,d", po::value(&dim_max)->default_value(1), diff --git a/src/Persistent_cohomology/example/rips_persistence_step_by_step.cpp b/src/Persistent_cohomology/example/rips_persistence_step_by_step.cpp index a503d983..6f37cf5c 100644 --- a/src/Persistent_cohomology/example/rips_persistence_step_by_step.cpp +++ b/src/Persistent_cohomology/example/rips_persistence_step_by_step.cpp @@ -112,7 +112,7 @@ void program_options(int argc, char * argv[] visible.add_options() ("help,h", "produce help message") ("output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog") + "Name of file in which the persistence diagram is written. Default print in standard output") ("max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.") diff --git a/src/Persistent_cohomology/example/rips_persistence_via_boundary_matrix.cpp b/src/Persistent_cohomology/example/rips_persistence_via_boundary_matrix.cpp index 8c5742aa..6b60f603 100644 --- a/src/Persistent_cohomology/example/rips_persistence_via_boundary_matrix.cpp +++ b/src/Persistent_cohomology/example/rips_persistence_via_boundary_matrix.cpp @@ -109,7 +109,7 @@ void program_options(int argc, char * argv[] visible.add_options() ("help,h", "produce help message") ("output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog") + "Name of file in which the persistence diagram is written. Default print in standard output") ("max-edge-length,r", po::value(&threshold)->default_value(0), "Maximal length of an edge for the Rips complex construction.") ("cpx-dimension,d", po::value(&dim_max)->default_value(1), diff --git a/src/Rips_complex/utilities/rips_correlation_matrix_persistence.cpp b/src/Rips_complex/utilities/rips_correlation_matrix_persistence.cpp index b473738e..72ddc797 100644 --- a/src/Rips_complex/utilities/rips_correlation_matrix_persistence.cpp +++ b/src/Rips_complex/utilities/rips_correlation_matrix_persistence.cpp @@ -118,7 +118,7 @@ void program_options(int argc, char* argv[], std::string& csv_matrix_file, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "min-edge-corelation,c", po::value(&correlation_min)->default_value(0), "Minimal corelation of an edge for the Rips complex construction.")( "cpx-dimension,d", po::value(&dim_max)->default_value(1), diff --git a/src/Rips_complex/utilities/rips_distance_matrix_persistence.cpp b/src/Rips_complex/utilities/rips_distance_matrix_persistence.cpp index 6306755d..77ad841a 100644 --- a/src/Rips_complex/utilities/rips_distance_matrix_persistence.cpp +++ b/src/Rips_complex/utilities/rips_distance_matrix_persistence.cpp @@ -79,7 +79,7 @@ void program_options(int argc, char* argv[], std::string& csv_matrix_file, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.")( diff --git a/src/Rips_complex/utilities/rips_persistence.cpp b/src/Rips_complex/utilities/rips_persistence.cpp index 9d7490b3..43194821 100644 --- a/src/Rips_complex/utilities/rips_persistence.cpp +++ b/src/Rips_complex/utilities/rips_persistence.cpp @@ -81,7 +81,7 @@ void program_options(int argc, char* argv[], std::string& off_file_points, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.")( diff --git a/src/Rips_complex/utilities/sparse_rips_persistence.cpp b/src/Rips_complex/utilities/sparse_rips_persistence.cpp index ac935b41..829c85e6 100644 --- a/src/Rips_complex/utilities/sparse_rips_persistence.cpp +++ b/src/Rips_complex/utilities/sparse_rips_persistence.cpp @@ -84,7 +84,7 @@ void program_options(int argc, char* argv[], std::string& off_file_points, std:: po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-edge-length,r", po::value(&threshold)->default_value(std::numeric_limits::infinity()), "Maximal length of an edge for the Rips complex construction.")( diff --git a/src/Witness_complex/utilities/strong_witness_persistence.cpp b/src/Witness_complex/utilities/strong_witness_persistence.cpp index 614de0d4..b2ecad82 100644 --- a/src/Witness_complex/utilities/strong_witness_persistence.cpp +++ b/src/Witness_complex/utilities/strong_witness_persistence.cpp @@ -108,7 +108,7 @@ void program_options(int argc, char* argv[], int& nbL, std::string& file_name, s visible.add_options()("help,h", "produce help message")("landmarks,l", po::value(&nbL), "Number of landmarks to choose from the point cloud.")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-sq-alpha,a", po::value(&max_squared_alpha)->default_value(default_alpha), "Maximal squared relaxation parameter.")( "field-charac,p", po::value(&p)->default_value(11), diff --git a/src/Witness_complex/utilities/weak_witness_persistence.cpp b/src/Witness_complex/utilities/weak_witness_persistence.cpp index 5ea31d6b..c7ead7de 100644 --- a/src/Witness_complex/utilities/weak_witness_persistence.cpp +++ b/src/Witness_complex/utilities/weak_witness_persistence.cpp @@ -108,7 +108,7 @@ void program_options(int argc, char* argv[], int& nbL, std::string& file_name, s visible.add_options()("help,h", "produce help message")("landmarks,l", po::value(&nbL), "Number of landmarks to choose from the point cloud.")( "output-file,o", po::value(&filediag)->default_value(std::string()), - "Name of file in which the persistence diagram is written. Default print in std::clog")( + "Name of file in which the persistence diagram is written. Default print in standard output")( "max-sq-alpha,a", po::value(&max_squared_alpha)->default_value(default_alpha), "Maximal squared relaxation parameter.")( "field-charac,p", po::value(&p)->default_value(11), diff --git a/src/Witness_complex/utilities/witnesscomplex.md b/src/Witness_complex/utilities/witnesscomplex.md index 3a3a7d83..e994e0b8 100644 --- a/src/Witness_complex/utilities/witnesscomplex.md +++ b/src/Witness_complex/utilities/witnesscomplex.md @@ -29,7 +29,7 @@ and `p` is the characteristic of the field *Z/pZ* used for homology coefficients * `-h [ --help ]` Produce help message * `-l [ --landmarks ]` Number of landmarks to choose from the point cloud. -* `-o [ --output-file ]` Name of file in which the persistence diagram is written. By default, print in std::clog. +* `-o [ --output-file ]` Name of file in which the persistence diagram is written. By default, print in standard output. * `-a [ --max-sq-alpha ]` (default = inf) Maximal squared relaxation parameter. * `-p [ --field-charac ]` (default = 11) Characteristic p of the coefficient field Z/pZ for computing homology. * `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals. @@ -60,7 +60,7 @@ and `p` is the characteristic of the field *Z/pZ* used for homology coefficients * `-h [ --help ]` Produce help message * `-l [ --landmarks ]` Number of landmarks to choose from the point cloud. -* `-o [ --output-file ]` Name of file in which the persistence diagram is written. By default, print in std::clog. +* `-o [ --output-file ]` Name of file in which the persistence diagram is written. By default, print in standard output. * `-a [ --max-sq-alpha ]` (default = inf) Maximal squared relaxation parameter. * `-p [ --field-charac ]` (default = 11) Characteristic p of the coefficient field Z/pZ for computing homology. * `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals. -- cgit v1.2.3 From 8227688cfb716bb881b836d245a869664bb10d36 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 19 Oct 2022 08:55:27 +0200 Subject: Update src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp Co-authored-by: Vincent Rouvreau <10407034+VincentRouvreau@users.noreply.github.com> --- src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp index 0085ae67..e7c261f1 100644 --- a/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp +++ b/src/Alpha_complex/test/Alpha_complex_dim3_unit_test.cpp @@ -9,7 +9,7 @@ */ #define BOOST_TEST_DYN_LINK -#define BOOST_TEST_MODULE "alpha_complex" +#define BOOST_TEST_MODULE "alpha_complex_dim3" #include #include -- cgit v1.2.3 From 475a20c0645b079aaf204fd47398b7bb278fb490 Mon Sep 17 00:00:00 2001 From: hschreiber Date: Wed, 19 Oct 2022 11:52:41 +0200 Subject: correction of indentation --- .../include/gudhi/Persistence_graph.h | 72 +++++++++++----------- .../test/bottleneck_unit_test.cpp | 62 +++++++++---------- 2 files changed, 67 insertions(+), 67 deletions(-) (limited to 'src') diff --git a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h index de989e1b..c1e10f8e 100644 --- a/src/Bottleneck_distance/include/gudhi/Persistence_graph.h +++ b/src/Bottleneck_distance/include/gudhi/Persistence_graph.h @@ -32,7 +32,7 @@ namespace persistence_diagram { * \ingroup bottleneck_distance */ class Persistence_graph { - public: +public: /** \internal \brief Constructor taking 2 PersistenceDiagrams (concept) as parameters. */ template Persistence_graph(const Persistence_diagram1& diag1, const Persistence_diagram2& diag2, double e); @@ -59,7 +59,7 @@ class Persistence_graph { /** \internal \brief Returns the corresponding internal point */ Internal_point get_v_point(int v_point_index) const; - private: +private: std::vector u; std::vector v; double b_alive; @@ -68,7 +68,7 @@ class Persistence_graph { template Persistence_graph::Persistence_graph(const Persistence_diagram1 &diag1, const Persistence_diagram2 &diag2, double e) - : u(), v(), b_alive(0.) { + : u(), v(), b_alive(0.) { std::vector u_alive; std::vector v_alive; std::vector u_nalive; @@ -79,28 +79,28 @@ Persistence_graph::Persistence_graph(const Persistence_diagram1 &diag1, double neginf = -inf; for (auto it = std::begin(diag1); it != std::end(diag1); ++it) { - if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ - if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) - u_inf++; - else if (std::get<0>(*it) == neginf) - u_nalive.push_back(std::get<1>(*it)); - else if (std::get<1>(*it) == inf) - u_alive.push_back(std::get<0>(*it)); - else if (std::get<1>(*it) - std::get<0>(*it) > e) - u.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), u.size())); - } + if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ + if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) + u_inf++; + else if (std::get<0>(*it) == neginf) + u_nalive.push_back(std::get<1>(*it)); + else if (std::get<1>(*it) == inf) + u_alive.push_back(std::get<0>(*it)); + else if (std::get<1>(*it) - std::get<0>(*it) > e) + u.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), u.size())); + } } for (auto it = std::begin(diag2); it != std::end(diag2); ++it) { - if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ - if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) - v_inf++; - else if (std::get<0>(*it) == neginf) - v_nalive.push_back(std::get<1>(*it)); - else if (std::get<1>(*it) == inf) - v_alive.push_back(std::get<0>(*it)); - else if (std::get<1>(*it) - std::get<0>(*it) > e) - v.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), v.size())); - } + if (std::get<0>(*it) != inf && std::get<1>(*it) != neginf){ + if (std::get<0>(*it) == neginf && std::get<1>(*it) == inf) + v_inf++; + else if (std::get<0>(*it) == neginf) + v_nalive.push_back(std::get<1>(*it)); + else if (std::get<1>(*it) == inf) + v_alive.push_back(std::get<0>(*it)); + else if (std::get<1>(*it) - std::get<0>(*it) > e) + v.push_back(Internal_point(std::get<0>(*it), std::get<1>(*it), v.size())); + } } if (u.size() < v.size()) swap(u, v); @@ -108,14 +108,14 @@ Persistence_graph::Persistence_graph(const Persistence_diagram1 &diag1, if (u_alive.size() != v_alive.size() || u_nalive.size() != v_nalive.size() || u_inf != v_inf) { b_alive = std::numeric_limits::infinity(); } else { - std::sort(u_alive.begin(), u_alive.end()); - std::sort(v_alive.begin(), v_alive.end()); - std::sort(u_nalive.begin(), u_nalive.end()); - std::sort(v_nalive.begin(), v_nalive.end()); + std::sort(u_alive.begin(), u_alive.end()); + std::sort(v_alive.begin(), v_alive.end()); + std::sort(u_nalive.begin(), u_nalive.end()); + std::sort(v_nalive.begin(), v_nalive.end()); for (auto it_u = u_alive.cbegin(), it_v = v_alive.cbegin(); it_u != u_alive.cend(); ++it_u, ++it_v) b_alive = (std::max)(b_alive, std::fabs(*it_u - *it_v)); - for (auto it_u = u_nalive.cbegin(), it_v = v_nalive.cbegin(); it_u != u_nalive.cend(); ++it_u, ++it_v) - b_alive = (std::max)(b_alive, std::fabs(*it_u - *it_v)); + for (auto it_u = u_nalive.cbegin(), it_v = v_nalive.cbegin(); it_u != u_nalive.cend(); ++it_u, ++it_v) + b_alive = (std::max)(b_alive, std::fabs(*it_u - *it_v)); } } @@ -129,17 +129,17 @@ inline bool Persistence_graph::on_the_v_diagonal(int v_point_index) const { inline int Persistence_graph::corresponding_point_in_u(int v_point_index) const { return on_the_v_diagonal(v_point_index) ? - v_point_index - static_cast (v.size()) : v_point_index + static_cast (u.size()); + v_point_index - static_cast (v.size()) : v_point_index + static_cast (u.size()); } inline int Persistence_graph::corresponding_point_in_v(int u_point_index) const { return on_the_u_diagonal(u_point_index) ? - u_point_index - static_cast (u.size()) : u_point_index + static_cast (v.size()); + u_point_index - static_cast (u.size()) : u_point_index + static_cast (v.size()); } inline double Persistence_graph::distance(int u_point_index, int v_point_index) const { if (on_the_u_diagonal(u_point_index) && on_the_v_diagonal(v_point_index)) - return 0.; + return 0.; Internal_point p_u = get_u_point(u_point_index); Internal_point p_v = get_v_point(v_point_index); return (std::max)(std::fabs(p_u.x() - p_v.x()), std::fabs(p_u.y() - p_v.y())); @@ -157,9 +157,9 @@ inline std::vector Persistence_graph::sorted_distances() const { std::vector distances; distances.push_back(0.); // for empty diagrams for (int u_point_index = 0; u_point_index < size(); ++u_point_index) { - distances.push_back(distance(u_point_index, corresponding_point_in_v(u_point_index))); - for (int v_point_index = 0; v_point_index < size(); ++v_point_index) - distances.push_back(distance(u_point_index, v_point_index)); + distances.push_back(distance(u_point_index, corresponding_point_in_v(u_point_index))); + for (int v_point_index = 0; v_point_index < size(); ++v_point_index) + distances.push_back(distance(u_point_index, v_point_index)); } #ifdef GUDHI_USE_TBB tbb::parallel_sort(distances.begin(), distances.end()); @@ -171,7 +171,7 @@ inline std::vector Persistence_graph::sorted_distances() const { inline Internal_point Persistence_graph::get_u_point(int u_point_index) const { if (!on_the_u_diagonal(u_point_index)) - return u.at(u_point_index); + return u.at(u_point_index); Internal_point projector = v.at(corresponding_point_in_v(u_point_index)); double m = (projector.x() + projector.y()) / 2.; return Internal_point(m, m, u_point_index); diff --git a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp index 38ed89a8..9872f20c 100644 --- a/src/Bottleneck_distance/test/bottleneck_unit_test.cpp +++ b/src/Bottleneck_distance/test/bottleneck_unit_test.cpp @@ -31,14 +31,14 @@ std::vector< std::pair > v1, v2; BOOST_AUTO_TEST_CASE(persistence_graph) { // Random construction for (int i = 0; i < n1; i++) { - double a = unif(re); - double b = unif(re); - v1.emplace_back(std::min(a, b), std::max(a, b)); + double a = unif(re); + double b = unif(re); + v1.emplace_back(std::min(a, b), std::max(a, b)); } for (int i = 0; i < n2; i++) { - double a = unif(re); - double b = unif(re); - v2.emplace_back(std::min(a, b), std::max(a, b)); + double a = unif(re); + double b = unif(re); + v2.emplace_back(std::min(a, b), std::max(a, b)); } Persistence_graph g(v1, v2, 0.); std::vector d(g.sorted_distances()); @@ -84,14 +84,14 @@ BOOST_AUTO_TEST_CASE(neighbors_finder) { Persistence_graph g(v1, v2, 0.); Neighbors_finder nf(g, 1.); for (int v_point_index = 1; v_point_index < ((n2 + n1)*9 / 10); v_point_index += 2) - nf.add(v_point_index); + nf.add(v_point_index); // int v_point_index_1 = nf.pull_near(n2 / 2); BOOST_CHECK((v_point_index_1 == -1) || (g.distance(n2 / 2, v_point_index_1) <= 1.)); std::vector l = nf.pull_all_near(n2 / 2); bool v = true; for (auto it = l.cbegin(); it != l.cend(); ++it) - v = v && (g.distance(n2 / 2, *it) > 1.); + v = v && (g.distance(n2 / 2, *it) > 1.); BOOST_CHECK(v); int v_point_index_2 = nf.pull_near(n2 / 2); BOOST_CHECK(v_point_index_2 == -1); @@ -101,7 +101,7 @@ BOOST_AUTO_TEST_CASE(layered_neighbors_finder) { Persistence_graph g(v1, v2, 0.); Layered_neighbors_finder lnf(g, 1.); for (int v_point_index = 1; v_point_index < ((n2 + n1)*9 / 10); v_point_index += 2) - lnf.add(v_point_index, v_point_index % 7); + lnf.add(v_point_index, v_point_index % 7); // int v_point_index_1 = lnf.pull_near(n2 / 2, 6); BOOST_CHECK((v_point_index_1 == -1) || (g.distance(n2 / 2, v_point_index_1) <= 1.)); @@ -119,7 +119,7 @@ BOOST_AUTO_TEST_CASE(graph_matching) { m1.set_r(0.); int e = 0; while (m1.multi_augment()) - ++e; + ++e; BOOST_CHECK(e > 0); BOOST_CHECK(e <= 2 * sqrt(2 * (n1 + n2))); Graph_matching m2 = m1; @@ -127,7 +127,7 @@ BOOST_AUTO_TEST_CASE(graph_matching) { m2.set_r(upper_bound); e = 0; while (m2.multi_augment()) - ++e; + ++e; BOOST_CHECK(e <= 2 * sqrt(2 * (n1 + n2))); BOOST_CHECK(m2.perfect()); BOOST_CHECK(!m1.perfect()); @@ -139,16 +139,16 @@ BOOST_AUTO_TEST_CASE(global) { std::default_random_engine re; std::vector< std::pair > v1, v2; for (int i = 0; i < n1; i++) { - double a = unif1(re); - double b = unif1(re); - double x = unif2(re); - double y = unif2(re); - v1.emplace_back(std::min(a, b), std::max(a, b)); - v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); - if (i % 5 == 0) - v1.emplace_back(std::min(a, b), std::min(a, b) + x); - if (i % 3 == 0) - v2.emplace_back(std::max(a, b), std::max(a, b) + y); + double a = unif1(re); + double b = unif1(re); + double x = unif2(re); + double y = unif2(re); + v1.emplace_back(std::min(a, b), std::max(a, b)); + v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); + if (i % 5 == 0) + v1.emplace_back(std::min(a, b), std::min(a, b) + x); + if (i % 3 == 0) + v2.emplace_back(std::max(a, b), std::max(a, b) + y); } BOOST_CHECK(bottleneck_distance(v1, v2, 0.) <= upper_bound / 100.); BOOST_CHECK(bottleneck_distance(v1, v2, upper_bound / 10000.) <= upper_bound / 100. + upper_bound / 10000.); @@ -166,16 +166,16 @@ BOOST_AUTO_TEST_CASE(neg_global) { std::default_random_engine re; std::vector< std::pair > v1, v2; for (int i = 0; i < n1; i++) { - double a = std::log(unif1(re)); - double b = std::log(unif1(re)); - double x = std::log(unif2(re)); - double y = std::log(unif2(re)); - v1.emplace_back(std::min(a, b), std::max(a, b)); - v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); - if (i % 5 == 0) - v1.emplace_back(std::min(a, b), std::min(a, b) + x); - if (i % 3 == 0) - v2.emplace_back(std::max(a, b), std::max(a, b) + y); + double a = std::log(unif1(re)); + double b = std::log(unif1(re)); + double x = std::log(unif2(re)); + double y = std::log(unif2(re)); + v1.emplace_back(std::min(a, b), std::max(a, b)); + v2.emplace_back(std::min(a, b) + std::min(x, y), std::max(a, b) + std::max(x, y)); + if (i % 5 == 0) + v1.emplace_back(std::min(a, b), std::min(a, b) + x); + if (i % 3 == 0) + v2.emplace_back(std::max(a, b), std::max(a, b) + y); } BOOST_CHECK(bottleneck_distance(v1, v2, 0.) <= upper_bound / 100.); BOOST_CHECK(bottleneck_distance(v1, v2, upper_bound / 10000.) <= upper_bound / 100. + upper_bound / 10000.); -- cgit v1.2.3 From 0a1d09caba097117e08c16d12bf682d611454c02 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Fri, 21 Oct 2022 16:54:09 +0200 Subject: Stores vertices on init_from_range and insert all of them in the Simplex_tree to avoid quadratic behavior --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 23 ++++++++++++++++++++--- 1 file changed, 20 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index aec8c1b1..d5d1a0f0 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -17,8 +17,7 @@ // to construct Alpha_complex from a OFF file of points #include -#include -#include // isnan, fmax +#include // isnan, fmax #include // for std::unique_ptr #include // for std::size_t @@ -45,6 +44,8 @@ #include // std::pair #include #include // for std::iota +#include // for std::sort +#include // for complex.insert_simplex_and_subfaces({...}, ...) // Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10 #if CGAL_VERSION_NR < 1041101000 @@ -146,6 +147,10 @@ class Alpha_complex { std::unique_ptr triangulation_; /** \brief Kernel for triangulation_ functions access.*/ A_kernel_d kernel_; + /** \brief Vertices to be inserted first by the create_complex method to avoid quadratic complexity. + * It isn't just [0, n) if some points have multiplicity (only one copy appears in the complex). + */ + std::vector vertices_; /** \brief Cache for geometric constructions: circumcenter and squared radius of a simplex.*/ std::vector cache_, old_cache_; @@ -279,6 +284,9 @@ class Alpha_complex { // structure to retrieve CGAL points from vertex handle - one vertex handle per point. // Needs to be constructed before as vertex handles arrives in no particular order. vertex_handle_to_iterator_.resize(point_cloud.size()); + // List of sorted unique vertices in the triangulation. We take advantage of the existing loop to construct it + // Vertices list avoids quadratic complexity with the Simplex_tree. We should not fill it up with Toplex_map e.g. + vertices_.reserve(triangulation_->number_of_vertices()); // Loop on triangulation vertices list for (CGAL_vertex_iterator vit = triangulation_->vertices_begin(); vit != triangulation_->vertices_end(); ++vit) { if (!triangulation_->is_infinite(*vit)) { @@ -286,8 +294,10 @@ class Alpha_complex { std::clog << "Vertex insertion - " << vit->data() << " -> " << vit->point() << std::endl; #endif // DEBUG_TRACES vertex_handle_to_iterator_[vit->data()] = vit; + vertices_.push_back(vit->data()); } } + std::sort(vertices_.begin(), vertices_.end()); // -------------------------------------------------------------------------------------------- } } @@ -384,12 +394,19 @@ class Alpha_complex { // -------------------------------------------------------------------------------------------- // Simplex_tree construction from loop on triangulation finite full cells list if (num_vertices() > 0) { + for (auto vertex : vertices_) { +#ifdef DEBUG_TRACES + std::clog << "SimplicialComplex insertion " << vertex << std::endl; +#endif // DEBUG_TRACES + complex.insert_simplex_and_subfaces({static_cast(vertex)}, std::numeric_limits::quiet_NaN()); + } + for (auto cit = triangulation_->finite_full_cells_begin(); cit != triangulation_->finite_full_cells_end(); ++cit) { Vector_vertex vertexVector; #ifdef DEBUG_TRACES - std::clog << "Simplex_tree insertion "; + std::clog << "SimplicialComplex insertion "; #endif // DEBUG_TRACES for (auto vit = cit->vertices_begin(); vit != cit->vertices_end(); ++vit) { if (*vit != nullptr) { -- cgit v1.2.3 From f4e90e943dd368b76f466c91075193124c191ef5 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Fri, 21 Oct 2022 17:16:01 +0200 Subject: Remove strange behavior with initializer list --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index d5d1a0f0..7b837ae0 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -45,7 +45,6 @@ #include #include // for std::iota #include // for std::sort -#include // for complex.insert_simplex_and_subfaces({...}, ...) // Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10 #if CGAL_VERSION_NR < 1041101000 @@ -394,11 +393,13 @@ class Alpha_complex { // -------------------------------------------------------------------------------------------- // Simplex_tree construction from loop on triangulation finite full cells list if (num_vertices() > 0) { + std::vector one_vertex(1); for (auto vertex : vertices_) { #ifdef DEBUG_TRACES std::clog << "SimplicialComplex insertion " << vertex << std::endl; #endif // DEBUG_TRACES - complex.insert_simplex_and_subfaces({static_cast(vertex)}, std::numeric_limits::quiet_NaN()); + one_vertex[0] = vertex; + complex.insert_simplex_and_subfaces(one_vertex, std::numeric_limits::quiet_NaN()); } for (auto cit = triangulation_->finite_full_cells_begin(); -- cgit v1.2.3 From 797191125ef4ef87d8b8e7689507ded1c5e20159 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 26 Oct 2022 18:03:56 +0200 Subject: Use std::optional to check C++17 support --- src/Tangential_complex/include/gudhi/Tangential_complex.h | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) (limited to 'src') diff --git a/src/Tangential_complex/include/gudhi/Tangential_complex.h b/src/Tangential_complex/include/gudhi/Tangential_complex.h index cc424810..56a24af0 100644 --- a/src/Tangential_complex/include/gudhi/Tangential_complex.h +++ b/src/Tangential_complex/include/gudhi/Tangential_complex.h @@ -36,7 +36,6 @@ #include #include // for EIGEN_VERSION_AT_LEAST -#include #include #include #include @@ -56,6 +55,7 @@ #include // for std::sqrt #include #include // for std::size_t +#include #ifdef GUDHI_USE_TBB #include @@ -994,7 +994,7 @@ class Tangential_complex { // circumspheres of the star of "center_vertex" // If th the m_max_squared_edge_length is set the maximal radius of the "star sphere" // is at most square root of m_max_squared_edge_length - boost::optional squared_star_sphere_radius_plus_margin = m_max_squared_edge_length; + std::optional squared_star_sphere_radius_plus_margin = m_max_squared_edge_length; // Insert points until we find a point which is outside "star sphere" for (auto nn_it = ins_range.begin(); nn_it != ins_range.end(); ++nn_it) { @@ -1036,7 +1036,7 @@ class Tangential_complex { // Let's recompute squared_star_sphere_radius_plus_margin if (triangulation.current_dimension() >= tangent_space_dim) { - squared_star_sphere_radius_plus_margin = boost::none; + squared_star_sphere_radius_plus_margin = std::nullopt; // Get the incident cells and look for the biggest circumsphere std::vector incident_cells; triangulation.incident_full_cells(center_vertex, std::back_inserter(incident_cells)); @@ -1044,7 +1044,7 @@ class Tangential_complex { cit != incident_cells.end(); ++cit) { Tr_full_cell_handle cell = *cit; if (triangulation.is_infinite(cell)) { - squared_star_sphere_radius_plus_margin = boost::none; + squared_star_sphere_radius_plus_margin = std::nullopt; break; } else { // Note that this uses the perturbed point since it uses @@ -2030,7 +2030,7 @@ class Tangential_complex { // and their center vertex Stars_container m_stars; std::vector m_squared_star_spheres_radii_incl_margin; - boost::optional m_max_squared_edge_length; + std::optional m_max_squared_edge_length; #ifdef GUDHI_TC_USE_ANOTHER_POINT_SET_FOR_TANGENT_SPACE_ESTIM Points m_points_for_tse; -- cgit v1.2.3 From ff0037d54f231d1cc52b40bba12d69c0a7788c51 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 26 Oct 2022 18:18:36 +0200 Subject: Update -std=c++14 for python --- src/cmake/modules/GUDHI_compilation_flags.cmake | 1 + src/python/CMakeLists.txt | 2 +- 2 files changed, 2 insertions(+), 1 deletion(-) (limited to 'src') diff --git a/src/cmake/modules/GUDHI_compilation_flags.cmake b/src/cmake/modules/GUDHI_compilation_flags.cmake index e2d3d872..b43ccf73 100644 --- a/src/cmake/modules/GUDHI_compilation_flags.cmake +++ b/src/cmake/modules/GUDHI_compilation_flags.cmake @@ -12,6 +12,7 @@ macro(add_cxx_compiler_flag _flag) endmacro() set (CMAKE_CXX_STANDARD 17) +# This number needs to be changed in python/CMakeLists.txt at the same time enable_testing() diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 5f323935..21602a13 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -131,7 +131,7 @@ if(PYTHONINTERP_FOUND) if(MSVC) set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'/fp:strict', ") else(MSVC) - set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-std=c++14', ") + set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-std=c++17', ") endif(MSVC) if(CMAKE_COMPILER_IS_GNUCXX) set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-frounding-math', ") -- cgit v1.2.3 From 20f25ed9db1fd1a9d781f3a4a0a976439a11d011 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 26 Oct 2022 18:31:19 +0200 Subject: For C++17, at least MacOS 10.14 --- src/python/CMakeLists.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 21602a13..704379e9 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -246,8 +246,8 @@ if(PYTHONINTERP_FOUND) # Specific for Mac if (${CMAKE_SYSTEM_NAME} MATCHES "Darwin") - set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-mmacosx-version-min=10.12', ") - set(GUDHI_PYTHON_EXTRA_LINK_ARGS "${GUDHI_PYTHON_EXTRA_LINK_ARGS}'-mmacosx-version-min=10.12', ") + set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-mmacosx-version-min=10.14', ") + set(GUDHI_PYTHON_EXTRA_LINK_ARGS "${GUDHI_PYTHON_EXTRA_LINK_ARGS}'-mmacosx-version-min=10.14', ") endif(${CMAKE_SYSTEM_NAME} MATCHES "Darwin") # Loop on INCLUDE_DIRECTORIES PROPERTY -- cgit v1.2.3 From b1639bdd54e5e34a8f801c147fcd7f64a829ec4c Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 26 Oct 2022 18:51:44 +0200 Subject: Visual Studio requires /std:c++17 --- src/python/CMakeLists.txt | 1 + 1 file changed, 1 insertion(+) (limited to 'src') diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 704379e9..8f8df138 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -129,6 +129,7 @@ if(PYTHONINTERP_FOUND) # Gudhi and CGAL compilation option if(MSVC) + set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'/std:c++17', ") set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'/fp:strict', ") else(MSVC) set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-std=c++17', ") -- cgit v1.2.3 From 26d7bcc518f3bdc9b0d8f854f2879ed9c219e440 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Thu, 3 Nov 2022 14:56:34 +0100 Subject: Translate n_jobs to workers for SciPy --- src/python/gudhi/point_cloud/knn.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) (limited to 'src') diff --git a/src/python/gudhi/point_cloud/knn.py b/src/python/gudhi/point_cloud/knn.py index de5844f9..7dc83817 100644 --- a/src/python/gudhi/point_cloud/knn.py +++ b/src/python/gudhi/point_cloud/knn.py @@ -314,7 +314,9 @@ class KNearestNeighbors: return None if self.params["implementation"] == "ckdtree": - qargs = {key: val for key, val in self.params.items() if key in {"p", "eps", "n_jobs"}} + qargs = {key: val for key, val in self.params.items() if key in {"p", "eps"}} + # SciPy renamed n_jobs to workers + qargs["workers"] = self.params.get("workers") or self.params.get("n_jobs") or 1 distances, neighbors = self.kdtree.query(X, k=self.k, **qargs) if k == 1: # SciPy decided to squeeze the last dimension for k=1 -- cgit v1.2.3 From 11dec62e555e446fb34a1f5ba3d0c9aba1a7e956 Mon Sep 17 00:00:00 2001 From: Vincent Rouvreau Date: Thu, 3 Nov 2022 17:49:33 +0100 Subject: code review: Use std::ptrdiff_t as internal vertex handle type when not yet known --- src/Alpha_complex/include/gudhi/Alpha_complex.h | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) (limited to 'src') diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex.h b/src/Alpha_complex/include/gudhi/Alpha_complex.h index 7b837ae0..a7372f19 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex.h @@ -101,13 +101,17 @@ template struct Is_Epeck_D> { static const bool val */ template, bool Weighted = false> class Alpha_complex { + private: + // Vertex_handle internal type (required by triangulation_ and vertices_). + using Internal_vertex_handle = std::ptrdiff_t; + public: /** \brief Geometric traits class that provides the geometric types and predicates needed by the triangulations.*/ using Geom_traits = std::conditional_t, Kernel>; // Add an int in TDS to save point index in the structure using TDS = CGAL::Triangulation_data_structure, + CGAL::Triangulation_vertex, CGAL::Triangulation_full_cell >; /** \brief A (Weighted or not) Delaunay triangulation of a set of points in \f$ \mathbb{R}^D\f$.*/ @@ -132,9 +136,6 @@ class Alpha_complex { // Vertex_iterator type from CGAL. using CGAL_vertex_iterator = typename Triangulation::Vertex_iterator; - // size_type type from CGAL. - using size_type = typename Triangulation::size_type; - // Structure to switch from simplex tree vertex handle to CGAL vertex iterator. using Vector_vertex_iterator = std::vector< CGAL_vertex_iterator >; @@ -149,7 +150,7 @@ class Alpha_complex { /** \brief Vertices to be inserted first by the create_complex method to avoid quadratic complexity. * It isn't just [0, n) if some points have multiplicity (only one copy appears in the complex). */ - std::vector vertices_; + std::vector vertices_; /** \brief Cache for geometric constructions: circumcenter and squared radius of a simplex.*/ std::vector cache_, old_cache_; @@ -261,11 +262,11 @@ class Alpha_complex { std::vector point_cloud(first, last); // Creates a vector {0, 1, ..., N-1} - std::vector indices(boost::counting_iterator(0), - boost::counting_iterator(point_cloud.size())); + std::vector indices(boost::counting_iterator(0), + boost::counting_iterator(point_cloud.size())); using Point_property_map = boost::iterator_property_map::iterator, - CGAL::Identity_property_map>; + CGAL::Identity_property_map>; using Search_traits_d = CGAL::Spatial_sort_traits_adapter_d; CGAL::spatial_sort(indices.begin(), indices.end(), Search_traits_d(std::begin(point_cloud))); -- cgit v1.2.3 From 95330fb658db4184813d8720c52da5e25bcdd38c Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 4 Nov 2022 11:01:43 +0100 Subject: Mention SciPy version in the doc --- src/python/doc/installation.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/python/doc/installation.rst b/src/python/doc/installation.rst index b704f778..276ac4e2 100644 --- a/src/python/doc/installation.rst +++ b/src/python/doc/installation.rst @@ -391,7 +391,7 @@ The :doc:`persistence graphical tools ` and mathematics, science, and engineering. :class:`~gudhi.point_cloud.knn.KNearestNeighbors` can use the Python package -`SciPy `_ as a backend if explicitly requested. +`SciPy `_ :math:`\geq` 1.6.0 as a backend if explicitly requested. TensorFlow ---------- -- cgit v1.2.3 From 454ef5370d5a188cd303c497ed473df9aac81708 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 4 Nov 2022 23:02:56 +0100 Subject: Use fetch_spiral_2d in the doc --- src/python/doc/clustering.rst | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/src/python/doc/clustering.rst b/src/python/doc/clustering.rst index c5a57d3c..62422682 100644 --- a/src/python/doc/clustering.rst +++ b/src/python/doc/clustering.rst @@ -17,9 +17,8 @@ As a by-product, we produce the persistence diagram of the merge tree of the ini :include-source: import gudhi - f = open(gudhi.__root_source_dir__ + '/data/points/spiral_2d.csv', 'r') - import numpy as np - data = np.loadtxt(f) + from gudhi.datasets.remote import fetch_spiral_2d + data = fetch_spiral_2d() import matplotlib.pyplot as plt plt.scatter(data[:,0],data[:,1],marker='.',s=1) plt.show() -- cgit v1.2.3 From 4246261962f5b3f37cf4e845ff192da2be023834 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 4 Nov 2022 23:39:46 +0100 Subject: Remove the deprecated std::iterator --- .../include/gudhi/Bitmap_cubical_complex.h | 16 ++++++++++++++-- .../include/gudhi/Bitmap_cubical_complex_base.h | 16 ++++++++++++++-- 2 files changed, 28 insertions(+), 4 deletions(-) (limited to 'src') diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h index 4a6af3a4..ae7e6ed3 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h @@ -241,10 +241,16 @@ class Bitmap_cubical_complex : public T { **/ class Filtration_simplex_range; - class Filtration_simplex_iterator : std::iterator { + class Filtration_simplex_iterator { // Iterator over all simplices of the complex in the order of the indexing scheme. // 'value_type' must be 'Simplex_handle'. public: + typedef std::input_iterator_tag iterator_category; + typedef Simplex_handle value_type; + typedef std::ptrdiff_t difference_type; + typedef value_type* pointer; + typedef value_type& reference; + Filtration_simplex_iterator(Bitmap_cubical_complex* b) : b(b), position(0) {} Filtration_simplex_iterator() : b(NULL), position(0) {} @@ -386,10 +392,16 @@ class Bitmap_cubical_complex : public T { **/ class Skeleton_simplex_range; - class Skeleton_simplex_iterator : std::iterator { + class Skeleton_simplex_iterator { // Iterator over all simplices of the complex in the order of the indexing scheme. // 'value_type' must be 'Simplex_handle'. public: + typedef std::input_iterator_tag iterator_category; + typedef Simplex_handle value_type; + typedef std::ptrdiff_t difference_type; + typedef value_type* pointer; + typedef value_type& reference; + Skeleton_simplex_iterator(Bitmap_cubical_complex* b, std::size_t d) : b(b), dimension(d) { if (globalDbg) { std::clog << "Skeleton_simplex_iterator ( Bitmap_cubical_complex* b , std::size_t d )\n"; diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h index bafe7981..253e4a54 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h @@ -251,8 +251,14 @@ class Bitmap_cubical_complex_base { * @brief Iterator through all cells in the complex (in order they appear in the structure -- i.e. * in lexicographical order). **/ - class All_cells_iterator : std::iterator { + class All_cells_iterator { public: + typedef std::input_iterator_tag iterator_category; + typedef T value_type; + typedef std::ptrdiff_t difference_type; + typedef value_type* pointer; + typedef value_type& reference; + All_cells_iterator() { this->counter = 0; } All_cells_iterator operator++() { @@ -355,8 +361,14 @@ class Bitmap_cubical_complex_base { * @brief Iterator through top dimensional cells of the complex. The cells appear in order they are stored * in the structure (i.e. in lexicographical order) **/ - class Top_dimensional_cells_iterator : std::iterator { + class Top_dimensional_cells_iterator { public: + typedef std::input_iterator_tag iterator_category; + typedef T value_type; + typedef std::ptrdiff_t difference_type; + typedef value_type* pointer; + typedef value_type& reference; + Top_dimensional_cells_iterator(Bitmap_cubical_complex_base& b) : b(b) { this->counter = std::vector(b.dimension()); // std::fill( this->counter.begin() , this->counter.end() , 0 ); -- cgit v1.2.3 From f5af745c4f671f477365115bedf17d167b8d07a5 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 4 Nov 2022 23:46:45 +0100 Subject: Fix bogus traits. Not great, but less broken. --- src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h | 4 ++-- .../include/gudhi/Bitmap_cubical_complex_base.h | 8 ++++---- 2 files changed, 6 insertions(+), 6 deletions(-) (limited to 'src') diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h index ae7e6ed3..29fabc6c 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex.h @@ -249,7 +249,7 @@ class Bitmap_cubical_complex : public T { typedef Simplex_handle value_type; typedef std::ptrdiff_t difference_type; typedef value_type* pointer; - typedef value_type& reference; + typedef value_type reference; Filtration_simplex_iterator(Bitmap_cubical_complex* b) : b(b), position(0) {} @@ -400,7 +400,7 @@ class Bitmap_cubical_complex : public T { typedef Simplex_handle value_type; typedef std::ptrdiff_t difference_type; typedef value_type* pointer; - typedef value_type& reference; + typedef value_type reference; Skeleton_simplex_iterator(Bitmap_cubical_complex* b, std::size_t d) : b(b), dimension(d) { if (globalDbg) { diff --git a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h index 253e4a54..2bf62f9b 100644 --- a/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h +++ b/src/Bitmap_cubical_complex/include/gudhi/Bitmap_cubical_complex_base.h @@ -254,10 +254,10 @@ class Bitmap_cubical_complex_base { class All_cells_iterator { public: typedef std::input_iterator_tag iterator_category; - typedef T value_type; + typedef std::size_t value_type; typedef std::ptrdiff_t difference_type; typedef value_type* pointer; - typedef value_type& reference; + typedef value_type reference; All_cells_iterator() { this->counter = 0; } @@ -364,10 +364,10 @@ class Bitmap_cubical_complex_base { class Top_dimensional_cells_iterator { public: typedef std::input_iterator_tag iterator_category; - typedef T value_type; + typedef std::size_t value_type; typedef std::ptrdiff_t difference_type; typedef value_type* pointer; - typedef value_type& reference; + typedef value_type reference; Top_dimensional_cells_iterator(Bitmap_cubical_complex_base& b) : b(b) { this->counter = std::vector(b.dimension()); -- cgit v1.2.3 From 8804aa1580c500ed927d65c25e8b78700725338e Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Tue, 8 Nov 2022 22:26:36 +0100 Subject: Clarify doc of RipsComplex.__init__ --- src/python/gudhi/rips_complex.pyx | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/rips_complex.pyx b/src/python/gudhi/rips_complex.pyx index c3470292..a0924cd6 100644 --- a/src/python/gudhi/rips_complex.pyx +++ b/src/python/gudhi/rips_complex.pyx @@ -45,20 +45,19 @@ cdef class RipsComplex: max_edge_length=float('inf'), sparse=None): """RipsComplex constructor. - :param max_edge_length: Rips value. - :type max_edge_length: float - :param points: A list of points in d-Dimension. - :type points: list of list of float + :type points: List[List[float]] Or :param distance_matrix: A distance matrix (full square or lower triangular). - :type points: list of list of float + :type distance_matrix: List[List[float]] And in both cases + :param max_edge_length: Rips value. + :type max_edge_length: float :param sparse: If this is not None, it switches to building a sparse Rips and represents the approximation parameter epsilon. :type sparse: float -- cgit v1.2.3 From 7c17408897a95a1f74626e8ff0ec8101ac4f92fd Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Tue, 8 Nov 2022 22:36:16 +0100 Subject: Reject positional arguments in RipsComplex.__init__ --- .github/next_release.md | 3 +++ src/python/gudhi/rips_complex.pyx | 4 ++-- src/python/test/test_simplex_generators.py | 2 +- 3 files changed, 6 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/.github/next_release.md b/.github/next_release.md index 81599b2c..d5fcef1c 100644 --- a/.github/next_release.md +++ b/.github/next_release.md @@ -9,6 +9,9 @@ Below is a list of changes made since GUDHI 3.6.0: - [Module](link) - ... +- [Rips complex](https://gudhi.inria.fr/python/latest/rips_complex_user.html) + - Construction now rejects positional arguments, you need to specify `points=X`. + - Installation - c++17 is the new minimal standard to compile the library. This implies Visual Studio minimal version is now 2017. diff --git a/src/python/gudhi/rips_complex.pyx b/src/python/gudhi/rips_complex.pyx index a0924cd6..d748f91e 100644 --- a/src/python/gudhi/rips_complex.pyx +++ b/src/python/gudhi/rips_complex.pyx @@ -41,7 +41,7 @@ cdef class RipsComplex: cdef Rips_complex_interface thisref # Fake constructor that does nothing but documenting the constructor - def __init__(self, points=None, distance_matrix=None, + def __init__(self, *, points=None, distance_matrix=None, max_edge_length=float('inf'), sparse=None): """RipsComplex constructor. @@ -64,7 +64,7 @@ cdef class RipsComplex: """ # The real cython constructor - def __cinit__(self, points=None, distance_matrix=None, + def __cinit__(self, *, points=None, distance_matrix=None, max_edge_length=float('inf'), sparse=None): if sparse is not None: if distance_matrix is not None: diff --git a/src/python/test/test_simplex_generators.py b/src/python/test/test_simplex_generators.py index 8a9b4844..c567d4c1 100755 --- a/src/python/test/test_simplex_generators.py +++ b/src/python/test/test_simplex_generators.py @@ -14,7 +14,7 @@ import numpy as np def test_flag_generators(): pts = np.array([[0, 0], [0, 1.01], [1, 0], [1.02, 1.03], [100, 0], [100, 3.01], [103, 0], [103.02, 3.03]]) - r = gudhi.RipsComplex(pts, max_edge_length=4) + r = gudhi.RipsComplex(points=pts, max_edge_length=4) st = r.create_simplex_tree(max_dimension=50) st.persistence() g = st.flag_persistence_generators() -- cgit v1.2.3 From ad1123c3c7cfddc1c15e9933b96af08ef3398b3c Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Wed, 9 Nov 2022 17:46:39 +0100 Subject: New write_points_to_off_file --- src/python/CMakeLists.txt | 4 +- src/python/doc/installation.rst | 4 +- src/python/doc/point_cloud.rst | 5 ++ src/python/gudhi/off_reader.pyx | 41 -------------- src/python/gudhi/off_utils.pyx | 57 ++++++++++++++++++++ src/python/test/test_subsampling.py | 103 ++++++++---------------------------- 6 files changed, 88 insertions(+), 126 deletions(-) delete mode 100644 src/python/gudhi/off_reader.pyx create mode 100644 src/python/gudhi/off_utils.pyx (limited to 'src') diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 8f8df138..35ddb778 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -53,7 +53,7 @@ if(PYTHONINTERP_FOUND) set(GUDHI_PYTHON_MODULES_EXTRA "${GUDHI_PYTHON_MODULES_EXTRA}'datasets', ") # Cython modules - set(GUDHI_PYTHON_MODULES "${GUDHI_PYTHON_MODULES}'off_reader', ") + set(GUDHI_PYTHON_MODULES "${GUDHI_PYTHON_MODULES}'off_utils', ") set(GUDHI_PYTHON_MODULES "${GUDHI_PYTHON_MODULES}'simplex_tree', ") set(GUDHI_PYTHON_MODULES "${GUDHI_PYTHON_MODULES}'rips_complex', ") set(GUDHI_PYTHON_MODULES "${GUDHI_PYTHON_MODULES}'cubical_complex', ") @@ -152,7 +152,7 @@ if(PYTHONINTERP_FOUND) set(GUDHI_PYTHON_EXTRA_COMPILE_ARGS "${GUDHI_PYTHON_EXTRA_COMPILE_ARGS}'-DCGAL_EIGEN3_ENABLED', ") endif (EIGEN3_FOUND) - set(GUDHI_CYTHON_MODULES "${GUDHI_CYTHON_MODULES}'off_reader', ") + set(GUDHI_CYTHON_MODULES "${GUDHI_CYTHON_MODULES}'off_utils', ") set(GUDHI_CYTHON_MODULES "${GUDHI_CYTHON_MODULES}'simplex_tree', ") set(GUDHI_CYTHON_MODULES "${GUDHI_CYTHON_MODULES}'rips_complex', ") set(GUDHI_CYTHON_MODULES "${GUDHI_CYTHON_MODULES}'cubical_complex', ") diff --git a/src/python/doc/installation.rst b/src/python/doc/installation.rst index 276ac4e2..5491542f 100644 --- a/src/python/doc/installation.rst +++ b/src/python/doc/installation.rst @@ -150,7 +150,7 @@ You shall have something like: Cython version 0.29.25 Numpy version 1.21.4 Boost version 1.77.0 - + Installed modules are: off_reader;simplex_tree;rips_complex;cubical_complex;periodic_cubical_complex; + + Installed modules are: off_utils;simplex_tree;rips_complex;cubical_complex;periodic_cubical_complex; persistence_graphical_tools;reader_utils;witness_complex;strong_witness_complex; + Missing modules are: bottleneck;nerve_gic;subsampling;tangential_complex;alpha_complex;euclidean_witness_complex; euclidean_strong_witness_complex; @@ -188,7 +188,7 @@ A complete configuration would be : GMPXX_LIBRARIES = /usr/lib/x86_64-linux-gnu/libgmpxx.so MPFR_LIBRARIES = /usr/lib/x86_64-linux-gnu/libmpfr.so TBB version 9107 found and used - + Installed modules are: bottleneck;off_reader;simplex_tree;rips_complex;cubical_complex;periodic_cubical_complex; + + Installed modules are: bottleneck;off_utils;simplex_tree;rips_complex;cubical_complex;periodic_cubical_complex; persistence_graphical_tools;reader_utils;witness_complex;strong_witness_complex;nerve_gic;subsampling; tangential_complex;alpha_complex;euclidean_witness_complex;euclidean_strong_witness_complex; + Missing modules are: diff --git a/src/python/doc/point_cloud.rst b/src/python/doc/point_cloud.rst index ffd8f85b..473b303f 100644 --- a/src/python/doc/point_cloud.rst +++ b/src/python/doc/point_cloud.rst @@ -13,6 +13,11 @@ File Readers .. autofunction:: gudhi.read_lower_triangular_matrix_from_csv_file +File Writers +------------ + +.. autofunction:: gudhi.write_points_to_off_file + Subsampling ----------- diff --git a/src/python/gudhi/off_reader.pyx b/src/python/gudhi/off_reader.pyx deleted file mode 100644 index a3200704..00000000 --- a/src/python/gudhi/off_reader.pyx +++ /dev/null @@ -1,41 +0,0 @@ -# This file is part of the Gudhi Library - https://gudhi.inria.fr/ - -# which is released under MIT. -# See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full -# license details. -# Author(s): Vincent Rouvreau -# -# Copyright (C) 2016 Inria -# -# Modification(s): -# - YYYY/MM Author: Description of the modification - -from __future__ import print_function -from cython cimport numeric -from libcpp.vector cimport vector -from libcpp.string cimport string -import errno -import os - -__author__ = "Vincent Rouvreau" -__copyright__ = "Copyright (C) 2016 Inria" -__license__ = "MIT" - -cdef extern from "Off_reader_interface.h" namespace "Gudhi": - vector[vector[double]] read_points_from_OFF_file(string off_file) - -def read_points_from_off_file(off_file=''): - """Read points from OFF file. - - :param off_file: An OFF file style name. - :type off_file: string - - :returns: The point set. - :rtype: List[List[float]] - """ - if off_file: - if os.path.isfile(off_file): - return read_points_from_OFF_file(off_file.encode('utf-8')) - else: - raise FileNotFoundError(errno.ENOENT, os.strerror(errno.ENOENT), - off_file) - diff --git a/src/python/gudhi/off_utils.pyx b/src/python/gudhi/off_utils.pyx new file mode 100644 index 00000000..155575d5 --- /dev/null +++ b/src/python/gudhi/off_utils.pyx @@ -0,0 +1,57 @@ +# This file is part of the Gudhi Library - https://gudhi.inria.fr/ - +# which is released under MIT. +# See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full +# license details. +# Author(s): Vincent Rouvreau +# +# Copyright (C) 2016 Inria +# +# Modification(s): +# - YYYY/MM Author: Description of the modification + +from __future__ import print_function +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.string cimport string +cimport cython +import errno +import os +import numpy as np + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Off_reader_interface.h" namespace "Gudhi": + vector[vector[double]] read_points_from_OFF_file(string off_file) + +def read_points_from_off_file(off_file=''): + """Read points from OFF file. + + :param off_file: An OFF file style name. + :type off_file: string + + :returns: The point set. + :rtype: List[List[float]] + """ + if off_file: + if os.path.isfile(off_file): + return read_points_from_OFF_file(off_file.encode('utf-8')) + else: + raise FileNotFoundError(errno.ENOENT, os.strerror(errno.ENOENT), + off_file) + +@cython.embedsignature(True) +def write_points_to_off_file(fname, points): + """Write points to an OFF file. + + A simple wrapper for `numpy.savetxt`. + + :param fname: Name of the OFF file. + :type fname: str or file handle + :param points: Point coordinates. + :type points: numpy array of shape (n, dim) + """ + points = np.array(points, copy=False) + assert len(points.shape) == 2 + np.savetxt(fname, points, header='nOFF\n{} {} 0 0'.format(points.shape[1], points.shape[0]), comments='') diff --git a/src/python/test/test_subsampling.py b/src/python/test/test_subsampling.py index 3431f372..c1cb4e3f 100755 --- a/src/python/test/test_subsampling.py +++ b/src/python/test/test_subsampling.py @@ -16,17 +16,9 @@ __license__ = "MIT" def test_write_off_file_for_tests(): - file = open("subsample.off", "w") - file.write("nOFF\n") - file.write("2 7 0 0\n") - file.write("1.0 1.0\n") - file.write("7.0 0.0\n") - file.write("4.0 6.0\n") - file.write("9.0 6.0\n") - file.write("0.0 14.0\n") - file.write("2.0 19.0\n") - file.write("9.0 17.0\n") - file.close() + gudhi.write_points_to_off_file( + "subsample.off", [[1.0, 1.0], [7.0, 0.0], [4.0, 6.0], [9.0, 6.0], [0.0, 14.0], [2.0, 19.0], [9.0, 17.0]] + ) def test_simple_choose_n_farthest_points_with_a_starting_point(): @@ -34,54 +26,29 @@ def test_simple_choose_n_farthest_points_with_a_starting_point(): i = 0 for point in point_set: # The iteration starts with the given starting point - sub_set = gudhi.choose_n_farthest_points( - points=point_set, nb_points=1, starting_point=i - ) + sub_set = gudhi.choose_n_farthest_points(points=point_set, nb_points=1, starting_point=i) assert sub_set[0] == point_set[i] i = i + 1 # The iteration finds then the farthest - sub_set = gudhi.choose_n_farthest_points( - points=point_set, nb_points=2, starting_point=1 - ) + sub_set = gudhi.choose_n_farthest_points(points=point_set, nb_points=2, starting_point=1) assert sub_set[1] == point_set[3] - sub_set = gudhi.choose_n_farthest_points( - points=point_set, nb_points=2, starting_point=3 - ) + sub_set = gudhi.choose_n_farthest_points(points=point_set, nb_points=2, starting_point=3) assert sub_set[1] == point_set[1] - sub_set = gudhi.choose_n_farthest_points( - points=point_set, nb_points=2, starting_point=0 - ) + sub_set = gudhi.choose_n_farthest_points(points=point_set, nb_points=2, starting_point=0) assert sub_set[1] == point_set[2] - sub_set = gudhi.choose_n_farthest_points( - points=point_set, nb_points=2, starting_point=2 - ) + sub_set = gudhi.choose_n_farthest_points(points=point_set, nb_points=2, starting_point=2) assert sub_set[1] == point_set[0] # Test the limits - assert ( - gudhi.choose_n_farthest_points(points=[], nb_points=0, starting_point=0) == [] - ) - assert ( - gudhi.choose_n_farthest_points(points=[], nb_points=1, starting_point=0) == [] - ) - assert ( - gudhi.choose_n_farthest_points(points=[], nb_points=0, starting_point=1) == [] - ) - assert ( - gudhi.choose_n_farthest_points(points=[], nb_points=1, starting_point=1) == [] - ) + assert gudhi.choose_n_farthest_points(points=[], nb_points=0, starting_point=0) == [] + assert gudhi.choose_n_farthest_points(points=[], nb_points=1, starting_point=0) == [] + assert gudhi.choose_n_farthest_points(points=[], nb_points=0, starting_point=1) == [] + assert gudhi.choose_n_farthest_points(points=[], nb_points=1, starting_point=1) == [] # From off file test for i in range(0, 7): - assert ( - len( - gudhi.choose_n_farthest_points( - off_file="subsample.off", nb_points=i, starting_point=i - ) - ) - == i - ) + assert len(gudhi.choose_n_farthest_points(off_file="subsample.off", nb_points=i, starting_point=i)) == i def test_simple_choose_n_farthest_points_randomed(): @@ -104,10 +71,7 @@ def test_simple_choose_n_farthest_points_randomed(): # From off file test for i in range(0, 7): - assert ( - len(gudhi.choose_n_farthest_points(off_file="subsample.off", nb_points=i)) - == i - ) + assert len(gudhi.choose_n_farthest_points(off_file="subsample.off", nb_points=i)) == i def test_simple_pick_n_random_points(): @@ -130,9 +94,7 @@ def test_simple_pick_n_random_points(): # From off file test for i in range(0, 7): - assert ( - len(gudhi.pick_n_random_points(off_file="subsample.off", nb_points=i)) == i - ) + assert len(gudhi.pick_n_random_points(off_file="subsample.off", nb_points=i)) == i def test_simple_sparsify_points(): @@ -152,31 +114,10 @@ def test_simple_sparsify_points(): ] assert gudhi.sparsify_point_set(points=point_set, min_squared_dist=2.001) == [[0, 1]] - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=0.0)) - == 7 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=30.0)) - == 5 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=40.1)) - == 4 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=89.9)) - == 3 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=100.0)) - == 2 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=324.9)) - == 2 - ) - assert ( - len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=325.01)) - == 1 - ) + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=0.0)) == 7 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=30.0)) == 5 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=40.1)) == 4 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=89.9)) == 3 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=100.0)) == 2 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=324.9)) == 2 + assert len(gudhi.sparsify_point_set(off_file="subsample.off", min_squared_dist=325.01)) == 1 -- cgit v1.2.3 From 4118bcb622c624130e768d9116a7e147a5e45c68 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 11 Nov 2022 18:00:42 +0100 Subject: Link to OFF format --- src/python/gudhi/off_utils.pyx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/python/gudhi/off_utils.pyx b/src/python/gudhi/off_utils.pyx index 155575d5..a8142791 100644 --- a/src/python/gudhi/off_utils.pyx +++ b/src/python/gudhi/off_utils.pyx @@ -26,7 +26,7 @@ cdef extern from "Off_reader_interface.h" namespace "Gudhi": vector[vector[double]] read_points_from_OFF_file(string off_file) def read_points_from_off_file(off_file=''): - """Read points from OFF file. + """Read points from an `OFF file `_. :param off_file: An OFF file style name. :type off_file: string @@ -43,7 +43,7 @@ def read_points_from_off_file(off_file=''): @cython.embedsignature(True) def write_points_to_off_file(fname, points): - """Write points to an OFF file. + """Write points to an `OFF file `_. A simple wrapper for `numpy.savetxt`. -- cgit v1.2.3 From 2d5039b7eeb16116ab859076aa0a93f092250d88 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 11 Nov 2022 18:59:42 +0100 Subject: Special case for writing 3d OFF --- src/python/CMakeLists.txt | 1 + src/python/gudhi/off_utils.pyx | 7 ++++++- src/python/test/test_off.py | 21 +++++++++++++++++++++ 3 files changed, 28 insertions(+), 1 deletion(-) create mode 100644 src/python/test/test_off.py (limited to 'src') diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 35ddb778..32ec13bd 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -546,6 +546,7 @@ if(PYTHONINTERP_FOUND) # Reader utils add_gudhi_py_test(test_reader_utils) + add_gudhi_py_test(test_off) # Wasserstein if(OT_FOUND) diff --git a/src/python/gudhi/off_utils.pyx b/src/python/gudhi/off_utils.pyx index a8142791..9276c7b0 100644 --- a/src/python/gudhi/off_utils.pyx +++ b/src/python/gudhi/off_utils.pyx @@ -54,4 +54,9 @@ def write_points_to_off_file(fname, points): """ points = np.array(points, copy=False) assert len(points.shape) == 2 - np.savetxt(fname, points, header='nOFF\n{} {} 0 0'.format(points.shape[1], points.shape[0]), comments='') + dim = points.shape[1] + if dim == 3: + head = 'OFF\n{} 0 0'.format(points.shape[0]) + else: + head = 'nOFF\n{} {} 0 0'.format(dim, points.shape[0]) + np.savetxt(fname, points, header=head, comments='') diff --git a/src/python/test/test_off.py b/src/python/test/test_off.py new file mode 100644 index 00000000..69bfa1f9 --- /dev/null +++ b/src/python/test/test_off.py @@ -0,0 +1,21 @@ +""" This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT. + See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. + Author(s): Marc Glisse + + Copyright (C) 2022 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +import gudhi as gd +import numpy as np +import pytest + + +def test_off_rw(): + for dim in range(2, 6): + X = np.random.rand(123, dim) + gd.write_points_to_off_file('rand.off', X) + Y = gd.read_points_from_off_file('rand.off') + assert Y == pytest.approx(X) -- cgit v1.2.3 From f4f930169fbf2db2707d69b334cfe5b941c64350 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Fri, 11 Nov 2022 19:03:54 +0100 Subject: black --- src/python/test/test_off.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/python/test/test_off.py b/src/python/test/test_off.py index 69bfa1f9..aea1941b 100644 --- a/src/python/test/test_off.py +++ b/src/python/test/test_off.py @@ -16,6 +16,6 @@ import pytest def test_off_rw(): for dim in range(2, 6): X = np.random.rand(123, dim) - gd.write_points_to_off_file('rand.off', X) - Y = gd.read_points_from_off_file('rand.off') + gd.write_points_to_off_file("rand.off", X) + Y = gd.read_points_from_off_file("rand.off") assert Y == pytest.approx(X) -- cgit v1.2.3 From 1e71120bdacb6f5ec4c90fb4c1365f76ecea7ff9 Mon Sep 17 00:00:00 2001 From: Marc Glisse Date: Sat, 12 Nov 2022 17:22:07 +0100 Subject: Handle QGLViewer >= 2.7.0 --- src/GudhUI/view/Viewer.cpp | 4 ++++ 1 file changed, 4 insertions(+) (limited to 'src') diff --git a/src/GudhUI/view/Viewer.cpp b/src/GudhUI/view/Viewer.cpp index 6b17c833..2c00f86f 100644 --- a/src/GudhUI/view/Viewer.cpp +++ b/src/GudhUI/view/Viewer.cpp @@ -31,7 +31,11 @@ void Viewer::set_bounding_box(const Point_3 & lower_left, const Point_3 & upper_ } void Viewer::update_GL() { +#if QGLVIEWER_VERSION >= 0x020700 + this->update(); +#else this->updateGL(); +#endif } void Viewer::init_scene() { -- cgit v1.2.3