CODEWITHSUNDEEP
matlabnumerical-integration
angryip
angryip

Reputation: 2290

Using trapz to find the area under a curve

I am attempting to use trapz to find the area under the curve.

My x data represents datetime, and my y data represents acceleration as a f(x). Accelerometer readings are in the form of SI.

Example data within x ( HH:mm:ss.SSS ) :

'01:06:48.330'
'01:06:48.352'
'01:06:48.374'

Example data within y ( accelerometer value * 9.81 ):

8.73750256159470
8.59619296907904
8.63520263017352

When I enter the following command (using the whole array of data ):

velocity = trapz(x,y);

I get back a duration array which looks like this:

velocity = 
    duration
    00:00:13

I am not entirely too sure I understand what 00:00:13 means. When I calculate for velocity, I'd expect to see something like 9.81 m/s , or 5m/s. Am I misusing the function or should I convert the duration array to a different object type?

Upvotes: 2

Views: 1268

Answers (2)

Wolfie
Wolfie

Reputation: 30047

The reason you expect m/s output from integrating acceleration is simply because you're doing a particular calculation involving (m/s^2)*s, i.e. y axis * x axis.

Let's use a simple example, where we first convert to seconds and then integrate.

x = datetime( {'00:00', '00:01', '00:02'}, 'inputformat', 'HH:mm' ); % Time
y = [1 2 4]; % Acceleration (m/s^2)

x_elapsed_seconds = datenum(x-min(x))*24*60*60; % Datenum returns days, convert to secs

z = trapz( x_elapsed_seconds, y ); % Velocity = 270 m/s

We can verify that, for this example, 270m/s is correct because there are simply 2 trapeziums in the calculation:

  1. Trapezium from 1m/s^2 to 2m/s^2 lasting 1min = 60secs: 60*(1+2)/2 = 60*1.5 = 90 m/s
  2. Trapezium from 2m/s^2 to 4m/s^2 lasting 1min = 60secs: 60*(2+4)/2 = 60*3 = 180 m/s

We sum the trapezium areas for the result: 90 + 180 = 270 as expected. This is why it's always best to use a simple example for verification before using real world data.

Upvotes: 1

Herman Wilén
Herman Wilén

Reputation: 195

Datetimes are usually tricky in MATLAB. I would convert it to seconds before doing any integrals. Here is an example solution. 

x=['01:06:48.330'; '01:06:48.352'; '01:06:48.374'];
datetime(x,'Format','HH:mm:ss.SSS');
secvec = datevec(x);
secvec = 3600*secvec(:,4)+60*secvec(:,5)+secvec(:,6); % times measured in seconds

y = [8.73750256159470; 8.59619296907904; 8.63520263017352]; % accelerations in m/s^2
trapz(secvec,y)
>> ans = 0.380216002428058 % Gained velocity, measured in m/s

Upvotes: 1

Related Questions