Delaunay 3D

a (very buggy) attempt to perform a  3D Delaunay tetrahedralization. computationally expensive, click to respawn.

by nicoptere 16 Aug 2012

Talk

bradsedito at 17 Aug 2012 03:55
Nicoptere!! You're the man! =D
nicoptere at 17 Aug 2012 13:37
thanks, it was a follow up to your last question ;) as far as I understood (from your reply on Makc's delaunay fork), you'd only need a distance function to calculate the distances of your planets. it's really not a good approach, as Makc mentioned, if you're willing to go Delaunay, you can check this Voronoi 3D: http://wonderfl.net/c/qilS
makc3d at 18 Aug 2012 04:02
with every edge colored the same and no faces my brain refuses to reconstruct it in 3D :S
nicoptere at 18 Aug 2012 16:13
haha! the worst is that even if you render each tetrahedron with a different color, your brain will still refuse to see it properly. as in the 2D version, lots of triangles alter the original shape and turn it into a visual clusterfuck. a while ago after studying the thing more indepth, I got pretty frustrated by the results (and gave up btw). check this for instance: http://www.flickr.com/photos/barradeau/sets/72157627806881267/ some renders are cool most are not "readable"
nicoptere at 18 Aug 2012 16:15
and an important step that is not implemented is the winding : here, the tetrahedra are not convex which is another reason for it not to work :)

Tags

3D delaunay geometry
Embed

/**
 * Copyright nicoptere ( http://wonderfl.net/user/nicoptere )
 * MIT License ( http://www.opensource.org/licenses/mit-license.php )
 * Downloaded from: http://wonderfl.net/c/4HZO
 */

package 
{
	import flash.display.Sprite;
	import flash.display.TriangleCulling;
	import flash.events.Event;
	import flash.events.MouseEvent;
	
	/**
	 * @author nicolas barradeau
	 * @website http://en.nicoptere.net/
	 */
	public class Delaunay3dTest extends Sprite 
	{
		private var mesh:Mesh;
		
		public function Delaunay3dTest() 
		{
			
			reset();	
			addEventListener( Event.ENTER_FRAME, oef );
			stage.addEventListener( MouseEvent.MOUSE_DOWN, onMouseDownHandler );
		}
		
		private function onMouseDownHandler(e:MouseEvent):void 
		{
			reset();
		}
		
		private function reset():void 
		{
			var size:Number = 400;
			var vertices:Vector.<Number> = new Vector.<Number>();
			for ( var i:int = 0; i < 10; i++ )
			{
				vertices.push( 	( Math.random() - .5 ) * size,
								( Math.random() - .5 ) * size,
								( Math.random() - .5 ) * size );
			}
			x = y = 465 * .5;;
			mesh = new Mesh( vertices );
		
			var d:Delaunay3D = new Delaunay3D();
			mesh.indices = d.compute( mesh.vertices );
		
		}
		
		private function oef(e:Event):void 
		{
			mesh.rotationY += 1;
			graphics.clear();
			graphics.lineStyle( 0, 0x0033DD, .8 );
			mesh.render( graphics, TriangleCulling.POSITIVE );
		}
		
	}

}

import flash.display.Graphics;
import flash.display.TriangleCulling;
import flash.events.EventDispatcher;
import flash.geom.Matrix3D;
import flash.geom.Utils3D;
import flash.geom.Vector3D;
class Delaunay3D
{
	static public var EPSILON:Number = Number.MIN_VALUE;
	static public var SUPER_TETRAHEDRON_RADIUS:Number = 1000000000;
	
	private var indices:Vector.<int>;
	public var tetras:Vector.<Tetra>;
	private var faces:Vector.<Face>;
	
	public function compute( vertices:Vector.<Number> ):Vector.<int>
	{
		if ( vertices.length < 4) return null;
		var nv:int = vertices.length;	
		
		//creates a super tetrahedraon
		var d:Number = SUPER_TETRAHEDRON_RADIUS;
		vertices.push( 	-0.798892958090	* d, -0.274518427439	* d, 0.535172563978		* d,
						0.337674633600	* d, -0.783974919468	* d, -0.520921458490	* d,
						0.717804165557	* d, 0.325344236568		* d, 0.615555282682		* d,
						-0.256585841067	* d, 0.733149110340		* d, -0.629806388170	* d);
	
		var sid0:int = ( vertices.length / 3 ) - 4;
		var sid1:int = ( vertices.length / 3 ) - 3;
		var sid2:int = ( vertices.length / 3 ) - 2;
		var sid3:int = ( vertices.length / 3 ) - 1;
		
		tetras = Vector.<Tetra>( [ new Tetra( vertices, sid0, sid1, sid2, sid3 ) ] );
		faces = new Vector.<Face>();
		
		var t:Tetra, f:Face, nf:Face;
		var i:int, j:int, k:int;
		var id0:int, id1:int, id2:int, id3:int
		for ( i = 0; i < nv; i+=3 )
		{
			//test if vertex is inside one of the circumspheres
			for ( j = 0; j < tetras.length; j++ )
			{
				t = tetras[ j ];
				
				if ( t.sphereContains( vertices[ i ], vertices[ i + 1 ], vertices[ i + 2 ] )	)
				{
					for each ( f in t.faces ) faces.push( f );
					tetras.splice( j, 1 );
					j--;
				}
			}
			//remoev duplicate faces
			for ( j = 0; j < faces.length; j++ )
			{
				f = faces[ j ];
				search : for ( k = j + 1; k < faces.length; k++ )
				{
					nf = faces[ k ];
					if ( f.equals( nf ) )
					{
						//Hmmmmm ... should be done but, no...
						//faces.splice( j, 1 );
						//j--;
						//if ( j < 0 ) break;
					
						faces.splice( k, 1 );
						k--;
						if ( k < 0 ) break;
					}
				}
			}
			for ( j = 0; j < faces.length; j++ )
			{
				tetras.push( new Tetra( vertices, i / 3, faces[ j ].id0, faces[ j ].id1, faces[ j ].id2 ) );
			}
			faces.length = 0;
			
		}
		
		//cleanup & export (remove the super tetrahedron and all tetrahedra attached to iits vertices )
		faces.length = 0;
		indices ||= new Vector.<int>();
		for each ( t in tetras ) 
		{
			var c:int = 0;
			if ( t.containsVertex( sid0 ) ) c++;
			if ( t.containsVertex( sid1 ) ) c++;
			if ( t.containsVertex( sid2 ) ) c++;
			if ( t.containsVertex( sid3 ) ) c++;
			if ( c == 0 )
			{
				for each ( f in t.faces) 
				{
					indices.push( f.id0, f.id1, f.id2 );
				}
			}
		}
		return indices
	}
}
/*************************************
 *
 * UTILS
 * 
*************************************/
class Object3D extends EventDispatcher
{
	static public const RADIANS:String = "radians";
	static public const DEGREES:String = "degrees";
	
	static private var DEG_TO_RAD:Number = Math.PI / 180;
	static private var RAD_TO_DEG:Number = 180 / Math.PI;
	
	protected var _transform:Matrix3D = new Matrix3D();		
	protected var _position:Vector3D = new Vector3D( 0,0,0 );
	protected var _rotation:Vector3D = new Vector3D( 0,0,0 );
	protected var _scale:Vector3D = new Vector3D( 1, 1, 1 );
	protected var geomDirty:Boolean = true;
	
	public function Object3D() {}
	public function get transform():Matrix3D 
	{
		if ( geomDirty )
		{
			_transform.recompose( Vector.<Vector3D>([ _position, _rotation, _scale ]) );
			geomDirty = false;
		}
		return _transform;
	}
	
	public function set transform(value:Matrix3D):void 
	{
		_transform = value;
		geomDirty = false;
	}
	
	
	public function get position():Vector3D {	return _position;	}
	public function set position(value:Vector3D):void 
	{
		_position = value;
		geomDirty = true;
	}
	
	public function get rotation():Vector3D {	return _rotation;	}
	public function set rotation(value:Vector3D):void 
	{
		_rotation = value;
		geomDirty = true;
	}
	
	public function set x(value:Number):void {		position.x = value;		geomDirty = true;		}
	public function set y(value:Number):void {		position.y = value;		geomDirty = true;		}
	public function set z(value:Number):void {		position.z = value;		geomDirty = true;		}
	
	public function get x():Number {		return position.x;		}
	public function get y():Number {		return position.y;		}
	public function get z():Number {		return position.z;		}
	
	public function set rotationX(value:Number):void {		rotation.x = value * DEG_TO_RAD;		geomDirty = true;		}
	public function set rotationY(value:Number):void {		rotation.y = value * DEG_TO_RAD;		geomDirty = true;		}
	public function set rotationZ(value:Number):void {		rotation.z = value * DEG_TO_RAD;		geomDirty = true;		}
	
	public function get rotationX():Number {		return rotation.x * RAD_TO_DEG;		}
	public function get rotationY():Number {		return rotation.y * RAD_TO_DEG;		}
	public function get rotationZ():Number {		return rotation.z * RAD_TO_DEG;		}
	
	public function set scale(value:Number):void {		_scale.x = _scale.y = _scale.z = value;		geomDirty = true;		}
	public function set scaleX(value:Number):void {		_scale.x = value;		geomDirty = true;		}
	public function set scaleY(value:Number):void {		_scale.y = value;		geomDirty = true;		}
	public function set scaleZ(value:Number):void {		_scale.z = value;		geomDirty = true;		}
	
	public function get scaleX():Number {		return _scale.x;		}
	public function get scaleY():Number {		return _scale.y;		}
	public function get scaleZ():Number {		return _scale.z;		}
	
	
	static public function set ANGLE_TYPE( value:String ):void
	{
		DEG_TO_RAD = value == DEGREES ? ( Math.PI / 180 ) : 1;
		RAD_TO_DEG = value == DEGREES ? ( 180 / Math.PI ) : 1;
	}
	
}
class Mesh extends Object3D 
{
	
	protected var _vertices:Vector.<Number> = new Vector.<Number>();
	protected var _indices:Vector.<int> = new Vector.<int>();
	protected var _uvts:Vector.<Number> = new Vector.<Number>();
	protected var projected:Vector.<Number> = new Vector.<Number>();
	
	public function Mesh( vertices:Vector.<Number> = null )
	{
		if ( vertices != null )
		{
			this.vertices = vertices;
		}
	}
	
	public function render( graphics:Graphics, culling:String = TriangleCulling.POSITIVE ):void
	{
	 
		if( geomDirty )compute();
		
		if ( indices != null )
		{
			
			graphics.drawTriangles( projected, indices, uvts, culling );
			graphics.lineStyle( 0, 0xCC2200 );
			for ( var i:int = 0; i < projected.length; i+=2 ) 
			{
				graphics.drawCircle( projected[ i ], projected[ i+1 ], 10 );
			}
		}
		
	}
	
	public function compute():void 
	{
		Utils3D.projectVectors( transform, vertices, projected, uvts );
	}
	
	public function get vertices():Vector.<Number> 
	{
		return _vertices;
	}
	
	public function set vertices(value:Vector.<Number>):void 
	{
		_vertices = value;
		geomDirty = true;
	}
	
	public function get indices():Vector.<int> 
	{
		return _indices;
	}
	
	public function set indices(value:Vector.<int>):void 
	{
		_indices = value;
		geomDirty = true;
	}
	
	public function get uvts():Vector.<Number> 
	{
		return _uvts;
	}
	
	public function set uvts(value:Vector.<Number>):void 
	{
		_uvts = value;
		geomDirty = true;
	}
	
}

internal class Tetra
{
	public var x:Number;
	public var y:Number;
	public var z:Number;
	public var radius:Number;
	
	public var id0:int;
	public var id1:int;
	public var id2:int;
	public var id3:int;
	public var faces:Vector.<Face>;
	
	public function Tetra( vertices:Vector.<Number>, id0:int, id1:int, id2:int, id3:int )
	{
		
		this.id0 = id0;
		this.id1 = id1;
		this.id2 = id2;
		this.id3 = id3;
		
		faces = Vector.<Face>( [ 
								new Face( id0, id1, id2 ),
								new Face( id0, id2, id3 ),
								new Face( id0, id3, id1 ),
								new Face( id1, id2, id3 )
								] );
		
		computeSphere( vertices, id0, id1, id2, id3 );
		
	}
	
	private function computeSphere( vertices:Vector.<Number>, id0:int, id1:int, id2:int, id3:int ):void
	{
		var ax:Number = vertices[ id0 ], ay:Number = vertices[ id0 + 1 ], az:Number = vertices[ id0 + 2 ];
		var bx:Number = vertices[ id1 ], by:Number = vertices[ id1 + 1 ], bz:Number = vertices[ id1 + 2 ];
		var cx:Number = vertices[ id2 ], cy:Number = vertices[ id2 + 1 ], cz:Number = vertices[ id2 + 2 ];
		var dx:Number = vertices[ id3 ], dy:Number = vertices[ id3 + 1 ], dz:Number = vertices[ id3 + 2 ];
		
		var xba:Number, yba:Number, zba:Number, xca:Number, yca:Number, zca:Number, xda:Number, yda:Number, zda:Number;
		var balength:Number, calength:Number, dalength:Number;
		var xcrosscd:Number, ycrosscd:Number, zcrosscd:Number;
		var xcrossdb:Number, ycrossdb:Number, zcrossdb:Number;
		var xcrossbc:Number, ycrossbc:Number, zcrossbc:Number;
		var denominator:Number;
		
		xba = bx - ax;
		yba = by - ay;
		zba = bz - az;
		
		xca = cx - ax;
		yca = cy - ay;
		zca = cz - az;
		
		xda = dx - ax;
		yda = dy - ay;
		zda = dz - az;
		
		balength = xba * xba + yba * yba + zba * zba;
		calength = xca * xca + yca * yca + zca * zca;
		dalength = xda * xda + yda * yda + zda * zda;
		
		xcrosscd = yca * zda - yda * zca;
		ycrosscd = zca * xda - zda * xca;
		zcrosscd = xca * yda - xda * yca;
		
		xcrossdb = yda * zba - yba * zda;
		ycrossdb = zda * xba - zba * xda;
		zcrossdb = xda * yba - xba * yda;
		
		xcrossbc = yba * zca - yca * zba;
		ycrossbc = zba * xca - zca * xba;
		zcrossbc = xba * yca - xca * yba;
		
		denominator = 0.5 / (xba * xcrosscd + yba * ycrosscd + zba * zcrosscd);
		
		//circumcenter
		x = ax + (balength * xcrosscd + calength * xcrossdb + dalength * xcrossbc) * denominator;
		y = ay + (balength * ycrosscd + calength * ycrossdb + dalength * ycrossbc) * denominator;
		z = az + (balength * zcrosscd + calength * zcrossdb + dalength * zcrossbc) * denominator;
		
		//circumRadius
		dx = x -= ax;
		dy = y -= ay;
		dz = z -= az;
		radius = dx * dx + dy * dy + dz * dz;
		
	}
	
	public function sphereContains( x:Number, y:Number, z:Number ):Boolean 
	{
		var dx:Number = this.x - x;
		var dy:Number = this.y - y;
		var dz:Number = this.z - z;
		return radius <= dx * dx + dy * dy + dz * dz;
	}
	
	public function containsVertex( id:int ):Boolean
	{
		return id == id0 || id == id1 || id == id2 || id == id3;
	}
	
}
class Face
{
	public var sum:Number;
	public var id0:int;
	public var id1:int;
	public var id2:int;
	
	public function Face( id0:int, id1:int, id2:int )
	{
		this.id0 = id0;
		this.id1 = id1;
		this.id2 = id2;
		sum = id0 + id1 + id2;
	}
	
	public function equals( face:Face ):Boolean
	{
		if ( sum == face.sum ) 
		{
			
			if ( id0 == face.id0 && id1 == face.id1 && id2 == face.id2 ) return true;
			if ( id0 == face.id0 && id2 == face.id1 && id1 == face.id2 ) return true;
			if ( id1 == face.id0 && id0 == face.id1 && id2 == face.id2 ) return true;
			if ( id1 == face.id0 && id2 == face.id1 && id0 == face.id2 ) return true;
			if ( id2 == face.id0 && id0 == face.id1 && id1 == face.id2 ) return true;
			if ( id2 == face.id0 && id1 == face.id1 && id0 == face.id2 ) return true;
			
		}		
		return false;
		
	}
}