test_hausdorff.py 6.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172
  1. import numpy as np
  2. from numpy.testing import (assert_allclose,
  3. assert_array_equal,
  4. assert_equal)
  5. import pytest
  6. from scipy.spatial.distance import directed_hausdorff
  7. from scipy.spatial import distance
  8. from scipy._lib._util import check_random_state
  9. class TestHausdorff:
  10. # Test various properties of the directed Hausdorff code.
  11. def setup_method(self):
  12. np.random.seed(1234)
  13. random_angles = np.random.random(100) * np.pi * 2
  14. random_columns = np.column_stack(
  15. (random_angles, random_angles, np.zeros(100)))
  16. random_columns[..., 0] = np.cos(random_columns[..., 0])
  17. random_columns[..., 1] = np.sin(random_columns[..., 1])
  18. random_columns_2 = np.column_stack(
  19. (random_angles, random_angles, np.zeros(100)))
  20. random_columns_2[1:, 0] = np.cos(random_columns_2[1:, 0]) * 2.0
  21. random_columns_2[1:, 1] = np.sin(random_columns_2[1:, 1]) * 2.0
  22. # move one point farther out so we don't have two perfect circles
  23. random_columns_2[0, 0] = np.cos(random_columns_2[0, 0]) * 3.3
  24. random_columns_2[0, 1] = np.sin(random_columns_2[0, 1]) * 3.3
  25. self.path_1 = random_columns
  26. self.path_2 = random_columns_2
  27. self.path_1_4d = np.insert(self.path_1, 3, 5, axis=1)
  28. self.path_2_4d = np.insert(self.path_2, 3, 27, axis=1)
  29. def test_symmetry(self):
  30. # Ensure that the directed (asymmetric) Hausdorff distance is
  31. # actually asymmetric
  32. forward = directed_hausdorff(self.path_1, self.path_2)[0]
  33. reverse = directed_hausdorff(self.path_2, self.path_1)[0]
  34. assert forward != reverse
  35. def test_brute_force_comparison_forward(self):
  36. # Ensure that the algorithm for directed_hausdorff gives the
  37. # same result as the simple / brute force approach in the
  38. # forward direction.
  39. actual = directed_hausdorff(self.path_1, self.path_2)[0]
  40. # brute force over rows:
  41. expected = max(np.amin(distance.cdist(self.path_1, self.path_2),
  42. axis=1))
  43. assert_allclose(actual, expected)
  44. def test_brute_force_comparison_reverse(self):
  45. # Ensure that the algorithm for directed_hausdorff gives the
  46. # same result as the simple / brute force approach in the
  47. # reverse direction.
  48. actual = directed_hausdorff(self.path_2, self.path_1)[0]
  49. # brute force over columns:
  50. expected = max(np.amin(distance.cdist(self.path_1, self.path_2),
  51. axis=0))
  52. assert_allclose(actual, expected)
  53. def test_degenerate_case(self):
  54. # The directed Hausdorff distance must be zero if both input
  55. # data arrays match.
  56. actual = directed_hausdorff(self.path_1, self.path_1)[0]
  57. assert_allclose(actual, 0.0)
  58. def test_2d_data_forward(self):
  59. # Ensure that 2D data is handled properly for a simple case
  60. # relative to brute force approach.
  61. actual = directed_hausdorff(self.path_1[..., :2],
  62. self.path_2[..., :2])[0]
  63. expected = max(np.amin(distance.cdist(self.path_1[..., :2],
  64. self.path_2[..., :2]),
  65. axis=1))
  66. assert_allclose(actual, expected)
  67. def test_4d_data_reverse(self):
  68. # Ensure that 4D data is handled properly for a simple case
  69. # relative to brute force approach.
  70. actual = directed_hausdorff(self.path_2_4d, self.path_1_4d)[0]
  71. # brute force over columns:
  72. expected = max(np.amin(distance.cdist(self.path_1_4d, self.path_2_4d),
  73. axis=0))
  74. assert_allclose(actual, expected)
  75. def test_indices(self):
  76. # Ensure that correct point indices are returned -- they should
  77. # correspond to the Hausdorff pair
  78. path_simple_1 = np.array([[-1,-12],[0,0], [1,1], [3,7], [1,2]])
  79. path_simple_2 = np.array([[0,0], [1,1], [4,100], [10,9]])
  80. actual = directed_hausdorff(path_simple_2, path_simple_1)[1:]
  81. expected = (2, 3)
  82. assert_array_equal(actual, expected)
  83. def test_random_state(self):
  84. # ensure that the global random state is not modified because
  85. # the directed Hausdorff algorithm uses randomization
  86. rs = check_random_state(None)
  87. old_global_state = rs.get_state()
  88. directed_hausdorff(self.path_1, self.path_2)
  89. rs2 = check_random_state(None)
  90. new_global_state = rs2.get_state()
  91. assert_equal(new_global_state, old_global_state)
  92. @pytest.mark.parametrize("seed", [None, 27870671])
  93. def test_random_state_None_int(self, seed):
  94. # check that seed values of None or int do not alter global
  95. # random state
  96. rs = check_random_state(None)
  97. old_global_state = rs.get_state()
  98. directed_hausdorff(self.path_1, self.path_2, seed)
  99. rs2 = check_random_state(None)
  100. new_global_state = rs2.get_state()
  101. assert_equal(new_global_state, old_global_state)
  102. def test_invalid_dimensions(self):
  103. # Ensure that a ValueError is raised when the number of columns
  104. # is not the same
  105. rng = np.random.default_rng(189048172503940875434364128139223470523)
  106. A = rng.random((3, 2))
  107. B = rng.random((3, 5))
  108. msg = r"need to have the same number of columns"
  109. with pytest.raises(ValueError, match=msg):
  110. directed_hausdorff(A, B)
  111. @pytest.mark.parametrize("A, B, seed, expected", [
  112. # the two cases from gh-11332
  113. ([(0,0)],
  114. [(0,1), (0,0)],
  115. 0,
  116. (0.0, 0, 1)),
  117. ([(0,0)],
  118. [(0,1), (0,0)],
  119. 1,
  120. (0.0, 0, 1)),
  121. # slightly more complex case
  122. ([(-5, 3), (0,0)],
  123. [(0,1), (0,0), (-5, 3)],
  124. 77098,
  125. # the maximum minimum distance will
  126. # be the last one found, but a unique
  127. # solution is not guaranteed more broadly
  128. (0.0, 1, 1)),
  129. ])
  130. def test_subsets(self, A, B, seed, expected):
  131. # verify fix for gh-11332
  132. actual = directed_hausdorff(u=A, v=B, seed=seed)
  133. # check distance
  134. assert_allclose(actual[0], expected[0])
  135. # check indices
  136. assert actual[1:] == expected[1:]
  137. @pytest.mark.xslow
  138. def test_massive_arr_overflow():
  139. # on 64-bit systems we should be able to
  140. # handle arrays that exceed the indexing
  141. # size of a 32-bit signed integer
  142. try:
  143. import psutil
  144. except ModuleNotFoundError:
  145. pytest.skip("psutil required to check available memory")
  146. if psutil.virtual_memory().available < 80*2**30:
  147. # Don't run the test if there is less than 80 gig of RAM available.
  148. pytest.skip('insufficient memory available to run this test')
  149. size = int(3e9)
  150. arr1 = np.zeros(shape=(size, 2))
  151. arr2 = np.zeros(shape=(3, 2))
  152. arr1[size - 1] = [5, 5]
  153. actual = directed_hausdorff(u=arr1, v=arr2)
  154. assert_allclose(actual[0], 7.0710678118654755)
  155. assert_allclose(actual[1], size - 1)