none
Need a free triangulation application or library... RRS feed

  • Question

  • I need a free triangulation Application or library which take some 3D points and segments as input and return me a triangulation. It will be better if it is a library.
     
    I need that just for a prove of concept. If that's concept worked I can write my own Triangulation application. But till that I need a simple triangulation application ASAP.
     
    Note:-
    >> Segments as input I need them conserved.
    >> 3D(x,y,z) triangulation will be better but for start I can some how work with 2D(x,y) triangulation as well. 
    >> For start I can try with only Point input and triangulation as output type also. If it takes Segment as input and conserve them, That will be great.

    Saurabh

    Sunday, March 10, 2013 10:11 AM

Answers

  • This procedure calculates the Delauney Triangulation of a set of geometry Points in SQL Server (from "Pro Spatial with SQL Server 2012"):

    public partial class Triangulation
      {
    
        [Microsoft.SqlServer.Server.SqlProcedure]
        public static void GeometryTriangulate3d(SqlGeometry MultiPoint)
        {
          // Retrieve the SRID
          int srid = (int)MultiPoint.STSrid;
    
          // Check valid input
          if (!(MultiPoint.STGeometryType() == "MULTIPOINT" && MultiPoint.STNumPoints() > 3))
          {
            throw new ArgumentException("Input must be a MultiPoint containing at least three points");
          }
    
          // Initialise a list of vertices
          List<SimplePoint3d> Vertices = new List<SimplePoint3d>();
          // Add all the original supplied points
          for (int i = 1; i <= MultiPoint.STNumPoints(); i++)
          {
            SqlGeometry p = MultiPoint.STPointN(i);
            SimplePoint3d Point = new SimplePoint3d((double)p.STX, (double)p.STY, p.HasZ ? (double)p.Z : 0);
            // MultiPoints can contain the same point twice, but this messes up Delauney
            if (!Vertices.Contains(Point))
            {
              Vertices.Add(Point);
            }
          }
          int numPoints = Vertices.Count;
          // Sort the list so that points sweep from left - right
          Vertices.Sort();
    
          // Calculate the "supertriangle" that encompasses the pointset
          SqlGeometry Envelope = MultiPoint.STEnvelope();
          // Width
          double dx = (double)(Envelope.STPointN(2).STX - Envelope.STPointN(1).STX);
          // Height 
          double dy = (double)(Envelope.STPointN(4).STY - Envelope.STPointN(1).STY);
          // Maximum dimension
          double dmax = (dx > dy) ? dx : dy;
          // Centre
          double avgx = (double)Envelope.STCentroid().STX;
          double avgy = (double)Envelope.STCentroid().STY;
          // Create the points at corners of the supertriangle
          SimplePoint3d a = new SimplePoint3d(avgx - 2 * dmax, avgy - dmax, 0);
          SimplePoint3d b = new SimplePoint3d(avgx + 2 * dmax, avgy - dmax, 0);
          SimplePoint3d c = new SimplePoint3d(avgx, avgy + 2 * dmax, 0);
    
          // Add the supertriangle vertices to the end of the vertex array
          Vertices.Add(a);
          Vertices.Add(b);
          Vertices.Add(c);
    
          double radius;
          SimplePoint3d circumcentre;
          CalculateCircumcircle3d(a, b, c, out circumcentre, out radius);
    
          // Create a triangle from the vertices
          SimpleTriangle3d SuperTriangle = new SimpleTriangle3d(numPoints, numPoints + 1, numPoints + 2, circumcentre, radius);
    
          // Add the supertriangle to the list of triangles
          List<SimpleTriangle3d> Triangles = new List<SimpleTriangle3d>();
          Triangles.Add(SuperTriangle);
    
          List<SimpleTriangle3d> CompletedTriangles = new List<SimpleTriangle3d>();
    
          // Loop through each point
          for (int i = 0; i < numPoints; i++)
          {
            // Initialise the edge buffer
            List<int[]> Edges = new List<int[]>();
    
            // Loop through each triangle
            for (int j = Triangles.Count - 1; j >= 0; j--)
            {
              // If the point lies within the circumcircle of this triangle
              if (Distance3d(Triangles[j].circumcentre, Vertices[i]) < Triangles[j].radius)
              {
                // Add the triangle edges to the edge buffer
                Edges.Add(new int[] {Triangles[j].a, Triangles[j].b});
                Edges.Add(new int[] { Triangles[j].b, Triangles[j].c });
                Edges.Add(new int[] { Triangles[j].c, Triangles[j].a });
                
                // Remove this triangle from the list
                Triangles.RemoveAt(j);
              }
    
              // If this triangle is complete
              else if (Vertices[i].x > Triangles[j].circumcentre.x + Triangles[j].radius)
              {
                {
                  CompletedTriangles.Add(Triangles[j]);
                }
                Triangles.RemoveAt(j);
              }
              
            }
    
            // Remove duplicate edges
            for (int j = Edges.Count - 1; j > 0; j--)
            {
              for (int k = j - 1; k >= 0; k--)
              {
                // Compare if this edge match in either direction
                if (Edges[j][0].Equals(Edges[k][1]) && Edges[j][1].Equals(Edges[k][0]))
                {
                  // Remove both duplicates
                  Edges.RemoveAt(j);
                  Edges.RemoveAt(k);
    
                 // We've removed an item from lower down the list than where j is now, so update j
                  j--;
                  break;
                }
              }
            }
    
            // Create new triangles for the current point
            for (int j = 0; j < Edges.Count; j++)
            {
              CalculateCircumcircle3d(Vertices[Edges[j][0]], Vertices[Edges[j][1]], Vertices[i], out circumcentre, out radius);
              SimpleTriangle3d T = new SimpleTriangle3d(Edges[j][0], Edges[j][1], i, circumcentre, radius);
              Triangles.Add(T);
            }
          }
    
          // We've finished triangulation. Move any remaining triangles onto the completed list
          CompletedTriangles.AddRange(Triangles);
    
          // Define the metadata of the results column
          SqlMetaData metadata = new SqlMetaData("Triangle", SqlDbType.Udt, typeof(SqlGeometry));
    
          // Create a record based on this metadata
          SqlDataRecord record = new SqlDataRecord(metadata);
    
          // Send the results back to the client
          SqlContext.Pipe.SendResultsStart(record);
          foreach (SimpleTriangle3d Tri in CompletedTriangles)
          {
            // Check that this is a triangle formed only from vertices in the original multipoint
            // i.e. not from the vertices of the supertriangle.
            if (Tri.a < numPoints && Tri.b < numPoints && Tri.c < numPoints)
            {
              SqlGeometry triangle = Triangle3dFromPoints(Vertices[Tri.a], Vertices[Tri.b], Vertices[Tri.c], srid);
              record.SetValue(0, triangle);
              SqlContext.Pipe.SendResultsRow(record);
            }
          }
          SqlContext.Pipe.SendResultsEnd();
        }
      }

    Used as:

    SELECT @MultiPoint;

    EXEC dbo.GeometryTriangulate @MultiPoint;


    twitter: @alastaira blog: http://alastaira.wordpress.com/ | Pro Spatial with SQL Server 2012

    Wednesday, March 13, 2013 8:37 PM
    Moderator
  • I wrote a blog post on 2D triangulation on the surface of spherical globes (maps) here: http://rbrundritt.wordpress.com/2008/10/14/triangulation/ which is a good starting point.

    http://rbrundritt.wordpress.com

    Wednesday, March 13, 2013 6:47 PM

All replies

  • Saurabh --
     
    I believe you have erroneously posted your question in a public user forum dedicated to Project Online, a feature of Office 365 SharePoint used for managing projects.  I would recommend you find a more relevant user forum for your question about triangulation software.  Hope this helps.
     

    Dale A. Howard [MVP]
    VP of Educational Services
    msProjectExperts
    http://www.msprojectexperts.com
    http://www.projectserverexperts.com
    "We write the books on Project Server"

    • Proposed as answer by Julie Sheets Monday, March 11, 2013 12:54 PM
    • Unproposed as answer by Ricky_Brundritt Wednesday, March 13, 2013 6:45 PM
    Monday, March 11, 2013 11:51 AM
  • I wrote a blog post on 2D triangulation on the surface of spherical globes (maps) here: http://rbrundritt.wordpress.com/2008/10/14/triangulation/ which is a good starting point.

    http://rbrundritt.wordpress.com

    Wednesday, March 13, 2013 6:47 PM
  • This procedure calculates the Delauney Triangulation of a set of geometry Points in SQL Server (from "Pro Spatial with SQL Server 2012"):

    public partial class Triangulation
      {
    
        [Microsoft.SqlServer.Server.SqlProcedure]
        public static void GeometryTriangulate3d(SqlGeometry MultiPoint)
        {
          // Retrieve the SRID
          int srid = (int)MultiPoint.STSrid;
    
          // Check valid input
          if (!(MultiPoint.STGeometryType() == "MULTIPOINT" && MultiPoint.STNumPoints() > 3))
          {
            throw new ArgumentException("Input must be a MultiPoint containing at least three points");
          }
    
          // Initialise a list of vertices
          List<SimplePoint3d> Vertices = new List<SimplePoint3d>();
          // Add all the original supplied points
          for (int i = 1; i <= MultiPoint.STNumPoints(); i++)
          {
            SqlGeometry p = MultiPoint.STPointN(i);
            SimplePoint3d Point = new SimplePoint3d((double)p.STX, (double)p.STY, p.HasZ ? (double)p.Z : 0);
            // MultiPoints can contain the same point twice, but this messes up Delauney
            if (!Vertices.Contains(Point))
            {
              Vertices.Add(Point);
            }
          }
          int numPoints = Vertices.Count;
          // Sort the list so that points sweep from left - right
          Vertices.Sort();
    
          // Calculate the "supertriangle" that encompasses the pointset
          SqlGeometry Envelope = MultiPoint.STEnvelope();
          // Width
          double dx = (double)(Envelope.STPointN(2).STX - Envelope.STPointN(1).STX);
          // Height 
          double dy = (double)(Envelope.STPointN(4).STY - Envelope.STPointN(1).STY);
          // Maximum dimension
          double dmax = (dx > dy) ? dx : dy;
          // Centre
          double avgx = (double)Envelope.STCentroid().STX;
          double avgy = (double)Envelope.STCentroid().STY;
          // Create the points at corners of the supertriangle
          SimplePoint3d a = new SimplePoint3d(avgx - 2 * dmax, avgy - dmax, 0);
          SimplePoint3d b = new SimplePoint3d(avgx + 2 * dmax, avgy - dmax, 0);
          SimplePoint3d c = new SimplePoint3d(avgx, avgy + 2 * dmax, 0);
    
          // Add the supertriangle vertices to the end of the vertex array
          Vertices.Add(a);
          Vertices.Add(b);
          Vertices.Add(c);
    
          double radius;
          SimplePoint3d circumcentre;
          CalculateCircumcircle3d(a, b, c, out circumcentre, out radius);
    
          // Create a triangle from the vertices
          SimpleTriangle3d SuperTriangle = new SimpleTriangle3d(numPoints, numPoints + 1, numPoints + 2, circumcentre, radius);
    
          // Add the supertriangle to the list of triangles
          List<SimpleTriangle3d> Triangles = new List<SimpleTriangle3d>();
          Triangles.Add(SuperTriangle);
    
          List<SimpleTriangle3d> CompletedTriangles = new List<SimpleTriangle3d>();
    
          // Loop through each point
          for (int i = 0; i < numPoints; i++)
          {
            // Initialise the edge buffer
            List<int[]> Edges = new List<int[]>();
    
            // Loop through each triangle
            for (int j = Triangles.Count - 1; j >= 0; j--)
            {
              // If the point lies within the circumcircle of this triangle
              if (Distance3d(Triangles[j].circumcentre, Vertices[i]) < Triangles[j].radius)
              {
                // Add the triangle edges to the edge buffer
                Edges.Add(new int[] {Triangles[j].a, Triangles[j].b});
                Edges.Add(new int[] { Triangles[j].b, Triangles[j].c });
                Edges.Add(new int[] { Triangles[j].c, Triangles[j].a });
                
                // Remove this triangle from the list
                Triangles.RemoveAt(j);
              }
    
              // If this triangle is complete
              else if (Vertices[i].x > Triangles[j].circumcentre.x + Triangles[j].radius)
              {
                {
                  CompletedTriangles.Add(Triangles[j]);
                }
                Triangles.RemoveAt(j);
              }
              
            }
    
            // Remove duplicate edges
            for (int j = Edges.Count - 1; j > 0; j--)
            {
              for (int k = j - 1; k >= 0; k--)
              {
                // Compare if this edge match in either direction
                if (Edges[j][0].Equals(Edges[k][1]) && Edges[j][1].Equals(Edges[k][0]))
                {
                  // Remove both duplicates
                  Edges.RemoveAt(j);
                  Edges.RemoveAt(k);
    
                 // We've removed an item from lower down the list than where j is now, so update j
                  j--;
                  break;
                }
              }
            }
    
            // Create new triangles for the current point
            for (int j = 0; j < Edges.Count; j++)
            {
              CalculateCircumcircle3d(Vertices[Edges[j][0]], Vertices[Edges[j][1]], Vertices[i], out circumcentre, out radius);
              SimpleTriangle3d T = new SimpleTriangle3d(Edges[j][0], Edges[j][1], i, circumcentre, radius);
              Triangles.Add(T);
            }
          }
    
          // We've finished triangulation. Move any remaining triangles onto the completed list
          CompletedTriangles.AddRange(Triangles);
    
          // Define the metadata of the results column
          SqlMetaData metadata = new SqlMetaData("Triangle", SqlDbType.Udt, typeof(SqlGeometry));
    
          // Create a record based on this metadata
          SqlDataRecord record = new SqlDataRecord(metadata);
    
          // Send the results back to the client
          SqlContext.Pipe.SendResultsStart(record);
          foreach (SimpleTriangle3d Tri in CompletedTriangles)
          {
            // Check that this is a triangle formed only from vertices in the original multipoint
            // i.e. not from the vertices of the supertriangle.
            if (Tri.a < numPoints && Tri.b < numPoints && Tri.c < numPoints)
            {
              SqlGeometry triangle = Triangle3dFromPoints(Vertices[Tri.a], Vertices[Tri.b], Vertices[Tri.c], srid);
              record.SetValue(0, triangle);
              SqlContext.Pipe.SendResultsRow(record);
            }
          }
          SqlContext.Pipe.SendResultsEnd();
        }
      }

    Used as:

    SELECT @MultiPoint;

    EXEC dbo.GeometryTriangulate @MultiPoint;


    twitter: @alastaira blog: http://alastaira.wordpress.com/ | Pro Spatial with SQL Server 2012

    Wednesday, March 13, 2013 8:37 PM
    Moderator
  • Thanks all for your Reply.

    I got a code of jonathan shewchuk on constrained delaunay triangulation in 2D  that helped me a lot.


    Saurabh

    • Marked as answer by Saurabh Saini Thursday, March 14, 2013 12:43 PM
    • Unmarked as answer by Ricky_Brundritt Thursday, March 14, 2013 2:12 PM
    • Marked as answer by Saurabh Saini Thursday, March 14, 2013 5:02 PM
    • Unmarked as answer by Saurabh Saini Thursday, March 14, 2013 5:04 PM
    Thursday, March 14, 2013 12:43 PM