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

Question

Sanufee 2024-2-20

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2084483-how-can-i-calculate-the-area-under-two-curves-that-intersect

评论： Star Strider 2024-2-24

采纳的回答： Star Strider

I have two curves that plot by x-y coordinates

How can I calculate this area?

3 个评论
显示 1更早的评论隐藏 1更早的评论

Sanufee 2024-2-21

Yess

Sanufee 2024-2-21

That's a two diffent cross section of the beach, So I want to know the volume that transiton

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2024-2-20

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2084483-how-can-i-calculate-the-area-under-two-curves-that-intersect#answer_1412878

编辑：Star Strider 2024-2-21

在 MATLAB Online 中打开

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 个评论
显示 25更早的评论隐藏 25更早的评论

Star Strider 2024-2-21

在 MATLAB Online 中打开

It would help to have your data.

My variables are all row vectors, so my guess is that yours are column vectors.

In that instance, one option would be:

xv = linspace(min([x1; x2]), max([x2; x2]), 100); % Common 'x' Vector

and another:

xv = linspace(min([x1 x2],[],'all'), max([x2 x2],[],'all'), 100); % Common 'x' Vector

At least one of those should work.

If you are still having problems with your data and my code, please provide your data as a text, spreadhseet, or .mat file.

Converting my initial vectors to column vectors works with the first option, and no other changes in the rest of my code —

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])))))

.

Star Strider 2024-2-21

在 MATLAB Online 中打开

I needed to tweak my code a bit to include the areas at both ends (defined by only one crossing). I could not understand what information you want from the spreadsheet files, so I used the supplied vectors.

Try this —

% Pre = readtable('Premonsoon.xlsx')

% Post = readtable('Postmonsoon.xlsx')

x1 = [0.8961 3.3569 5.8782 8.3384 10.8515 13.3112 15.8254 18.3049 20.8188 23.3092 25.771 28.2613 29.7442 32.2171 33.1798 35.9226].';

y1 = [-0.3838 -0.5309 -0.762 -0.9384 -1.1361 -1.2114 -1.4222 -1.2912 -1.3113 -1.5192 -1.8071 -2.0278 -2.1672 -2.4575 -2.5841 -2.9109].';

x2 = [2.0252 4.4448 6.8848 9.4151 12.0184 14.3916 17.01 19.3367 21.8709 24.3984 26.8193 29.3898 31.8715 34.2408 35.8614].';

y2 = [-0.3586 -0.5645 -0.75 -1.0408 -1.0724 -1.1438 -1.164 -1.5574 -1.8426 -1.7756 -1.7499 -1.9891 -2.2325 -2.3882 -2.5566].';

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'

x1a = x1;

y1a = y1;

x2a = x2;

y2a = y2;

x = xv;

y1 = y1v;

y2 = y2v;

L = numel(xv);

idx = find(diff(sign(y1-y2))); % indices Of Intersections

idx = [1 idx L]; % Include Both Ends

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

if isnan(xi(1)) % The 'NaN' Values Result Of There Are No Intersections (At Each End), So Test For Them & Make Appropriate Substitutions

xi(1) = x(1);

elseif isnan(xi(2))

xi(2) = x(end);

end

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

xc(k) = mean(xv); % Used For 'text' Loations

yc(k) = mean([y1v y2v]); % Used For 'text' Loations

end

A % Segment Areas

A = 1×5

-0.5386 0.0191 -1.1436 2.0686 -2.0872

figure

hp(1) = plot(x1a, y1a, '.b', 'MarkerFaceColor','b', 'DisplayName','Premonsoon');

hold on

plot(x, y1, '-b')

hp(2) = plot(x2a, y2a, '.r', 'MarkerFaceColor','r', 'DisplayName','Postmonsoon');

plot(x, y2, '-r')

hold off

grid

axis([0 45 -3.5 0])

legend(hp,'Location','best')

for k = 1:numel(A)

text(xc(k), yc(k), sprintf('\\leftarrow Area = %.3f', A(k)))

end

text(7, -2.25, sprintf('Increasing = %.3f\nDecreasing = %.3f',abs(sum(A([1:3 5]))),abs(A(4))))

I did not use patch to colour the individual areas, although I can if that is important.

Star Strider 2024-2-21

在 MATLAB Online 中打开

Postmonsoon.xlsx

It is best to use the scatteredInterpolant function to calcualte the surface from the data. (I did it two ways here with the same scatteredInterpolant function, first with the original unsorted vectors to do the interpolation to calculate the surface, and second with sorted vectors. The interpolation method choices are 'linear', 'nearest', and 'natural'. I experimented with them, and like 'natural' best in this instance.) Choose the appropriate colormap and view to depict it as you want. I also used the surfc function to draw the surfaces, since the contour plot lets me see details that I would otherwise not be able to. Use it or surf to drtaw the surfaces.

Post = readtable('Postmonsoon.xlsx')

Post = 529×4 table

Point E N Elevation _____ ______ _______ _________ NaN 2.0642 -0.0163 -0.1885 NaN 0 0 0 1 4.4896 -0.0456 -0.3672 2 4.4866 -0.0455 -0.3671 3 6.9964 -0.0443 -0.4542 4 9.4778 -0.0158 -0.7015 5 12.005 -0.0133 -0.8593 6 14.603 -0.0125 -1.0918 7 16.983 -0.0117 -1.4529 8 19.3 -0.011 -1.6625 9 21.899 -0.0101 -1.7106 10 23.91 -0.0145 -1.6612 11 26.93 -0.0143 -1.8578 12 29.437 -0.014 -2.123 13 34.331 -0.0138 -2.564 14 36.5 -0.049 -2.8317

E = Post.E(3:end);

N = Post.N(3:end);

Elev = Post.Elevation(3:end);

F = scatteredInterpolant(E,N,Elev,'natural')

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

F =

scatteredInterpolant with properties: Points: [526×2 double] Values: [526×1 double] Method: 'natural' ExtrapolationMethod: 'linear'

[Em,Nm] = ndgrid(E, N);

Elevm = F(Em,Nm);

figure

surfc(Em, Nm, Elevm, 'EdgeColor','none')

colormap(turbo)

xlabel('E')

ylabel('N')

zlabel('Elevation')

title('Unsorted Vectors')

view(135, 30)

[Em,Nm] = ndgrid(sort(E), sort(N));

Elevm = F(Em,Nm);

figure

surfc(Em, Nm, Elevm, 'EdgeColor','interp')

colormap(turbo)

xlabel('E')

ylabel('N')

zlabel('Elevation')

title('Sorted Vectors')

view(135, 30)

.

Star Strider 2024-2-22

在 MATLAB Online 中打开

Postmonsoon.xlsx

The Z = -2 plane is straightforward. You can either use surf or patch. The'FaceAlpha' and 'EdgeAlpha' name-value pairs control the transparency of the surface. I set them both to 0.5 here, increasing that value makes them less transparent, and decreasing it makes them more transparent. The allowed range is from 0 to 1.

I am not certain what you want with respect to expanding the E axis. One option is to set the limits to be larger with

xlim([-10 50])

since the ‘x’ axis is the ‘E’ axis, however it would not be a good idea to expand the surface to match it. That would require extrapolation, and that is never appropriate because it is impossible to know what the surface (or any such function) actually does outside the region that originally defines it. Interpolating is appropriate, however extrapolating never is appropriate.

Post = readtable('Postmonsoon.xlsx')

Post = 529x4 table

Point E N Elevation _____ ______ _______ _________ NaN 2.0642 -0.0163 -0.1885 NaN 0 0 0 1 4.4896 -0.0456 -0.3672 2 4.4866 -0.0455 -0.3671 3 6.9964 -0.0443 -0.4542 4 9.4778 -0.0158 -0.7015 5 12.005 -0.0133 -0.8593 6 14.603 -0.0125 -1.0918 7 16.983 -0.0117 -1.4529 8 19.3 -0.011 -1.6625 9 21.899 -0.0101 -1.7106 10 23.91 -0.0145 -1.6612 11 26.93 -0.0143 -1.8578 12 29.437 -0.014 -2.123 13 34.331 -0.0138 -2.564 14 36.5 -0.049 -2.8317

E = Post.E(3:end);

N = Post.N(3:end);

Elev = Post.Elevation(3:end);

F = scatteredInterpolant(E,N,Elev,'natural');

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

[Em,Nm] = ndgrid(sort(E), sort(N));

Elevm = F(Em,Nm);

figure

surfc(Em, Nm, Elevm, 'EdgeColor','interp')

hold on

surf(Em, Nm, -2*ones(size(Em)), 'FaceAlpha',0.5, 'EdgeAlpha',0.5, 'EdgeColor','none')

hold off

colormap(turbo)

xlabel('E')

ylabel('N')

zlabel('Elevation')

title('Sorted Vectors')

view(135, 30)

xlim([-10 50])

.

Sanufee 2024-2-22

编辑：Sanufee 2024-2-22

在 MATLAB Online 中打开

Now, I want to plot the 2 surface in same fig for comparison. So i cant change different color of surface (Colormap). Can you help please?

I tried it by change in fiqure manualy but It's not color that i want. I want to change to different colormap to see the different clearly.

Then that's better if I can show the legend also.

Thankss

Here's my code

Pre = readtable('Premonsoon.xlsx')
Post = readtable('Postmonsoon.xlsx')
E1 = Pre.E(3:end);
N1 = Pre.N(3:end);
Elev1 = Pre.Elev(3:end);
E2 = Post.E(3:end);
N2 = Post.N(3:end);
Elev2 = Post.Elevation(3:end);
F1 = scatteredInterpolant(E1,N1,Elev1,'natural')
[Em1,Nm1] = ndgrid(sort(E1), sort(N1));
Elevm1 = F1(Em1,Nm1);
F2 = scatteredInterpolant(E2,N2,Elev2,'natural')
[Em2,Nm2] = ndgrid(sort(E2), sort(N2));
Elevm2 = F2(Em2,Nm2);
% Plot
hold on
% Pre-monsoon surface
surfc(Em1, Nm1, Elevm1, 'EdgeColor','interp')
colormap("turbo") % Change colormap for pre-monsoon surface
title('Pre-monsoon')
view(140, 30)
hold off
hold on
% Post-monsoon surface
surfc(Em2, Nm2, Elevm2, 'EdgeColor','interp')
colormap("parula") % Change colormap for post-monsoon surface
title('Post-monsoon')
view(140, 30)
hold off
% Legend
h = colorbar;
ylabel(h,'Elevation')
hold on
[X,Y] = meshgrid(min(Em1(:)):max(Em1(:)), min(Nm1(:)):max(Nm1(:)));
Z = -2 * ones(size(X)); % Set Z values to -2
surf(X, Y, Z, 'FaceAlpha', 0.3, 'FaceColor', '#4DBEEE', 'EdgeColor', 'none');
hold off

Star Strider 2024-2-22

在 MATLAB Online 中打开

That is apparently not possible. The colormap has to apply to an axis, not a specific surf plot.

One option is to set each surf plot to a different colour, howefvr the best approach may be to interpolate the larger values (‘Premonsoon’) to the size of the smaller values (‘Postmonsoon’) so they are both the same size (this requires a third scatteredInterpolant call) and then subtract them, as I did at the end (using a slightly different approach to plotting it) —

Pre = readtable('Premonsoon.xlsx')

Pre = 606×4 table

No E N Elev __ ______ ______ _______ 1 2.0028 0.0003 -0.247 2 4.4752 0.1113 -0.3718 3 6.9702 0.1733 -0.4969 4 9.471 0.2356 -0.6752 5 11.946 0.2971 -0.9908 6 14.44 0.359 -1.1287 7 16.893 0.4018 -1.3546 8 19.425 0.4804 -1.4187 9 21.903 0.5552 -1.4402 10 24.411 0.6455 -1.5291 11 26.882 0.673 -1.7164 12 29.372 0.7274 -1.9758 13 31.862 0.7887 -2.2727 14 34.376 0.8114 -2.646 15 36.864 0.836 -2.9736 16 1.7236 9.9816 -0.2711

Post = readtable('Postmonsoon.xlsx')

Post = 531×6 table

Var1 Var2 Var3 Var4 Var5 Var6 ____ ______ _______ _______ ____________________ ____ NaN 0 0 0 {'10 Feb, Saturday'} NaN 1 4.4896 -0.0456 -0.3672 {0×0 char } NaN 2 4.4866 -0.0455 -0.3671 {0×0 char } NaN 3 6.9964 -0.0443 -0.4542 {0×0 char } NaN 4 9.4778 -0.0158 -0.7015 {0×0 char } NaN 5 12.005 -0.0133 -0.8593 {0×0 char } NaN 6 14.603 -0.0125 -1.0918 {0×0 char } NaN 7 16.983 -0.0117 -1.4529 {0×0 char } NaN 8 19.3 -0.011 -1.6625 {0×0 char } NaN 9 21.899 -0.0101 -1.7106 {0×0 char } NaN 10 23.91 -0.0145 -1.6612 {0×0 char } NaN 11 26.93 -0.0143 -1.8578 {0×0 char } NaN 12 29.437 -0.014 -2.123 {0×0 char } NaN 13 34.331 -0.0138 -2.564 {0×0 char } NaN 14 36.5 -0.049 -2.8317 {0×0 char } NaN 15 2.1108 9.9106 -0.2643 {0×0 char } NaN

E1 = Pre.E(3:end);

N1 = Pre.N(3:end);

Elev1 = Pre.Elev(3:end);

% E2 = Post.E(3:end);

% N2 = Post.N(3:end);

% Elev2 = Post.Elevation(3:end);

E2 = Post{3:end,2}; % No Appropriate Variable Names In 'Post'

N2 = Post{3:end,3};

Elev2 = Post{3:end,4};

F1 = scatteredInterpolant(E1,N1,Elev1,'natural');

[Em1,Nm1] = ndgrid(sort(E1), sort(N1));

Elevm1 = F1(Em1,Nm1);

F2 = scatteredInterpolant(E2,N2,Elev2,'natural');

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

[Em2,Nm2] = ndgrid(sort(E2), sort(N2));

Elevm2 = F2(Em2,Nm2);

% Plot

hold on

% Pre-monsoon surface

surfc(Em1, Nm1, Elevm1, 'EdgeColor','interp', 'FaceColor','r')

% colormap(turbo) % Change colormap for pre-monsoon surface

title('Pre-monsoon')

view(140, 30)

hold off

hold on

% Post-monsoon surface

surfc(Em2, Nm2, Elevm2, 'EdgeColor','interp', 'FaceColor','g');

% colormap(parula) % Change colormap for post-monsoon surface

title('Post-monsoon')

view(140, 30)

hold off

% Legend

h = colorbar;

ylabel(h,'Elevation')

hold on

[X,Y] = meshgrid(min(Em1(:)):max(Em1(:)), min(Nm1(:)):max(Nm1(:)));

Z = -2 * ones(size(X)); % Set Z values to -2

surf(X, Y, Z, 'FaceAlpha', 0.3, 'FaceColor', '#4DBEEE', 'EdgeColor', 'none');

hold off

size(Elevm1)

ans = 1×2

604 604

size(Elevm2)

ans = 1×2

529 529

F3 = scatteredInterpolant(Em1(:),Nm1(:),Elevm1(:));

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

Elevm1i = F3(Em2,Nm2); % 'Elev1m' Interpolated To The dimensions of 'Elev2m'

figure

hsl = surfl(Em2,Nm2,Elevm2-Elevm1i, 'light')

hsl =

1×2 graphics array: Surface Light

hsl(1).EdgeColor = 'none';

% surfc(Em2,Nm2,Elevm2-Elevm1i, 'EdgeColor','interp')

colormap(turbo)

Ax = gca;

Ax.YDir = 'reverse';

Ax.XDir = 'reverse';

title('Post-Pre Difference: Elevm2-Elevm1 (Uses Interpolated Elevm1)')

zlim([-2 1])

figure

contourf(Em2,Nm2,Elevm2-Elevm1i, (-2:0.25:1), 'ShowText','on')

Ax = gca;

colormap(turbo)

colorbar

Ax.YDir = 'reverse';

Ax.XDir = 'reverse';

title('Post-Pre Difference: Elevm2-Elevm1 (Uses Interpolated Elevm1)')

.

Sanufee 2024-2-24

It's so useful for me! Thanks so much. maybe I could give my graduate certificate to you. 🥹

Star Strider 2024-2-24

As always, my pleasure!

(I already have a few of my own! I am happy to be able to help you.)

请先登录，再进行评论。

Answer 2

John D'Errico 2024-2-20

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2084483-how-can-i-calculate-the-area-under-two-curves-that-intersect#answer_1412873

编辑：John D'Errico 2024-2-21

在 MATLAB Online 中打开

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.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Sanufee 2024-2-21

Thanks!! I will try it

请先登录，再进行评论。

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

3 个评论
显示 1更早的评论隐藏 1更早的评论

采纳的回答

27 个评论
显示 25更早的评论隐藏 25更早的评论

更多回答（1 个）

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

另请参阅

类别

标签

产品

Community Treasure Hunt

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

3 个评论 显示 1更早的评论隐藏 1更早的评论

采纳的回答

27 个评论 显示 25更早的评论隐藏 25更早的评论

更多回答（1 个）

1 个评论 显示 -1更早的评论隐藏 -1更早的评论

另请参阅

类别

标签

产品

Community Treasure Hunt

3 个评论
显示 1更早的评论隐藏 1更早的评论

27 个评论
显示 25更早的评论隐藏 25更早的评论

1 个评论
显示 -1更早的评论隐藏 -1更早的评论