title: How to average all coordinates within a given distance in a vectorized way
body: <p>I did find a way to calculate the center coordinate of a cluster of points. However, my method is quite slow when the number of initial coordinates is increased (I have about 100 000 coordinates). </p>
  
  <p>The bottleneck is the for-loop in the code. I tried to remove it by using np.apply_along_axis, but discovered that this is nothing more than a hidden python-loop. </p>
  
  <p>Is it possible to detect and average out various sized clusters of too close points in a vectorized way?</p>
  
  <pre><code>import numpy as np
  from scipy.spatial import cKDTree
  np.random.seed(7)
  max_distance=1
  
  #Create random points
  points = np.array([[1,1],[1,2],[2,1],[3,3],[3,4],[5,5],[8,8],[10,10],[8,6],[6,5]])
  
  #Create trees and detect the points and neighbours which needs to be fused
  tree = cKDTree(points)
  rows_to_fuse = np.array(list(tree.query_pairs(r=max_distance))).astype('uint64')
  
  #Split the points and neighbours into two groups
  points_to_fuse = points[rows_to_fuse[:,0], :2]
  neighbours = points[rows_to_fuse[:,1], :2]
  
  #get unique points_to_fuse
  nonduplicate_points = np.ascontiguousarray(points_to_fuse)
  unique_points = np.unique(nonduplicate_points.view([('', nonduplicate_points.dtype)]\
                                                   *nonduplicate_points.shape[1]))
  unique_points = unique_points.view(nonduplicate_points.dtype).reshape(\
                                            (unique_points.shape[0],\
                                             nonduplicate_points.shape[1]))
  #Empty array to store fused points
  fused_points = np.empty((len(unique_points), 2))
  
  ####BOTTLENECK LOOP####
  for i, point in enumerate(unique_points):
      #Detect all locations where a unique point occurs
      locs=np.where(np.logical_and((points_to_fuse[:,0] == point[0]), (points_to_fuse[:,1]==point[1])))
      #Select all neighbours on these locations take the average
      fused_points[i,:] = (np.average(np.hstack((point[0],neighbours[locs,0][0]))),np.average(np.hstack((point[1],neighbours[locs,1][0]))))
  
  #Get original points that didn't need to be fused
  points_without_fuse = np.delete(points, np.unique(rows_to_fuse.reshape((1, -1))), axis=0)
  
  #Stack result
  points = np.row_stack((points_without_fuse, fused_points))
  </code></pre>
  
  <p><strong>Expected output</strong></p>
  
  <pre><code>&gt;&gt;&gt; points
  array([[  8.        ,   8.        ],
         [ 10.        ,  10.        ],
         [  8.        ,   6.        ],
         [  1.33333333,   1.33333333],
         [  3.        ,   3.5       ],
         [  5.5       ,   5.        ]])
  </code></pre>
  
  <p><strong>EDIT 1: Example of 1 loop with desired result</strong></p>
  
  <p><em>Step 1</em>: Create variables for the loop</p>
  
  <pre><code>#outside loop
  points_to_fuse = np.array([[100,100],[101,101],[100,100]])
  neighbours = np.array([[103,105],[109,701],[99,100]])
  unique_points = np.array([[100,100],[101,101]])
  
  #inside loop
  point = np.array([100,100])
  i = 0
  </code></pre>
  
  <p><em>Step 2</em>: Detect all locations where a unique point occurs in the points_to_fuse array</p>
  
  <pre><code>locs=np.where(np.logical_and((points_to_fuse[:,0] == point[0]), (points_to_fuse[:,1]==point[1])))
  &gt;&gt;&gt; (array([0, 2], dtype=int64),)
  </code></pre>
  
  <p><em>Step 3</em>: Create an array of the point and the neighbouring points at these locations and calculate the average</p>
  
  <pre><code>array_of_points = np.column_stack((np.hstack((point[0],neighbours[locs,0][0])),np.hstack((point[1],neighbours[locs,1][0]))))
  &gt;&gt;&gt; array([[100, 100],
             [103, 105],
             [ 99, 100]])
  fused_points[i, :] = np.average(array_of_points, 0)
  &gt;&gt;&gt; array([ 100.66666667,  101.66666667])
  </code></pre>
  
  <p><em>Loop output after a complete run</em>:</p>
  
  <pre><code>&gt;&gt;&gt; print(fused_points)
  &gt;&gt;&gt; array([[ 100.66666667,  101.66666667],
             [ 105.        ,  401.        ]])
  </code></pre>