How to calculate area between curve and horizontal line?

14 vues (au cours des 30 derniers jours)
Hamzah Mahmood
Hamzah Mahmood le 21 Août 2020
Commenté : esat gulhan le 23 Août 2020
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;

Réponse acceptée

esat gulhan
esat gulhan le 21 Août 2020
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 commentaires
Hamzah Mahmood
Hamzah Mahmood le 22 Août 2020
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 le 23 Août 2020
Sorry, i dont understand your question. If you ask this question in new stage, other users can see and answer...

Connectez-vous pour commenter.

Plus de réponses (2)

esat gulhan
esat gulhan le 21 Août 2020
Modifié(e) : esat gulhan le 21 Août 2020
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 commentaire
Hamzah Mahmood
Hamzah Mahmood le 21 Août 2020
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.

Connectez-vous pour commenter.


esat gulhan
esat gulhan le 22 Août 2020
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

Catégories

En savoir plus sur Graph and Network Algorithms dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by