Geometric Ojbect

  • 핵심 Objects는 Points, Lines, Polygons
  • Shapely로 할 수 있는 일
  • Line과 Polygon은 Point의 Collection으로 만들 수 있음
  • Input geometries의 area, length, bounds 등을 계산
  • Union, Difference, Distance 등의 geometric operation 수행 가능
  • Intersects, Touches, Crosses, Within 등의 spatial queries도 가능

Geometric Objects

  • Point : (x, y) 또는 (x, y, z)
  • LineString
  • Polygon : 최소 3개의 coordinate tuples로 구성된 폴리곤
  • MultiPoint : coordinate tuples의 리스트로 구성된 포인트들
  • MultiLineString : sequences 같이 리스트로 된 것들
  • MultiPolygon : 폴리곤의 리스트로 구성

Point


In [3]:
from shapely.geometry import Point, LineString, Polygon

# Create Point geometric object(s) with coordinates
point1 = Point(2.2, 4.2)
point2 = Point(7.2, -25.1)
point3 = Point(9.26, -2.456)
point3D = Point(9.26, -2.456, 0.57)

# What is the type of the point?
point_type = type(point1)

In [5]:
print(point1)
print(point3D)
print(type(point1))
print(type(point3D))


POINT (2.2 4.2)
POINT Z (9.26 -2.456 0.57)
<class 'shapely.geometry.point.Point'>
<class 'shapely.geometry.point.Point'>
  • shapely의 Points는 GEOS C++ 라이브러리 기반

Point 속성과 함수


In [6]:
point_coords = point1.coords
print(point_coords)
print(type(point_coords))


<shapely.coords.CoordinateSequence object at 0x116679518>
<class 'shapely.coords.CoordinateSequence'>

In [7]:
# Get x and y coordinates
xy = point_coords.xy

# Get only x coordinates of Point1
x = point1.x

# Whatabout y coordinate?
y = point1.y

# Print out
print("xy variable:\n", xy, "\n")
print("x variable:\n", x, "\n")
print("y variable:\n", y)


xy variable:
 (array('d', [2.2]), array('d', [4.2])) 

x variable:
 2.2 

y variable:
 4.2
  • .xy는 x, y를 tuple로 리턴하고 .x는 x 리턴
  • .distance로 거리를 구할 수 있음
  • 이 거리는 좌표계에 기반한 거리임

In [8]:
point_dist = point1.distance(point2)

In [9]:
print("Distance between the points is {0:.2f} decimal degrees".format(point_dist))


Distance between the points is 29.72 decimal degrees

LineString


In [11]:
# Create a LineString from our Point objects
line = LineString([point1, point2, point3])

# It is also possible to use coordinate tuples having the same outcome
line2 = LineString([(2.2, 4.2), (7.2, -25.1), (9.26, -2.456)])

# Print the results
print("line variable: \n", line, "\n")
print("line2 variable: \n", line2, "\n")
print("type of the line: \n", type(line))


line variable: 
 LINESTRING (2.2 4.2, 7.2 -25.1, 9.26 -2.456) 

line2 variable: 
 LINESTRING (2.2 4.2, 7.2 -25.1, 9.26 -2.456) 

type of the line: 
 <class 'shapely.geometry.linestring.LineString'>
  • line은 multiple 좌표 쌍을 가지고 있음

LineString 속성과 함수

  • .xy로 x, y를 리턴

In [13]:
lxy = line.xy

print(lxy)


(array('d', [2.2, 7.2, 9.26]), array('d', [4.2, -25.1, -2.456]))

In [14]:
# Extract x coordinates
line_x = lxy[0]

# Extract y coordinates straight from the LineObject by referring to a array at index 1
line_y = line.xy[1]

print('line_x:\n', line_x, '\n')

print('line_y:\n', line_y)


line_x:
 array('d', [2.2, 7.2, 9.26]) 

line_y:
 array('d', [4.2, -25.1, -2.456])
  • length나 중심점을 추출할 수 있음

In [15]:
# Get the lenght of the line
l_length = line.length

# Get the centroid of the line
l_centroid = line.centroid

# What type is the centroid?
centroid_type = type(l_centroid)

# Print the outputs
print("Length of our line: {0:.2f}".format(l_length))
print("Centroid of our line: ", l_centroid)
print("Type of the centroid:", centroid_type)


Length of our line: 52.46
Centroid of our line:  POINT (6.229961354035622 -11.89241115757239)
Type of the centroid: <class 'shapely.geometry.point.Point'>

Polyon


In [16]:
# Create a Polygon from the coordinates
poly = Polygon([(2.2, 4.2), (7.2, -25.1), (9.26, -2.456)])

# We can also use our previously created Point objects (same outcome)
# --> notice that Polygon object requires x,y coordinates as input
poly2 = Polygon([[p.x, p.y] for p in [point1, point2, point3]])

# Geometry type can be accessed as a String
poly_type = poly.geom_type

# Using the Python's type function gives the type in a different format
poly_type2 = type(poly)

# Let's see how our Polygon looks like
print('poly:', poly)
print('poly2:', poly2)
print("Geometry type as text:", poly_type)
print("Geometry how Python shows it:", poly_type2)


poly: POLYGON ((2.2 4.2, 7.2 -25.1, 9.26 -2.456, 2.2 4.2))
poly2: POLYGON ((2.2 4.2, 7.2 -25.1, 9.26 -2.456, 2.2 4.2))
Geometry type as text: Polygon
Geometry how Python shows it: <class 'shapely.geometry.polygon.Polygon'>

In [24]:
help(Polygon)


Help on class Polygon in module shapely.geometry.polygon:

class Polygon(shapely.geometry.base.BaseGeometry)
 |  Polygon(shell=None, holes=None)
 |  
 |  A two-dimensional figure bounded by a linear ring
 |  
 |  A polygon has a non-zero area. It may have one or more negative-space
 |  "holes" which are also bounded by linear rings. If any rings cross each
 |  other, the feature is invalid and operations on it may fail.
 |  
 |  Attributes
 |  ----------
 |  exterior : LinearRing
 |      The ring which bounds the positive space of the polygon.
 |  interiors : sequence
 |      A sequence of rings which bound all existing holes.
 |  
 |  Method resolution order:
 |      Polygon
 |      shapely.geometry.base.BaseGeometry
 |      builtins.object
 |  
 |  Methods defined here:
 |  
 |  __eq__(self, other)
 |      Return self==value.
 |  
 |  __init__(self, shell=None, holes=None)
 |      Parameters
 |      ----------
 |      shell : sequence
 |          A sequence of (x, y [,z]) numeric coordinate pairs or triples
 |      holes : sequence
 |          A sequence of objects which satisfy the same requirements as the
 |          shell parameters above
 |      
 |      Example
 |      -------
 |      Create a square polygon with no holes
 |      
 |        >>> coords = ((0., 0.), (0., 1.), (1., 1.), (1., 0.), (0., 0.))
 |        >>> polygon = Polygon(coords)
 |        >>> polygon.area
 |        1.0
 |  
 |  __ne__(self, other)
 |      Return self!=value.
 |  
 |  svg(self, scale_factor=1.0, fill_color=None)
 |      Returns SVG path element for the Polygon geometry.
 |      
 |      Parameters
 |      ==========
 |      scale_factor : float
 |          Multiplication factor for the SVG stroke-width.  Default is 1.
 |      fill_color : str, optional
 |          Hex string for fill color. Default is to use "#66cc99" if
 |          geometry is valid, and "#ff3333" if invalid.
 |  
 |  ----------------------------------------------------------------------
 |  Class methods defined here:
 |  
 |  from_bounds(xmin, ymin, xmax, ymax) from builtins.type
 |      Construct a `Polygon()` from spatial bounds.
 |  
 |  ----------------------------------------------------------------------
 |  Data descriptors defined here:
 |  
 |  __array_interface__
 |      Provide the Numpy array protocol.
 |  
 |  __geo_interface__
 |      Dictionary representation of the geometry
 |  
 |  coords
 |      Access to geometry's coordinates (CoordinateSequence)
 |  
 |  ctypes
 |      Return ctypes buffer
 |  
 |  exterior
 |  
 |  interiors
 |  
 |  ----------------------------------------------------------------------
 |  Data and other attributes defined here:
 |  
 |  __hash__ = None
 |  
 |  ----------------------------------------------------------------------
 |  Methods inherited from shapely.geometry.base.BaseGeometry:
 |  
 |  __and__(self, other)
 |  
 |  __del__(self)
 |  
 |  __or__(self, other)
 |  
 |  __reduce__(self)
 |      Helper for pickle.
 |  
 |  __setstate__(self, state)
 |  
 |  __str__(self)
 |      Return str(self).
 |  
 |  __sub__(self, other)
 |  
 |  __xor__(self, other)
 |  
 |  almost_equals(self, other, decimal=6)
 |      Returns True if geometries are equal at all coordinates to a
 |      specified decimal place
 |      
 |      Refers to approximate coordinate equality, which requires coordinates be
 |      approximately equal and in the same order for all components of a geometry.
 |  
 |  buffer(self, distance, resolution=16, quadsegs=None, cap_style=1, join_style=1, mitre_limit=5.0)
 |      Returns a geometry with an envelope at a distance from the object's
 |      envelope
 |      
 |      A negative distance has a "shrink" effect. A zero distance may be used
 |      to "tidy" a polygon. The resolution of the buffer around each vertex of
 |      the object increases by increasing the resolution keyword parameter
 |      or second positional parameter. Note: the use of a `quadsegs` parameter
 |      is deprecated and will be gone from the next major release.
 |      
 |      The styles of caps are: CAP_STYLE.round (1), CAP_STYLE.flat (2), and
 |      CAP_STYLE.square (3).
 |      
 |      The styles of joins between offset segments are: JOIN_STYLE.round (1),
 |      JOIN_STYLE.mitre (2), and JOIN_STYLE.bevel (3).
 |      
 |      The mitre limit ratio is used for very sharp corners. The mitre ratio
 |      is the ratio of the distance from the corner to the end of the mitred
 |      offset corner. When two line segments meet at a sharp angle, a miter
 |      join will extend the original geometry. To prevent unreasonable
 |      geometry, the mitre limit allows controlling the maximum length of the
 |      join corner. Corners with a ratio which exceed the limit will be
 |      beveled.
 |      
 |      Example:
 |      
 |        >>> from shapely.wkt import loads
 |        >>> g = loads('POINT (0.0 0.0)')
 |        >>> g.buffer(1.0).area        # 16-gon approx of a unit radius circle
 |        3.1365484905459389
 |        >>> g.buffer(1.0, 128).area   # 128-gon approximation
 |        3.1415138011443009
 |        >>> g.buffer(1.0, 3).area     # triangle approximation
 |        3.0
 |        >>> list(g.buffer(1.0, cap_style='square').exterior.coords)
 |        [(1.0, 1.0), (1.0, -1.0), (-1.0, -1.0), (-1.0, 1.0), (1.0, 1.0)]
 |        >>> g.buffer(1.0, cap_style='square').area
 |        4.0
 |  
 |  contains(self, other)
 |      Returns True if the geometry contains the other, else False
 |  
 |  covers(self, other)
 |      Returns True if the geometry covers the other, else False
 |  
 |  crosses(self, other)
 |      Returns True if the geometries cross, else False
 |  
 |  difference(self, other)
 |      Returns the difference of the geometries
 |  
 |  disjoint(self, other)
 |      Returns True if geometries are disjoint, else False
 |  
 |  distance(self, other)
 |      Unitless distance to other geometry (float)
 |  
 |  empty(self, val=140590495520064)
 |  
 |  equals(self, other)
 |      Returns True if geometries are equal, else False
 |      
 |      Refers to point-set equality (or topological equality), and is equivalent to
 |      (self.within(other) & self.contains(other))
 |  
 |  equals_exact(self, other, tolerance)
 |      Returns True if geometries are equal to within a specified
 |      tolerance
 |      
 |      Refers to coordinate equality, which requires coordinates to be equal 
 |      and in the same order for all components of a geometry
 |  
 |  geometryType(self)
 |  
 |  hausdorff_distance(self, other)
 |      Unitless hausdorff distance to other geometry (float)
 |  
 |  interpolate(self, distance, normalized=False)
 |      Return a point at the specified distance along a linear geometry
 |      
 |      If the normalized arg is True, the distance will be interpreted as a
 |      fraction of the geometry's length.
 |  
 |  intersection(self, other)
 |      Returns the intersection of the geometries
 |  
 |  intersects(self, other)
 |      Returns True if geometries intersect, else False
 |  
 |  overlaps(self, other)
 |      Returns True if geometries overlap, else False
 |  
 |  project(self, other, normalized=False)
 |      Returns the distance along this geometry to a point nearest the
 |      specified point
 |      
 |      If the normalized arg is True, return the distance normalized to the
 |      length of the linear geometry.
 |  
 |  relate(self, other)
 |      Returns the DE-9IM intersection matrix for the two geometries
 |      (string)
 |  
 |  relate_pattern(self, other, pattern)
 |      Returns True if the DE-9IM string code for the relationship between
 |      the geometries satisfies the pattern, else False
 |  
 |  representative_point(self)
 |      Returns a point guaranteed to be within the object, cheaply.
 |  
 |  simplify(self, tolerance, preserve_topology=True)
 |      Returns a simplified geometry produced by the Douglas-Peucker
 |      algorithm
 |      
 |      Coordinates of the simplified geometry will be no more than the
 |      tolerance distance from the original. Unless the topology preserving
 |      option is used, the algorithm may produce self-intersecting or
 |      otherwise invalid geometries.
 |  
 |  symmetric_difference(self, other)
 |      Returns the symmetric difference of the geometries
 |      (Shapely geometry)
 |  
 |  to_wkb(self)
 |  
 |  to_wkt(self)
 |  
 |  touches(self, other)
 |      Returns True if geometries touch, else False
 |  
 |  union(self, other)
 |      Returns the union of the geometries (Shapely geometry)
 |  
 |  within(self, other)
 |      Returns True if geometry is within the other, else False
 |  
 |  ----------------------------------------------------------------------
 |  Data descriptors inherited from shapely.geometry.base.BaseGeometry:
 |  
 |  __dict__
 |      dictionary for instance variables (if defined)
 |  
 |  __weakref__
 |      list of weak references to the object (if defined)
 |  
 |  area
 |      Unitless area of the geometry (float)
 |  
 |  array_interface_base
 |  
 |  boundary
 |      Returns a lower dimension geometry that bounds the object
 |      
 |      The boundary of a polygon is a line, the boundary of a line is a
 |      collection of points. The boundary of a point is an empty (null)
 |      collection.
 |  
 |  bounds
 |      Returns minimum bounding region (minx, miny, maxx, maxy)
 |  
 |  centroid
 |      Returns the geometric center of the object
 |  
 |  convex_hull
 |      Imagine an elastic band stretched around the geometry: that's a
 |      convex hull, more or less
 |      
 |      The convex hull of a three member multipoint, for example, is a
 |      triangular polygon.
 |  
 |  envelope
 |      A figure that envelopes the geometry
 |  
 |  geom_type
 |      Name of the geometry's type, such as 'Point'
 |  
 |  has_z
 |      True if the geometry's coordinate sequence(s) have z values (are
 |      3-dimensional)
 |  
 |  is_closed
 |      True if the geometry is closed, else False
 |      
 |      Applicable only to 1-D geometries.
 |  
 |  is_empty
 |      True if the set of points in this geometry is empty, else False
 |  
 |  is_ring
 |      True if the geometry is a closed ring, else False
 |  
 |  is_simple
 |      True if the geometry is simple, meaning that any self-intersections
 |      are only at boundary points, else False
 |  
 |  is_valid
 |      True if the geometry is valid (definition depends on sub-class),
 |      else False
 |  
 |  length
 |      Unitless length of the geometry (float)
 |  
 |  minimum_rotated_rectangle
 |      Returns the general minimum bounding rectangle of
 |      the geometry. Can possibly be rotated. If the convex hull
 |      of the object is a degenerate (line or point) this same degenerate
 |      is returned.
 |  
 |  type
 |  
 |  wkb
 |      WKB representation of the geometry
 |  
 |  wkb_hex
 |      WKB hex representation of the geometry
 |  
 |  wkt
 |      WKT representation of the geometry
 |  
 |  xy
 |      Separate arrays of X and Y coordinate values
 |  
 |  ----------------------------------------------------------------------
 |  Data and other attributes inherited from shapely.geometry.base.BaseGeometry:
 |  
 |  __geom__ = 140590495520064
 |  
 |  __p__ = None
 |  
 |  impl = <GEOSImpl object: GEOS C API version (1, 8, 0)>

  • Polygon은 안에 구멍이 있을수도 있음

In [25]:
# First we define our exterior
world_exterior = [(-180, 90), (-180, -90), (180, -90), (180, 90)]

# Let's create a single big hole where we leave ten decimal degrees at the boundaries of the world
# Notice: there could be multiple holes, thus we need to provide a list of holes
hole = [[(-170, 80), (-170, -80), (170, -80), (170, 80)]]

# World without a hole
world = Polygon(shell=world_exterior)

# Now we can construct our Polygon with the hole inside
world_has_a_hole = Polygon(shell=world_exterior, holes=hole)

In [26]:
print('world:', world)
print('world_has_a_hole:', world_has_a_hole)
print('type:', type(world_has_a_hole))


world: POLYGON ((-180 90, -180 -90, 180 -90, 180 90, -180 90))
world_has_a_hole: POLYGON ((-180 90, -180 -90, 180 -90, 180 90, -180 90), (-170 80, -170 -80, 170 -80, 170 80, -170 80))
type: <class 'shapely.geometry.polygon.Polygon'>

Polygon 속성과 함수

  • .area, .centroid, .bounds, .exteror, .length

In [27]:
# Get the centroid of the Polygon
world_centroid = world.centroid

# Get the area of the Polygon
world_area = world.area

# Get the bounds of the Polygon (i.e. bounding box)
world_bbox = world.bounds

# Get the exterior of the Polygon
world_ext = world.exterior

# Get the length of the exterior
world_ext_length = world_ext.length

# Print the outputs
print("Poly centroid: ", world_centroid)
print("Poly Area: ", world_area)
print("Poly Bounding Box: ", world_bbox)
print("Poly Exterior: ", world_ext)
print("Poly Exterior Length: ", world_ext_length)


Poly centroid:  POINT (-0 -0)
Poly Area:  64800.0
Poly Bounding Box:  (-180.0, -90.0, 180.0, 90.0)
Poly Exterior:  LINEARRING (-180 90, -180 -90, 180 -90, 180 90, -180 90)
Poly Exterior Length:  1080.0

Geometry Collections

  • Y모양의 선 형상(도로) 또는 다중 다각형(섬)은 MultiLineString 또는 MultiPolygon을 적절하게 사용해 표현할 수 있음

In [28]:
from shapely.geometry import MultiPoint, MultiLineString, MultiPolygon, box

# Create a MultiPoint object of our points 1,2 and 3
multi_point = MultiPoint([point1, point2, point3])

# It is also possible to pass coordinate tuples inside
multi_point2 = MultiPoint([(2.2, 4.2), (7.2, -25.1), (9.26, -2.456)])

# We can also create a MultiLineString with two lines
line1 = LineString([point1, point2])
line2 = LineString([point2, point3])
multi_line = MultiLineString([line1, line2])

# MultiPolygon can be done in a similar manner
# Let's divide our world into western and eastern hemispheres with a hole on the western hemisphere
# --------------------------------------------------------------------------------------------------

# Let's create the exterior of the western part of the world
west_exterior = [(-180, 90), (-180, -90), (0, -90), (0, 90)]

# Let's create a hole --> remember there can be multiple holes, thus we need to have a list of hole(s).
# Here we have just one.
west_hole = [[(-170, 80), (-170, -80), (-10, -80), (-10, 80)]]

# Create the Polygon
west_poly = Polygon(shell=west_exterior, holes=west_hole)

# Let's create the Polygon of our Eastern hemisphere polygon using bounding box
# For bounding box we need to specify the lower-left corner coordinates and upper-right coordinates
min_x, min_y = 0, -90
max_x, max_y = 180, 90

# Create the polygon using box() function
east_poly_box = box(minx=min_x, miny=min_y, maxx=max_x, maxy=max_y)

# Let's create our MultiPolygon. We can pass multiple Polygon -objects into our MultiPolygon as a list
multi_poly = MultiPolygon([west_poly, east_poly_box])

# Print outputs
print("MultiPoint:", multi_point)
print("MultiLine: ", multi_line)
print("Bounding box: ", east_poly_box)
print("MultiPoly: ", multi_poly)


MultiPoint: MULTIPOINT (2.2 4.2, 7.2 -25.1, 9.26 -2.456)
MultiLine:  MULTILINESTRING ((2.2 4.2, 7.2 -25.1), (7.2 -25.1, 9.26 -2.456))
Bounding box:  POLYGON ((180 -90, 180 90, 0 90, 0 -90, 180 -90))
MultiPoly:  MULTIPOLYGON (((-180 90, -180 -90, 0 -90, 0 90, -180 90), (-170 80, -170 -80, -10 -80, -10 80, -170 80)), ((180 -90, 180 90, 0 90, 0 -90, 180 -90)))

Geometry collection의 속성과 함수

  • .convex_hull : 다양한 점들을 감싸는 볼록 껍질을 구함
    • (알고리즘 이해하기 위해)
    • 외곽이 되는 점 : 상하좌우 최 외각에 있는 점
    • 선택된 점에서 모든 점 사이의 각도를 구함
    • 각도를 기준으로 정렬
    • 0, 1을 stack에 넣고 그 다음 점이 반 시계 방향이면 넣고 아니면 pop
    • 시간 복잡도 : O(nlogn)
    • 참고 슬라이드
  • .is_valid : polygon이나 선들이 서로 교차하는지

In [29]:
# Convex Hull of our MultiPoint --> https://en.wikipedia.org/wiki/Convex_hull
convex = multi_point.convex_hull

# How many lines do we have inside our MultiLineString?
lines_count = len(multi_line)

# Let's calculate the area of our MultiPolygon
multi_poly_area = multi_poly.area

# We can also access different items inside our geometry collections. We can e.g. access a single polygon from
# our MultiPolygon -object by referring to the index

# Let's calculate the area of our Western hemisphere (with a hole) which is at index 0
west_area = multi_poly[0].area

# We can check if we have a "valid" MultiPolygon. MultiPolygon is thought as valid if the individual polygons
# does notintersect with each other. Here, because the polygons have a common 0-meridian, we should NOT have
# a valid polygon. This can be really useful information when trying to find topological errors from your data
valid = multi_poly.is_valid

# Print outputs
print("Convex hull of the points: ", convex)
print("Number of lines in MultiLineString:", lines_count)
print("Area of our MultiPolygon:", multi_poly_area)
print("Area of our Western Hemisphere polygon:", west_area)
print("Is polygon valid?: ", valid)


Convex hull of the points:  POLYGON ((7.2 -25.1, 2.2 4.2, 9.26 -2.456, 7.2 -25.1))
Number of lines in MultiLineString: 2
Area of our MultiPolygon: 39200.0
Area of our Western Hemisphere polygon: 6800.0
Is polygon valid?:  False

In [ ]: