[python] 지도에 대한민국 행정 단위 경계 그리기

- 개요

저번 포스트에서 ASOS 관측소 위치를 표시했습니다.

한반도 위의 ASOS 관측소 위치를 보면 대충 어느 지역에 있다는 감은 오겠지만 정확히 어디에 있는지 알 수 없습니다.

이 때 행정구역 경계를 그려주면 ASOS 관측소 위치를 정확히 알 수 있습니다.

또한, 서울 대기질 관측소 위치를 그린다고 치면 서울의 25개구 경계 자료를 같이 그려줘야 대기질 관측소 위치를 알기 쉽습니다.

이번 포스트에서는 행정구역 경계를 그려보겠습니다.

 

- 대한민국 행정 단위 경계 자료

행정구역 경계 자료는 shp 파일의 형태로 이 자료를 그림에 얹어주면 됩니다.

대한민국 행정구역 자료는 아래의 블로그에서 다운로드 받아주세요.

http://www.gisdeveloper.co.kr/?p=2332

대한민국 최신 행정구역(SHP) 다운로드 – GIS Developer

 

다운로드를 받아서 압축을 풀면 dbf, shp, shx 파일이 있습니다.

세 파일 모두 경계를 그리는데 필요하니 같은 경로에 둡시다.

저는 시군구 파일을 다운로드 받았습니다.

 

- 행정구역 그리기

필요 라이브러리: matplotlib, cartopy, geopandas

1. 한반도 전체 시군구 그리기

import cartopy.crs as ccrs
import geopandas as gpd
import cartopy.crs as ccrs
import matplotlib.pyplot as plt

fname = 'sig.shp'
# cp949로 인코딩되어있습니다.
gdf = gpd.read_file(fname, encoding='cp949')

"""
set_crs로 현재 파일의 좌표계를 설정합니다.
블로그 설명을 보시면 우리가 받은 파일의 좌표계는 EPSG:5179입니다.
to_crs를 이용해 EPSG:4326 좌표계로 변환합니다.
EPSG:4326 좌표계는 우리가 흔히 아는 위도, 경도 좌표계입니다.
"""
gdf.geometry = gdf.geometry.set_crs("EPSG:5179")
gdf.geometry = gdf.geometry.to_crs("EPSG:4326")

plt.figure(figsize=(10, 8))
# EPSG:4326 좌표계는 PlateCarree 좌표계
ax = plt.axes(projection=ccrs.PlateCarree()) 
gdf.plot(ax=ax, lw=1, color='white', edgecolor='black')
plt.show()

2. 서울 시군구 그리기

import cartopy.crs as ccrs
import geopandas as gpd
import cartopy.crs as ccrs
import matplotlib.pyplot as plt

fname = 'sig.shp'
gdf = gpd.read_file(fname, encoding='cp949')
gdf.geometry = gdf.geometry.set_crs("EPSG:5179")
gdf.geometry = gdf.geometry.to_crs("EPSG:4326")

"""
gdf 파일에는 SIG_CD: 지역 코드 정보가 있습니다.
지역 코드의 앞부분 숫자 두개는 도시의 지역 코드입니다.
서울의 지역 코드는 11이므로 SIG_CD의 앞에 두 자리가 11인 행을 선택해야 합니다
"""
gdf_seoul = gdf[gdf['SIG_CD'].str.startswith('11')]

plt.figure(figsize=(10, 8))
ax = plt.axes(projection=ccrs.PlateCarree())
gdf_seoul.plot(ax=ax, lw=1, color='white', edgecolor='black')
plt.show()

3. geopandas 쓰지 않고 행정구역 경계 그리기

앞선 코드에서는 모두 geopandas를 이용해서 그림을 그렸습니다.

cartopy의 함수를 이용하면 geopandas 없이도 행정구역 경계를 그릴 수 있지만 코드가 꽤 느립니다.

geopandas를 썼을 땐 5초만에 그려지던게 geopandas를 안 쓰니까 1분이 넘게 걸립니다.

시군구보다 더 큰 규모의 행정구역은 쓸만할지도 모르겠네요.

import matplotlib.pyplot as plt
import cartopy.io.shapereader as shpreader
from cartopy.feature import ShapelyFeature
import cartopy.crs as ccrs

fname = 'sig.shp'

# crrs.epsg(5179) = EPSG:5179 좌표계로 읽어야 함
shape_feature = ShapelyFeature(shpreader.Reader(fname).geometries(),
                               ccrs.epsg(5179), facecolor='none')

plt.figure(figsize=(10, 8))
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_extent([124, 132, 33, 38], crs=ccrs.PlateCarree())
ax.add_feature(shape_feature)
plt.show()