Title: Lesson 7: Spatial Data
Course: Basic scripting and programming > Section: Course Content

Using Python and Shapely for the manipulation of spatial objects.

Objectives:

- Spatial processing using Python and Shapely
- Coordinate transformations using PyProj

                (note: PyProj is not needed for the exercises below)

Entrance requirements:

- Lessons 1-6

Required materials:

- Computer with Python
- Special document on “Spatial Data Handling with Python”
- The data files deaths.txt and pumps.txt
- Documents can be downloaded to Resilio Sync using the following link: https://link.resilio.com/#f=2_Basic_Scripting_Programming&sz=73E5&t=1&s=X6D7CSK7YN4GKZMQ4R5GKMRWR23UZHMH&i=COWFEDPMR7N3VHWF4E5CKNP4Z3PF4OCNJ&v=2.6&a=2

Learning outcomes:

At the end of the lesson the student:

- has installed the Python modules Shapely and PyProj
- knows how spatial data is represented
- knows about spatial objects such as Points, LineStrings, Polygons, Multipoints
- knows about spatial reference systems and projections
- knows about coordinate transformations
- is able to apply all of the above

Learning activities:

- Read chapter 1 and 2 of the  document on “Spatial Data Handling with Python.”
Chapter 3 is optional.
- Duplicate the Python code of the John Snow example and study it!
- Does the problematic pump change if you change the size of the buffer?
- Determine the death count for every pump.

Exercises

- Install Shapely and PyProj on your system. See:

                https://pypi.python.org/pypi/Shapely

                https://shapely.readthedocs.io/

•         Do the following two exercises

Exercise 1 – Shapely objects

This exercise presumes that you have at least skimmed through the Shapely User Manual. Keep it close at hand, because you will need it during the exercise!

In previous exercises we've created our own class for Points. The development of a full set of geometric objects and their relationships is a tedious one. Fortunately, people have done this before us so we don't have to reinvent the wheel. One of the packages for working with geometric objects in Python is called Shapely. First we'll use shapely for Creating some simple spatial objects. After that, we'll use shapely for some more complex calculations.

Suppose we want to know whether the rectangle with corners at (0, 0) and (2, 2) contains the points at (1, 1) and (3, 1). The first problem we encounter is that shapely does not have a representation for rectangles. However, rectangles are polygons too, and we can represent rectangles by polygons. In our case the rectangle can be represented by a polygon defined by the four corner points of the rectangle: (0, 0), (2, 0), (2, 2) and (0, 2). Note that the points must be entered in the correct order, or we might end up with a twisted polygon. The points of a polygon must be entered as a sequence of tuples, where the tuples represents the points of the polygon. Polygons are supposed to be closed so there is no need to enter the starting point as the end point as well. Point objects can be created very much in the same way we've seen before.

Here's the code to check whether the rectangle contains the two points. Run it! And feel free to experiment. Note: this code will only run if you've properly installed shapely.

from shapely.geometry import Point, Polygon

r = Polygon(((0, 0), (2, 0), (2, 2), (0, 2)))

print(r.contains(Point(1, 1)))

print(r.contains(Point(3, 1)))

Task 1

Two polygons are defined. Polygon one has the following points: (0, 0), (5, -3), (5, 0) and (0, 3). Polygon two has the following points: (6, 5), (1, 2) and (6, 2).

Determine, using shapely, whether the two polygons overlap.

from shapely.geometry import Polygon

...

Task 2

Using shapely, determine the distance between the line running between the points (8, 5) and (0, 7) and the point (6, 4).

from shapely.geometry import Point, LineString

...

Task 3

Two lines are defined. One line runs from (8, 0) to (0, 7). The second line runs from (2, 3) to (8, 8).

Using shapely, determine the intersection between the two lines.

Print the coordinates of the intersection.

from shapely.geometry import LineString

...

Task 4

Determine the intersection of the two polygons used in the first task of this exercise.

- What is the type of the object returned by the intersection method?
- How many points does it contain?
- Print the coordinates of the intersection.

from shapely.geometry import Polygon

...

 Exercise 2 - advanced geometry

The code below defines a set of 100 points and 1 triangle. Use these for the following tasks.

Tasks

- Convert the points and the triangle to shapely objects.
- Determine how many points fall within the triangle.
- Print the coordinates of these points.
- Determine how many of all the points are within a distance of 0.1 of the triangle.
- Print the coordinates of these points.
- Calculate the area of that part of the triangle that is not inside the bounding box of the set of the 100 points.

Use the IDLE editor for developing your code.

from shapely.geometry import MultiPoint, Polygon

pointlist = [

(0.35, 0.28), (0.94, 0.70), (0.98, 0.59), (0.96, 0.07), (0.05, 0.57),

(0.48, 0.59), (0.33, 0.91), (0.20, 0.85), (0.32, 0.87), (0.13, 0.73),

(0.91, 0.90), (0.93, 0.48), (0.29, 0.33), (0.44, 0.06), (0.35, 0.19),

(0.73, 0.05), (0.08, 0.20), (0.58, 0.32), (0.61, 0.45), (0.65, 0.30),

(0.97, 0.34), (0.39, 0.15), (0.30, 0.83), (0.56, 0.55), (0.14, 0.37),

(0.09, 0.10), (0.97, 0.22), (0.42, 0.82), (0.06, 0.94), (0.82, 0.39),

(0.46, 0.09), (0.01, 0.81), (0.33, 0.23), (0.04, 0.02), (0.60, 0.16),

(0.68, 0.21), (0.61, 0.40), (0.53, 0.92), (0.96, 0.50), (0.14, 0.72),

(0.96, 0.60), (0.53, 0.29), (0.50, 0.16), (0.92, 0.06), (0.85, 0.81),

(0.24, 0.67), (0.05, 0.67), (0.57, 0.79), (0.36, 0.63), (0.63, 0.47),

(0.62, 0.65), (0.96, 0.47), (0.41, 0.48), (0.43, 0.57), (0.04, 0.18),

(0.90, 0.90), (0.71, 0.12), (0.94, 0.54), (0.30, 0.64), (0.48, 0.71),

(0.10, 0.33), (0.20, 0.11), (0.61, 0.19), (0.65, 0.34), (0.36, 0.02),

(0.74, 0.62), (0.13, 0.89), (0.09, 0.50), (0.79, 0.50), (0.17, 0.96),

(0.98, 0.77), (0.26, 0.68), (0.56, 0.39), (0.68, 0.71), (0.65, 0.69),

(0.26, 0.60), (0.89, 0.03), (0.96, 0.94), (0.74, 0.57), (0.54, 0.48),

(0.38, 0.65), (0.53, 0.85), (0.82, 0.68), (0.16, 0.12), (0.16, 0.15),

(0.31, 0.72), (0.13, 0.36), (0.12, 0.88), (0.57, 0.36), (0.13, 0.07),

(0.88, 0.34), (0.53, 0.67), (0.86, 0.16), (0.16, 0.09), (0.24, 0.06),

(0.00, 0.94), (0.53, 0.79), (0.34, 0.68), (0.11, 0.11), (0.49, 0.73)]

triangle = [(0.52, 0.99), (0.33, 0.76), (0.32, 0.22)]

# add your code here!

...

Check the figure below for an overview of the various geometric objects defined in the code above.

SUBMIT all your answers in a PDF document!

Planned time frame for this lesson:

6 hours