Algorithm/approach for regridding 3D coordinates

92 Views Asked by At

I have a question concerning surfaces of 3-dimensional structures, specifically if there is a way to interpolate or ‘regrid’ points/coordinates on a surface generated from the structure. I ask here on Math instead of SO, because it is not so much a question of programming, but instead a question of how to approach the problem mathematically. I know how to interpolate data in 1D (e.g. splines) or regrid them using binning (group data points into fewer points), but my usual methods fail with this 3D case. Importantly, there may be duplicates in my list of coordinates, the original list of coordinates is not sorted nor regularly spaced so there may be areas with high density of points.

For an example of what the structure look like with heaving down-sampled surface points (1.5 A away form surface):

Screenshot of structure (green), with the surface I am interested in regridding (blue). Only every 50th point in my list of coordinates for the blue surface has been plotted

The fundamental difference between what I usually do and what I need to do here is that I usually have a function f(x), which I need to have one another x-axis with different spacing between the points. In the present case I have coordinates (x,y,z) which I want to re-group in into a (much) smaller number of points, and ensure that they are regularly spaced – from e.g. 2 million points down to say 4000 points. By regularly spaced I mean that the Euclidian distance between neighboring points are the same down to some tolerance e.g. 5%.

I have considered using 3D interpolation routines such as Scipys Regular grid interpolator ( https://docs.scipy.org/doc/scipy-0.16.0/reference/generated/scipy.interpolate.RegularGridInterpolator.html ), but it is for interpolating a function of (x,y,z) not just the coordinates.

I have also considered using some sort of clustering method e.g. KMeans clustering, but I do not know what the impact of dublicates will have on that, and the method is very slow. From my initial small test, the clusters are also not so uniformly spread, so the final grid will also not be uniform.

Can someone suggest an algorithm/approach for solving this? Any suggestions would be highly appreciated. If someone would like to test something,

I have made this minimal example using python which will produce a grid which is not uniform.

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def generate_sphere(n, skew=0):
  """
  Returns list of coordinates on a sphere using the Golden-
  Section Spiral algorithm.
  """
  points = []
  inc = np.pi * (3 - np.sqrt(5))
  offset = 2 / float(n+skew)
  for k in range(n):
    y = k * offset - 1 + (offset / 2)
    r = np.sqrt(1 - y*y)
    phi = k * inc
    points.append([np.cos(phi)*r, y, np.sin(phi)*r])
  return np.array(points)



#%%
offset = 0.95
sph1 = generate_sphere(120)
sph2 = generate_sphere(120,25)
sph3 = generate_sphere(120)+offset

struc = np.vstack((sph1,sph2))
norms1 = np.linalg.norm(struc-np.mean(sph3,axis=0), axis=1)
struc = struc[norms1 > offset,:]
norms2 = np.linalg.norm(sph3,axis=1)
sph3 = sph3[norms2 > offset,:]
struc = np.vstack((struc,sph3))

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(struc[:,0],struc[:,1],struc[:,2])

EDIT: To answer Siddharth Bhat question, I am not interested in having my new data points/coordinates reflect the original distribution of coordinates. In fact, since my original distribution may contain duplicates and areas of higher density, it is not useful for me. I can round off the floating points down to some precision (say one decimal e.g. 13.1) and get rid of the duplicates, but the packing of coordinates will still be uneven. For the case you mention (one third of the points at x=0, one third of the points at x=1, and the final one third at x=100000) for 100000 points would like to regrid them into an evenly spaced list of e.g. 100: 0., 1010.1, 2020.2, 3030.3 ... 100000.

The initial distribution of one third of each of the values is not useful for me, and it is an artifact of my algorithm for generating the surface.

Apologies for being unclear :)