How to evaluate area between 2 curve

1 view (last 30 days)
Here is my code.I try to evaluate surface area between of these 2 curves, Can anyone help me with this issue please.
unitcellsize = 1; % Cell size [mm]
length_x = 1; % Samplesize in x direction [mm]
length_y = 1; % Samplesize in y direction [mm]
iso = 0.5;
% Define the function
syms x y z
f = @(x,y,z) cos((2.*pi.*x)./unitcellsize)+cos((2.*pi.*y)./unitcellsize)+cos((2.*pi.*z)./unitcellsize);
% Set the z-value for the slice
z = 0; % Adjust this value to the desired z-slice
N = 100;
a = linspace(-length_x/2,length_x/2,N); % Adjust the resolution as needed
b = linspace(-length_y/2,length_y/2,N); % Adjust the resolution as needed
c = linspace(-0,0,N);
[X,Y] = meshgrid(a,b);
Z = meshgrid(c);
% Evaluate the function at the current slice
V = f(X,Y,Z);
% Create a contour plot to find the enclosed area
figure
contour(a,b,V,[-iso iso]);
[contour_matrix,h] = contour(a,b,V,[-iso,iso]); % Adjust iso values as needed
axis equal
xlabel('X')
ylabel('Y')
xticks(-0.5:0.1:0.5)
yticks(-0.5:0.1:0.5)
% Extract the contour coordinates
contour_x = contour_matrix(1,2:end);
contour_y = contour_matrix(2,2:end);
plot(contour_x,contour_y)
% Calculate the area of the enclosed region using the shoelace formula
enclosed_area = polyarea(contour_x,contour_y);
% Find total area of slice from sample dimensions
total_slice_area = length_x.*length_y;
% Calculate the area fraction for the slice
area_fraction = enclosed_area/total_slice_area;
Area = trapz(abs(contour_x),abs(contour_y))
Area = 2.8742
fprintf('The total surface area is %f\n',enclosed_area);
The total surface area is 14.621238
fprintf('The area fraction is %f\n',area_fraction);
The area fraction is 14.621238

Answers (1)

Harsha Vardhan
Harsha Vardhan on 22 Sep 2023
Hi,
I understand that you want to calculate the area between 2 contour lines.
‘Contour’ matrix returned by the ‘contour’ function defines contour lines in its columns.
In the following line of your code,
[contour_matrix,h] = contour(a,b,V,[-iso,iso]);
The ‘contour_matrix’ returned by the ‘contour’ function is of the 2x416 dimension.
Contour matrix, returned as a two-row matrix is of the following form
The columns of the matrix define the contour lines. Each contour line starts with a column containing Z and N values where:
  • Zi — The height of the ith contour line
  • Ni — The number of vertices in the ith contour line
  • (xij, yij) — The coordinates of the vertices for the ith contour line, where j ranges from 1 to Ni
So, the Zi, Ni column pair may appear not only in the first column but also in other columns.
Kindly refer to the following link for more information on ‘contour’ function - https://www.mathworks.com/help/matlab/ref/contour.html#mw_27cd6c94-d861-4e0a-837c-0a19f2574186
The contour_matrix returned by the following line of your code is as follows:
[contour_matrix,h] = contour(a,b,V,[-iso,iso]);
So, to extract the actual coordinates of the contour lines, all the columns like column 1, column 36, column 71 etc., which contain the height and number of vertices of ith contour line must be excluded.
Therefore, to extract the contour coordinates I would recommend changing the following lines of code:
% Extract the contour coordinates
contour_x = contour_matrix(1,2:end);
contour_y = contour_matrix(2,2:end);
I provided correct code below to replace the above lines of code.
% Extract the contour coordinates
currentColumnWithHeight=1;
contour_x = [];
contour_y = [];
while currentColumnWithHeight <= length(contour_matrix)
nextColumnWithHeight=currentColumnWithHeight+contour_matrix(2,currentColumnWithHeight)+1;
x = contour_matrix(1,currentColumnWithHeight+1:nextColumnWithHeight-1);
y = contour_matrix(2,currentColumnWithHeight+1:nextColumnWithHeight-1);
contour_x = [contour_x ,x];
contour_y = [contour_y ,y];
currentColumnWithHeight = nextColumnWithHeight;
end
Hope it helps in resolving your query!

Categories

Find more on Contour Plots in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!