球坐标中的插值曲线

Interpolating curve in spherical coordinates

我有点,属于3D中的一些非参数曲线space。 我需要在球坐标系中工作,并在给定的 azimuthelevation 值中插入给定的曲线。

你能帮我看看这种插值算法吗?

P.S. 我不确定,从我的问题想法 1D 曲线绘制在 1D space 很清楚。因此,以防万一下面的代码创建了这样的参数。

x=linspace(-2*pi,2*pi,10^3);
y=sin(x);
z=sinh(x);

更新:

感谢@jodag 非常好的 post 关于上采样。 slerp 完美!尽管如此,我的问题是不同的,或者我不明白如何使用上采样来解决它。 假设我将我的点定义为数组(只是一些用于演示的虚拟数据):

x_given = rand(9, 1) - 0.5; 
y_given = rand(9, 1) - 0.5; 
z_given = rand(9, 1) - 0.5; 
x_given(end+1) = x_given(1); % should work for closed curves
y_given(end+1) = y_given(1); % should work for closed curves
z_given(end+1) = x_given(1); % should work for closed curves

我想找到给定 elevation_qazimuth_q 的所有点 (x_i, y_i, z_i):

elevation_q = pi./4;
azimuth_q = pi./7;
[x_i, y_i ,z_i] = interpolateCurve(x_given, y_given, y_given, elevation_q, azimuth_q);

那是我真正需要的。我明白,我们不需要从笛卡尔坐标系来回跳到球坐标系。但即使我将 x_given, y_given, y_given 替换为 elevation_given, azimuth_given, r_given,我也不知道如何编写 interpolateCurve() 函数。

这是一个在更高维度中对曲线进行线性插值的示例 space。请注意,我引入了一个参数化变量 t 对应于点的索引。

function [xq,yq,zq] = interp_sph(x,y,z,upsample)
    [az,el,r] = cart2sph(x,y,z);
    t = 1:numel(x);
    tq = 1:(1/upsample):numel(x);
    azq = interp1(t,az,tq);
    elq = interp1(t,el,tq);
    rq = interp1(t,r,tq);
    [xq,yq,zq] = sph2cart(azq,elq,rq);
end

从问题中不清楚这是否是您真正要查找的内容。球坐标中的直接插值肯定存在问题。例如,如果两个点彼此靠近,但一个点的方位角 = -179 度,另一个点的方位角 = 179,则插值将产生跨越 [-179,179] 的样本,而不是采用最短的 2 度路径。您可能感兴趣的是 slerp,它给出了更好的插值结果。

我写了下面的slerp例子来演示。请注意,slerp 要求定义 1/sin(theta),其中 theta 是函数中两点之间的角度。如果 theta == pi 那么我们有问题,slerp 将无法工作,所以如果您的数据接近原点,请小心。

function [xq,yq,zq] = interp_slerp(x,y,z,upsample)
    r = sqrt(x.^2+y.^2+z.^2);
    vq = zeros(3,(upsample-1)*(numel(x)-1)+1);
    for idx = 1:numel(x)-1
        t = 0:1;
        tq = linspace(0,1,upsample);
        % compute interpolated direction
        dq = slerp([x(idx),y(idx),z(idx)],[x(idx+1),y(idx+1),z(idx+1)],tq);
        % linearly interpolate radius
        rq = interp1(t,r(idx:(idx+1)),tq);
        % get interpolated path
        idxq = (idx-1)*upsample + 1;
        vq(:,idxq:(idxq+(upsample-1))) = bsxfun(@times,rq,dq);
    end
    xq = vq(1,:);
    yq = vq(2,:);
    zq = vq(3,:);
end

function pq = slerp(p0,p1,t)
    p0 = p0(:) / norm(p0);
    p1 = p1(:) / norm(p1);
    theta = acos( dot(p0,p1) );
    if(0==theta)
        pq = a;
    elseif(pi==theta)
        error('Angle between points cannot be exactly pi.');
    else
        pq = bsxfun(@times,(sin((1.0-t)*theta)/sin(theta)),p0) + bsxfun(@times,(sin(t*theta)/sin(theta)),p1);
    end
end

为了演示 interp_sph 的问题,请考虑以下示例。

% create test path
az = pi+sin(linspace(pi/2,3*pi/2,10));
el = cos(linspace(-pi,pi,10));
r = 1 + 1/5*sin(6*az).*sin(5*el);
[x,y,z] = sph2cart(az,el,r);

% test both spherical interpolation and slerp
upsample = 20;
[xq_sph,yq_sph,zq_sph] = interp_sph(x,y,z,upsample);
[xq_slerp,yq_slerp,zq_slerp] = interp_slerp(x,y,z,upsample);

% plot results
plot3(x,y,z,'LineWidth',2); hold on;
plot3(xq_sph,yq_sph,zq_sph,'LineWidth',2);
plot3(xq_slerp,yq_slerp,zq_slerp,'LineWidth',2);
axis('vis3d');
axis([-1,1,-1,1,-1,1]);
grid('on');
legend('Cartisian', 'Spherical', 'SLERP');

这产生了下面的情节。请注意,当方位角从正变为负时,球面插值会走很长一段路,但 slerp 会走预期路线。