Borrar filtros
Borrar filtros

Plot the integral of a discrete function

2 visualizaciones (últimos 30 días)
Antonio
Antonio el 14 de Mayo de 2014
Editada: Antonio el 14 de Mayo de 2014
Good morning, I have to plot the integral of a discrete function (t,df/dt) defined in a file *.txt. I want to plot (t,f). f is an angle. I've tried with the command trapz but it gives only the numeric value of the area. Thanks

Respuestas (2)

Youssef  Khmou
Youssef Khmou el 14 de Mayo de 2014
Trapz function returns a scalar value, numeric primitive can be calculated using the following function :
function itg=integral(f,dx,smooth);
% INTEGRAL ANS=INTEGRAL(F,DX)
% This function computes the integral of F(X)DX where the integrand
% is specified at discrete points F spaced DX apart (F is a vector,
% DX is a scalar). Simpsons Rule is used, so that the error
% is O(dx^5*F4). (F4 is the 4th derivative of F).
%
% If F is a matrix, then the integration is done for each column.
%
% If F is really spiky, then INTEGRAL(F,DX,'smooth') may
% provide a better looking result (the result is smoothed
% with a 3 point triangular filter).
%
% Author: RP (WHOI) 15/Aug/92
[N,M]=size(f);
if (N==1 | M==1),
N=max(size(f));
itg=zeros(size(f));
itg(1)=0; % first element
itg(2)=(5*f(1)+8*f(2)-f(3))*dx/12; % Parabolic approx to second
itg(3:N)=(f(1:N-2)+4*f(2:N-1)+f(3:N))*dx/3; % Simpsons rule for 2-segment
% intervals
itg(1:2:N)=cumsum(itg(1:2:N)); % Sum up 2-seg integrals
itg(2:2:N)=cumsum(itg(2:2:N));
if (nargin>2), % ... apply smoothing
itg(2:N-1)=(itg(1:N-2)+2*itg(2:N-1)+itg(3:N))/4;
itg(N)= (itg(N-1)+itg(N))/2;
end;
else
itg=zeros(size(f));
itg(1,:)=zeros(1,M);
itg(2,:)=(5*f(1,:)+8*f(2,:)-f(3,:))*dx/12;
itg(3:N,:)=(f(1:N-2,:)+4*f(2:N-1,:)+f(3:N,:))*dx/3;
itg(1:2:N,:)=cumsum(itg(1:2:N,:)); % Sum up 2-seg integrals
itg(2:2:N,:)=cumsum(itg(2:2:N,:));
if (nargin>2), % ... apply smoothing
itg(2:N-1,:)=(itg(1:N-2,:)+2*itg(2:N-1,:)+itg(3:N,:))/4;
itg(N,:)= (itg(N-1,:)+itg(N,:))/2;
end;
end;
Example :
t=0:0.1:10;
y=sin(t);
z=integral(y,0.1);
figure; plot(t,y,t,z);

Antonio
Antonio el 14 de Mayo de 2014
Editada: Antonio el 14 de Mayo de 2014
Thanks for the answer.
The code gives me the error
In an assignment A(I) = B, the number of elements in B and I must be the same.
could you apply your code to the real case of study I've attached here?
I need to plot f.
load data_hw1.txt
dx = data_hw1(:,1)
F = data_hw1(:,3) % F=df/dx
ps my homework suggests the trapezoidal integration.
  2 comentarios
Youssef  Khmou
Youssef Khmou el 14 de Mayo de 2014
the sampling rate is :
DX=dx(300)-dx(299);
FF=integral(F,DX);
figure; plot(dx,F,dx,FF),legend(' F','\int F');
Antonio
Antonio el 14 de Mayo de 2014
Editada: Antonio el 14 de Mayo de 2014
It gives me error :( I've tried with this code (I need to plot angle in the interval 0deg:360deg. is it correct?
load data_hw1.txt
dx = data_hw1(:,1)
f = data_hw1(:,3);
FF=cumsum(f)
a=FF./10;
for n = 1:1:length(dx)
if a(n,1) < -360
FF(n)=(a(n,1)+720)
else
FF(n)=(a(n,1)+360)
end
end
figure(1)
plot(dx,f,dx,FF)

Iniciar sesión para comentar.

Categorías

Más información sobre Numerical Integration and Differentiation en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by