Clear Filters
Clear Filters

Find the number of regions and calculate the area of each region

5 views (last 30 days)
I used Gauss Random Field equation to create the bi-phase region, I want to find the number of regions separated by the equation and calculate the area of every single region of them. can anyone help? I attached the code and the image of the bi-phase field.
clear; format short e;
rng('shuffle');
%% Calculate
export = 1; % Whether or not to save data for ABAQUS
N = 1000; % number of waves to sum
mesh_size = 301; % coarseness of pattern
max_phi = 2*pi; % maximum phase shift
pct_stiff = 0.57; % fraction of material representing hard phase
lambda = 0.320; % characteristic wavelength (mm)
sample_size = 3.5 ; % characteristic length of volume (mm)
% Generate physical coordinates
x = linspace(-sample_size/2,sample_size/2,mesh_size);
y = linspace(-sample_size/2,sample_size/2,mesh_size);
% Create 3D coordinate matrices
[X, Y] = meshgrid(x, y);
phase = zeros(size(X));
for n=1:N
% Generate random phase shift
rand_phase = max_phi * rand(1,1);
% Generate random vector
rand_theta = 2 * pi * rand(1,1);
rand_vec = [cos(rand_theta) sin(rand_theta)];
% Dot direction vector with coordinate matrices
vec = X*rand_vec(1) + Y*rand_vec(2);
% Update the continuos phase field with new wave function
phase = phase + sqrt(2/N)*cos(2*pi/lambda*vec + rand_phase);
end
%% Threshold
% Allocate space for new binary matrix
biphase = phase;
% Find average of all phase data
val_avg = mean(biphase,'all');
% Find standard deviation of all phase data
val_std = std(biphase,0,'all');
% Find the value at which the population splits into desired proportions
split = val_avg - sqrt(2) * erfinv((pct_stiff-0.5) * 2) * val_std;
% Apply threshold
biphase(phase>=split) = 1;
biphase(phase<split) = 0;
contourf(X,-Y,biphase)

Answers (1)

Hassaan
Hassaan on 15 Apr 2024
clear;
format short e;
rng('shuffle');
%% Parameters
export = 1; % Whether or not to save data for ABAQUS
N = 1000; % number of waves to sum
mesh_size = 301; % coarseness of pattern
max_phi = 2 * pi; % maximum phase shift
pct_stiff = 0.57; % fraction of material representing hard phase
lambda = 0.320; % characteristic wavelength (mm)
sample_size = 3.5; % characteristic length of volume (mm)
%% Generate physical coordinates
x = linspace(-sample_size / 2, sample_size / 2, mesh_size);
y = linspace(-sample_size / 2, sample_size / 2, mesh_size);
[X, Y] = meshgrid(x, y);
phase = zeros(size(X));
%% Generate Gaussian Random Field
for n = 1:N
rand_phase = max_phi * rand(1,1);
rand_theta = 2 * pi * rand(1,1);
rand_vec = [cos(rand_theta) sin(rand_theta)];
vec = X * rand_vec(1) + Y * rand_vec(2);
phase = phase + sqrt(2/N) * cos(2 * pi / lambda * vec + rand_phase);
end
%% Threshold to create a bi-phase region
val_avg = mean(phase, 'all');
val_std = std(phase, 0, 'all');
split = val_avg - sqrt(2) * erfinv((pct_stiff - 0.5) * 2) * val_std;
biphase = double(phase >= split);
%% Find connected components
CC = bwconncomp(biphase);
numRegions = CC.NumObjects;
disp(['Number of Regions: ', num2str(numRegions)]);
%% Calculate the area of each region
pixelArea = (sample_size / mesh_size)^2; % Area represented by each pixel
areas = zeros(numRegions, 1);
for k = 1:numRegions
numPixels = numel(CC.PixelIdxList{k});
areas(k) = numPixels * pixelArea;
end
disp('Area of each region (mm^2):');
disp(areas);
%% Visualization
labelMatrix = labelmatrix(CC);
figure;
contourf(X, -Y, labelMatrix); % Using negative Y for correct orientation
colormap(jet(numRegions)); % Apply a colormap with a distinct color per region
colorbar;
title('Contour Plot of Regions');
  1 Comment
Huan
Huan on 15 Apr 2024
Thank you very much for your answer. however, it seems it didn't find all the regions separated by the equation.

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!