python - matlab cumtrapz 和 scipy.integrate cumtrapz 的不同结果
问题描述
我将一些matlab代码翻译成python代码并调试这两个代码我从调用cumtrapz
函数得到不同的结果,我还验证了两者的输入数据是相似的。这是代码:
Python代码
from numpy import zeros, ceil, array, mean, ptp, abs, sqrt, power
from scipy.integrate import cumtrapz
def step_length_vector(ics_y, fcs_y, acc_y, l, sf):
step_length_m1 = zeros(int(ceil(len(ics_y)/2))-1)
for i in range(0, len(ics_y)-2, 2):
av = acc_y[int(ics_y[i]):int(ics_y[i+2])+1]
t = array(range(1, int((ics_y[i+2]-ics_y[i])+2)))/sf
hvel = cumtrapz(t, av - mean(av), initial=0)
h = cumtrapz(t, hvel - mean(hvel), initial=0)
hend = ptp(h)
sl = 6.8*(sqrt(abs(2*l*hend - hend**2)))
step_length_m1[int(ceil(i/2))] = sl
return step_length_m1
Matlab代码
function [StepLengthM1] = StepLengthVector(ICsY,FCsY,ACCY,l,sf)
StepLengthM1 = zeros(1,ceil(length(ICsY)/2)-1);
for i= 1:2:length(ICsY)-2
av = ACCY(ICsY(i):ICsY(i+2));
t = (1:(ICsY(i+2)-ICsY(i))+1)/sf;
hvel = cumtrapz(t,av-mean(av));
h = cumtrapz(t,hvel-mean(hvel));
hend = peak2peak(h);
sl = 6.8*(sqrt(abs(2*l*hend - hend.^2)));
StepLengthM1(ceil(i/2)) = sl;
end
end
两个代码的hvel
变量不同。也许我用错了,scipy
cumtrapz
因为我认为接收的initial
价值是0
。在这两种情况下,输入ics_y(ICsy)
、fcs_y(FCsY)
、 acc_y(ACCY)
都是一维数组,并且l
和sf
是标量。
谢谢!!!
解决方案
(如果这个问题cumtrapz
是cumtrapz
关于对应的 matlab 函数的完全相同的副本。)
问题是,当您同时给出x
和y
值时,它们在 matlab/octave 中给出的顺序是x, y
,但在 SciPy 版本中,它是y, x
.
例如,
octave:11> t = [0 1 1.5 4 4.5 6]
t =
0.00000 1.00000 1.50000 4.00000 4.50000 6.00000
octave:12> y = [1 2 3 -2 0 1]
y =
1 2 3 -2 0 1
octave:13> cumtrapz(t, y)
ans =
0.00000 1.50000 2.75000 4.00000 3.50000 4.25000
要获得相同的结果scipy.integrate.cumtrapz
:
In [22]: from scipy.integrate import cumtrapz
In [23]: t = np.array([0, 1, 1.5, 4, 4.5, 6])
In [24]: y = np.array([1, 2, 3, -2, 0, 1])
In [25]: cumtrapz(y, t, initial=0)
Out[25]: array([0. , 1.5 , 2.75, 4. , 3.5 , 4.25])
推荐阅读
- xcode - 在 xcode 中添加另一个团队
- c# - 将列表复制到剪贴板以将富文本粘贴到 Outlook
- aws-lambda - (AWS) 通过 Lambda 获取 RDS 状态
- java - 应用网络流
- javascript - Mongoose + Mongodb User.update 不工作
- qt5 - 如何启动 SystemC 模拟的 Qt 前端
- php - 在 Elasticsearch 中实施分析器后没有结果
- javascript - 在Web应用程序上按下按钮时下载文件?
- firebase - React Native:com.android.build.api.transform.TransformException:生成主 dex 列表时出错
- java - 在 Intellij 中删除和重新导入项目