使用径向基函数对球上的函数进行插值 [英] Using Radial Basis Functions to Interpolate a Function on a Sphere

问题描述

首先，有一点背景知识:

我正在使用球谐函数作为球体表面(如此图像中的前球体)上的函数的示例:

我制作了其中一个球体，并根据其表面上各点的谐波函数的值对其进行了着色.我首先针对大量点进行此操作，因此我的功能非常准确.我称它为我的fine领域.

现在我有了我的fine球体，我在球体上获得的点相对较少.这些是我想从训练数据中插值的点，我称它们为interp点.这是我的interp点，用它们的值着色，并绘制在我的fine球面上.

现在，该项目的目标是使用这些interp点来训练解决方案

您可以更改标准距离:

def euclidean_norm(x1, x2):
    return np.sqrt( ((x1 - x2)**2).sum(axis=0) )

按球面距离(例如，参见以下问题SciPy Radial Basis Function to interpolate my function on the sphere. I was able to do this using:

# Train the interpolation using interp coordinates
rbf = Rbf(interp.phi, interp.theta, harmonic13_coarse)
# The result of the interpolation on fine coordinates
interp_values = rbf(fine.phi, fine.theta)

Which produced this interpolation, plotted on the sphere:

Hopefully, through this last image, you can see my problem. Notice the line running through the interpolation? This is because the interpolation data has a boundary. The boundary is because I trained the radial basis function using spherical coordinates (boundaries at [0,pi] and [0,2pi]).

rbf = Rbf(interp.phi, interp.theta, harmonic13_coarse)

My goal, and why I'm posting this problem, is to interpolate my function on the surface of the sphere using the x,y,z Cartesian coordinates of the data on the sphere. This way, since spheres don't have boundaries, I won't have this boundary error like I do in spherical coordinates. However, I just can't figure out how to do this.

I've tried simply giving the Rbf function the x,y,z coordinates and the value of the function.

rbf=Rbf(interp.x, interp.y, interp.z, harmonic13_coarse)
interp_values=rbf(fine.x,fine.y,fine.z)

But NumPy throws me a Singular Matrix Error

numpy.linalg.linalg.LinAlgError: singular matrix

Is there any way for me to give Rbf my data sites in Cartesian coordinates, with the function values at each site and have it behave like it does with spherical coordinates but without that boundaries? From the Rbf documentation, there is the attribute norm for defining a different distance norm, could I have to use a spherical distance to get this to work?

I'm pretty much stumped on this. Let me know if you have any ideas for interpolating my function on a sphere without the boundaries of spherical coordinates.

Here is my code in full:

import matplotlib.pyplot as plt
from matplotlib import cm, colors
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
from scipy import special
from scipy.interpolate import Rbf
from collections import namedtuple
from mayavi import mlab

# Nice aliases
pi = np.pi
cos = np.cos
sin = np.sin

# Creating a sphere in Cartesian and Sphereical
# Saves coordinates as named tuples


def coordinates(r, n):
    phi, theta = np.mgrid[0:pi:n, 0:2 * pi:n]
    Coor = namedtuple('Coor', 'r phi theta x y z')
    r = r
    x = r * sin(phi) * cos(theta)
    y = r * sin(phi) * sin(theta)
    z = r * cos(phi)
    return Coor(r, phi, theta, x, y, z)

# Creating a sphere
# fine is coordinates on a fine grid
# interp is coordinates on coarse grid for training interpolation
fine = coordinates(1, 100j)
interp = coordinates(1, 5j)


# Defining finection to colour sphere
# Here we are using a spherical harmonic
def harmonic(m, n, theta, phi):
    return special.sph_harm(m, n, theta, phi).real
norm = colors.Normalize()

# One example of the harmonic function, for testing
harmonic13_fine = harmonic(1, 3, fine.theta, fine.phi)
harmonic13_coarse = harmonic(1, 3, interp.theta, interp.phi)


# Train the interpolation using interp coordinates
rbf = Rbf(interp.phi, interp.theta, harmonic13_coarse)
# The result of the interpolation on fine coordinates
interp_values = rbf(fine.phi, fine.theta)


rbf=Rbf(interp.x, interp.y, interp.z, harmonic13_coarse)
interp_values=rbf(fine.x,fine.y,fine.z)

#Figure of harmoinc function on sphere in fine cordinates
#Points3d showing interpolation training points coloured to their value
mlab.figure()
vmax, vmin = np.max(harmonic13_fine), np.min(harmonic13_fine)
mlab.mesh(fine.x, fine.y, fine.z, scalars=harmonic13_fine, vmax=vmax, vmin=vmin)
mlab.points3d(interp.x, interp.y, interp.z, harmonic13_coarse,
              scale_factor=0.1, scale_mode='none', vmax=vmax, vmin=vmin)


#Figure showing results of rbf interpolation
mlab.figure()
vmax, vmin = np.max(harmonic13_fine), np.min(harmonic13_fine)
mlab.mesh(fine.x, fine.y, fine.z, scalars=interp_values)
# mlab.points3d(interp.x, interp.y, interp.z, scalars, scale_factor=0.1, scale_mode='none',vmax=vmax, vmin=vmin)

mlab.show()

解决方案

You can change the standard distance:

def euclidean_norm(x1, x2):
    return np.sqrt( ((x1 - x2)**2).sum(axis=0) )

by the sphere distance (see, for instance, this question Haversine Formula in Python (Bearing and Distance between two GPS points)).

这篇关于使用径向基函数对球上的函数进行插值的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！