Vision¶

Inspire a community of professionals skilled in spatiotemporal analysis for community health and safety.

Mission¶

Prototype an algorithm that places service drop poles near customers.

Owner¶

Roy Hyunjin Han

Context¶

Placing service drop poles is the first step to connecting a village. It is also possible that we can abstract this algorithm to similar problem definitions.

Timeframe¶

20170912-1715 - 20170913-1715: 1 day estimated

20170912-1715 - 20170919-2215: 1 week actual

Objectives¶

+ Consider different ways to implement algorithm
+ Implement at least one version

Log¶

20170912-2130 - 20170912-2200: 30 minutes

While walking home, I thought of a potential solution.

Let's take each point (or geometry) that represents a customer. Then we'll draw a buffer around the geometry. That buffer represents the locus of possible locations of a service drop pole for that customer.
Then let's take the intersection of all the buffers. If there are multiple buffers, there there are likely to be multiple areas of intersection. For each area of intersection, we'll count how many times it intersects with a buffer, which indicates the number of customers a service drop pole in that area will be able to reach.
Then we'll start with the intersection area that has the most number of overlapping buffers and put as many service drop poles in it until it satisfies its customers.
Then we'll go through the intersection areas in order, from the most number of overlapping buffers to the least until we've satisfied them all.
Each time we satisfy a customer, we'll remove it from the unsatisfied customers list.
Finally, we'll have the unsatisfied buffers, which we'll keep.
The end result is a set of geometries, preserving a degree of freedom on where exactly to place each service drop pole for future steps.

We should refine these thoughts in more precise pseudocode.

+ Draft potential solution in words

Start with a set of geometries that represent each customer.
Draw a buffer around each geometry that represents the possible locations of a service drop pole for that customer.
Compute the intersection of all buffers.
Split the buffer intersection into separate geometries which we will call buffer intersections.
For each buffer intersection, loop through the buffers and count the number of times the buffer intersection is contained in a buffer.
Sort the buffer intersections in order of most number of overlapping buffers.
For each buffer intersection, assign the minimum number of service drop poles that will cover the customers for that buffer intersection. We'll call these social service drop poles.
Remove each satisfied buffer from the list of unsatisfied buffers.
Assign a service drop pole to each unsatified buffer, which we will call lonely service drop poles.
Count the number of service drop poles.
Output the buffer intersections and unsatisfied buffers, which is our solution.

I think our metric for now is simply the number of service drop poles combined with the sum of distances from customer to service drop pole for each customer.

score = A * sum(customer_distances) - B * pole_count

where A is $/length and B is $/pole

+ Draft potential solution in pseudocode
+ Think of different ways to implement algorithm
+ Propose a metric

For the toy problem, we'll start with points, though we should confirm that we can get a buffer around any geometry.

In [ ]:

from shapely.geometry import Point
point = Point(0, 0)
point

In [ ]:

point.buffer(1)

In [ ]:

from shapely.geometry import LineString
line = LineString([(0, 0), (1, 1)])
line

In [ ]:

line.buffer(1)

In [ ]:

from shapely.geometry import Polygon
polygon = Polygon([(0, 0), (1, 1), (1, 0)])
polygon

In [ ]:

polygon.buffer(1)

+ Review use of buffers

In [ ]:

# Start with a set of geometries that represent each customer
from shapely.geometry import MultiPolygon, Polygon
customers = MultiPolygon([
    Polygon([(0, 0), (1, 1), (1, 0)]),
    Polygon([(0, 0.5), (0, 1.5), (1, 1.5)]),
    Polygon([(2.5, 0.5), (2.5, 1), (3, 1), (3, 0.5)]),   
])
customers

In [ ]:

# Draw a buffer around each geometry that represents
# possible locations of a service drop pole for that customer
maximum_drop_length = 0.5
customers_buffered = customers.buffer(maximum_drop_length)
customers_buffered

I think we want to keep each customer separate.

In [ ]:

buffered_customers = [x.buffer(maximum_drop_length) for x in customers]
MultiPolygon(buffered_customers)

In [ ]:

buffered_customers[0].intersection(buffered_customers[1])

Now we are trying to compute the intersection of all buffers. There seems to be one possible solution here: https://gis.stackexchange.com/questions/187402/how-to-find-the-intersection-areas-of-overlapping-buffer-zones-in-single-shapefi/187499#187499

In [ ]:

from shapely.ops import unary_union, polygonize
pieces = list(polygonize(unary_union(buffered_customers)))
print(len(pieces))  # Expect 4 pieces
MultiPolygon(pieces)  # This is not what we want

In [ ]:

# Compute the intersection of all buffers
from shapely.geometry import LineString
from shapely.ops import unary_union, polygonize
rings = [LineString(list(x.exterior.coords)) for x in buffered_customers]
pieces = list(polygonize(unary_union(rings)))
print(len(pieces))  # Expect 4 pieces
MultiPolygon(pieces)

+ Start with a set of geometries that represent each customer
+ Draw a buffer around each geometry
+ Compute the intersection of all buffers

In [ ]:

buffered_customers[0]

In [ ]:

buffered_customers[1]

In [ ]:

pieces[1]

In [ ]:

pieces[2]

In [ ]:

index = 2
print([
    buffered_customers[0].contains(pieces[index]),
    buffered_customers[0].overlaps(pieces[index]),
    buffered_customers[0].intersects(pieces[index]),
    buffered_customers[0].touches(pieces[index]),
])
print([
    buffered_customers[1].contains(pieces[index]),
    buffered_customers[1].overlaps(pieces[index]),
    buffered_customers[1].intersects(pieces[index]),
    buffered_customers[1].touches(pieces[index]),
])

In [ ]:

# Isolate buffer intersections
intersections = []
buffers = buffered_customers
for piece in pieces:
    membership_count = 0
    for buffer in buffers:
        if buffer.contains(piece):
            membership_count += 1
    if membership_count > 1:
        piece.membership_count = membership_count
        intersections.append(piece)
MultiPolygon(intersections)

20170919-2100 - 20170919-2200: 60 minutes

I was trying to isolate the buffer intersections last time. Maybe we can get a random point inside each intersection (such as the centroid) and use that as the test point.

In [ ]:

piece.centroid

In [ ]:

# Compute the intersection of all buffers
from shapely.geometry import LineString
from shapely.ops import unary_union, polygonize
rings = [LineString(list(x.exterior.coords)) for x in buffered_customers]
pieces = list(polygonize(unary_union(rings)))

intersections = []
buffers = buffered_customers
for piece in pieces:
    membership_count = 0
    for buffer in buffers:
        if buffer.contains(piece.centroid):
            membership_count += 1
    if membership_count > 1:
        piece.membership_count = membership_count
        intersections.append(piece)
MultiPolygon(intersections)

We were able to isolate the areas of intersection.

+ Split the buffer intersection into separate geometries which we will call buffer intersections.
+ For each buffer intersection, loop through the buffers and count the number of times the buffer intersection is contained in a buffer.

In [ ]:

for intersection in intersections:
    print(intersection.membership_count)

In [ ]:

sorted_intersections = sorted(intersections, key=lambda x: -x.membership_count)
sorted_intersections

Now we will assign service drop poles for each intersection area.

+ Sort the buffer intersections in order of most number of overlapping buffers.

The number of service drop poles is membership_count divided by maximum_customer_count_per_drop_pole, rounded up to nearest integer.

+ For each buffer intersection, assign the minimum number of service drop poles that will cover the customers for that buffer intersection. We'll call these social service drop poles.

In [ ]:

maximum_customer_count_per_drop_pole = 4
membership_count = 10
from math import ceil
ceil(membership_count / maximum_customer_count_per_drop_pole)

In [ ]:

# Compute the intersection of all buffers
from shapely.geometry import LineString
from shapely.ops import unary_union, polygonize
rings = [LineString(list(x.exterior.coords)) for x in buffered_customers]
pieces = list(polygonize(unary_union(rings)))

import copy
intersections = []
buffers = buffered_customers
lonely_buffers = copy.copy(buffers)
for piece in pieces:
    piece.buffers = []
    for buffer in buffers:
        if buffer.contains(piece.centroid):
            piece.buffers.append(buffer)
    if len(piece.buffers) > 1:
        for buffer in piece.buffers:
            lonely_buffers.remove(buffer)
        intersections.append(piece)
MultiPolygon(lonely_buffers)

In [ ]:

drop_pole_count = 0
drop_pole_count += len(lonely_buffers)
drop_pole_count

The next step is to process the isolated customers.

+ Remove each satisfied buffer from the list of unsatisfied buffers.
+ Assign a service drop pole to each unsatified buffer, which we will call lonely service drop poles.
+ Count the number of service drop poles.
+ Output the buffer intersections and unsatisfied buffers, which is our solution.

In [ ]:

# Start with a set of geometries that represent each customer
from shapely.geometry import MultiPolygon, Polygon
customers = MultiPolygon([
    Polygon([(0, 0), (1, 1), (1, 0)]),
    Polygon([(0, 0.5), (0, 1.5), (1, 1.5)]),
    Polygon([(2.5, 0.5), (2.5, 1), (3, 1), (3, 0.5)]),   
])
customers

In [ ]:

maximum_drop_length = 0.5

# Compute buffers
buffers = [x.buffer(maximum_drop_length) for x in customers]

# Compute the intersection of all buffers
from shapely.geometry import LineString
from shapely.ops import unary_union, polygonize
rings = [LineString(list(x.exterior.coords)) for x in buffers]
pieces = list(polygonize(unary_union(rings)))

import copy
intersections = []
lonely_buffers = copy.copy(buffers)
for piece in pieces:
    piece.buffers = []
    for buffer in buffers:
        if buffer.contains(piece.centroid):
            piece.buffers.append(buffer)
    if len(piece.buffers) > 1:
        for buffer in piece.buffers:
            lonely_buffers.remove(buffer)
        intersections.append(piece)

from math import ceil
maximum_customer_count_per_drop_pole = 4
drop_pole_count = 0
for intersection in intersections:
    drop_pole_count += ceil(len(intersection.buffers) / maximum_customer_count_per_drop_pole)
drop_pole_count += len(lonely_buffers)

print('drop_pole_count = %s' % drop_pole_count)

In [ ]:

MultiPolygon(intersections)

In [ ]:

MultiPolygon(lonely_buffers)

+ Combine the steps into a single algorithm

20170919-2200 - 20170919-2215: 15 minutes

In [ ]:

maximum_customer_count_per_drop_pole = 4
maximum_drop_length = 0.5
drop_pole_count = 0

from copy import copy
from math import ceil
from shapely.geometry import LineString
from shapely.ops import unary_union, polygonize

def get_disjoint_polygons(overlapping_polygons):
    rings = [LineString(list(
        x.exterior.coords)) for x in overlapping_polygons]
    return list(polygonize(unary_union(rings)))

# Compute buffers
buffers = [x.buffer(maximum_drop_length) for x in customers]
lonely_buffers = copy(buffers)

# Identify social service drop poles
intersections = []
for polygon in get_disjoint_polygons(buffers):
    polygon.buffers = []
    for buffer in buffers:
        if buffer.contains(polygon.centroid):
            polygon.buffers.append(buffer)
    if len(polygon.buffers) > 1:
        for buffer in polygon.buffers:
            lonely_buffers.remove(buffer)
        intersections.append(polygon)

# Count service drop poles        
for intersection in intersections:
    drop_pole_count += ceil(len(
        intersection.buffers) / maximum_customer_count_per_drop_pole)
drop_pole_count += len(lonely_buffers)
drop_pole_count

Okay, now we have one version.

+ Clean code
+ Implement at least one version

I think we can call service drop poles as drop poles for short.

20171116-0100 - 20171116-0130: 30 minutes

Deril noted that the service drop pole algorithm is not solving the triangle problem correctly. It seems like I forgot the step of sorting the intersections by connection_count.

In [13]:

from copy import copy
from shapely.geometry import LineString
from shapely.ops import polygonize, unary_union

def get_source_polygons(target_geometries, maximum_distance):
    # Convert target_geometries into target_polygons using maximum_distance
    target_polygons = [x.buffer(maximum_distance) for x in target_geometries]
    # Identify overlapping areas
    sliced_polygons = get_disjoint_polygons(target_polygons)
    for polygon in sliced_polygons:
        candidate_polygons = []
        for target_polygon in target_polygons:
            if target_polygon.contains(polygon.centroid):
                candidate_polygons.append(target_polygon)
        polygon.candidate_polygons = candidate_polygons
    # Sort overlapping areas by overlap count
    sorted_polygons = sorted(
        sliced_polygons, key=lambda x: -len(x.candidate_polygons))
    # Assign target_polygons to each sorted_polygon
    social_polygons, lonely_polygons = [], []
    for polygon in sorted_polygons:
        connected_polygons = []
        for candidate_polygon in polygon.candidate_polygons:
            try:
                target_polygons.remove(candidate_polygon)
            except ValueError:
                continue
            connected_polygons.append(candidate_polygon)
        connection_count = len(connected_polygons)
        if connection_count > 1:
            social_polygons.append(polygon)
        elif connection_count == 1:
            lonely_polygons.append(polygon)
        polygon.connected_polygons = connected_polygons
    return social_polygons + lonely_polygons

def get_disjoint_polygons(overlapping_polygons):
    'Split overlapping polygons into disjoint polygons'
    rings = [LineString(list(
        x.exterior.coords)) for x in overlapping_polygons]
    return list(polygonize(unary_union(rings)))

In [19]:

from shapely.geometry import GeometryCollection, Point

points = [
    Point(-1, 0),
    Point(+1, 0),
    Point(0, +2),
]
source_polygons = get_source_polygons(points, 1.5)
GeometryCollection(points + source_polygons)

Out[19]:

<shapely.geometry.collection.GeometryCollection at 0x7fa0d7373320>

Pay Notebook Creator: Roy Hyunjin Han	250
Set Container: Numerical CPU with TINY Memory for 10 Minutes	0
Total	0

Prototype Solar Mini Grid Layout (1st Implementation)

Vision¶

Mission¶

Owner¶

Context¶

Timeframe¶

Objectives¶

Log¶