ST_Buffer와 공간 분석 — 500m 버퍼가 397m로 줄어든다
건물 12,912개가 사라졌다
"서울 소방서에서 500m 안에 있는 건물이 몇 개인가." 문장이 짧으니 쿼리도 짧습니다. 소방서마다 반경 500m 원을 그리고, 그 원에 걸치는 건물을 세면 되죠.
ST_Buffer(way, 500)한 번이면 됩니다. 답은 23,695개가 나왔습니다.그런데 같은 질문의 답이 36,607개이기도 합니다. 같은 데이터베이스, 같은 소방서, 같은 숫자
500. 건물이 어디로 옮겨간 것도 아닌데 12,912개가 첫 번째 셈에서 통째로 빠져 있었어요. 에러는 없었습니다. 경고도 없었고요.
오늘은 이 과제 하나를 처음부터 끝까지 풀어 보겠습니다. 두 숫자가 왜 다른지 찾아내고, 이어서 "서울 땅의 몇 %가 소방서에서 500m 안인가" 까지 답을 냅니다.
그 과정에서 필요할 때마다 함수가 하나씩 불려 나옵니다 — 반경을 그리는 ST_Buffer, 겹친 걸 합치는 ST_Union, 잘라내는 ST_Intersection·ST_Difference, 도형에서 수를 뽑는 ST_Centroid·ST_Area. 지난 편이 있는 도형에 관계를 묻는 일이었다면, 오늘은 없는 도형을 만드는 쪽입니다.
측정 환경: PostGIS 3.5.2, GEOS 3.9.0, PROJ 7.2.1, PostgreSQL 17, 측정일 2026-09-09.
오늘 쓸 데이터
앞선 글들에서 적재해 둔 테이블을 그대로 씁니다.
| 테이블 | 내용 | SRID | 어디서 넣었나 |
|---|---|---|---|
planet_osm_point · planet_osm_polygon | 서울 OSM 점과 면 | 3857 | osm2pgsql 기본 출력 |
sgg | 통계청 2012년 시군구 경계 251개 (서울 25개 구 포함) | 5174 | shp2pgsql |
convenience | 서울 편의점 7,100개 | 4326 | osm2pgsql flex 출력 |
없어도 오늘 글은 읽힙니다 — 숫자는 전부 본문에 적어 두었어요. 직접 만들어 보실 거라면 각 링크에 명령이 그대로 있고, PostGIS가 깔린 데이터베이스 자체가 없다면 PostGIS 시작하기 편에 도커로 띄우는 방법이 있습니다.
소방서를 모으는 데 이미 함정이 하나 있다
과제의 주인공인 소방서부터 꺼내야 하는데, 여기서부터 걸립니다. OSM에서 amenity=fire_station은 점으로도 면으로도 찍혀 있거든요.
SELECT (SELECT count(*) FROM planet_osm_point WHERE amenity='fire_station') AS 점,
(SELECT count(*) FROM planet_osm_polygon WHERE amenity='fire_station') AS 면;
점 | 면
----+----
35 | 86
점 35개, 면 86개. 점 테이블만 보면 서울에 소방서가 35개뿐인 걸로 나옵니다. 누군가는 위치에 점 하나를 찍었고, 누군가는 건물 외곽선을 그린 뒤 거기에 태그를 붙인 것뿐인데요. 지난 편에서 편의점 7,100개를 모을 때 닫힌 Way의 중심점을 뽑았던 것과 같은 상황입니다.
그러니 둘을 합쳐야 합니다. 면은 ST_Centroid로 점을 만들어 붙이고요.
CREATE VIEW fire AS
SELECT osm_id, name, way FROM planet_osm_point
WHERE amenity = 'fire_station'
UNION ALL
SELECT osm_id, name, ST_Centroid(way) FROM planet_osm_polygon
WHERE amenity = 'fire_station';
SELECT count(*) AS 소방서 FROM fire;
소방서
--------
121
121개. 점만 봤다면 35개로 시작할 뻔했습니다.
💡 여기서 쓴
ST_Centroid는 나중에 다시 붙잡을 겁니다. 면을 점으로 바꾸는 게 안전한 일인지 아직 확인하지 않았거든요. 5절에서 다시 보겠습니다.
1. 반경을 그린다 — ST_Buffer
ST_Buffer는 도형을 주면 그만큼 부풀린 새 도형을 돌려줍니다. 점에 쓰면 원이 되죠. 과제 문장을 그대로 SQL로 옮기면 이렇게 됩니다.
SELECT count(DISTINCT b.osm_id) AS 건물
FROM planet_osm_polygon b
JOIN fire f ON ST_Intersects(b.way, ST_Buffer(f.way, 500))
WHERE b.building IS NOT NULL;
건물
-------
23695
23,695개. 읽기에 이상한 데가 없습니다. 그런데 여기서 물어야 할 게 하나 있어요.
500은 무엇의 500인가
ST_Buffer(f.way, 500)에서 500에는 단위가 없습니다. 미터라고 쓴 적이 없어요. 그냥 숫자 500입니다.
PostGIS는 이 숫자를 좌표계의 단위 그대로 씁니다. PostGIS 시작하기 편에서 ST_Distance가 4326 좌표에 대해 0.08을 돌려줬던 게 같은 이유였죠 — 4326의 단위는 도(degree)니까요.
그럼 planet_osm_polygon의 단위는 뭘까요. 지난 편에서 확인했듯 osm2pgsql이 우리에게 묻지 않고 붙여준 좌표계, EPSG:3857 웹 메르카토르입니다. 3857의 단위는 미터입니다.
그런데 그 미터는 줄자로 재는 미터가 아닙니다.
메르카토르 도법은 지구를 원통에 펼치면서 위도가 높아질수록 가로로 늘립니다. 지도 투영법 편에서 본 "그린란드가 아프리카만 해 보이는" 그 왜곡이에요.
그러니 늘려 그린 지도 위에서 눈금 500칸을 가도, 실제 땅에서는 500m를 못 갑니다. 지도가 부풀어 있으니 눈금 한 칸이 실제 1m보다 짧은 거리를 나타내는 거죠.
얼마나 짧은지 재봅시다. 소방서 한 곳을 골라 3857 좌표에서 정확히 500만큼 옆으로 옮긴 다음, 그 두 지점이 실제 땅 위에서 몇 미터 떨어져 있는지 묻는 겁니다. 아래 ::geography가 그 "현장에서 줄자로 재기"에 해당해요 — 평면 눈금이 아니라 지구 표면을 따라 재 줍니다.
SELECT round(ST_Y(ST_Transform(way, 4326))::numeric, 4) AS 위도,
round(ST_Distance(ST_Transform(way, 4326)::geography,
ST_Transform(ST_Translate(way, 500, 0), 4326)::geography)::numeric, 2) AS 실제_m,
round((cos(radians(ST_Y(ST_Transform(way, 4326)))) * 500)::numeric, 2) AS "cosφ×500"
FROM fire WHERE name = '세검정 소방서';
위도 | 실제_m | cosφ×500
---------+--------+----------
37.6058 | 396.61 | 396.11
396.61m. 500을 줬는데 실제로는 397m짜리 원을 그리고 있었습니다.
왜 397m일까요. 지도가 부풀어 있는 만큼 나눠지기 때문입니다.
메르카토르가 위도에 따라 지도를 옆으로 늘린다고 했죠. 그 늘어난 배율이 1 / cos φ 입니다. 세검정 소방서의 위도 37.6058°를 넣으면 이렇게 돼요.
cos(37.6058°) = 0.79223
확대 배율 = 1 / 0.79223 = 1.2623배
서울에서 3857 지도는 실제 땅보다 1.26배 크게 그려져 있습니다. 지도가 1.26배 부풀어 있으니 그 위의 눈금 500칸은 실제 거리로 환산하면 그만큼 나눈 값이 됩니다.
500 ÷ 1.2623 = 396.11 m
앞 표의 셋째 칸(cosφ×500)이 바로 이 값이고, 실측한 396.61m와 0.5m 차이로 맞습니다. 남은 0.5m는 3857이 지구를 완전한 구로 놓고 계산하는 반면 ::geography는 타원체로 재기 때문이에요.
같은 계산을 소방서 121곳 전부에 돌리면 최소 396.1m, 평균 396.9m, 최대 397.6m. 서울은 남북 폭이 좁아 위도 차가 크지 않으니 어디서 재도 397m 언저리예요. "반경 500m"라고 쓴 쿼리가 실제로는 397m 반경을 재고 있었습니다.
정말 397m 원인지 확인한다
추론이 맞는지 검산할 방법이 있습니다. 3857에서 그린 그 버퍼가 진짜로 397m짜리라면, 눈금 1칸이 진짜 1m인 좌표계에서 반경 396.61m로 검색한 결과와 개수가 같아야 합니다.
그런 좌표계가 필요한데, 3857로는 안 됩니다. 3857 좌표를 3857로 재면 당연히 500이 나오거든요 — 자기 눈금으로 자기를 재는 셈이라 의심하던 숫자가 그대로 돌아옵니다. 그래서 바깥의 자를 가져옵니다. 한국이라면 EPSG:5179예요(한국이 EPSG:5179를 쓰는 이유 편에서 다뤘습니다).
이 자가 믿을 만한지부터 확인하죠. 3857에 했던 것과 똑같은 실험을 5179에 해봅니다 — 500만큼 옆으로 옮기고 실제 거리를 재는 겁니다.
3857 에서 500 이동 → 실제 396.61 m (103m 부족)
5179 에서 500 이동 → 실제 500.19 m (0.19m 오차)
5179는 500을 주면 실제로 500m를 갑니다. 한국 지역에 맞춰 만든 좌표계라 이 정도 범위에서는 눈금이 실제 거리와 거의 일치해요. 이제 이 자로 3857을 재봅니다.
WITH fs AS (SELECT way AS g3857, ST_Transform(way, 5179) AS g FROM fire WHERE name = '세검정 소방서')
SELECT (SELECT count(*) FROM planet_osm_polygon b, fs
WHERE b.building IS NOT NULL AND ST_Intersects(b.way, ST_Buffer(fs.g3857, 500))) AS "3857_버퍼500",
(SELECT count(*) FROM planet_osm_polygon b, fs
WHERE b.building IS NOT NULL AND ST_DWithin(ST_Transform(b.way,5179), fs.g, 396.61)) AS "5179_반경396.61",
(SELECT count(*) FROM planet_osm_polygon b, fs
WHERE b.building IS NOT NULL AND ST_DWithin(ST_Transform(b.way,5179), fs.g, 500)) AS "5179_반경500";
3857_버퍼500 | 5179_반경396.61 | 5179_반경500
--------------+-----------------+--------------
10 | 10 | 34
10개와 10개. 3857에서 "500"으로 그린 버퍼는 미터 좌표계의 396.61m 반경과 정확히 같은 답을 냅니다. 우리가 원한 500m 반경은 같은 소방서에서 34개였고요.
PostGIS는 SRID를 계산에 쓰지 않는다
여기서 한 발 더 들어가 봅시다. 방금 일이 "3857이라서" 생긴 걸까요? 아닙니다. 좌표계만 바꿔 같은 버퍼를 세 번 만들어 보면 압니다.
SELECT ST_Area(ST_Buffer(ST_SetSRID(ST_MakePoint(127, 37.5), 4326), 500)) AS "4326",
ST_Area(ST_Buffer(ST_SetSRID(ST_MakePoint(127, 37.5), 5174), 500)) AS "5174",
ST_Area(ST_Buffer(ST_SetSRID(ST_MakePoint(127, 37.5), 0), 500)) AS "SRID없음";
4326 | 5174 | SRID없음
-------------------+-------------------+-------------------
780361.2880645131 | 780361.2880645131 | 780361.2880645131
소수점 끝자리까지 같습니다. 위경도로 주든, 미터 좌표계로 주든, 아예 안 주든 결과가 똑같아요. ST_Buffer는 SRID를 꼬리표로만 들고 다닐 뿐 계산에 쓰지 않습니다. 반지름 500짜리 32각형 넓이를 세 번 똑같이 구한 거고, 그 숫자를 제곱미터로 읽을지 제곱도로 읽을지는 사람 몫입니다.
경계상자를 찍어 보면 계산이 얼마나 단순한지 더 분명해집니다.
SELECT ST_AsText(ST_Envelope(ST_Buffer('SRID=4326;POINT(127 37.5)'::geometry, 500))) AS 경계상자,
ST_IsValid(ST_Buffer('SRID=4326;POINT(127 37.5)'::geometry, 500)) AS 유효한가;
경계상자 | 유효한가
--------------------------------------------------------------------+----------
POLYGON((-373 -462.5,-373 537.5,627 537.5,627 -462.5,-373 -462.5)) | t
중심에서 그냥 500을 더하고 뺐습니다 — 127 − 500 = −373, 127 + 500 = 627. 경도가 ±180도, 위도가 ±90도까지뿐이라는 건 고려 대상이 아니었어요. 그런데도 ST_IsValid는 t 입니다. 유효성은 "링이 닫혔나, 자기교차가 없나" 만 보지 "이 좌표가 지구에 실재하나" 는 안 보거든요.
그래서 3857이 제일 위험하다
PostGIS의 침묵은 어느 좌표계에서나 똑같습니다. 달라지는 건 그 침묵이 얼마나 큰 사고로 이어지느냐뿐이에요.
| 좌표계 | 눈금 1칸의 실제 크기 | 500을 주면 | 결과 |
|---|---|---|---|
| 4326 | 위도 1도 ≈ 110.6 km | 반지름 55,287 km | 지구 둘레(40,075 km)를 넘는다 |
| 3857 (서울) | ≈ 0.79 m | 반지름 397 m | 그럴듯하다 |
4326 쪽은 실제로 지구를 통째로 삼킵니다. ST_Contains로 확인해 보면 그 버퍼 하나가 경도 −180180, 위도 −9090 전체를 포함한다고 나와요. 이런 실수는 곧바로 들킵니다 — 건물이 전부 매치되고, 쿼리는 인덱스를 못 타 하염없이 돕니다.
진짜 위험한 건 그럴듯한 쪽입니다. 397m는 500m와 충분히 비슷해서 결과를 봐도 이상한 데가 없어요. 우리가 밟은 게 정확히 그 자리였고, 건물 12,720개가 그렇게 조용히 빠졌습니다.
고쳐서 다시 센다
해법은 이미 나왔습니다. 검산에 쓴 그 5179로 아예 옮겨놓고 재는 것.
다만 쿼리마다 ST_Transform을 부르면 인덱스를 못 씁니다. 인덱스는 원래 좌표에 걸려 있는데 계산은 변환된 좌표로 하니까요. 그래서 변환한 결과를 테이블로 만들어 두고 인덱스를 새로 겁니다.
CREATE TABLE fire_m AS SELECT osm_id, name, ST_Transform(way, 5179) AS geom FROM fire;
CREATE INDEX ON fire_m USING gist(geom);
CREATE TABLE bld_m AS
SELECT osm_id, building, ST_Transform(way, 5179) AS geom
FROM planet_osm_polygon WHERE building IS NOT NULL;
CREATE INDEX ON bld_m USING gist(geom);
ANALYZE fire_m; ANALYZE bld_m;
SELECT (SELECT count(*) FROM fire_m) AS 소방서, (SELECT count(*) FROM bld_m) AS 건물;
소방서 | 건물
--------+--------
121 | 240092
건물 240,092개가 미터 좌표계로 옮겨졌습니다. 이제 같은 질문을 다시 던집니다.
SELECT count(DISTINCT b.osm_id) AS 건물
FROM bld_m b JOIN fire_m f ON ST_Intersects(b.geom, ST_Buffer(f.geom, 500));
건물
-------
36415
23,695 → 36,415. 건물 12,720개가 돌아왔습니다.
다른 데를 뒤진 게 아니라, 좁게 뒤진 것이다
여기서 한 번 짚고 갈 게 있습니다. 좌표계를 바꿔서 답이 달라졌다고 하면 "두 쿼리가 애초에 서로 다른 자리를 조회한 것 아닌가" 싶어집니다. 그럴듯한데, 아닙니다.
ST_Transform은 같은 땅 위 같은 자리를 다른 숫자로 다시 적을 뿐입니다. 소방서는 두 쿼리에서 물리적으로 같은 지점에 있고, 건물도 마찬가지예요. 조회 위치는 한 번도 안 움직였습니다. 움직인 건 그 지점을 감싼 원의 반지름뿐이에요 — 397m에서 500m로.
둘은 구별할 수 있습니다. 조회 자리가 옮겨갔다면 새로 들어온 건물이 있는 만큼 빠져나간 건물도 있어야 하죠. 반면 같은 중심에서 원만 커진 거라면 나가는 건 하나도 없고, 새로 들어온 것들은 전부 두 반지름 사이 도넛에 모여 있어야 합니다.
WITH a AS (SELECT DISTINCT b.osm_id FROM planet_osm_polygon b JOIN fire f
ON ST_Intersects(b.way, ST_Buffer(f.way,500)) WHERE b.building IS NOT NULL),
b5179 AS (SELECT DISTINCT b.osm_id FROM bld_m b JOIN fire_m f
ON ST_Intersects(b.geom, ST_Buffer(f.geom,500)))
SELECT (SELECT count(*) FROM (SELECT osm_id FROM b5179 EXCEPT SELECT osm_id FROM a) x) AS 늘어난건물,
(SELECT count(*) FROM (SELECT osm_id FROM a EXCEPT SELECT osm_id FROM b5179) x) AS 빠진건물;
늘어난건물 | 빠진건물
------------+----------
12720 | 0
빠진 건물이 0개. 교환이 아니라 한 방향으로만 늘었습니다. 새로 잡힌 것들이 어디 있었는지도 재봅니다.
WITH a AS (SELECT DISTINCT b.osm_id FROM planet_osm_polygon b JOIN fire f
ON ST_Intersects(b.way, ST_Buffer(f.way,500)) WHERE b.building IS NOT NULL)
SELECT round(min(d)::numeric,1) AS 최단_m, round(max(d)::numeric,1) AS 최장_m
FROM (SELECT b.osm_id, min(ST_Distance(b.geom, f.geom)) AS d
FROM bld_m b JOIN fire_m f ON ST_Intersects(b.geom, ST_Buffer(f.geom,500))
WHERE b.osm_id NOT IN (SELECT osm_id FROM a) GROUP BY b.osm_id) t;
최단_m | 최장_m
--------+--------
393.1 | 499.6
393.1m에서 499.6m 사이. 12,720개가 전부 그 좁은 띠에 몰려 있습니다. 소방서 바로 옆 건물은 한 채도 없어요 — 그것들은 이미 첫 번째 셈에서 잡혔으니까요.
💡 하한이 397m가 아니라 393.1m인 게 눈에 걸리실 겁니다. 397m 원 안이면 진작 잡혔어야 하니까요. 실제로 그 최단 건물을 꺼내 보면 3857 눈금으로 498.08 — 500 안쪽인데도 버퍼가 놓쳤습니다. 반지름과는 다른 두 번째 원인이 있다는 뜻인데, 바로 다음 절의 주제입니다.
건물은 처음부터 거기 있었고, 조회한 자리도 같았습니다. 그물이 작았을 뿐이에요.
2. 그런데 그건 원이 아니었다
숫자가 하나 더 남아 있습니다. 앞에서 세검정 소방서 한 곳을 ST_DWithin으로 쟀을 때 34개가 나왔었죠. 방금 전체 집계에는 ST_Buffer를 썼고요. 같은 소방서에 둘을 나란히 걸어 봅니다.
WITH fs AS (SELECT geom AS g FROM fire_m WHERE osm_id = 2279744031)
SELECT (SELECT count(*) FROM bld_m b, fs WHERE ST_Intersects(b.geom, ST_Buffer(fs.g, 500))) AS 버퍼,
(SELECT count(*) FROM bld_m b, fs WHERE ST_DWithin(b.geom, fs.g, 500)) AS dwithin;
버퍼 | dwithin
------+---------
33 | 34
같은 좌표계, 같은 500인데 답이 다릅니다. 이유는 ST_Buffer가 돌려준 도형을 열어 보면 나옵니다.
SELECT ST_NPoints(ST_Buffer(geom, 500)) AS 점개수,
round(ST_Area(ST_Buffer(geom, 500))::numeric, 2) AS 버퍼넓이,
round((pi() * 500 * 500)::numeric, 2) AS 진짜원,
round((ST_Area(ST_Buffer(geom,500)) / (pi()*500*500) * 100)::numeric, 3) AS 비율
FROM fire_m WHERE osm_id = 2279744031;
점개수 | 버퍼넓이 | 진짜원 | 비율
--------+-----------+-----------+--------
33 | 780361.29 | 785398.16 | 99.359
점이 33개뿐입니다. 원이 아니라 32각형이에요.
당연하다면 당연합니다. PostGIS의 geometry에는 원이라는 타입이 없습니다. 저장할 수 있는 건 좌표 목록뿐이니, "원"을 만들라고 하면 원에 가까운 다각형을 그려줄 수밖에 없어요.
💡 33인데 왜 32각형인가요?
ST_NPoints가 하나를 더해 준 게 아니라, 도형 자체에 점이 33개 들어 있습니다. 폴리곤의 링은 반드시 닫혀 있어야 해서ST_Buffer가 만들 때부터 첫 점을 마지막에 한 번 더 저장하거든요.quad_segs = 1버퍼를 그대로 출력해 보면 보입니다.textPOLYGON((952911.01 1956402.16, 952411.01 1955902.16, 951911.01 1956402.16, 952411.01 1956902.16, 952911.01 1956402.16)) └─ 첫 점과 마지막 점이 같다. 모서리는 4개, 저장된 점은 5개그래서 링이 늘면 그만큼 같이 늡니다 — 구멍 하나 뚫린 사각형은 링이 둘이라
ST_NPoints가 10이에요.
그래서 넓이가 모자랍니다. 반지름 500m 원의 진짜 넓이는 π × 500² = 785,398.16 m²인데 버퍼는 780,361.29 m², 99.359% 입니다. 다각형이 원 안쪽에 내접하니까 5,037 m²가 부족한 거죠.
몇 각형으로 그릴지는 세 번째 인자로 정합니다. quad_segs — 사분원 하나를 몇 개의 선분으로 나눌지입니다. 기본값이 8이라 사분원마다 8개, 한 바퀴에 32개가 되는 거예요.
SELECT q AS quad_segs,
ST_NPoints(ST_Buffer(geom, 500, q)) AS 점개수,
round(ST_Area(ST_Buffer(geom, 500, q))::numeric, 1) AS 넓이,
round((ST_Area(ST_Buffer(geom,500,q)) / (pi()*500*500) * 100)::numeric, 3) AS 원_대비
FROM fire_m, unnest(ARRAY[1,2,4,8,16,64]) q
WHERE osm_id = 2279744031 ORDER BY q;
quad_segs | 점개수 | 넓이 | 원_대비
-----------+--------+----------+---------
1 | 5 | 500000.0 | 63.662
2 | 9 | 707106.8 | 90.032
4 | 17 | 765366.9 | 97.450
8 | 33 | 780361.3 | 99.359
16 | 65 | 784137.1 | 99.839
64 | 257 | 785319.3 | 99.990
맨 윗줄이 재미있습니다. quad_segs = 1이면 점 5개, 즉 모서리 4개 — 정사각형입니다. 넓이 비율 **63.662%**는 2/π고요. 원에 내접하는 정사각형의 넓이가 정확히 그만큼이거든요.
그 0.641%에 유치원이 하나 있었다
넓이가 0.641% 모자란 게 실제로 문제가 되냐면, 됩니다. 32각형은 원 안쪽에 내접하니 원 안이지만 다각형 밖인 초승달 조각이 변마다 하나씩, 서른두 군데 생기거든요. 거기 뭔가 놓여 있으면 조용히 빠집니다.
세검정 소방서에서 빠진 그 하나를 찾아봤습니다.
WITH fs AS (SELECT geom AS g FROM fire_m WHERE osm_id = 2279744031)
SELECT b.osm_id, b.building,
round(ST_Distance(b.geom, fs.g)::numeric, 2) AS 거리_m,
ST_Intersects(b.geom, ST_Buffer(fs.g, 500)) AS 기본버퍼,
ST_Intersects(b.geom, ST_Buffer(fs.g, 500, 64)) AS 정밀버퍼
FROM bld_m b, fs
WHERE ST_DWithin(b.geom, fs.g, 500) AND NOT ST_Intersects(b.geom, ST_Buffer(fs.g, 500));
osm_id | building | 거리_m | 기본버퍼 | 정밀버퍼
------------+--------------+--------+----------+----------
1280181582 | kindergarten | 499.93 | f | t
유치원입니다. 소방서에서 499.93m — 500m 선 안쪽으로 7cm 들어와 있어요. 기본 버퍼로는 빠지고, quad_segs를 64로 올리면 잡힙니다.
그럴 거면 도형을 만들 이유가 없다
여기서 방향을 바꿔야 합니다. 우리가 묻고 싶은 건 "이 건물이 소방서에서 500m 안인가" 였어요. 그건 거리에 대한 질문이지 도형에 대한 질문이 아닙니다. 원을 만든 건 그 질문을 표현하는 한 가지 방법이었을 뿐인데, 그 과정에서 32각형이라는 불필요한 근사가 끼어든 거죠.
거리를 직접 물으면 근사가 없습니다. ST_DWithin(A, B, r) — "A와 B의 거리가 r 이하인가".
SELECT count(DISTINCT b.osm_id) AS 건물
FROM bld_m b JOIN fire_m f ON ST_DWithin(b.geom, f.geom, 500);
건물
-------
36607
36,607개. 이게 과제의 답입니다.
정말 같은 질문의 두 표현인지 확인해 봅시다. 버퍼의 quad_segs를 올려 원에 가깝게 만들면 ST_DWithin으로 수렴해야 합니다.
SELECT (SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire_m f ON ST_Intersects(b.geom, ST_Buffer(f.geom,500))) AS "버퍼_기본8",
(SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire_m f ON ST_Intersects(b.geom, ST_Buffer(f.geom,500,64))) AS "버퍼_64",
(SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire_m f ON ST_Intersects(b.geom, ST_Buffer(f.geom,500,256))) AS "버퍼_256",
(SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire_m f ON ST_DWithin(b.geom, f.geom, 500)) AS "ST_DWithin";
버퍼_기본8 | 버퍼_64 | 버퍼_256 | ST_DWithin
------------+---------+----------+------------
36415 | 36605 | 36607 | 36607
36,415 → 36,605 → 36,607 = 36,607. quad_segs를 256까지 올리면 ST_DWithin과 한 개도 다르지 않습니다. 둘은 같은 질문이었고, 차이는 전부 다각형 근사에서 나왔던 거예요.
그러니 반경 안인지만 묻는다면 ST_DWithin 입니다. ST_Buffer는 원 자체가 결과물일 때 — 지도에 그리거나 다른 면과 겹쳐야 할 때 꺼내는 거고요. 3절부터가 정확히 그 경우입니다.
💡
ST_DWithin을 쓸 이유가 하나 더 있는데, 속도입니다. 같은 답을 내면서 훨씬 빠르거든요. 다만 "인덱스가 있는데도 왜 어떤 함수는 안 빨라지나"는 그 자체로 한 편 분량이라 Phase 7에서 따로 다룹니다. 오늘은 정확성만 봤습니다.
3. 겹친 것을 합친다 — ST_Union
이제 과제의 뒷부분입니다. "서울 땅의 몇 %가 소방서에서 500m 안인가." 건물 개수가 아니라 면적 비율을 묻는 거예요.
넓이를 구하려면 원을 도형으로 들고 있어야 하니 여기서는 ST_Buffer가 제대로 쓰입니다. 그런데 그냥 더하면 안 됩니다. 소방서 121곳의 원들이 서로 겹치거든요.
겹침을 무시하면 어떻게 되는지는 이미 겪었습니다. 앞의 쿼리들에 붙어 있던 DISTINCT가 그거예요. 빼 보면 차이가 드러납니다.
SELECT count(*) AS 단순합,
count(DISTINCT b.osm_id) AS 중복제거
FROM bld_m b JOIN fire_m f ON ST_DWithin(b.geom, f.geom, 500);
단순합 | 중복제거
--------+----------
39943 | 36607
건물 3,336개가 두 번 이상 세어졌습니다. 소방서 두 곳에서 동시에 500m 안인 건물들이죠.
넓이도 똑같습니다. 버퍼 121개의 넓이를 그냥 더하면 겹친 부분을 여러 번 더하게 돼요. ST_Union이 그걸 하나의 도형으로 녹여 줍니다.
SELECT round((sum(ST_Area(ST_Buffer(geom,500))) / 1e6)::numeric, 3) AS 낱개합_km2,
round((ST_Area(ST_Union(ST_Buffer(geom,500))) / 1e6)::numeric, 3) AS 합친뒤_km2,
ST_NumGeometries(ST_Union(ST_Buffer(geom,500))) AS 덩어리
FROM fire_m;
낱개합_km2 | 합친뒤_km2 | 덩어리
------------+------------+--------
94.424 | 87.313 | 106
94.424 km²가 87.313 km²로. 7.1 km²가 중복이었습니다.
셋째 칸도 봐 두세요. 원 121개를 합쳤더니 덩어리가 106개입니다. 15번의 합쳐짐이 일어났다는 뜻이에요 — 서로 1km 안에 있는 소방서들이 그만큼 있었던 겁니다. ST_Union은 떨어진 것은 떨어진 채로 두고 MultiPolygon으로 묶습니다.
지난 편의 1.02 m²가 여기서 사라진다
ST_Union이 겹침을 녹인다는 걸 더 확실히 보여주는 데이터가 있습니다. 지난 편에서 서울 25개 구 중 21쌍이 서로 미세하게 겹쳐 있다는 걸 발견했었죠. 경계선을 디지털로 그리는 과정에서 생긴 틈, 슬리버입니다. 21쌍의 겹친 넓이를 전부 더하면 1.02 m² 였고요.
그럼 25개 구를 ST_Union으로 합치면 그 1.02 m²가 사라져야 합니다. 겹친 부분을 한 번만 세게 되니까요.
SELECT round(sum(ST_Area(geom))::numeric, 2) AS 낱개합,
round(ST_Area(ST_Union(geom))::numeric, 2) AS 합친뒤,
round((sum(ST_Area(geom)) - ST_Area(ST_Union(geom)))::numeric, 2) AS 사라진넓이
FROM sgg WHERE code LIKE '11%';
낱개합 | 합친뒤 | 사라진넓이
--------------+--------------+------------
606993569.91 | 606993568.89 | 1.02
정확히 1.02 m². 지난 편은 ST_Overlaps로 겹친 쌍을 찾아 ST_Intersection의 넓이를 더했고, 오늘은 ST_Union으로 합쳐 낱개 합과의 차이를 봤습니다. 완전히 다른 두 경로가 같은 숫자에 도착했어요.
⚠️ 이 절의
sgg계산은 원래 좌표계인 5174에서 했습니다. 소방서 커버리지 쪽은 5179고요. 두 좌표계는 서울 넓이를 606.99 km²와 606.5 km²로 조금 다르게 냅니다(0.08% 차이). 같은 표에 두 숫자를 나란히 놓으면 안 됩니다.
4. 겹치고 깎아낸다 — ST_Intersection · ST_Difference
이제 답을 낼 수 있습니다. 커버리지 도형과 서울 도형을 겹치면 되죠.
ST_Intersection(A, B)— 둘 다에 속하는 부분ST_Difference(A, B)— A에서 B를 뺀 부분
WITH cover AS (SELECT ST_Union(ST_Buffer(geom, 500)) AS c FROM fire_m),
seoul AS (SELECT ST_Union(ST_Transform(geom, 5179)) AS s FROM sgg WHERE code LIKE '11%')
SELECT round((ST_Area(s) / 1e6)::numeric, 1) AS 서울,
round((ST_Area(ST_Intersection(s, c)) / 1e6)::numeric, 1) AS 덮인곳,
round((ST_Area(ST_Difference(s, c)) / 1e6)::numeric, 1) AS 안덮인곳,
round((ST_Area(ST_Intersection(s, c)) / ST_Area(s) * 100)::numeric, 1) AS 커버율
FROM seoul, cover;
서울 | 덮인곳 | 안덮인곳 | 커버율
-------+--------+----------+--------
606.5 | 63.6 | 543.0 | 10.5
서울의 10.5%. 과제의 두 번째 답입니다.
검산도 됩니다. 63.6 + 543.0 = 606.6 — 서울 전체 606.5와 반올림 오차 안에서 맞아요. ST_Intersection과 ST_Difference를 더하면 원본입니다.
앞에서 잰 커버리지 전체는 87.313 km²였는데 덮인 곳은 63.6 km²죠. 나머지 23.7 km²는 서울 밖입니다 — 지난 편에서 자른 bbox가 서울보다 넓어 경기도 소방서까지 들어와 있거든요(121곳 중 서울 안은 89곳).
구별로 쪼개면 더 쓸 만해집니다.
WITH cover AS (SELECT ST_Union(ST_Buffer(geom, 500)) AS c FROM fire_m)
SELECT s.name,
round((ST_Area(ST_Intersection(ST_Transform(s.geom, 5179), c.c))
/ ST_Area(ST_Transform(s.geom, 5179)) * 100)::numeric, 1) AS 커버율
FROM sgg s, cover c
WHERE s.code LIKE '11%'
ORDER BY 2 DESC;
name | 커버율
----------+--------
중구 | 23.5
양천구 | 19.6
동대문구 | 18.9
중랑구 | 16.6
마포구 | 16.3
…
노원구 | 4.9
강북구 | 3.8
관악구 | 3.0
중구 23.5%, 관악구 3.0%. 8배 가까이 차이가 납니다.
⚠️ 이 숫자를 소방 서비스의 품질로 읽으면 안 됩니다. 소방서 위치는 OSM 자원봉사자가 그린 것이라 누락이 있을 수 있고, 무엇보다 500m 직선거리는 출동 시간이 아닙니다. 실제로 중요한 건 도로를 따라 몇 분인지고, 그건 직선 반경으로 답할 수 없어요. 도로망을 그래프로 놓고 도달권을 구하는 이야기는 Phase 8의 주제입니다.
5. 도형에서 수를 뽑는다 — ST_Centroid · ST_Area
과제는 풀렸습니다. 그런데 1절에서 미뤄 둔 게 하나 있었죠. 소방서 면 86개를 ST_Centroid로 점을 만들어 붙인 것 — 그게 안전했을까요?
ST_Centroid는 도형의 무게중심을 돌려줍니다. 균질한 판으로 오려내 손가락에 올렸을 때 균형이 잡히는 지점이죠. 문제는 그 지점이 도형 밖일 수 있다는 겁니다. 도넛이나 초승달을 떠올리면 바로 보여요.
시군구 251개로 확인해 봅시다.
SELECT count(*) FILTER (WHERE NOT ST_Contains(geom, ST_Centroid(geom))) AS centroid_밖,
count(*) FILTER (WHERE NOT ST_Contains(geom, ST_PointOnSurface(geom))) AS onsurface_밖,
count(*) AS 전체
FROM sgg;
centroid_밖 | onsurface_밖 | 전체
-------------+--------------+------
8 | 0 | 251
251개 중 8개가 자기 중심점을 자기 안에 담지 못합니다. 누군지 봅시다.
SELECT name, ST_NumGeometries(geom) AS 조각,
round(ST_Distance(geom, ST_Centroid(geom))::numeric, 0) AS 밖으로_m
FROM sgg WHERE NOT ST_Contains(geom, ST_Centroid(geom)) ORDER BY 3 DESC;
name | 조각 | 밖으로_m
--------------+------+----------
옹진군 | 129 | 9201
안산시단원구 | 16 | 5224
청원군 | 1 | 3028
신안군 | 310 | 1920
여수시 | 184 | 1460
완도군 | 137 | 1019
사천시 | 32 | 97
남해군 | 49 | 30
옹진군의 중심점은 자기 땅에서 9.2km 떨어진 바다 위에 있습니다. 조각이 129개 — 서해에 흩어진 섬들이니 그 무게중심이 섬 사이 바다에 찍히는 게 당연하죠. 신안군은 조각이 310개고요.
그런데 청원군은 조각이 1개입니다. 섬이 아닌데도 3km 밖으로 나가요. 청주시를 도넛처럼 감싸고 있는 모양이라 그렇습니다. 가운데 구멍이 청주시고, 무게중심은 그 구멍 한복판에 찍히는 거죠.
지도에 구 이름을 중심점 자리에 찍는 흔한 작업을 생각해 보면, 여덟 개의 라벨이 남의 땅이나 바다에 떨어집니다.
해법은 이름부터 그렇게 생긴 함수입니다. ST_PointOnSurface — "표면 위의 점". 무게중심은 포기하는 대신 반드시 도형 안에 있는 점을 돌려주고, 앞 결과에서 봤듯 251개 전부에서 성공합니다. 라벨이나 대표점이 필요한 거라면 이쪽입니다.
그래서 소방서는 어땠나
1절로 돌아갑니다. 소방서 면 86개에 ST_Centroid를 썼는데, 그중 밖으로 나간 게 있었을까요?
SELECT osm_id, name, building,
round(ST_Distance(way, ST_Centroid(way))::numeric, 1) AS 밖으로_m
FROM planet_osm_polygon
WHERE amenity='fire_station' AND NOT ST_Contains(way, ST_Centroid(way));
osm_id | name | building | 밖으로_m
-----------+------------+----------+----------
542869466 | 은평소방서 | yes | 5.1
278615353 | 서초소방서 | yes | 0.2
두 곳이 나갔습니다. 은평소방서는 자기 건물에서 5.1m 떨어진 곳에 중심점이 찍혔어요 — ㄱ자나 ㄷ자로 꺾인 건물이면 이렇게 됩니다.
최종 집계가 달라지는지도 재봤습니다.
WITH fire2 AS (
SELECT ST_Transform(way,5179) AS geom FROM planet_osm_point WHERE amenity='fire_station'
UNION ALL
SELECT ST_Transform(ST_PointOnSurface(way),5179) FROM planet_osm_polygon WHERE amenity='fire_station')
SELECT (SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire_m f ON ST_DWithin(b.geom,f.geom,500)) AS centroid,
(SELECT count(DISTINCT b.osm_id) FROM bld_m b JOIN fire2 f ON ST_DWithin(b.geom,f.geom,500)) AS pointonsurface;
centroid | pointonsurface
----------+----------------
36607 | 36620
13개 차이. 5m짜리 어긋남 두 개가 최종 답을 13건 움직였습니다. 12,912개에 비하면 작지만, 1절에서 아무 생각 없이 고른 함수 하나가 여기까지 흘러왔다는 게 요점이에요.
ST_Area도 좌표계를 탄다
ST_Area는 넓이를 돌려줍니다. 이미 3·4절에서 계속 썼죠. 그런데 이 함수도 ST_Buffer와 똑같은 이유로 좌표계를 탑니다. 종로구 하나를 네 가지로 재봤습니다.
SELECT round(ST_Area(geom)::numeric, 1) AS "5174_m2",
ST_Area(ST_Transform(geom, 4326)) AS "4326_제곱도",
round(ST_Area(ST_Transform(geom, 4326)::geography)::numeric, 1) AS "geography_m2",
round(ST_Area(ST_Transform(geom, 3857))::numeric, 1) AS "3857_m2"
FROM sgg WHERE code = '11010';
5174_m2 | 4326_제곱도 | geography_m2 | 3857_m2
------------+----------------------+--------------+------------
23892610.4 | 0.002437515019652436 | 23892213.1 | 38121750.1
믿을 만한 건 첫 칸입니다 — 미터 좌표계에서 잰 23.9 km². ::geography가 타원체 위에서 잰 값도 여기서 0.002%밖에 안 벗어나고요. 4326은 넓이가 아니라 제곱도라 아예 읽을 수 없는 숫자입니다.
문제는 마지막 칸이에요. 3857이 종로구를 38.1 km²로 봅니다. 실제의 1.6배입니다.
1절과 같은 원인인데, 제곱으로 나타납니다. 거리가 cos φ 배로 줄었으니 넓이는 1/cos² φ 배로 붑니다. 종로구 중심의 위도가 37.59°이니 1 / cos²(37.59°) = 1.5928이고, 실측 배율은 38,121,750 / 23,892,610 = 1.5955. 0.2% 안에서 맞습니다 — 종로구가 남북으로 폭이 있어 위도마다 배율이 조금씩 다르니 중심 한 점으로 계산한 값과 완전히 같을 수는 없어요.
그러니 규칙은 하나로 정리됩니다. 거리든 넓이든, 재기 전에 미터 좌표계로 옮긴다.
6. 거꾸로 물어본다 — 가장 가까운 것
마지막으로 근접 질문의 다른 형태 하나. 지금까지는 "반경 안에 무엇이 있나" 였는데, 실무에서 그만큼 자주 나오는 게 "가장 가까운 것이 무엇이고 얼마나 머나" 입니다. 반경을 미리 정할 수 없을 때죠.
PostGIS에는 이걸 위한 연산자가 있습니다. <-> — 두 도형 사이 거리를 돌려주는데, ORDER BY에 쓰면 인덱스가 가까운 순서대로 꺼내 줍니다. 전부 재서 정렬하는 게 아니라요. 이런 검색을 KNN(K-Nearest Neighbor)이라고 부릅니다.
편의점 7,100개마다 가장 가까운 소방서를 찾아봤습니다.
CREATE TABLE conv_m AS SELECT osm_id, name, ST_Transform(geom, 5179) AS geom FROM convenience;
CREATE INDEX ON conv_m USING gist(geom);
ANALYZE conv_m;
WITH d AS (
SELECT c.osm_id,
(SELECT c.geom <-> f.geom FROM fire_m f ORDER BY c.geom <-> f.geom LIMIT 1) AS dist
FROM conv_m c)
SELECT count(*) AS 편의점, round(avg(dist)::numeric,0) AS 평균_m,
round(min(dist)::numeric,0) AS 최소_m, round(max(dist)::numeric,0) AS 최대_m,
count(*) FILTER (WHERE dist <= 500) AS "500m_안"
FROM d;
편의점 | 평균_m | 최소_m | 최대_m | 500m_안
--------+--------+--------+--------+---------
7100 | 1037 | 14 | 4844 | 1083
평균 1,037m, 최대 4,844m. 500m 안에 있는 건 1,083개로 전체의 15%고요.
⚠️ 인덱스를 안 걸면 이 쿼리는 안 끝납니다. 실제로 처음엔 임시 뷰로 돌렸다가 2분을 넘겨 중단했어요.
<->가 인덱스를 타야 "가까운 순으로 몇 개만" 꺼내는데, 인덱스가 없으면 편의점 7,100개마다 소방서 121개를 전부 재고 정렬합니다. 위 코드에서CREATE INDEX가 그래서 붙어 있습니다.
정리 — 한 문장으로
PostGIS는 좌표를 숫자로만 볼 뿐 단위를 모르고, 원을 저장할 수 없어 다각형으로 근사한다. 그래서 "반경 500m"라고 쓴 쿼리가 3857에서는 397m 32각형이 되고, 아무 에러 없이 건물 12,912개를 빠뜨린 답을 돌려준다.
같은 질문에 대한 세 개의 답이 어디서 갈렸는지 되짚으면 이렇습니다.
| 쿼리 | 건물 | 무엇이 틀렸나 |
|---|---|---|
ST_Buffer(way, 500) — 3857 | 23,695 | 500이 미터가 아니었다 (실제 397m) |
ST_Buffer(geom, 500) — 5179 | 36,415 | 원이 아니라 32각형이었다 |
ST_DWithin(geom, geom, 500) | 36,607 | — |
오늘 쓴 함수를 묻고 싶은 것에서 거꾸로 찾아가면 이렇게 됩니다.
| 묻고 싶은 것 | 쓸 함수 | 조심할 것 |
|---|---|---|
| 반경 N 안에 있는가 | ST_DWithin | 좌표계 단위 |
| 그 반경을 도형으로 갖고 싶다 | ST_Buffer | 다각형 근사(quad_segs) |
| 겹치는 도형들을 하나로 | ST_Union | 안 하면 중복 집계 |
| 두 도형이 겹치는 부분 | ST_Intersection | — |
| A에서 B를 뺀 부분 | ST_Difference | — |
| 넓이 | ST_Area | 좌표계 단위(제곱으로 틀린다) |
| 면을 대표하는 점 | ST_PointOnSurface | ST_Centroid는 밖으로 나갈 수 있다 |
| 가장 가까운 것 | <-> + ORDER BY | 인덱스 없으면 안 끝난다 |
오늘 이 표에 한 줄이 두 번 등장했습니다. 좌표계 단위요. ST_Buffer에서 한 번, ST_Area에서 제곱으로 또 한 번. 그리고 그때마다 우회로가 ST_Transform이었죠.
그런데 gis-21부터 오늘까지, 거리를 잴 때마다 슬쩍 등장한 또 다른 방법이 있었습니다. ::geography 캐스팅이요. 오늘도 3857 왜곡을 측정할 때 썼고, ST_Area 표에도 한 줄 있었는데 한 번도 설명하지 않았습니다.
다음 편은 그 빚을 갚습니다. geometry vs geography — 같은 두 점의 거리가 타입에 따라 8,770.9m와 8,767.6m로 갈리는데 그 3.3m가 어디서 오는지, 그리고 좌표계 변환이 필요 없는 그 편함의 대가로 무엇을 포기하는지 봅니다.
참고
- https://postgis.net/docs/ST_Buffer.html
- https://postgis.net/docs/ST_DWithin.html
- https://postgis.net/docs/ST_Union.html
- https://postgis.net/docs/ST_Intersection.html
- https://postgis.net/docs/ST_Difference.html
- https://postgis.net/docs/ST_Centroid.html
- https://postgis.net/docs/ST_PointOnSurface.html
- https://postgis.net/docs/ST_Area.html
- https://postgis.net/docs/geometry_distance_knn.html
- https://postgis.net/workshops/postgis-intro/geometry_returning.html
- https://libgeos.org/doxygen/classgeos_1_1operation_1_1buffer_1_1BufferParameters.html
- https://proj.org/en/stable/operations/projections/webmerc.html
- https://epsg.io/3857
- https://www.postgresql.org/docs/17/gist.html