update of point distance

changelog.rst

@@ -6,6 +6,14 @@ Pyinterpolate is the Python library for **geostatistics** and **spatial statisti
 Changes by date
 ===============
 
+2024-02-
+--------
+
+**version 0.5.1** (*pre production development*)
+
+* (refactoring) `point_distance` and `select_values_in_range` functions,
+*
+
 2023-09-16
 ----------
 

developer/dev_checks/debug_sandbox/check_cent_pk/check.py

@@ -170,7 +170,7 @@ def centroid_poisson_kriging(semivariogram_model: TheoreticalVariogram,
 
     # Distances between known blocks
     coordinates = kriging_data[:, :values_column_index]
-    block_distances = calc_point_to_point_distance(coordinates).flatten()
+    block_distances = point_distance(coordinates, coordinates).flatten()
     known_blocks_semivars = semivariogram_model.predict(block_distances)
     predicted = np.array(known_blocks_semivars.reshape(n, n))
     predicted = sem_to_cov(predicted, sill)

docs/source/api/distance/distance.rst

@@ -1,14 +1,19 @@
 Distance calculations
 =====================
 
-.. autofunction:: pyinterpolate.calc_point_to_point_distance
+.. autofunction:: pyinterpolate.point_distance
    :noindex:
 
-.....
-
 .. autofunction:: pyinterpolate.calc_block_to_block_distance
    :noindex:
 
+deprecated
+..........
+
+.. autofunction:: pyinterpolate.calc_point_to_point_distance
+   :noindex:
+
+
 Gridding
 ========
 

pyinterpolate/__init__.py

@@ -1,5 +1,5 @@
 # Distance
-from pyinterpolate.distance import calc_point_to_point_distance, calc_block_to_block_distance
+from pyinterpolate.distance import calc_point_to_point_distance, calc_block_to_block_distance, point_distance
 # from pyinterpolate.distance import aggregate_cluster, ClusterDetector
 from pyinterpolate.distance import create_grid, points_to_grid
 
@@ -42,4 +42,4 @@
 from pyinterpolate.viz import interpolate_raster
 
 
-__version__ = "0.5.0.post1"
+__version__ = "0.5.1"

pyinterpolate/distance/__init__.py

@@ -1,3 +1,4 @@
 from .distance import calc_point_to_point_distance, calc_block_to_block_distance
 # from .clusters import ClusterDetector, aggregate_cluster
-from .gridding import create_grid, points_to_grid
+from .gridding import create_grid, points_to_grid
+from .point import point_distance

pyinterpolate/distance/distance.py

@@ -4,6 +4,8 @@
 Authors
 -------
 1. Szymon Moliński | @SimonMolinsky
+
+TODO: the file will be removed in version 1.0
 """
 
 from typing import Dict, Union, Iterable
@@ -18,6 +20,8 @@
 from pyinterpolate.processing.preprocessing.blocks import PointSupport
 from pyinterpolate.processing.transform.transform import point_support_to_dict, block_dataframe_to_dict
 
+from deprecation import deprecated
+
 
 def _calc_b2b_dist_from_array(blocks: np.ndarray) -> Dict:
     """Function calculates distances between blocks.
@@ -305,6 +309,12 @@ def calc_angles(points_b: Iterable, point_a: Iterable = None, origin: Iterable =
     return angles
 
 
+@deprecated(
+    deprecated_in='0.5.1',
+    removed_in='1.0',
+    current_version='0.5.1',
+    details="Use `point_distance()` instead"
+)
 def calc_point_to_point_distance(points_a, points_b=None):
     """Function calculates distances between two group of points of a single group to itself.
 

pyinterpolate/distance/point.py

@@ -0,0 +1,89 @@
+"""
+Distance calculation functions.
+
+Authors
+-------
+1. Szymon Moliński | @SimonMolinsky
+"""
+import numpy as np
+from numpy.typing import ArrayLike
+from scipy.spatial.distance import cdist
+
+
+# noinspection PyTypeChecker
+def point_distance(points: ArrayLike,
+                   other: ArrayLike,
+                   metrics: str = 'euclidean') -> np.ndarray:
+    """
+    Calculates the euclidean distance from one group of points to another
+    group of points.
+
+    Parameters
+    ----------
+    points : array
+        Spatial coordinates.
+
+    other : array
+        Other array with spatial coordinates.
+
+    metrics : str, default = 'euclidean'
+        Metrics used to calculate distance.
+        See ``scipy.spatial.distance.cdist`` for more details.
+
+    Returns
+    -------
+    distances : array
+        Distances matrix. Row index = ``points`` point index, and column
+        index = ``other`` point index.
+
+    Notes
+    -----
+    The function creates array of size MxN, where M = number of ``points``
+    and N = number of ``other``. Very big coordinates array may cause
+    memory errors.
+
+    Examples
+    --------
+    >> points = [(0, 0), (0, 1), (0, 2)]
+    >> other = [(2, 2), (3, 3)]
+    >> distances = point_distance(points=points, other=other)
+    >> print(distances)
+    [[2.82842712 4.24264069]
+     [2.23606798 3.60555128]
+     [2.         3.16227766]]
+    """
+
+    distances = cdist(points, other, metrics)
+    return distances
+
+
+def select_values_in_range(data: np.ndarray,
+                           current_lag: float,
+                           previous_lag: float):
+    """
+    Function selects set of values which are greater than
+    (lag - step_size size) and smaller or equal to (lag).
+
+    Parameters
+    ----------
+    data : numpy array
+           Distances between points.
+
+    current_lag : float
+
+    previous_lag : float
+
+    Returns
+    -------
+    : numpy array
+        Mask with distances within a specified radius.
+    """
+
+    # Check conditions
+    condition_matrix = np.logical_and(
+            np.greater(data, previous_lag),
+            np.less_equal(data, current_lag))
+
+    # Find positions
+    position_matrix = np.where(condition_matrix)
+    return position_matrix

pyinterpolate/idw/idw.py

@@ -72,7 +72,7 @@ def inverse_distance_weighting(known_points: np.ndarray,
         unknown_location = unknown_location[np.newaxis, ...]
 
     # Calculate distances
-    distances = calc_point_to_point_distance(unknown_location, known_points[:, :-1])
+    distances = point_distance(unknown_location, known_points[:, :-1])
 
     # Check if any distance is equal to 0 - then return this value
     if not np.all(distances[0]):

pyinterpolate/kriging/models/block/centroid_based_poisson_kriging.py

@@ -135,7 +135,7 @@ def centroid_poisson_kriging(semivariogram_model: TheoreticalVariogram,
 
     # Distances between known blocks
     coordinates = kriging_data[:, :values_column_index]
-    block_distances = calc_point_to_point_distance(coordinates).flatten()
+    block_distances = point_distance(coordinates, coordinates).flatten()
     known_blocks_semivars = semivariogram_model.predict(block_distances)
     predicted = np.array(known_blocks_semivars.reshape(n, n))
     predicted = sem_to_cov(predicted, sill)

pyinterpolate/kriging/utils/process.py

@@ -77,7 +77,7 @@ def get_predictions(theoretical_model: TheoreticalVariogram,
     k = theoretical_model.predict(unknown_distances)
     k = k.T
 
-    dists = calc_point_to_point_distance(prepared_data[:, :2])
+    dists = point_distance(prepared_data[:, :2], prepared_data[:, :2])
 
     predicted_weights = theoretical_model.predict(dists.ravel())
     predicted = np.array(predicted_weights.reshape(n, n))

pyinterpolate/processing/select_values.py

@@ -368,7 +368,7 @@ def get_distances_within_unknown(point_support: np.ndarray):
         [[value1, value-n, distance between points 1-n], ..., [value-n, value1, distance between points n-1]]
     """
 
-    distances = calc_point_to_point_distance(point_support[:, :-1])
+    distances = point_distance(point_support[:, :-1], point_support[:, :-1])
     fdistances = distances.flatten()
 
     values = []
@@ -437,7 +437,7 @@ def prepare_pk_known_areas(point_support_dict: Dict,
             coordinates_b = ps_b[:, :-1]
             values_b = ps_b[:, -1]
             if bid_a != bid_b:
-                distances = calc_point_to_point_distance(coordinates_a, coordinates_b)
+                distances = point_distance(coordinates_a, coordinates_b)
             else:
                 distances = np.zeros(len(values_a) * len(values_b))
             fdistances = distances.flatten()
@@ -553,7 +553,7 @@ def select_kriging_data_from_direction(unknown_position: Iterable,
     r = np.array([unknown_position])
 
     known_pos = data_array[:, :-1]
-    dists = calc_point_to_point_distance(r, known_pos)
+    dists = point_distance(r, known_pos)
     angles = calc_angles(known_pos, origin=unknown_position)
     angle_diffs = calculate_angular_distance(angles, direction)
 
@@ -606,7 +606,7 @@ def select_kriging_data(unknown_position: Iterable,
     r = np.array([unknown_position])
 
     known_pos = data_array[:, :-1]
-    dists = calc_point_to_point_distance(r, known_pos)
+    dists = point_distance(r, known_pos)
 
     # Prepare data for kriging
     neighbors_and_dists = np.c_[data_array, dists.T]
@@ -753,7 +753,7 @@ def select_poisson_kriging_data(u_block_centroid: np.ndarray,
         point_s = k_point_support_dict[idx]
 
         # Distances between points
-        distances = calc_point_to_point_distance(u_point_support[:, :-1],
+        distances = point_distance(u_point_support[:, :-1],
                                                  point_s[:, :-1])
         fdistances = distances.flatten()
         ldist = len(fdistances)
@@ -972,7 +972,7 @@ def select_centroid_poisson_kriging_data(u_block_centroid: np.ndarray,
         dists = _calculate_weighted_distances(k_point_support_dict, u_index, u_point_support)
     else:
         # Calc from centroids
-        dists = calc_point_to_point_distance(k_centroids[:, :-1], [u_coordinates])
+        dists = point_distance(k_centroids[:, :-1], [u_coordinates])
 
     if direction is not None:
         angles = calc_angles(k_centroids[:, :-1], origin=u_coordinates)

pyinterpolate/variogram/empirical/cloud.py

@@ -44,7 +44,7 @@ def omnidirectional_point_cloud(input_array: np.array,
     variogram_cloud : dict
                       {Lag: array of semivariances within a given lag}
     """
-    distances = calc_point_to_point_distance(input_array[:, :-1])
+    distances = point_distance(input_array[:, :-1], input_array[:, :-1])
     lags = np.arange(step_size, max_range, step_size)
     variogram_cloud = OrderedDict()
 

pyinterpolate/variogram/empirical/covariance.py

@@ -60,7 +60,7 @@ def omnidirectional_covariogram(points: np.array, lags: np.array, step_size: flo
     """
 
     covariances_and_lags = list()
-    distances = calc_point_to_point_distance(points[:, :-1])
+    distances = point_distance(points[:, :-1], points[:, :-1])
 
     for h in lags:
         distances_in_range = select_values_in_range(distances, h, step_size)

pyinterpolate/variogram/empirical/semivariance.py

@@ -42,7 +42,7 @@ def omnidirectional_semivariogram(points: np.array, lags: np.array, step_size: f
 
     semivariances_and_lags = list()
     pts = points[:, :-1]
-    distances = calc_point_to_point_distance(pts)
+    distances = point_distance(pts, pts)
 
     for h in lags:
         distances_in_range = select_values_in_range(distances, h, step_size)

pyinterpolate/variogram/regularization/block/block_to_block_semivariance.py

@@ -42,7 +42,7 @@ def block_pair_semivariance(block_a: Collection,
                                   block a to every point in block b).
     """
 
-    distances_between_points = calc_point_to_point_distance(block_a, block_b).flatten()
+    distances_between_points = point_distance(block_a, block_b).flatten()
 
     predictions = semivariogram_model.predict(distances_between_points)
 

pyinterpolate/variogram/regularization/block/inblock_semivariance.py

@@ -34,7 +34,8 @@ def inblock_semivariance(points_of_block: np.ndarray, variogram_model: Theoretic
     number_of_points_within_block = len(points_of_block)  # P
     p = number_of_points_within_block * number_of_points_within_block  # P^2
 
-    distances_between_points = calc_point_to_point_distance(points_of_block[:, :-1])  # Matrix of size PxP
+    distances_between_points = point_distance(points_of_block[:, :-1],
+                                              points_of_block[:, :-1])  # Matrix of size PxP
     flattened_distances = distances_between_points.flatten()
     semivariances = variogram_model.predict(flattened_distances)
 

pyinterpolate/viz/raster.py

@@ -130,7 +130,7 @@ def interpolate_raster(data,
     # Calculate semivariance if not provided
 
     if semivariogram_model is None:
-        distances = calc_point_to_point_distance(data[:, :-1])
+        distances = point_distance(data[:, :-1], data[:, :-1])
 
         maximum_range = np.max(distances)
         number_of_divisions = 100

requirements.txt

@@ -15,4 +15,5 @@ prettytable
 pandas
 dask
 pyarrow
-pylibtiff==0.5.1
+pylibtiff==0.5.1
+deprecation==2.1.0

setup.cfg

@@ -51,6 +51,7 @@ install_requires =
     tqdm==4.64.0; python_version<='3.9'
     pyarrow
     pylibtiff==0.5.1
+    deprecation==2.1.0
 
 [options.packages.find]
 exclude =