mean surface distance and residual mean square distance of 2 borders

hi every body. I have 2 borders of 2 surfaces called S1 and S2. I need to compute the surface distance and after that the mean surface distance and residual mean square distance from that.

A complete description of the concept and the code in python is presented in https://mlnotebook.github.io/post/surface-distance-function/

Is there any matlab code for that?

Thanks in advance

Answer 1

[Image Analyst]

You can simply use sqrt() and min():

x1 = S1(:, 1);
y1 = S1(:, 2);
x2 = S2(:, 1);
y2 = S2(:, 2);
for k = 1 : length(S1)
    % Get the distance from the kth point of S1 to all points in S2.
    distances = sqrt((x1(k) - x2) .^ 2 + (y1(k) - y2) .^ 2);
    % Find the min of those distances to find how how close point k came to curve 2
    minDistance(k) = min(distances);
end

If you have 3-D (x,y,z) coordinates, the extension to 3-D should be obvious.

[talayeh ghodsi]

thanks for your answer

and could you please help me to calculate The mean surface distance (MSD) and Residual mean Square distance (RMS) from that? the formula is attached below

[Image Analyst]

It's pretty simple so you probably already have it by now. Did you just simply try two loops with the reference changed in each? You probably did. So did you get this:

x1 = S1(:, 1);
y1 = S1(:, 2);
x2 = S2(:, 1);
y2 = S2(:, 2);
for k = 1 : length(S1)
    % Get the distance from the kth point of S1 to all points in S2.
    distances = sqrt((x1(k) - x2) .^ 2 + (y1(k) - y2) .^ 2);
    % Find the min of those distances to find how how close point k came to curve 2
    minDistance1(k) = min(distances);
end
for k = 1 : length(S2)
    % Get the distance from the kth point of S2 to all points in S1.
    distances = sqrt((x2(k) - x1) .^ 2 + (y2(k) - y1) .^ 2);
    % Find the min of those distances to find how how close point k came to curve 1
    minDistance2(k) = min(distances);
end
numPairs = numel(x1) + numel(x2)
MSD = (sum(minDistance1) + sum(minDistance2)) / numPairs
RMS = sqrt((sum(minDistance1.^2) + sum(minDistance2.^2)) / numPairs)

Does that agree with your code?

[talayeh ghodsi]

Thanks for your answer

the code returns 0 for MSD and RMS is 0

[talayeh ghodsi]

I have used the Sobel filter in order to extract the border of my image as attached below, but the contents in workspace are all zero. I also attach my images

S1=imread('ECHO-mask.jpg');
S2=rgb2gray(imread('CTBWslice.jpg'));
S2=imresize(S2,[512 512]);
S1 = edge(S1,'sobel');
S2 = edge(S2,'sobel');
% BW1=imcomplement(BW1);
% BW2=imcomplement(BW2);
C = rgb2gray(imfuse(S1,S2));
imshow(C)
%%%% S1 and S2 are two contours
x1 = S1(:, 1);
y1 = S1(:, 2);
x2 = S2(:, 1);
y2 = S2(:, 2);
for k = 1 : length(S1)
    % Get the distance from the kth point of S1 to all points in S2.
    distances = sqrt((x1(k) - x2) .^ 2 + (y1(k) - y2) .^ 2);
    % Find the min of those distances to find how how close point k came to curve 2
    minDistance1(k) = min(distances);
end
for k = 1 : length(S2)
    % Get the distance from the kth point of S2 to all points in S1.
    distances = sqrt((x2(k) - x1) .^ 2 + (y2(k) - y1) .^ 2);
    % Find the min of those distances to find how how close point k came to curve 1
    minDistance2(k) = min(distances);
end
numPairs = numel(x1) + numel(x2)
MSD = (sum(minDistance1) + sum(minDistance2)) / numPairs
RMS = sqrt((sum(minDistance1.^2) + sum(minDistance2.^2)) / numPairs)

[Image Analyst]

You didn't get the boundaries correctly. The boundaries are not an edge image, you need to call bwboundaries().

[Image Analyst]

talayeh, try this:

clc;    % Clear the command window.
fprintf('Beginning to run %s.m.\n', mfilename);
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;

% Read in and display the images.
image1=imread('ECHO-mask.jpg');
image1 = imbinarize(image1);
% Fill holes
image1 = imfill(image1, 'holes');
subplot(1, 2, 1);
imshow(image1, []);
title('ECHO-mask.jpg');
image2=rgb2gray(imread('CTBWslice.jpg'));
image2 = imbinarize(image2);
% Fill holes
image2 = imfill(image2, 'holes');
subplot(1, 2, 2);
imshow(image2, []);
title('CTBWslice.jpg');

% Get the boundaries
boundaries1 = bwboundaries(image1);
boundaries2 = bwboundaries(image2);
for k1 = 1 : length(boundaries1)
	b1 = boundaries1{k1};
	S1 = fliplr(b1); % Flip so it's (x,y), not (row, column).
	for k2 = 1 : length(boundaries2)
		b2 = boundaries2{k2};
		S2 = fliplr(b2); % Flip so it's (x,y), not (row, column).
		
		%%%% S1 and S2 are two contours
		x1 = S1(:, 1);
		y1 = S1(:, 2);
		x2 = S2(:, 1);
		y2 = S2(:, 2);
		for k = 1 : length(S1)
			% Get the distance from the kth point of S1 to all points in S2.
			distances = sqrt((x1(k) - x2) .^ 2 + (y1(k) - y2) .^ 2);
			% Find the min of those distances to find how how close point k came to curve 2
			minDistance1(k) = min(distances);
		end
		for k = 1 : length(S2)
			% Get the distance from the kth point of S2 to all points in S1.
			distances = sqrt((x2(k) - x1) .^ 2 + (y2(k) - y1) .^ 2);
			% Find the min of those distances to find how how close point k came to curve 1
			minDistance2(k) = min(distances);
		end
		
		% Compute the metrics we're interested in.
		numPairs = numel(x1) + numel(x2);
		MSD = (sum(minDistance1) + sum(minDistance2)) / numPairs;
		RMS = sqrt((sum(minDistance1.^2) + sum(minDistance2.^2)) / numPairs);
		fprintf('Comparing image1 boundary #%d with image2 boundary #%d,\n  we get MSD = %f, and RMS = %f.\n',...
			k1, k2, MSD, RMS);
	end
end
fprintf('Done running %s.m.\n', mfilename);

In the command window, you'll see:

Comparing image1 boundary #1 with image2 boundary #1,

we get MSD = 5.545358, and RMS = 7.403608.

Comparing image1 boundary #1 with image2 boundary #2,

we get MSD = 60.859496, and RMS = 75.365100.

Comparing image1 boundary #2 with image2 boundary #1,

we get MSD = 57.429919, and RMS = 74.059752.

Comparing image1 boundary #2 with image2 boundary #2,

we get MSD = 3.204780, and RMS = 8.873175.