Manhattan Voronoi

Some experiments on Manhattan Voronoi. I had to test 3 strategies in order to have descent quality and speed. I used Delaunay Triangulation, then Brep intersections (with pyramids) then the last using RTree and CreateBooleanIntersection

Manahattan Voronoi







22 Likes

may i ask what this was/is aimed for?

1 Like

It is another type of Voronoi, named Manhattan voronoi because it is good at creating city type plots.

https://www.shadertoy.com/view/ltV3R3

2 Likes

It is just a different distance used to calculate this Voronoi Tessellation, here the Manhattan distance between 2 points in plane XY P1(x1, y1) , P2(x2, y2)
d(P1, P2) = |x1-x2|+|y1-y2|

Classical distance is
d(P1, P2) = (|x1-x2|²+|y1-y2|²)^0.5

For the use it is a Tessellation, mainly for design.

2 Likes

Might add it’s called Manhattan as the distance is calculated in an orthogonal manner, so it’s like going from one point to another in Manhattan (where the city plan is grid-like)

1 Like

Nice! Would be interesting to see this on meshes, somehow applying it for planarization similar to this : https://github.com/formateng/giraffe

Wouldn’t be that hard to use them for generative floor plans also.

@Petras_Vestartas and @felipegutierrezduque i don’t see many uses as it is just a slightly different tessellations from classical Voronoi. Meanwhile The Manhattan and Chebychev distance are already in Grasshopper.

Voronoi => intersection of cones


Manhattan => Intersection of pyramids


Chebychev => Intersection of pyramids


Equilateral => Intersection of cubes (rotated 45° then Atan(1/sqrt(2))


19 Likes

I’m trying to understand too and I’m not completely clear

It’s possible to align every “pyramid” to a local plane?

As it is a 3d object the answer is yes. I included rotation and non uniform scale to make anisotropy. But as the result must be projected in a plane I don’t add rotation in X or Y.

Dont know why but all these make me feel satisfied…

Hello Laurent,
I didn’t know Manhattan distance Voronoi can be created using pyramide intersections, thank you for the hint!

Is Delaunay triangulation method geometrically correct or does it only produce similar looking results? Would you share your script with me or point me in a direction to help me create my own?

Cheers
Jonas

Hi. @laurent_delrieu I been wondering if you have an example file for the Manhattan Voronoi.
I believe it can be use for splitting building floor plans, and to get a better “equal” area partitions.

Hi Laurent,

Do you use only components to get the results or make use of codes to achieve it? Could you give some tips to achieve the Manhattan or Chebychev form, that would help a lot. I just feel it’s some kind of inspiration for the next cyberpunk project.

Best,
Anti-Vis!on

Hello
I use code and Clipper. I never ended nicely the tool. If you want to redo the algorithm. This one work.
voronoi_ManhattanClipper.gh (12.4 KB)


To have the good dll
Right click on component “Manage Assemblies”
image

Go to Component Folder


Use the dll named ClipperTools.dll

using System;
using System.Collections;
using System.Collections.Generic;

using Rhino;
using Rhino.Geometry;

using Grasshopper;
using Grasshopper.Kernel;
using Grasshopper.Kernel.Data;
using Grasshopper.Kernel.Types;

using System.Linq;
using System.Windows.Forms;
using ClipperLib;


/// <summary>
/// This class will be instantiated on demand by the Script component.
/// </summary>
public class Script_Instance : GH_ScriptInstance
{
#region Utility functions
  /// <summary>Print a String to the [Out] Parameter of the Script component.</summary>
  /// <param name="text">String to print.</param>
  private void Print(string text) { /* Implementation hidden. */ }
  /// <summary>Print a formatted String to the [Out] Parameter of the Script component.</summary>
  /// <param name="format">String format.</param>
  /// <param name="args">Formatting parameters.</param>
  private void Print(string format, params object[] args) { /* Implementation hidden. */ }
  /// <summary>Print useful information about an object instance to the [Out] Parameter of the Script component. </summary>
  /// <param name="obj">Object instance to parse.</param>
  private void Reflect(object obj) { /* Implementation hidden. */ }
  /// <summary>Print the signatures of all the overloads of a specific method to the [Out] Parameter of the Script component. </summary>
  /// <param name="obj">Object instance to parse.</param>
  private void Reflect(object obj, string method_name) { /* Implementation hidden. */ }
#endregion

#region Members
  /// <summary>Gets the current Rhino document.</summary>
  private readonly RhinoDoc RhinoDocument;
  /// <summary>Gets the Grasshopper document that owns this script.</summary>
  private readonly GH_Document GrasshopperDocument;
  /// <summary>Gets the Grasshopper script component that owns this script.</summary>
  private readonly IGH_Component Component;
  /// <summary>
  /// Gets the current iteration count. The first call to RunScript() is associated with Iteration==0.
  /// Any subsequent call within the same solution will increment the Iteration count.
  /// </summary>
  private readonly int Iteration;
#endregion

  /// <summary>
  /// This procedure contains the user code. Input parameters are provided as regular arguments,
  /// Output parameters as ref arguments. You don't have to assign output parameters,
  /// they will have a default value.
  /// </summary>
  private void RunScript(List<Point3d> lst_points, List<double> lst_weights, bool automaticSize, double size, bool closestNumber, int number, double scale, ref object A, ref object B)
  {

    A = Manhattan(lst_points, lst_weights, automaticSize, size, closestNumber, number, scale);

    //B = BooleanIntersection(lst_polylines, scale);

  }

  // <Custom additional code> 
  public static List<IntPoint> GeoPoints3dToIntPoints(Polyline poyline, double scale)
  {
    List<IntPoint> lst_IntPoints = new List<IntPoint>();
    foreach (Point3d point in poyline)
    {
      IntPoint ptInt = new IntPoint(point.X * scale, point.Y * scale);
      lst_IntPoints.Add(ptInt);
    }
    return lst_IntPoints;
  }


  public static List<IntPoint> BooleanIntersection(List<IntPoint> pl_1, List<IntPoint> pl_2)
  {
    Clipper cp = new Clipper();

    cp.AddPath(pl_1, PolyType.ptSubject, true);
    cp.AddPath(pl_2, PolyType.ptClip, true);
    List<List<IntPoint>> lst_polygons = new List<List<IntPoint>>();


    List<IntPoint> output = pl_1;
    if (cp.Execute(ClipType.ctIntersection, lst_polygons, PolyFillType.pftNonZero, PolyFillType.pftNonZero))
    {
      if (lst_polygons.Count > 0)
      {
        output = lst_polygons[0];
      }
    }
    return output;
  }

  public static Polyline BooleanIntersection(Polyline pl_1, Polyline pl_2, double arg_scale)
  {
    Clipper cp = new Clipper();

    cp.AddPath(GeoPoints3dToIntPoints(pl_1, arg_scale), PolyType.ptSubject, true);
    cp.AddPath(GeoPoints3dToIntPoints(pl_2, arg_scale), PolyType.ptClip, true);
    List<List<IntPoint>> lst_polygons = new List<List<IntPoint>>();
    cp.Execute(ClipType.ctIntersection, lst_polygons, PolyFillType.pftNonZero, PolyFillType.pftNonZero);

    List<Polyline> lst_polylines_out = new List<Polyline>();

    foreach (List<IntPoint> lst_intPoints in lst_polygons)
    {
      //On ajoute chaque liste de points du polygone fermé dans une liste (
      lst_polylines_out.Add(new Polyline(ClipperClosedPolygonToRealGeopPoint3d(lst_intPoints, arg_scale)));
    }
    return lst_polylines_out[0];
  }

  public static Polyline BooleanIntersection(List<Polyline> lst_polylines, double arg_scale)
  {
    List<IntPoint> lst_pl = GetClipperClosedPolygon(lst_polylines[0].ToList(), arg_scale);

    for (int i = 1; i < lst_polylines.Count; i++)
    {
      lst_pl = BooleanIntersection(lst_pl, GetClipperClosedPolygon(lst_polylines[i].ToList(), arg_scale));
    }
    return new Polyline(ClipperClosedPolygonToReal(lst_pl, arg_scale));
  }

  public static List<Point3d> ClipperClosedPolygonToRealGeopPoint3d(List<IntPoint> lst_intPoint, double arg_scale)
  {
    List<Point3d> lst_outerPoints = new List<Point3d>();
    foreach (IntPoint point in lst_intPoint)
    {
      lst_outerPoints.Add(new Point3d((double) point.X / arg_scale, (double) point.Y / arg_scale, 0.0));
    }
    lst_outerPoints.Add(lst_outerPoints[0]);//On ferme le polygone
    return lst_outerPoints;
  }

  public List<Point3d> PolylineCircle(Point3d arg_point, double arg_radius)
  {
    List<Point3d> output = new  List<Point3d>();
    double step = Math.PI / 10;
    for (double angle = 0; angle < Math.PI * 2.0; angle += step)
    {
      Point3d pt = new Point3d(arg_point.X + arg_radius * Math.Cos(angle), arg_point.Y + arg_radius * Math.Sin(angle), 0.0);
      output.Add(pt);
    }
    return output;
  }
  public static List<Point3d> ClipperClosedPolygonToReal(List<IntPoint> lst_intPoint, double arg_scale)
  {
    List<Point3d> lst_outerPoints = new List<Point3d>();
    foreach (IntPoint point in lst_intPoint)
    {
      lst_outerPoints.Add(new Point3d((double) point.X / arg_scale, (double) point.Y / arg_scale, 0.0));
    }
    lst_outerPoints.Add(lst_outerPoints[0]);//On ferme le polygone
    return lst_outerPoints;
  }
  public static List<IntPoint> GetClipperClosedPolygon(List<Point3d> lst_points, double arg_scale)
  {
    List<IntPoint> output = new List<IntPoint>();
    foreach (Point3d point in lst_points)
    {
      output.Add(new IntPoint(point.X * arg_scale, point.Y * arg_scale));
    }
    return output;
  }
  /// <summary>
  ///
  /// </summary>
  /// <param name="arg_lst_points"></param>
  /// <param name="arg_automaticSize"></param>
  /// <param name="arg_size"></param>
  /// <param name="arg_closestNumber"></param>
  /// <param name="arg_maxCount"></param>
  /// <param name="scale"></param>
  /// <returns></returns>
  public List<Polyline> Manhattan(List<Point3d> arg_lst_points, List<double> arg_lst_weights, bool arg_automaticSize, double arg_size, bool arg_closestNumber, int arg_maxCount, double scale)
  {
    //Bounding box to calculate the size of the bounding box
    BoundingBox bb = new BoundingBox(arg_lst_points);
    if (arg_automaticSize)
    {
      arg_size = bb.Diagonal.Length / 2;
    }
    IEnumerable<int[]> ie_tab_connectivity;

    //Use of RTree closest point by distance
    if (arg_closestNumber)
    {
      ie_tab_connectivity = RTree.Point3dKNeighbors(arg_lst_points, arg_lst_points, Math.Min(arg_lst_points.Count, arg_maxCount));
    }
      //Use of RTree number of closest point
    else
    {
      ie_tab_connectivity = RTree.Point3dClosestPoints(arg_lst_points, arg_lst_points, arg_size / 2);
    }
    List<int[]> lst_tab_connectivity = ie_tab_connectivity.ToList();
    List<Polyline> lst_outputCurves = new List<Polyline>();
    for (int vertexIndex = 0; vertexIndex < arg_lst_points.Count; vertexIndex++)
    {
      int[] tab_connectivity = lst_tab_connectivity[vertexIndex];
      List<Polyline> lst_curves = new List<Polyline>();
      //For all points aroung the point except itself
      foreach(int vertexIndexCutter in  tab_connectivity)
      {
        if (vertexIndex != vertexIndexCutter)
        {
          //Add closed curves
          lst_curves.Add(ManhattanPolylineBetweenPoints(arg_lst_points[vertexIndex], arg_lst_weights[vertexIndex], arg_lst_points[vertexIndexCutter], arg_lst_weights[vertexIndexCutter], arg_size));
        }
      }
      //Calculate the intersections
      if (lst_curves != null && lst_curves.Count > 0)
      {
        lst_outputCurves.Add(BooleanIntersection(lst_curves, scale));
      }
    }
    return lst_outputCurves;
  }

  /// <summary>
  /// Method to generate Manhattan Voronoi in plane XY
  /// </summary>
  /// <param name="arg_lst_points">List of points in XY plane</param>
  /// <param name="arg_automaticSize">If true, automatic will be used</param>
  /// <param name="arg_size">If arg_automaticSize is set to false this parameter is used to make the cells</param>
  /// <param name="arg_closestNumber">If true, the method cloest number is used in RTree</param>
  /// <param name="arg_maxCount">Number of closest points used</param>
  /// <returns></returns>
  public List<Curve> Manhattan(List < Point3d > arg_lst_points, List < double > arg_lst_weights, bool arg_automaticSize, double arg_size, bool arg_closestNumber, int arg_maxCount)
  {

    //Bounding box to calculate the size of the bounding box
    BoundingBox bb = new BoundingBox(arg_lst_points);
    if (arg_automaticSize)
    {
      arg_size = bb.Diagonal.Length / 2;
    }
    double tolerance = 1e-4;
    IEnumerable<int[]> ie_tab_connectivity;

    //Use of RTree closest point by distance
    if (arg_closestNumber)
    {
      ie_tab_connectivity = RTree.Point3dKNeighbors(arg_lst_points, arg_lst_points, Math.Min(arg_lst_points.Count, arg_maxCount));
    }
      //Use of RTree number of closest point
    else
    {

      ie_tab_connectivity = RTree.Point3dClosestPoints(arg_lst_points, arg_lst_points, arg_size / 2);
    }
    List<int[]> lst_tab_connectivity = ie_tab_connectivity.ToList();
    List<Curve> lst_outputCurves = new List<Curve>();
    for (int vertexIndex = 0; vertexIndex < arg_lst_points.Count; vertexIndex++)
    {
      int[] tab_connectivity = lst_tab_connectivity[vertexIndex];
      List<Curve> lst_curves = new List<Curve>();
      //For all points aroung the point except itself
      foreach(int vertexIndexCutter in  tab_connectivity)
      {
        if (vertexIndex != vertexIndexCutter)
        {
          //Add closed curves
          lst_curves.Add(ManhattanCurveBetweenPoints(arg_lst_points[vertexIndex], arg_lst_weights[vertexIndex], arg_lst_points[vertexIndexCutter], arg_lst_weights[vertexIndexCutter], arg_size));
        }
      }
      //Calculate the intersections
      if (lst_curves != null && lst_curves.Count > 0)
      {
        Curve intersection = lst_curves[0];
        for (int i = 1; i < lst_curves.Count; i++)
        {
          Curve[] curves = Curve.CreateBooleanIntersection(intersection, lst_curves[i], tolerance);

          if (curves != null && curves.Length > 0)
          {
            intersection = curves[0];
          }
        }
        lst_outputCurves.Add(intersection);
      }
    }

    return lst_outputCurves;
  }

  public Polyline ManhattanPolylineBetweenPoints(Point3d arg_point1, double weight1, Point3d arg_point2, double weight2, double arg_length)
  {
    List<Point3d> lst_points = new List<Point3d> ();

    Line outputLine = Line.Unset;

    Point3d middle = (arg_point1 * weight1 + arg_point2 * weight2) / (weight1 + weight2);
    Point3d manhattanPoint0 = Point3d.Unset;
    Point3d manhattanPoint1 = middle;
    Point3d manhattanPoint2 = middle;
    Point3d manhattanPoint3 = Point3d.Unset;

    double dx = arg_point2.X - arg_point1.X;
    double dy = arg_point2.Y - arg_point1.Y;

    //Middle curve calculation
    if (Math.Abs(dx * dy) > 1e-12)
    {
      double minWidth = Math.Min(Math.Abs(dx), Math.Abs(dy));
      Vector3d direction;

      if (dx * dy > 0)
      {
        direction = new Vector3d(1, -1, 0);
      }
      else
      {
        direction = new Vector3d(1, 1, 0);
      }
      manhattanPoint1 = middle + direction * minWidth / 2;
      manhattanPoint2 = middle - direction * minWidth / 2;
      Vector3d dir = manhattanPoint2 - manhattanPoint1;

      //Just diagonal line
      if (Math.Abs(Math.Abs(dx) - Math.Abs(dy)) < 1e-12)
      {
        dir.Unitize();
      }
      else
      {
        if (Math.Abs(dx) < Math.Abs(dy))
        {
          dir.Y = 0;
          dir.Unitize();

        }
        else
        {
          dir.X = 0;
          dir.Unitize();
        }

      }
      manhattanPoint0 = manhattanPoint1 - dir * arg_length;
      manhattanPoint3 = manhattanPoint2 + dir * arg_length;
    }
    else
    {

      if (Math.Abs(dx) > 1e-6)
      {
        manhattanPoint0 = middle - Vector3d.YAxis * arg_length;
        manhattanPoint3 = middle + Vector3d.YAxis * arg_length;
      }
      else
      {
        manhattanPoint0 = middle - Vector3d.XAxis * arg_length;
        manhattanPoint3 = middle + Vector3d.XAxis * arg_length;
      }
    }

    List<Point3d> lst_points_pl = new List<Point3d>();
    lst_points_pl.Add(manhattanPoint0);
    lst_points_pl.Add(manhattanPoint1);
    lst_points_pl.Add(manhattanPoint2);
    lst_points_pl.Add(manhattanPoint3);

    return GetClosedPolyline(lst_points_pl, arg_point1, arg_length * 2);
  }

  /// <summary>
  /// Make a cell between 2 points,
  /// </summary>
  /// <param name="arg_point1">the cell will contain this point</param>
  /// <param name="arg_point2"></param>
  /// <param name="arg_length">Length of 2 laterals lines</param>
  /// <param name="arg_angle">Angle of rotation</param>
  /// <returns>A closed curve</returns>
  public Curve ManhattanCurveBetweenPoints(Point3d arg_point1, double weight1, Point3d arg_point2, double weight2, double arg_length)
  {
    List<Curve> lst_curves = new List<Curve> ();

    Line outputLine = Line.Unset;

    Point3d middle = (arg_point1 * weight1 + arg_point2 * weight2) / (weight1 + weight2);

    double dx = arg_point2.X - arg_point1.X;
    double dy = arg_point2.Y - arg_point1.Y;

    //Middle curve calculation
    if (Math.Abs(dx * dy) > 1e-12)
    {
      double minWidth = Math.Min(Math.Abs(dx), Math.Abs(dy));
      Vector3d direction;

      if (dx * dy > 0)
      {
        direction = new Vector3d(1, -1, 0);
      }
      else
      {
        direction = new Vector3d(1, 1, 0);
      }
      Point3d manhattanPoint1 = middle + direction * minWidth / 2;
      Point3d manhattanPoint2 = middle - direction * minWidth / 2;
      outputLine = new Line(manhattanPoint1, manhattanPoint2);
      lst_curves.Add(outputLine.ToNurbsCurve());

      Vector3d dir = manhattanPoint2 - manhattanPoint1;

      //Just diagonal line
      if (Math.Abs(Math.Abs(dx) - Math.Abs(dy)) < 1e-12)
      {
        dir.Unitize();
      }
      else
      {
        if (Math.Abs(dx) < Math.Abs(dy))
        {
          dir.Y = 0;
          dir.Unitize();

        }
        else
        {
          dir.X = 0;
          dir.Unitize();
        }

      }
      Line line1 = new Line(manhattanPoint2, manhattanPoint2 + dir * arg_length);
      Line line2 = new Line(manhattanPoint1, manhattanPoint1 - dir * arg_length);
      lst_curves.Add(line1.ToNurbsCurve());
      lst_curves.Add(line2.ToNurbsCurve());
    }
    else
    {

      if (Math.Abs(dx) > 1e-6)
      {
        Line line1 = new Line(middle - Vector3d.YAxis * arg_length, middle + Vector3d.YAxis * arg_length);
        lst_curves.Add(line1.ToNurbsCurve());
      }
      else
      {
        Line line1 = new Line(middle - Vector3d.XAxis * arg_length, middle + Vector3d.XAxis * arg_length);
        lst_curves.Add(line1.ToNurbsCurve());
      }

    }

    Curve[] tab_curves = Curve.JoinCurves(lst_curves);
    Curve cc = null;
    if (tab_curves != null && tab_curves.Length > 0)
    {
      cc = GetClosedCurve(tab_curves[0], arg_point1, arg_length * 2);
    }

    return cc;
  }

  public Polyline GetClosedPolyline(List < Point3d > lst_points, Point3d p, double l)
  {
    List<Point3d> lst_pointsNew = new List<Point3d>();

    Vector3d direction = p - (lst_points[1] + lst_points[2]) / 2;
    direction.Unitize();
    lst_pointsNew.Add(lst_points[0] + direction * l);
    lst_pointsNew.Add(lst_points[0]);
    lst_pointsNew.Add(lst_points[1]);
    lst_pointsNew.Add(lst_points[2]);
    lst_pointsNew.Add(lst_points[3]);
    lst_pointsNew.Add(lst_points[3] + direction * l);
    lst_pointsNew.Add(lst_points[0] + direction * l);

    return new Polyline(lst_pointsNew);
  }

  public Curve GetClosedCurve(Curve c, Point3d p, double l)
  {
    Vector3d direction = p - c.PointAtNormalizedLength(0.5);
    direction.Unitize();
    Point3d p1 = c.PointAtEnd;
    Point3d p2 = c.PointAtEnd + direction * l;

    Point3d p4 = c.PointAtStart;
    Point3d p3 = c.PointAtStart + direction * l;

    Line l1 = new Line(p1, p2);
    Line l2 = new Line(p2, p3);
    Line l3 = new Line(p3, p4);
    List<Curve> curves = new List<Curve>();
    curves.Add(c);
    curves.Add(l1.ToNurbsCurve());
    curves.Add(l2.ToNurbsCurve());
    curves.Add(l3.ToNurbsCurve());

    Curve[] cs = Curve.JoinCurves(curves);


    return cs[0];

  }


  // </Custom additional code> 
}
2 Likes

Thank you so much for your patience, I’ve been playing these components for two days at my home and office, but never got the final result. :sweat_smile: :sweat_smile:

It is not easy, you will see that there are some bugs in my code, Clipper is very robust but there is still place to improvements in my code. You will see many parameters they are just useful to make the code work.
here some bug in a Unicorn with


Bug because there are some duplicates !!! Far more better now

If you want simple things use Mesh Pyramids and colors. I am quite sure this could lead to very robust but less precise way of making Voronoi.

2 Likes