How can i calculate the area under two curves that intersect?

14 visualizzazioni (ultimi 30 giorni)
I have two curves that plot by x-y coordinates
How can I calculate this area?
  3 Commenti
Sanufee
Sanufee il 21 Feb 2024
That's a two diffent cross section of the beach, So I want to know the volume that transiton

Accedi per commentare.

Risposta accettata

Star Strider
Star Strider il 20 Feb 2024
Modificato: Star Strider il 21 Feb 2024
Integrate them using trapz then do what you want with the results (keep them as they are, add them, add their absolute values, etc.) —
x = linspace(0, 40, 120);
y1 = sin(2*pi*x/30);
y2 = cos(2*pi*x/30);
L = numel(x);
idx = find(diff(sign(y1-y2))); % indices Of Intersections
% q = x(idx)
for k = 1:numel(idx)-1 % Integrate Between Indices
idxrng = max(1,idx(k)) : min(L,idx(k+1)); % Index Range
A(k) = trapz(x(idxrng), y1(idxrng)) - trapz(x(idxrng), y2(idxrng));
end
A
A = 1×2
13.4886 -13.4862
for k = 1:numel(idx)-1 % Interpolate Then Integrate
idxrng1 = max(1,idx(k)-1) : min(L,idx(k)+1); % Index Range (First Inmtersection)
idxrng2 = max(1,idx(k+1)-1) : min(L,idx(k+1)+1); % Index Range (Second Inmtersection)
xi(1) = interp1(y1(idxrng1)-y2(idxrng1), x(idxrng1), 0); % Find First 'x' Intersection
xi(2) = interp1(y1(idxrng2)-y2(idxrng2), x(idxrng2), 0); % Find Second 'x' Intersection
xv = x((x>=xi(1)) & (x<=xi(2))); % 'x' Vector Between Them
y1v = interp1(x, y1, xv); % Corresponding 'y1' Vector
y2v = interp1(x, y2, xv); % Corresponding 'y2' Vector
A(k) = trapz(xv, y1v) - trapz(xv, y2v); % Integrate
end
A
A = 1×2
13.4771 -13.4956
figure
plot(x, y1)
hold on
plot(x, y2)
hold off
grid
Here are two options. I suggest using the second one (interpolate first), since that seems to match your data more closely.
EDIT — (21 Feb 2024 02:35)
Guessing the data —
x1 = [2 4 7 9 10 12 14 17 19 22 24 27 29 32 34 37];
y1 = [-0.25 -0.4 -0.45 -0.6 -1.0 -1.2 -1.4 -1.45 -1.45 -1.5 -1.51 -1.6 -1.9 -2.3 -2.7 -2.9];
x2 = [2 4 7 9 10 12 13 17 19 22 23 27 29 34 36];
y2 = [-0.26 -0.4 -0.50 -0.7 -0.8 -1.1 -1.5 -1.7 -1.7 -1.7 -1.7 -1.8 -2.1 -2.6 -2.8];
% x1y1 = [x1; y1]
% x2y2 = [x2; y2]
xv = linspace(min([x1 x2]), max([x2 x2]), 100); % Common 'x' Vector
y1v = interp1(x1, y1, xv, 'pchip'); % Interpolate Using 'pchip'
y2v = interp1(x2, y2, xv, 'pchip'); % Interpolate Using 'pchip'
x = xv;
y1 = y1v;
y2 = y2v;
L = numel(xv);
idx = find(diff(sign(y1-y2))); % indices Of Intersections
for k = 1:numel(idx)-1 % Interpolate Then Integrate
idxrng1 = max(1,idx(k)-1) : min(L,idx(k)+1); % Index Range (First Inmtersection)
idxrng2 = max(1,idx(k+1)-1) : min(L,idx(k+1)+1); % Index Range (Second Inmtersection)
xi(1) = interp1(y1(idxrng1)-y2(idxrng1), x(idxrng1), 0); % Find First 'x' Intersection
xi(2) = interp1(y1(idxrng2)-y2(idxrng2), x(idxrng2), 0); % Find Second 'x' Intersection
xv = x((x>=xi(1)) & (x<=xi(2))); % 'x' Vector Between Them
y1v = interp1(x, y1, xv); % Corresponding 'y1' Vector
y2v = interp1(x, y2, xv); % Corresponding 'y2' Vector
A(k) = trapz(xv, y1v) - trapz(xv, y2v); % Integrate
end
figure
plot(x, y1, '-b', 'MarkerFaceColor','b')
hold on
plot(x, y2, '-r', 'MarkerFaceColor','r')
hold off
grid
text(10.5, -1, sprintf('\\leftarrow Area = %.3f', A(3)))
text(26, -1.7, sprintf('\\leftarrow Area = %.3f', A(4)))
text(20, -0.5, sprintf('Total Area = %.3f', sum(abs(A([3 4])))))
.
  27 Commenti
Sanufee
Sanufee il 24 Feb 2024
It's so useful for me! Thanks so much. maybe I could give my graduate certificate to you. 🥹
Star Strider
Star Strider il 24 Feb 2024
As always, my pleasure!
(I already have a few of my own! I am happy to be able to help you.)

Accedi per commentare.

Più risposte (1)

John D'Errico
John D'Errico il 20 Feb 2024
Modificato: John D'Errico il 21 Feb 2024
1. Subtract the two curves.
2. Take the absolute value of the difference.
3. Compute the integral.
Note that you may need the curves in a functional form if you want to do this accurately, so you might need an interpolation tool to do that. This is often in the form of a spline.
X = linspace(0,10,11);
Y1 = rand(1,11);
Y2 = rand(1,11);
f1 = spline(X,Y1);
f2 = spline(X,Y2);
fnplt(f1,'r')
hold on
fnplt(f2,'b')
plot(X,Y1,'ko',X,Y2,'kx')
hold off
AreaBetween = integral(@(x) abs(fnval(f1,x) - fnval(f2,x)),0,10)
AreaBetween = 3.2668
Note the improtance of using integral here, and of the use of functions in the form of splines to interpolate. The problem is, if you do not, then when the curves cross, you will get an error, and that error could be significant. For example, we could just use trapz directly on the data.
trapz(X,abs(Y1 - Y2))
ans = 3.5494
You should see there is a very significant difference.
plot(X,abs(Y1 - Y2),'o-')
hold on
fplot(@(x) abs(fnval(f1,x) - fnval(f2,x)),[0,10])
hold off
Where those curves actually crossed, the difference function should go all the way to zero. That is what the red curve does. In fact, the red curve is what you actually want to integrate.
But if the curves cross between two of the data points, the result will be the blue line. And that will totally mess up the trapezoidal integration, since it will see only the blue curve. You can see here, the difference is approximately a 10% error in the area estimate. The error from trapz will usually be an overestimate of the area.

Categorie

Scopri di più su Contour Plots in Help Center e File Exchange

Prodotti

Community Treasure Hunt

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

Start Hunting!

Translated by