Homec4science
Diffusion ShoulderDB_ImportSegmentMeasureAnalyze (aterrier-plotUpdate) Browse

R8218/ShoulderCase/@Plane/Plane.m
997424434cadaterrier-plotUpdate

Last Change
Raw File

Plane.m
View Options

classdef Plane < handle
% A plane is decribed by a point and a normal.
%
% Plane can define the desired plane's normal and point
% by fitting an array of points.
% Plane has a method to project points on the current plane.
properties
normal
point
fitPerformance
end
methods
function obj = Plane()
obj.normal = [];
obj.point = [];
obj.fitPerformance.points = [];
obj.fitPerformance.residuals = [];
obj.fitPerformance.RMSE = [];
obj.fitPerformance.R2 = [];
end
function plot(obj,varargin)
% Plot a circle in the plane. The circle's radius can be given in argument
radius = 100;
if nargin == 2
assert(isnum)
radius = varargin{1};
end
circleInPlane = Circle(obj.point,radius,obj.normal,36);
XYZ = circleInPlane.points;
patch(XYZ(:,1),XYZ(:,2),XYZ(:,3),'black','FaceAlpha',0.3); % Transparent plane
end
function fit(obj,points)
% Fit current plane to given points
[coeff,score,roots] = pca(points);
normal = cross(coeff(:,1),coeff(:,2));
normal = normal/norm(normal);
meanPoint = mean(points,1);
estimatedPoints = repmat(meanPoint,size(points,1),1) +...
score(:,1:2)*coeff(:,1:2)';
residuals = points - estimatedPoints;
error = vecnorm(residuals,2,2);
RMSE = norm(error)/sqrt(size(points,1));
sumSquareError = sum(error.^2);
total = vecnorm(points - meanPoint,2,2)';
sumSquareTotal = sum(total.^2);
R2 = 1-(sumSquareError/sumSquareTotal); %http://en.wikipedia.org/wiki/Coefficient_of_determination
obj.normal = normal';
obj.point = meanPoint;
obj.setFitPerformance(points,residuals,RMSE,R2);
end
function setFitPerformance(obj,points,residuals,RMSE,R2)
obj.fitPerformance.points = points;
obj.fitPerformance.residuals = residuals;
obj.fitPerformance.RMSE = RMSE;
obj.fitPerformance.R2 = R2;
end
function setPointOnPlane(obj,point)
obj.point = point;
end
function projectedPoints = projectOnPlane(obj,points)
N2 = obj.normal.'*obj.normal;
projectedPoints = points*(eye(3)-N2)+repmat(obj.point*N2,size(points,1),1);
end
end
end
c4science · Help