将 MATLAB interp 1 转换为 Python interp1d

Question

我正在将 MATLAB 代码转换为 Python 代码。

代码使用了MATLAB中的函数interp1。我发现 scipy 函数 interp1d 应该是我所追求的，但我不确定。你能告诉我我实现的代码是否正确吗？我的Python版本是3.4.1，MATLAB版本是R2013a。不过代码已经在2010年左右实现了。

MATLAB:

S_T = [0.0, 2.181716948, 4.363766232, 6.546480392, 8.730192373, ...
      10.91523573, 13.10194482, 15.29065504, 17.48170299, 19.67542671, ...
      21.87216588, 24.07226205, 26.27605882, 28.48390208; ...
     1.0, 1.000382662968538, 1.0020234819906781, 1.0040560245904753, ...
       1.0055690037530718, 1.0046180687475195, 1.000824223678225, ...
       0.9954866694014762, 0.9891408937764872, 0.9822543350571298, ...
       0.97480163751874, 0.9666158376141503, 0.9571711322843011, ...
       0.9460998105962408; ...
     1.0, 0.9992731388936672, 0.9995093132493109, 0.9997021748479805, ...
       0.9982835412406582, 0.9926319477117723, 0.9833685776596993, ...
       0.9730725288209638, 0.9626092685176822, 0.9525234896714959, ...
       0.9426698515488858, 0.9326788630704709, 0.9218100196936996, ...
       0.9095717918978693];

S = transpose(S_T);

dist = 0.00137;
old = 15.61;
ll = 125;
ref = 250;
start = 225;
high = 7500;
low = 2;

U = zeros(low,low,high);

for ii=1:high
    g0= start-ref*dist*ii;
    g1= g0+ll;
    if(g0 <=0.0 && g1 >= 0.0)
        temp= old/2*(1-cos(2*pi*g0/ll));
        for jj=1:low
            U(jj,jj,ii)= temp;
        end
    end
end

for ii=1:low
 S_mod(ii,1,:)=interp1(S(:,1),S(:,ii+1),U(ii,ii,:),'linear');
end

Python:

import numpy
import os
from scipy import interpolate

S = [[0.0, 2.181716948, 4.363766232, 6.546480392, 8.730192373, 10.91523573, 13.10194482, 15.29065504, \
      17.48170299, 19.67542671, 21.87216588, 24.07226205, 26.27605882, 28.48390208], \
     [1.0, 1.000382662968538, 1.0020234819906781, 1.0040560245904753, 1.0055690037530718, 1.0046180687475195, \
      1.000824223678225, 0.9954866694014762, 0.9891408937764872, 0.9822543350571298, 0.97480163751874, \
      0.9666158376141503, 0.9571711322843011, 0.9460998105962408], \
     [1.0, 0.9992731388936672, 0.9995093132493109, 0.9997021748479805, 0.9982835412406582, 0.9926319477117723, \
      0.9833685776596993, 0.9730725288209638, 0.9626092685176822, 0.9525234896714959, 0.9426698515488858, \
      0.9326788630704709, 0.9218100196936996, 0.9095717918978693]]

dist = 0.00137
old = 15.61
ll = 125
ref = 250
start = 225
high = 7500
low = 2

U = [numpy.zeros( [low, low] ) for _ in range(high)]

for ii in range(high):
    g0 = start - ref * dist * (ii+1)   
    g1 = g0 + ll
    if g0 <=0.0 and g1 >= 0.0:
        for jj in range(low):
            U[ii][jj,jj] = old / 2 * (1 - numpy.cos( 2 * numpy.pi * g0 / ll) )

S_mod = []

for jj in range(high):
    temp = []
    for ii in range(low):
        temp.append(interpolate.interp1d( S[0], S[ii+1], U[jj][ii,ii]))

    S_mod.append(temp)

Answer 1

有一个 python 库可让您通过包装器使用 MATLAB 函数：mlabwrap。如果您不需要更改函数本身的代码，这可以节省您一些时间。

Answer 2

我不知道scipy，但我可以告诉你 MATLAB 中的 interp1 调用在做什么：

http://www.mathworks.com/help/matlab/ref/interp1.html

您使用的语法：

vq = interp1(x,v,xq,method)

"Vector x contains the sample points, and v contains the corresponding values, v(x). Vector xq contains the coordinates of the query points."

因此，在您的代码中，S(:,1) 包含定义网格的采样点，S(:,ii+1) 包含一维函数的采样值，而 U(ii) ,ii,:) 包含您要在其中进行插值以在网格中的已知值之间查找新函数值的查询点。您正在使用线性插值。

一维插值是一个定义非常明确的操作，interp1 是该操作的一个相对简单的接口。你到底不明白什么？你清楚什么是插值吗？

本质上，你有一个离散定义的函数 f[x]，interp1 的第一个参数是 x，第二个参数是 f[x]，第三个参数是任意定义的查询点 Xq找到新的函数值 f[Xq]。由于这些值未知，因此您必须使用插值方法来近似 f[Xq]。 'linear' 表示您将使用最接近 Xq 的两个已知采样邻居（左右邻居）的线性加权平均值。

Answer 3

好的，我已经解决了我自己的问题（感谢 Alex 对 MATLAB interp1 的解释！）。

python interp1d 本身没有查询点，而是创建一个函数，然后您可以使用该函数来获取新的数据点。因此，它应该是：

        f = interpolate.interp1d( S[0], S[ii+1])
        temp.append(f(U[jj][ii,ii]))

将 MATLAB interp 1 转换为 Python interp1d

Converting MATLAB's interp1 to Python interp1d

python

matlab