Code covered by the BSD License  

Highlights from
octree - partitioning 3D points into spatial subvolumes

4.875

4.9 | 8 ratings Rate this file 74 Downloads (last 30 days) File Size: 16.3 KB File ID: #40732
image thumbnail

octree - partitioning 3D points into spatial subvolumes

by

 

11 Mar 2013 (Updated )

OcTree recursively splits a large set of points into smaller subvolumes. A QuadTree but in 3D.

| Watch this File

File Information
Description

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.
http://en.wikipedia.org/wiki/Octree
   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)

  Please post comments to this FEX page for this entry if you find any bugs or have any feature requests. There's plenty of room for improved efficiency and enhancement.

Acknowledgements

This file inspired 3 D Hrtf Interpolation.

Required Products MATLAB
MATLAB release MATLAB 8.1 (R2013a)
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (22)
31 Oct 2014 Peter Morovic

Nice work, thanks for sharing!

30 Oct 2014 Sven

Hi Jonathan,

I haven't implemented specifically a nearest neighbour routine. You could use OcTree to find which bin your query location belongs to, and then find the nearest neighbour within that bin. This will be faster than querying all points in a large point cloud, however this alone won't guarantee correctness (since a point in an adjacent bin may be closer) unless the distance to the found nearest neighbour is less than the distance from your query point to any edges of its bin.

Did that help you out? If your point cloud isn't massive (or even if it is), John D'Errico's ipdm is pretty efficient, and there are some nearest neighbour kdtree implementations that may be helpful too.

10 Oct 2014 Jonathan

Hello Sven,
I would like to use your Octree code in order to find the N closest point of an arbitrary position in my domain. I was using the quedtree code on the Matlab File Exchange but now that I am going to 3D, I would like to use your code for this for a fast localisation of neighbour points of a specific coordinates. Is it something already possible with your code? I found the OT.query function but this is not exactly what I want to do. Any idea of I can make it from your code? Thank you very much anyway for your great code. Best regards,Jonathan

04 Oct 2014 Zhouxin Xi  
30 Sep 2014 Sven

Hi Christophe, for your second question, you can force the decomposition to go to a specific level with the optional parameters:
OT = OcTree(pts,'binCapacity',inf,'maxDepth',3,'maxSize',0);

The OT will still contain the higher levels above that, but you can easily get only the desired level bins with a mask such as:
OT.BinDepths==2

For the neighbours that's an interesting question. It will be easy to find the 3 "siblings" of a given bin (they will be the same level and have the same parent bin No. as that bin). It might need a bit more work to find all neighbours... what exactly would you define as a neighbour? Can only "leaf-node" bins be neighbours or do you include any level's bin that shares a face with your target a neighbour?

29 Sep 2014 Christophe Langrenne

Hi,
An other question : Is it possible to have only one level tree (for exemple, only bins of level 3) ?

29 Sep 2014 Christophe Langrenne

Hi Sven,
I have to find the neighbors of a given box (for a Fast Multipole Method application). Do you have an easy way to do this ?

Thanks and "bravo" for your great job.

05 Sep 2014 Sven

@Jing: Yep, thanks for finding that bug. Fixed now.

05 Sep 2014 Sven

@Wu Jun: It's strange but QuadTree decomposition in images is a different beast to QuadTree decomposition of 2d points, but they both have the same name. I'm afraid it's the same situation for OcTree decomposition (of image volumes vs 3d points). OcTree.m is targeted at points rather than image volume data.

05 Sep 2014 Sven

@Greg: (sorry for the late reply) This is best done *after* plotting the bins. You can just find which bins don't contain a point and then delete or deactivate that bin's plot handle:

% Run example 1 and then:
binHasPts = ismember(1:OT.BinCount, OT.PointBins)
set(boxH(~binHasPts),'Visible','off')

04 Sep 2014 Jing

A bug in the code???
should "this.BinParents(binNo)+1" be "this.BinDepths(binNo)+1" ??

% Prevent dividing beyond the maximum depth

if this.BinParents(binNo)+1 >= this.Properties.maxDepth
continue;
end

19 Dec 2013 Wu Jun

Hi Sven, Thanks for your great job. I want to use Octree method to divides the 3D images volume into little 3D cube just like Quadree method (matlab qtdecomp()function). Do you have any plan to add this feature to your code ? Regards

08 Nov 2013 Greg

Hi Sven,

Thank you for your great job.
I am an absolute beginner in Matlab but I managed to run your code and got quite good results. I am wondering if it is possible to delete the empty bins (i.e. bins of the higher levels)and keep only the child bins of the lowest level before plotting the final OcTree.

Kind regards
Greg

17 Apr 2013 Siyi Deng  
04 Apr 2013 Moein

@Sven Thank you very much.

02 Apr 2013 Sven

@Moein: Thanks for your feedback and request, I'm looking at implementing the concept of element connectivity so that you can partition based on collections of points (elements) as well as the points themselves. Note that it seems you have an immediate need (2:1 balancing) for a particular application (FE meshing). I would suggest that if your problem is urgent then you might consider either adding to the class or inheriting and making a new class specific to your needs yourself. I will work on element connectivity first, and node balancing will be some time away, so my schedule won't necessarily fit in with yours.

01 Apr 2013 Moein

Hi Sven
I need element connectivity in your code , How can I get element connectivity ? Thanks for your help

29 Mar 2013 Moein

Also, in order to balance your mesh, you need information about elements connectivity. Because if you have second or higher level hanging nodes in one edge, you should subdivide all the elements which connected to that edge.you can find some algorithms for finding neighbors in octree mesh in Journals.

29 Mar 2013 Moein

@Sven Thanks for your reply. In octree decomposition, there will arise “hanging” nodes, i.e., nodes that lie on the edge of an undivided box, rather than at one of its corners. If, during the subdivision, a second hanging node is placed on the edge of a box, that box should be itself subdivided. In this way, the final mesh does not have any “second-level” hanging nodes.
Also I think it would be better to define the binary mask to your code showing which elements of a set are to be subdivided. It gives more freedom to users. The suitable outputs for Finite element method are the positions of the all nodes , The connectivity array which shows the number of 8 points of each cubic element.(for octree case).

28 Mar 2013 Sven

@Moein: Thanks for the feedback. Yes, it doesn't yet consider 2:1 balancing. I have a few plans for new features, starting with including the concept of elements (ie, triangles or quads or hexas... hopefully I can do N-dim elements).
I hadn't yet thought through node balancing, do you have a specific example you can provide and the type of output you'd want? For example are you simply looking for adjusting a current decomposition (swapping neigbours) based on particular criteria? If so, is that criteria based solely on points or do you need information about "element connectivity" too?

28 Mar 2013 Moein

Hi , Thanks for your great job. I tested this code and I noticed that the code doesn't consider 2:1 rule for neighboring mesh elements. ( I mean the final tree is not balanced). Is it true? If yes, Do you have any plan to add this feature to your code ? (Its very important in FEM mesh generation). Regards

13 Mar 2013 Anton Semechko

superb job

Updates
18 Mar 2013

Added ability to shrink bins to minimally cover their child points/bins

05 Sep 2014

Thanks to Jing for a BinDepths bugfix

Contact us