classdef OcTree < handle
% OcTree point decomposition in 3D
% OcTree is used to create a tree data structure of bins containing 3D
% points. Each bin may be recursively decomposed into 8 child bins.
%
% OT = OcTree(PTS) creates an OcTree from an N-by-3 matrix of point
% coordinates.
%
% OT = OcTree(...,'PropertyName',VALUE,...) takes any of the following
% property values:
%
% binCapacity - Maximum number of points a bin may contain. If more
% points exist, the bin will be recursively subdivided.
% Defaults to ceil(numPts/10).
% maxDepth - Maximum number of times a bin may be subdivided.
% Defaults to INF.
% maxSize - Maximum size of a bin edge. If any dimension of a bin
% exceeds maxSize, it will be recursively subdivided.
% Defaults to INF.
% minSize - Minimum size of a bin edge. Subdivision will stop after
% any dimension of a bin gets smaller than minSize.
% Defaults to 1000*eps.
% style - Either 'equal' (default) or 'weighted'. 'equal'
% subdivision splits bins at their central coordinate
% (ie, one bin subdivides into 8 equally sized bins).
% 'weighted' subdivision divides bins based on the mean
% of all points they contain. Weighted subdivision is
% slightly slower than equal subdivision for a large
% number of points, but it can produce a more efficient
% decomposition with fewer subdivisions.
%
% Example 1: Decompose 200 random points into bins of 20 points or less,
% then display each bin with its points in a separate colour.
% pts = (rand(200,3)-0.5).^2;
% OT = OcTree(pts,'binCapacity',20);
% figure
% boxH = OT.plot;
% cols = lines(OT.BinCount);
% doplot3 = @(p,varargin)plot3(p(:,1),p(:,2),p(:,3),varargin{:});
% for i = 1:OT.BinCount
% set(boxH(i),'Color',cols(i,:),'LineWidth', 1+OT.BinDepths(i))
% doplot3(pts(OT.PointBins==i,:),'.','Color',cols(i,:))
% end
% axis image, view(3)
%
% Example 2: Decompose 200 random points into bins of 10 points or less,
% shrunk to minimallly encompass their points, then display.
% pts = rand(200,3);
% OT = OcTree(pts,'binCapacity',10,'style','weighted');
% OT.shrink
% figure
% boxH = OT.plot;
% cols = lines(OT.BinCount);
% doplot3 = @(p,varargin)plot3(p(:,1),p(:,2),p(:,3),varargin{:});
% for i = 1:OT.BinCount
% set(boxH(i),'Color',cols(i,:),'LineWidth', 1+OT.BinDepths(i))
% doplot3(pts(OT.PointBins==i,:),'.','Color',cols(i,:))
% end
% axis image, view(3)
%
%
% OcTree methods:
% shrink - Shrink each bin to tightly encompass its children
% query - Ask which bins a new set of points belong to.
% plot, plot3 - Plots bin bounding boxes to the current axes.
%
% OcTree properties:
% Points - The coordinate of points in the decomposition.
% PointBins - Indices of the bin that each point belongs to.
% BinCount - Total number of bins created.
% BinBoundaries - BinCount-by-6 [MIN MAX] coordinates of bin edges.
% BinDepths - The # of subdivisions to reach each bin.
% BinParents - Indices of the bin that each bin belongs to.
% Properties - Name/Val pairs used for creation (see help above)
%
% See also qtdecomp.
% Created by Sven Holcombe.
% 1.0 - 2013-03 Initial release
% 1.1 - 2013-03 Added shrinking bins and allocate/deallocate space
%
% Please post comments to the FEX page for this entry if you have any
% bugs or feature requests.
properties
Points;
PointBins;
BinCount;
BinBoundaries;
BinDepths;
BinParents = zeros(0,1);
Properties;
end
methods
function this = OcTree(pts,varargin)
% This is the OcTree header line
validateattributes(pts,{'numeric'},...
{'real','finite','nonnan','ncols', 3},...
mfilename,'PTS')
% Initialise a single bin surrounding all given points
numPts = size(pts,1);
this.BinBoundaries = [min(pts,[],1) max(pts,[],1)];
this.Points = pts;
this.PointBins = ones(numPts,1);
this.BinDepths = 0;
this.BinParents(1) = 0;
this.BinCount = 1;
% Allow custom setting of Properties
IP = inputParser;
IP.addParamValue('binCapacity',ceil(numPts)/10);
IP.addParamValue('maxDepth',inf);
IP.addParamValue('maxSize',inf);
IP.addParamValue('minSize',1000 * eps);
IP.addParamValue('style','equal');
IP.parse(varargin{:});
this.Properties = IP.Results;
% Return on empty or trivial bins
if numPts<2, return; end
% Start dividing!
this.preallocateSpace;
this.divide(1);
this.deallocateSpace;
end
% MATLAB performs better if arrays that grow are initialised,
% rather than grown during a loop. These two functions do just that
% before and after the identification of new beens.
function preallocateSpace(this)
numPts = size(this.Points,1);
numBins = numPts;
if isfinite(this.Properties.binCapacity)
numBins = ceil(2*numPts/this.Properties.binCapacity);
end
this.BinDepths(numBins) = 0;
this.BinParents(numBins) = 0;
this.BinBoundaries(numBins,1) = 0;
end
function deallocateSpace(this)
this.BinDepths(this.BinCount+1:end) = [];
this.BinParents(this.BinCount+1:end) = [];
this.BinBoundaries(this.BinCount+1:end,:) = [];
end
function divide(this, startingBins)
% Loop over each bin we will consider for division
for i = 1:length(startingBins)
binNo = startingBins(i);
% Prevent dividing beyond the maximum depth
if this.BinParents(binNo)+1 >= this.Properties.maxDepth
continue;
end
% Prevent dividing beyond a minimum size
thisBounds = this.BinBoundaries(binNo,:);
binEdgeSize = diff(thisBounds([1:3;4:6]));
minEdgeSize = min(binEdgeSize);
maxEdgeSize = max(binEdgeSize);
if minEdgeSize < this.Properties.minSize
continue;
end
% There are two conditions under which we should divide
% this bin. 1: It's bigger than maxSize. 2: It contains
% more points than binCapacity.
oldCount = this.BinCount;
if nnz(this.PointBins==binNo) > this.Properties.binCapacity
this.divideBin(binNo);
this.divide(oldCount+1:this.BinCount);
continue;
end
if maxEdgeSize>this.Properties.maxSize
this.divideBin(binNo);
this.divide(oldCount+1:this.BinCount);
continue;
end
end
end
function divideBin(this,binNo)
% Gather the new points (a bit more efficient to copy once
评论0