# Fortune's Voronoi algorithm implemented in C#

By , 21 Apr 2013

 fortunevoronoi.zip FortuneVoronoi FortuneVoronoi.csproj.user ```using System; using System.Collections; using BenTools.Data; namespace BenTools.Mathematics { public class VoronoiGraph { public HashSet Vertizes = new HashSet(); public HashSet Edges = new HashSet(); } public class VoronoiEdge { public Vector RightData, LeftData; public Vector VVertexA = Fortune.VVUnkown, VVertexB = Fortune.VVUnkown; public void AddVertex(Vector V) { if(VVertexA==Fortune.VVUnkown) VVertexA = V; else if(VVertexB==Fortune.VVUnkown) VVertexB = V; else throw new Exception("Tried to add third vertex!"); } } // VoronoiVertex or VoronoiDataPoint are represented as Vector internal abstract class VNode { private VNode _Parent = null; private VNode _Left = null, _Right = null; public VNode Left { get{return _Left;} set { _Left = value; value.Parent = this; } } public VNode Right { get{return _Right;} set { _Right = value; value.Parent = this; } } public VNode Parent { get{return _Parent;} set{_Parent = value;} } public void Replace(VNode ChildOld, VNode ChildNew) { if(Left==ChildOld) Left = ChildNew; else if(Right==ChildOld) Right = ChildNew; else throw new Exception("Child not found!"); ChildOld.Parent = null; } public static VDataNode FirstDataNode(VNode Root) { VNode C = Root; while(C.Left!=null) C = C.Left; return (VDataNode)C; } public static VDataNode LeftDataNode(VDataNode Current) { VNode C = Current; //1. Up do { if(C.Parent==null) return null; if(C.Parent.Left == C) { C = C.Parent; continue; } else { C = C.Parent; break; } }while(true); //2. One Left C = C.Left; //3. Down while(C.Right!=null) C = C.Right; return (VDataNode)C; // Cast statt 'as' damit eine Exception kommt } public static VDataNode RightDataNode(VDataNode Current) { VNode C = Current; //1. Up do { if(C.Parent==null) return null; if(C.Parent.Right == C) { C = C.Parent; continue; } else { C = C.Parent; break; } }while(true); //2. One Right C = C.Right; //3. Down while(C.Left!=null) C = C.Left; return (VDataNode)C; // Cast statt 'as' damit eine Exception kommt } public static VEdgeNode EdgeToRightDataNode(VDataNode Current) { VNode C = Current; //1. Up do { if(C.Parent==null) throw new Exception("No Left Leaf found!"); if(C.Parent.Right == C) { C = C.Parent; continue; } else { C = C.Parent; break; } }while(true); return (VEdgeNode)C; } public static VDataNode FindDataNode(VNode Root, double ys, double x) { VNode C = Root; do { if(C is VDataNode) return (VDataNode)C; if(((VEdgeNode)C).Cut(ys,x)<0) C = C.Left; else C = C.Right; }while(true); } /// /// Will return the new root (unchanged except in start-up) /// public static VNode ProcessDataEvent(VDataEvent e, VNode Root, VoronoiGraph VG, double ys, out VDataNode[] CircleCheckList) { if(Root==null) { Root = new VDataNode(e.DataPoint); CircleCheckList = new VDataNode[] {(VDataNode)Root}; return Root; } //1. Find the node to be replaced VNode C = VNode.FindDataNode(Root, ys, e.DataPoint[0]); //2. Create the subtree (ONE Edge, but two VEdgeNodes) VoronoiEdge VE = new VoronoiEdge(); VE.LeftData = ((VDataNode)C).DataPoint; VE.RightData = e.DataPoint; VE.VVertexA = Fortune.VVUnkown; VE.VVertexB = Fortune.VVUnkown; VG.Edges.Add(VE); VNode SubRoot; if(Math.Abs(VE.LeftData[1]-VE.RightData[1])<1e-10) { if(VE.LeftData[0]=ys) return VC; return null; } } internal class VDataNode : VNode { public VDataNode(Vector DP) { this.DataPoint = DP; } public Vector DataPoint; } internal class VEdgeNode : VNode { public VEdgeNode(VoronoiEdge E, bool Flipped) { this.Edge = E; this.Flipped = Flipped; } public VoronoiEdge Edge; public bool Flipped; public double Cut(double ys, double x) { if(!Flipped) return Math.Round(x-Fortune.ParabolicCut(Edge.LeftData[0], Edge.LeftData[1], Edge.RightData[0], Edge.RightData[1], ys),10); return Math.Round(x-Fortune.ParabolicCut(Edge.RightData[0], Edge.RightData[1], Edge.LeftData[0], Edge.LeftData[1], ys),10); } } internal abstract class VEvent : IComparable { public abstract double Y {get;} public abstract double X {get;} #region IComparable Members public int CompareTo(object obj) { if(!(obj is VEvent)) throw new ArgumentException("obj not VEvent!"); int i = Y.CompareTo(((VEvent)obj).Y); if(i!=0) return i; return X.CompareTo(((VEvent)obj).X); } #endregion } internal class VDataEvent : VEvent { public Vector DataPoint; public VDataEvent(Vector DP) { this.DataPoint = DP; } public override double Y { get { return DataPoint[1]; } } public override double X { get { return DataPoint[0]; } } } internal class VCircleEvent : VEvent { public VDataNode NodeN, NodeL, NodeR; public Vector Center; public override double Y { get { return Math.Round(Center[1]+MathTools.Dist(NodeN.DataPoint[0],NodeN.DataPoint[1],Center[0],Center[1]),10); } } public override double X { get { return Center[0]; } } public bool Valid = true; } public abstract class Fortune { public static readonly Vector VVInfinite = new Vector(double.PositiveInfinity, double.PositiveInfinity); public static readonly Vector VVUnkown = new Vector(double.NaN, double.NaN); internal static double ParabolicCut(double x1, double y1, double x2, double y2, double ys) { // y1=-y1; // y2=-y2; // ys=-ys; // if(Math.Abs(x1-x2)<1e-10 && Math.Abs(y1-y2)<1e-10) { // if(y1>y2) // return double.PositiveInfinity; // if(y1xs2) { double h = xs1; xs1=xs2; xs2=h; } if(y1>=y2) return xs2; return xs1; } internal static Vector CircumCircleCenter(Vector A, Vector B, Vector C) { if(A==B || B==C || A==C) throw new Exception("Need three different points!"); double tx = (A[0] + C[0])/2; double ty = (A[1] + C[1])/2; double vx = (B[0] + C[0])/2; double vy = (B[1] + C[1])/2; double ux,uy,wx,wy; if(A[0] == C[0]) { ux = 1; uy = 0; } else { ux = (C[1] - A[1])/(A[0] - C[0]); uy = 1; } if(B[0] == C[0]) { wx = -1; wy = 0; } else { wx = (B[1] - C[1])/(B[0] - C[0]); wy = -1; } double alpha = (wy*(vx-tx)-wx*(vy - ty))/(ux*wy-wx*uy); return new Vector(tx+alpha*ux,ty+alpha*uy); } public static VoronoiGraph ComputeVoronoiGraph(IEnumerable Datapoints) { BinaryPriorityQueue PQ = new BinaryPriorityQueue(); Hashtable CurrentCircles = new Hashtable(); VoronoiGraph VG = new VoronoiGraph(); VNode RootNode = null; foreach(Vector V in Datapoints) { PQ.Push(new VDataEvent(V)); } while(PQ.Count>0) { VEvent VE = PQ.Pop() as VEvent; VDataNode[] CircleCheckList; if(VE is VDataEvent) { RootNode = VNode.ProcessDataEvent(VE as VDataEvent,RootNode,VG,VE.Y,out CircleCheckList); } else if(VE is VCircleEvent) { CurrentCircles.Remove(((VCircleEvent)VE).NodeN); if(!((VCircleEvent)VE).Valid) continue; RootNode = VNode.ProcessCircleEvent(VE as VCircleEvent,RootNode,VG,VE.Y,out CircleCheckList); } else throw new Exception("Got event of type "+VE.GetType().ToString()+"!"); foreach(VDataNode VD in CircleCheckList) { if(CurrentCircles.ContainsKey(VD)) { ((VCircleEvent)CurrentCircles[VD]).Valid=false; CurrentCircles.Remove(VD); } VCircleEvent VCE = VNode.CircleCheckDataNode(VD,VE.Y); if(VCE!=null) { PQ.Push(VCE); CurrentCircles[VD]=VCE; } } if(VE is VDataEvent) { Vector DP = ((VDataEvent)VE).DataPoint; foreach(VCircleEvent VCE in CurrentCircles.Values) { if(MathTools.Dist(DP[0],DP[1],VCE.Center[0],VCE.Center[1])1e-10) VCE.Valid = false; } } } return VG; } } } ```

If a file you wish to view isn't highlighted, and is a text file (not binary), please let us know and we'll add colourisation support for it.