How can I calculate the area between two curves?
140 次查看(过去 30 天)
显示 更早的评论
Hello everyone. I am trying to calculate the area between two curves. How can I do this having only the x and y coordinates for each curve? Thanks in advance.
1 个评论
Muhammad Usman
2019-2-2
Please Let me know if this is of any help...
area = 0;
for i=1:size(x1,2)-1
wid = ( sqrt((x1(i+1) - x1(i) )^2 + (y1(i+1) - y1(i) )^2) + sqrt((x2(i+1) - x2(i))^2 + (y2(i+1) - y2(i))^2) )/2;
height = ( sqrt((y2(i) - y1(i))^2 + (x2(i) - x1(i))^2) + sqrt( ( y2(i+1) - y1(i+1) )^2 + (x2(i+1) - x1(i+1))^2 ))/2
area(i) = wid*height;
end
sum(area)
采纳的回答
Star Strider
2016-9-9
Use the trapz (or cumtrapz if their x-coordinates are the same and you want to subtact them element-wise) function to integrate each one, and then subtract the integral of one (calculated by trapz) from the other.
7 个评论
Vanika Sharma
2021-5-20
I know the value of area of intersection and data points for one of the curve. Is there any way to get the second curve or points for that?
更多回答(2 个)
John D'Errico
2016-9-9
编辑:John D'Errico
2016-9-9
If the two curves are defined over the same support, so identical x values for each curve, then the simple answer is to use trapz on the difference if the y-values. So...
trapz(x,y1-y2)
This is equivalent to computing trapz for each curve, then subtracting the two results, since integration is a linear operator. So, IF they are defined on the same interval in x...
trapz(x1,y1) - trapz(x2-y2)
If the two curves cross at one or more points, then, assuming you want to compute the area of each section as being a positive area, summing them up, then you COULD use trapz on the absolute difference.
trapz(x,abs(y1-y2))
the problem is this fails to be fully accurate because at the intersection point, if that falls between two of your sampled data points, the abs introduces an error.
So a better scheme would be to identify the crossing points using linear interpolation. Insert new values into each curve at those locations, then trapz will have higher accuracy. Doug Schwarz has a nice tool to locate those intersections on the file exchange, here:
If the two curves are not defined over the same support, then you will best use some interpolation method to put the two curves onto the same support, then apply one of the schemes I describe above.
4 个评论
John D'Errico
2016-9-10
编辑:John D'Errico
2016-9-10
Oh well, I knew that by not being complete there in show how to do the interpolation, I'd be forced to add that part anyway. :)
It is easy enough. intersections does the interpolation for you.
I'll keep the points with the same support foreach function, so x will be used for both curves.
x = linspace(0,2*pi,10);
y1 = sin(x);
y2 = cos(x);
plot(x1,y1,'r-o',x2,y2,'b-o')
[xnew,ynew] = intersections(x,y1,x,y2)
Warning: NARGCHK will be removed in a future release. Use NARGINCHK or NARGOUTCHK instead.
> In intersections (line 91)
xnew =
0.790220733000934
3.92363040454986
ynew =
0.687902741912839
-0.667001165195401
(It looks like Doug needs to fix intersections to remove that warning message.) Those are the new points we need to add to each curve.
x = [x,xnew'];
y1 = [y1,ynew'];
y2 = [y2,ynew'];
xold = x; % save for a later plot
[x,tags] = sort(x);
y1 = y1(tags);
y2 = y2(tags);
plot(x,y1,'r-o',x,y2,'b-o')
Now, I happen to know that the integral of abs(sin(x)-cos(x)), over the limits [0,2*pi] is 4*sqrt(2).
First, compare the two absolute differences.
plot(xold,abs(sin(xold) - cos(xold)),'r--',x,abs(y1-y2),'b:')
As you can see, just taking the abs of the difference WITHOUT inserting the points of intersection of the curves, I will get a poor result.
trapz(x,abs(y1-y2))
ans =
5.43086292018145
Test that result:
syms u
int(abs(sin(u) - cos(u)),u,[0,2*pi])
ans =
4*2^(1/2)
vpa(ans)
ans =
5.6568542494923801952067548968388
IF the two functions did not live on the same support, then use interp1 to put them on the same support.
Zeeshan Salam
2019-11-4
Can anyon solve this question?
How i find area between two curves?
clear;
clc;
%define x for both curves
x_ref = 0:0.01:1;
x = 0:0.01:1;
%define the reference curve and the other curve
c_ref = sin(x_ref*2*pi);
y = 0.5+0.5*sin(x*2*pi);
%plot both curves
figure(1)
plot (x_ref,c_ref,'r');
hold on
plot(x,y,'b');
axis([min(x_ref)-1,max(x_ref)+1,min(c_ref)-1,max((c_ref))+1]);
legend('c_{ref}','y')
%the relative x&y positions
x_rel = abs(x-min(x)/(max(x)-min(x)));
y_rel = abs((y-min(y))/(max(y)-min(y)));
%% Parameter before shifting
%corner points (in p1 matrix)
p1(:,1) = [min(x); max(x); min(x); max(x)]; % x colomn
p1(:,2) = [min(y); min(y); max(y); max(y)]; % y Spalte
for i = 1:size(p1,1)
text(p1(i,1),p1(i,2),[' p',num2str(i)])
end
%plot corner Points as circles
plot(p1(:,1),p1(:,2),'ko')
% hold off
%finding Xup, Xlow, Yright, Yleft
x_up_temp = (p1(4,1)-p1(3,1))*x_rel;
x_up1 = x_up_temp;
x_low_temp = (p1(2,1)-p1(1,1))*x_rel;
x_low1 = x_low_temp;
y_right_temp = (p1(4,2)-p1(2,2))*y_rel;
y_right1 = y_right_temp;
y_left_temp = (p1(3,2)-p1(1,2))*y_rel;
y_left1 = y_left_temp;
%% Parameter after shifting the first corner
%p1 = before shifting
%p2 = first point shifting
for i=1:3
p = [i-2 0; 0 0; 0 0; 0 0];
for j=1:3
p2 = p1 + p + [0 j-2; 0 0; 0 0; 0 0];
%finding new Xup, Xlow, Yright, Yleft
f1 = p2(4,1)-p2(3,1);
x_up_temp = f1*x_rel;
x_up2 = x_up_temp;
f2 = p2(2,1)-p2(1,1);
x_low_temp = f2*x_rel;
x_low2 = x_low_temp;
f3 = p2(4,2)-p2(2,2);
y_right_temp = f3*y_rel;
y_right2 = y_right_temp;
f4 = p2(3,2)-p2(1,2);
y_left_temp = f4*y_rel;
y_left2 = y_left_temp;
%calculating delta_x & delta_y
x_up_d = x_up2 - x_up1;
x_low_d = x_low2 - x_low1;
delta_x = x_up_d - x_low_d;
delta_y = (y_right2 - y_right1) - (y_left2 - y_left1);
%% the new x&y values are
xn = x + delta_x .* y_rel;
yn = y + delta_y .* x_rel;
cn = [xn; yn];% should drow this as a new curve
% calculate the Area
%P2(i;j) = flaeche();
end
end
%find the minimum area and apply the shifting
%a = find(flaeche == min(min(flaeche)));
%P2{a}
%% p3, p4, p5 would be the same
1 个评论
Steven Lord
2019-11-4
Rather than posting this as an answer on a more than three year old question, you should ask this as a separate question of its own.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Surface and Mesh Plots 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!