How to calculate area between curve and horizontal line?

18 ビュー (過去 30 日間)
Hamzah Mahmood
Hamzah Mahmood 2020 年 8 月 21 日
コメント済み: esat gulhan 2020 年 8 月 23 日
Apologies if this has already been answered in another form elsewhere on the forums, however I couldn't find any answers that helped me directly with my code. I currently have code that maps the profile and area of a channel. I wish to have a horizontal y value across my channel width to simulate the fluctuating water level, the example line is shown in the image as the green line. So for example I would want to calculate and shade the area underneath the greenline only.
With the function would I be able to loop the same function for a range of values for the horizontal values, to provide an output I could then plot? A vector form would be ideal if possible.
I've already tried to use the area and fill commands unsuccesfully, or maybe I was not using the correct format for them.
I've also attached the relevant code below, any help would be much appreciated.
x = [ 0 1 2 3 4 5 6 7 8 9 10 11 12];
y = -[0 -0.5 -0.8 -0.8 -1 -1.1 -1.2 -1.2 -1.4 -1.2 -1.1 -1 0]; % y values in image were produced by rand function
subplot(2, 1, 1);
area(x,-y)
grid on;
hold on;
for k = 1 : length(x)
xline(x(k), 'Color', 'k');
end
yt = yticks;
for k = 1 : length(yt)
yline(yt(k), 'Color', 'k');
end
plot (x, -y, 'r.-', 'LineWidth', 3, 'MarkerSize', 30)
hold on
yline(-0.5,'g',11);
xlabel('Chainage (m)', 'FontSize', 11);
ylabel('Depth (m)', 'FontSize', 11);
An = (trapz(x,y)); % Overall area.
% Compute areas within each pair of x locations.
midpoints = (y(1:end-1) + y(2:end))/2;
deltax = diff(x);
Ani = deltax .* midpoints;
subplot(2, 1, 2);
bar(x(1:end-1)-x(1)+ 0.5, Ani);
title('Areas of Slices', 'FontSize', 14);
xlabel('x', 'FontSize', 11);
ylabel('Area (m^2)', 'FontSize', 11);
grid on;

採用された回答

esat gulhan
esat gulhan 2020 年 8 月 21 日
clc;clear;
x = [ 0 1 2 3 4 5 6 7 8 9 10 11 12];
y = [0 -0.5 -0.8 -0.8 -1 -1.1 -1.2 -1.2 -1.4 -1.2 -1.1 -1 0];
h0= -0.5; %the axis in the graph
h=ones(size(x))*h0;
Int=pchip(x,y);
xx=linspace(x(1),x(end),500);
yy=y-h0;
yyy=pchip(x,yy,xx);
gt=yyy<0;
theArea =-trapz(xx(gt), (yyy(gt)))
hold off
area(xx, min(ppval(Int,xx), h(1)), h(1), 'EdgeColor', 'none', 'FaceColor', 'b'),grid,
hold on
plot(xx,ppval(Int,xx),'k',xx,h0,'Linewidth',2)
set(gca, 'XLim', [x(1) x(end)], 'YLim', [min(y) 0]);
title(strcat('BlueArea=', num2str(theArea)))
xlabel('Chainage (m)', 'FontSize', 11);
ylabel('Depth (m)', 'FontSize', 11);
You can change h0 which is horizontal line. when h0=-1
when h0=-0.5
when h0=0
I RESHAPED YOUR X AND Y DATAS 500 INTERVALS WITH PCHIP INTERPOLATION, IF NOT THE DATA IS NOT ENOUGH FOR GOOD SOLUTION.
  4 件のコメント
Hamzah Mahmood
Hamzah Mahmood 2020 年 8 月 22 日
I'll look into that, cheers.
Another thing has come up, how would I be able to then keep these values as an array, im trying this however it isn't working.
clc;
clear;
x = [ 0 1 2 3 4 5 6 7 8 9 10 11 12];
y = [0 -0.5 -0.8 -0.8 -1 -1.1 -1.2 -1.2 -1.4 -1.2 -1.1 -1 0];
for stage= [-2 -1.8 -1.6 -1.4 -1.2 -1 -0.8 -0.6 -0.4 -0.2 0] %the axis in the graph
stage_dim= length(h)
h=ones(size(x))*stage;
Int=pchip(x,y);
xx=linspace(x(1),x(end),500);
yy=y-stage;
yyy=pchip(x,yy,xx);
gt=yyy<0;
Area =zeros(stage_dim);
for i = 1:stage_dim
Area(i) =-trapz(xx(gt), (yyy(gt)))
end
hold off
area(xx, min(ppval(Int,xx), h(1)), h(1), 'EdgeColor', 'none', 'FaceColor', 'b'),grid,
hold on
plot(xx,ppval(Int,xx),'k',xx,stage,'Linewidth',2)
set(gca, 'XLim', [x(1) x(end)], 'YLim', [min(y) 0]);
title(strcat('Area=', num2str(Area)))
xlabel('Chainage (m)', 'FontSize', 11);
ylabel('Depth (m)', 'FontSize', 11);
end
esat gulhan
esat gulhan 2020 年 8 月 23 日
Sorry, i dont understand your question. If you ask this question in new stage, other users can see and answer...

サインインしてコメントする。

その他の回答 (2 件)

esat gulhan
esat gulhan 2020 年 8 月 21 日
編集済み: esat gulhan 2020 年 8 月 21 日
x = [ 0 1 2 3 4 5 6 7 8 9 10 11 12]
y = -[0 -0.5 -0.8 -0.8 -1 -1.1 -1.2 -1.2 -1.4 -1.2 -1.1 -1 0]
Int=pchip(x,y)
plot(x,ppval(Int,x))
MIItemp = (fnint(Int))
Area = fnval(MIItemp,12)%Area between 0 12
If you want to find area between 0-5 , Area = fnval(MIItemp,5)
There are other ways but this is the easiest way .
Note:If your matlab is latetest version you can write makima(x,y),instead of phcip(x,y). If your lines are strict use makima
  1 件のコメント
Hamzah Mahmood
Hamzah Mahmood 2020 年 8 月 21 日
Hi Esat,
I think you may have misinterpered which axes I'm looking to cacluate the area from. The -0.5 value is a line that shifts up and down the y-axis rather than along the x-axis like so (the shade blue area).
Would your code allow me to calculate the area n the verticle directon if moodified slightly( I tried for a box with y = ones(1,13) however this didnt work). I'm relatively new to Matlab, so your code has helped introduce me to some new commands I hadn't come across before too.
Thanks.

サインインしてコメントする。


esat gulhan
esat gulhan 2020 年 8 月 22 日
Did you see my second answer. I think it helps.
If you add pause to the code it will seem like an animation.
clc;clear;
x = [ 0 1 2 3 4 5 6 7 8 9 10 11 12];
y = [0 -0.5 -0.8 -0.8 -1 -1.1 -1.2 -1.2 -1.4 -1.2 -1.1 -1 0];
for h0= 0:-0.04:-2; %the axis in the graph
h=ones(size(x))*h0;Int=pchip(x,y);xx=linspace(x(1),x(end),500);
yy=y-h0;yyy=pchip(x,yy,xx);gt=yyy<0;theArea =-trapz(xx(gt), (yyy(gt)))
hold off
area(xx, min(ppval(Int,xx), h(1)), h(1), 'EdgeColor', 'none', 'FaceColor', 'b'),grid,
hold on
plot(xx,ppval(Int,xx),'k',xx,h0,'Linewidth',2)
set(gca, 'XLim', [x(1) x(end)], 'YLim', [min(y) 0]);
title(strcat('BlueArea=', num2str(theArea)))
xlabel('Chainage (m)', 'FontSize', 11);ylabel('Depth (m)', 'FontSize', 11);pause(0.00001)
end

カテゴリ

Help Center および File ExchangeVisual Exploration についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by