MATLAB 中具有任意数量交叉点的曲线之间的区域
Area between curves with arbitary number of intersections in MATLAB
我想计算 matlab 中具有不同 x 范围的两条曲线之间的面积。有两种情况。当一条曲线(图中的红色)高于另一条曲线时,我希望该区域为正,当该曲线低于蓝色曲线时,该区域为负。我意识到这通常是一个简单的问题,通常会找到两条曲线之间的所有交点,但在我的例子中,交点的数量可能非常大。我的困难是,根据所讨论的数据集,两条曲线将相交任意次数或根本不相交。
示例:
n=4
x = n*pi*rand(100,1);
x=[0;x;pi];
x=sort(x);
y=.5+exp(-.5*x).*sin(x);
x2 = n*pi*rand(1000,1);
x2=[0;x2;pi];
x2=sort(x2);
y2=.5+exp(-.5*x2).*0.6.*sin(x2/1.2);
figure,hold on
plot(x,y,'r-')
plot(x2,y2,'b-')
Area1 = trapz(x, y);
Area2 = trapz(x2, y2);
有什么建议?
您计算交叉点之间面积的方法是合理的。您可以使用 this library to calculate intersection points. Direct Link 压缩文件。
假设您有交叉点,您可以使用 trapz
函数遍历该交叉点以减去两个间隔之间的面积。如果曲线 1 高于曲线 2,唯一棘手的部分就是面积符号。为此,您可以设置一个阈值 epsilon
并在 x2 - epsilon
处计算 y1
和 y2
其中 (x1, x2)
是 X
交点。 y1 > y2
意味着你必须做 Area 1 - Area2
.
综上所述,代码如下所示:
%Find all intersections
[xIntersect,yIntersect] = intersections(x,y,x2,y2,1);
%Number of intersections
numIntersections = size(xIntersect, 1);
%Threshold to calculate sign of area
epsilon = 0.1;
curveArea = 0; %Initial value
lastX1Index = 1; % End X1 Index processes in last iteration
lastX2Index = 1; % End X2 Index processes in last iteration
%Loop over intersections processing pair of points hence until
%numIntersections - 1
for i = 1:numIntersections-1
% startX = xIntersect(i); %Start of interval
% startY = yIntersect(i);
endX = xIntersect(i+1); % Next Intersection point
% endY = yIntersect(i+1);
xMinus = endX - epsilon; % get a X value before intersection
%Sign of area is positive if y1Minus > y2Minus
y1Minus = .5+exp(-.5*xMinus).*sin(xMinus); %Convert this to function1
y2Minus = .5+exp(-.5*xMinus).*0.6.*sin(xMinus/1.2); %Convert this to function 2
x1Index = find(x < xIntersect(i+1), 1, 'last' ); % Find last X1 element before intersection
x2Index = find(x2 < xIntersect(i+1), 1, 'last' ); % Find last X2 element before intersection
%Calculate area for interval
Area1 = trapz(x(lastX1Index:x1Index), y(lastX1Index:x1Index));
Area2 = trapz(x2(lastX2Index:x2Index), y2(lastX2Index:x2Index));
%Store last X1, X2 Index for next run
lastX1Index = x1Index+1;
lastX2Index = x2Index+1;
%Save curveArea
if y1Minus > y2Minus
curveArea = curveArea + (Area1 - Area2);
else
curveArea = curveArea + (Area2 - Area1);
end
end
fprintf('curveArea= %f \n',curveArea);
我想计算 matlab 中具有不同 x 范围的两条曲线之间的面积。有两种情况。当一条曲线(图中的红色)高于另一条曲线时,我希望该区域为正,当该曲线低于蓝色曲线时,该区域为负。我意识到这通常是一个简单的问题,通常会找到两条曲线之间的所有交点,但在我的例子中,交点的数量可能非常大。我的困难是,根据所讨论的数据集,两条曲线将相交任意次数或根本不相交。
示例:
n=4
x = n*pi*rand(100,1);
x=[0;x;pi];
x=sort(x);
y=.5+exp(-.5*x).*sin(x);
x2 = n*pi*rand(1000,1);
x2=[0;x2;pi];
x2=sort(x2);
y2=.5+exp(-.5*x2).*0.6.*sin(x2/1.2);
figure,hold on
plot(x,y,'r-')
plot(x2,y2,'b-')
Area1 = trapz(x, y);
Area2 = trapz(x2, y2);
有什么建议
您计算交叉点之间面积的方法是合理的。您可以使用 this library to calculate intersection points. Direct Link 压缩文件。
假设您有交叉点,您可以使用 trapz
函数遍历该交叉点以减去两个间隔之间的面积。如果曲线 1 高于曲线 2,唯一棘手的部分就是面积符号。为此,您可以设置一个阈值 epsilon
并在 x2 - epsilon
处计算 y1
和 y2
其中 (x1, x2)
是 X
交点。 y1 > y2
意味着你必须做 Area 1 - Area2
.
综上所述,代码如下所示:
%Find all intersections
[xIntersect,yIntersect] = intersections(x,y,x2,y2,1);
%Number of intersections
numIntersections = size(xIntersect, 1);
%Threshold to calculate sign of area
epsilon = 0.1;
curveArea = 0; %Initial value
lastX1Index = 1; % End X1 Index processes in last iteration
lastX2Index = 1; % End X2 Index processes in last iteration
%Loop over intersections processing pair of points hence until
%numIntersections - 1
for i = 1:numIntersections-1
% startX = xIntersect(i); %Start of interval
% startY = yIntersect(i);
endX = xIntersect(i+1); % Next Intersection point
% endY = yIntersect(i+1);
xMinus = endX - epsilon; % get a X value before intersection
%Sign of area is positive if y1Minus > y2Minus
y1Minus = .5+exp(-.5*xMinus).*sin(xMinus); %Convert this to function1
y2Minus = .5+exp(-.5*xMinus).*0.6.*sin(xMinus/1.2); %Convert this to function 2
x1Index = find(x < xIntersect(i+1), 1, 'last' ); % Find last X1 element before intersection
x2Index = find(x2 < xIntersect(i+1), 1, 'last' ); % Find last X2 element before intersection
%Calculate area for interval
Area1 = trapz(x(lastX1Index:x1Index), y(lastX1Index:x1Index));
Area2 = trapz(x2(lastX2Index:x2Index), y2(lastX2Index:x2Index));
%Store last X1, X2 Index for next run
lastX1Index = x1Index+1;
lastX2Index = x2Index+1;
%Save curveArea
if y1Minus > y2Minus
curveArea = curveArea + (Area1 - Area2);
else
curveArea = curveArea + (Area2 - Area1);
end
end
fprintf('curveArea= %f \n',curveArea);