CODEWITHSUNDEEP
pythonmatlabnumpyscipy
Carlos Bermudez
Carlos Bermudez

Reputation: 73

Different results with matlab cumtrapz and scipy.integrate cumtrapz

Im translating some matlab code to python code and debugging both codes i get a different result from a call to the cumtrapz function, i also verified that the input data of both is similar. This are the codes:

Python Code

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 Code

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

The hvel variable is different for both codes. Maybe im using wrong the scipy cumtrapz because i asume that the initial value that recives is 0. In both cases the inputs ics_y(ICsy), fcs_y(FCsY), acc_y(ACCY) are one dimensional arrays and l and sf are scalars.

Thanks!!!

Upvotes: 2

Views: 1878

Answers (1)

Warren Weckesser
Warren Weckesser

Reputation: 114831

(If this question is about cumtrapz, you should simplify your tests to just a single call to cumtrapz with the same input arrays in matlab and Python. Also, be sure you read the matlab and SciPy documentation of each function carefully. The SciPy functions are typically not exact duplicates of the corresponding matlab function.)

The problem is that when you give both x and y values, the order in which they are given in matlab/octave is x, y, but in the SciPy version, it is y, x.

For example,

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

To get the same result with 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])

Upvotes: 5

Related Questions