Note: I could not find any code in C#. I found a few examples for the clustering in PHP and the Convex Hull code was in C++. I had to translate it to C#.
My latest project required the use of Google Maps API to fetch and display a KML file. “KML is a file format used to display geographic data in an earth browser, such as Google Earth, Google Maps, and Google Maps for mobile” [More on KML] The whole point of my project was to display “sites” on the map. Each site represents a customer’s location. The problem is that there could be any number of sites up to 50,000 and more. That many markers on a map is just worthless.
There are plug-ins that you can use to auto cluster markers on a map, but they all run client side. Client side trying to cluster 50,000 markers isn’t going to work. I decided to use the GGeoXml() function to request a KML file on the server. The KML file is generated when the user makes a request to the server. After the data comes back from the database it is processed, clustered and KML is generated for consumption by Google Maps.
The first step is to cluster the sites by distance. This gives us our grouping that we could display on the map, but since I can’t have any more than 100 markers on the map at any one time, I have to take it a step further. The seconds step is to take the new set of points (after clustering) and generate a polygon from them. The end result is a polygon that represent each cluster.
To start, we need a point class. Most examples you see will use integers but in the real world, longitude and latitude are doubles. We define out Point class as
Point.cs
public class Point
{
public Point(double x, double y)
{
this.X = x;
this.Y = y;
}
public double X { get; set; }
public double Y { get; set; }
public override string ToString()
{
return String.Format("{0} {1}", this.X, this.Y);
}
}
Next we define our ClusterPoint class that will hold the point where the cluster marker will be and a list of all points within that cluster. We do this so we can generate the polygon later on.
ClusterPoint.cs
public class ClusterPoint
{
public Point Point { get; set; }
public List<Point> PointsInCluster { get; set; }
public int TotalPointsInCluster { get; set; }
}
To do the actual clustering, we use two algorithms. The first calculates a distance value based on the zoom level and then loops through each point in the list. For each point, it gets removed from the list since it’s already being processed and then an inner loop compares the remaining points to the current point by calculating the distance between each point. If the calculated distance is within the threshold it gets grouped and removed from the list.
GeoClustering.cs
public class GeoClustering
{
public List<GeoCluster> Cluster(List<Point> points, int Zoom)
{
return Cluster(points, Zoom, false);
}
public List<ClusterPoint> Cluster(List<Point> points, int Zoom, bool applyConvexHull)
{
List<ClusterPoint> clusterMarkers = new List<ClusterPoint>();
double distance = (10000000 >> Zoom) / 100000.00;
while (points.Count > 0)
{
Point marker = points.Last();
points.Remove(marker);
ClusterPoint cp = new ClusterPoint() { Point = marker, PointsInCluster = new List<Point>() { marker }, TotalPointsInCluster = 0 };
bool group = false;
for (int i = 0; i < points.Count; i++)
{
double pixels = Math.Abs(marker.X - points[i].X) + Math.Abs(marker.Y - points[i].Y);
if (pixels <= distance)
{
group = true;
cp.PointsInCluster.Add(points[i]);
points.RemoveAt(i);
i--;
cp.TotalPointsInCluster++;
}
}
if (group)
clusterMarkers.Add(cp);
}
if (applyConvexHull)
foreach (ClusterPoint c in clusterMarkers)
c.PointsInCluster = ProcessConvexHull(c.PointsInCluster);
return clusterMarkers;
}
public List<Point> ProcessConvexHull(List<Point> points)
{
ConvexHull ch = new ConvexHull();
return ch.ConvexHullPoints(points);
}
}
You’ll notice that we have options for processing Convex Hull. Convex Hull is the algorithm that we use in the next step to generate the polygon from a set of points. [More on Convex Hull]
How does it work? The process starts by taking a set of points and partitions them. Partitioning orders the points by X then gets the first and last point on the list (leftmost and rightmost points) and both are removed from the list. The remaining points are looped and the direction of each point is calculated based on the left and right points extracted earlier. Based on the direction value, the point is put into either the upper or lower partition (an array).
After partitioning we build the “hull”. This process used the upper and lower partitions we created previously. Each partition is processed separately. The points in the partition are looped and 1 by 1 they are added and compared. A calculation is done to determine if point B is on the inside or outside of a line between point A & C, if it is on the inside then it gets removed. If it’s on the outside then it’s left in. Point A,B & C are referenced as Points[Len-2], Points[Len] & Points[Len-1] respectively. There is a great animation demonstrating this process here. After the upper and lower partitions are processed (giving us 2 sets of points that we could generate polygons from), we then combine them. Notice that the UpperHull array is reversed before the merge. If the points were not reversed you would en up with an incorrect polygon being drawn.
ConvexHull.cs
public class ConvexHull
{
private Point _left;
private Point _right;
private List<Point> _upper;
private List<Point> _lower;
private List<Point> _points;
public ConvexHull() { }
/// <summary>
/// Takes a bulk set of points and returns points of outer edge (for drawing polygon)
/// </summary>
/// <param name="points">Expects Point(double x, double y) x = Lat, y = Lng</param>
/// <returns></returns>
public List<Point> ConvexHullPoints(List<Point> points)
{
if (points == null)
throw new ArgumentNullException("points");
if (points.Count < 3)
return points;
_points = points;
Partition();
return BuildHull();
}
private void Partition()
{
_points = _points.OrderBy(c => c.X).ToList();
_left = _points.First();
_right = _points.Last();
_points.Remove(_left);
_points.Remove(_right);
_upper = new List<Point>();
_lower = new List<Point>();
for (int i = 0; i < _points.Count; i++)
{
int dir = Direction(_left, _right, _points[i]);
if (dir < 0)
_upper.Add(_points[i]);
else
_lower.Add(_points[i]);
}
}
private int Direction(Point p0, Point p1, Point p2)
{
double res = (((p0.X - p1.X) * 1000) * ((p2.Y - p1.Y) * 1000)) - (((p2.X - p1.X) * 1000) * ((p0.Y - p1.Y) * 1000));
return Convert.ToInt32(res);
}
private List<Point> BuildHull()
{
List<Point> LowerHull;
List<Point> UpperHull;
BuildHalfHull(_lower, out LowerHull, 1);
BuildHalfHull(_upper, out UpperHull, -1);
UpperHull.Reverse();
LowerHull.AddRange(UpperHull);
return LowerHull;
}
private void BuildHalfHull(List<Point> input, out List<Point> output, int factor)
{
output = new List<Point>();
input.Add(_right);
output.Add(_left);
while (input.Count != 0)
{
output.Add(input.First());
input.Remove(input.First());
while (output.Count >= 3)
{
int end = output.Count - 1;
int res = (factor * Direction(output[end - 2], output[end], output[end - 1]));
if (res <= 0)
{
output.RemoveAt(end - 1);
}
else
break;
}
}
}
}
Putting it all together
var MyPoints = (from s in databaseResults select new Point(double.Parse(s.Lat), double.Parse(s.Lng))).OrderBy(c=>c.X).ThenBy(c=>c.Y).ToList();
GeoClustering gc = new GeoClustering();
var t = gc.Cluster(MyPoints, Zoom, true);
XDocument kml = new XDocument();
XElement doc = new XElement(XName.Get("Document"));
kml.Add(new XElement(XName.Get("kml", "http://earth.google.com/kml/2.2")));
foreach (ClusterPoint s in t)
{
XElement placemark = new XElement(XName.Get("Placemark"));
XElement name = new XElement(XName.Get("name"));
XElement desc = new XElement(XName.Get("description"));
XElement polygon = new XElement(XName.Get("Polygon"));
XElement altitudeMode = new XElement(XName.Get("altitudeMode"));
XElement outerBoundaryIs = new XElement(XName.Get("outerBoundaryIs"));
XElement LinearRing = new XElement(XName.Get("LinearRing"));
XElement coords = new XElement(XName.Get("coordinates"));
XElement style = new XElement(XName.Get("Style"));
XElement polystyle = new XElement(XName.Get("PolyStyle"));
XElement polycolor = new XElement(XName.Get("color"));
polycolor.Value = "8814B446";
polystyle.Add(polycolor);
style.Add(polystyle);
name.Value = "Group";
desc.Value = string.Format("There are {0} points in this cluster", s.TotalPointsInCluster);
StringBuilder sb = new StringBuilder();
foreach (Point p in s.PointsInCluster)
{
sb.AppendFormat("{0},{1},100 ", p.X, p.Y);
}
coords.Value = sb.ToString();
LinearRing.Add(coords);
outerBoundaryIs.Add(LinearRing);
polygon.Add(altitudeMode);
polygon.Add(outerBoundaryIs);
placemark.Add(name);
placemark.Add(desc);
placemark.Add(polygon);
placemark.Add(style);
doc.Add(placemark);
}
kml.Root.Add(doc);
return kml;
Microsoft MVP, MCTS, PostSharp MVP,
