GIS&GPS

매핑정보가 있는 .points 파일과 온맵이미지로 좌표 정보가 삽입된 geotiff 파일을 생성할

차례입니다.

파이썬으로 gdal 명렁어를 생성하는 코드를 짰습니다.

onmap2tif.py

OSGeo4W 명렁어창에 python onmap2tif.py 를 입력 후 엔터키를 누르면

D:\50k 폴더에 38815.bat 파일이 생성되었습니다.

OSGeo4W 명렁어창에 d:\50k\38815.bat 을 입력 후 엔터키를 누르면  

gdal 명령어가 실행되고 잠시 후 38815.tif 가 만들어졌습니다.

배치파일이 만들어서 하나하나 실행하려면 번거로우니

배치파일을 불러서 실행하는 배치파일을 만들었습니다.

1. D:\50k\bat 폴더를 만들고 모든 배치파일을 만든 폴더로 이동

2. 도스창에 dir /s/b d:\50k\bat > filelist50k.bat 을 입력 후 엔터

 -> d:\50k\bat 폴더안의 파일리스트가 담긴 filelist50k.bat 를 생성합니다.

     filelist5k.bat 는 cmd 창의 현재폴더에 생성됩니다.

3. filelist50k.bat 를 편집기로 열어서 첫째줄에 @echo off 을 삽입하고

파일리스트 앞에 call 을 삽입.

@eocho off
call d:\50k\bat\34504.bat
call d:\50k\bat\34505.bat
call d:\50k\bat\34506.bat
call d:\50k\bat\34508.bat

............. 이하 생략 ............

4. OSGeo4W 에서 filelist50k.bat 를 실행시키면 1:50,000 축척의 온맵 geotiff 파일이

생성됩니다.

작업에 사용했던 파일입니다.

250k_bat.zip

50k_bat.zip

QGIS 는 오픈소스 GIS 프로그램입니다.

전에 소개했던 OSGeo4W 를 설치했다면 이미 깔려있습니다.

http://gisngps.tistory.com/2

qgis 플러그인 중 GDAL 지오레퍼런서에 대해 예제를 통해서 살펴보겠습니다.

시작 > 모든프로그램 > OSGeo4W > qgis desktop 을 실행시키고

메뉴> 플러그인 > 플러그인 설치 및 관리... 를 누릅니다.

설치됨 항목에서 GDAL 지오레퍼런서를 체크해서 플러그인을 활성화시키면

메뉴 > 래스터 > 지오레퍼런서 메뉴가 생깁니다.

지오레퍼런서 메뉴 > 설정 > 지오레퍼런스 설정 에서 "ID 표시" 와 "좌표보기"를 체크.

1. 래스터 열기 아이콘을 눌러서 매칭할 이미지를 열면 좌표계 선택 창이 뜹니다.

2. 필터창에 4326을 입력하면 

3. 세계의 좌표계 란에 WGS84(EPSG:4326)이 나옵니다. 선택하고 확인을 누릅니다.

1.포인트 추가 아이콘을 누르고

2. 좌표를 알고 있는 지점을 마우스 좌클릭하면

3. 지도 좌표 입력창이 뜨는데 x, y 값(경도, 위도)을 써줍니다.

아래가 좌표점을 하나 입력한 후 모습입니다.

지도 아래쪽 GCP 테이블에 포인트에 대한 값이 나옵니다.

id, srcX, srcY, defX, defY,..(포인트ID, 이미지X 픽셀값, 이미지Y픽셀값, 경도값, 위도값) 

계속해서 총 16개의 포인트를 추가한 화면입니다.

이미지 픽셀값의 원점은 좌상단입니다.

꽤 번거로운 작업입니다. 

그래서 GCP 테이블을 .points 파일로 저장하고 불러오는 기능이 있습니다.

아래는 예제에 쓰인 1:50,000 축척 지형도(33602_추자)에 해당하는 points 파일입니다.

지형도 도엽코드 규칙을 이용하면 각 지형도에 해당하는 points 파일을 따로 만들어서

좀더 쉽게 지오레퍼런싱을 할 수 있습니다.

mapX,mapY,pixelX,pixelY,enable

126.25,33.75,543.68037159763355248,-8381.70678816568033653,1

126.33333333299999879,33.75,2671.43259171597628665,-8396.54119289940899762,1

126.41666666700000121,33.75,4799.14044497041413706,-8409.49250177514841198,1

126.5,33.75,6926.38575147928986553,-8420.59362366863933858,1

126.25,33.83333333329999704,562.22660828402388233,-5834.97752426035458484,1

126.33333333299999879,33.83333333329999704,2687.62172781065100935,-5847.46628639053233201,1

126.41666666700000121,33.83333333329999704,4813.47939408284037199,-5860.41759526627174637,1

126.5,33.83333333329999704,6938.8745136094676127,-5873.36890414201207022,1

126.25,33.91666666670000296,581.19102485207122299,-3287.75280473372777124,1

126.33333333299999879,33.91666666670000296,2704.27341065088785399,-3300.24156686390506366,1

126.41666666700000121,33.91666666670000296,4828.28088994082827412,-3313.65542248520705471,1

126.5,33.91666666670000296,6951.36327573964445037,-3326.60673136094646907,1

126.25,34,599.23034792899420609,-740.06553846153838094,1

126.33333333299999879,34,2720.92509349112424388,-754.86703431952651044,1

126.41666666700000121,34,4842.15729230769193236,-768.28088994082827412,1

126.5,34,6963.85203786982219754,-778.45691834319518421,1

이미지에 좌표를 부여할 차례입니다.

georeference, map matching, mapping 등으로 불리죠.

맵매칭 방법도 여러가지가 있지만 OSGeo4W 라는 오픈소스 GIS 프로그램을 사용했습니다.

이미지와 매칭포인트를 알고 있다면 OSGeo4W shell 으로 매핑작업을 하면 됩니다.
아래가 시험삼아 5만 축척 도엽코드 36905 온맵으로 작업한 겁니다.

gdal_translate -of GTiff -gcp 671.84 8494.73 127 36.5 -gcp 2728.18 8493.45 127.08333 36.5 -gcp 4784.65 8491.31 127.16667 36.5 -gcp 6841.04 8486.74 127.25 36.5 -gcp 671.413 5946.92 127 36.58333 -gcp 2726.04 5944.79 127.08333 36.58333 -gcp 4779.95 5941.91 127.16667 36.58333 -gcp 6834.21 5939.35 127.25 36.58333 -gcp 671.413 3400.01 127 36.66667 -gcp 2723.96 3397.45 127.08333 36.66667 -gcp 4775.87 3394.65 127.16667 36.66667 -gcp 6827.79 3392.08 127.25 36.66667 -gcp 672.156 853.006 127 36.75 -gcp 2721.4 851.724 127.08333 36.75 -gcp 4771.17 849.586 127.16667 36.75 -gcp 6820.81 845.151 127.25 36.75 "D:/georef_test/50k/36905.png" "D:/georef_test/50k/tt36905.png"
gdalwarp -t_srs EPSG:4326 -r cubic -tps -co COMPRESS=LZW  "D:/georef_test/50k/tt36905.png" "D:/georef_test/50k/36905_10p.tif"
gdal_translate -projwin 127.0 36.75 127.25 36.5 -of GTiff D:/georef_test/50k/36905_10p.tif D:/georef_test/50k/36905_10pc.tif

세가지 명령어로 되어있습니다.

gdal 로 시작되는 gdal 유틸리티를 이용한 겁니다. (아래 사이트 참조)

http://www.gdal.org/gdal_utilities.html

1. gdal_translate -of GTiff -gcp xxx xxx xxx xxx ...............

-> 이미지에 그라운드 컨트롤 포인트 (Ground Control Point) 를 부여해줍니다.

    -gcp 이후에 나오는 숫자 4개는 이미지픽셀x, 이미지픽셀y, 경도, 위도 순서입니다.

2. gdalwarp -t_srs EPSG:4326 .....

-> 매핑한 이미지를 goetiff 파일로 변환합니다.

3. gdal_translate -prowin xxx xxx xxx xxx ....

-> 지도에서 필요한 부분만 크롭

   -prowin 이후의 4개의 숫자는 좌측경도, 상단위도, 우측경도, 하단위도 순서입니다.

온맵 도엽마다 매핑에 필요한 gcp(Ground Control Point) 와 크롭할 좌표를 알고 있다면

gdal 명령어로 간단하게 매칭할 수 있습니다.

하지만, gcp 에서 1,2번째 숫자인 이미지펙셀값을 계산을 통해 알아낼 수 있는 방법을

못찾았습니다. 다행히 3,4번째 숫자인 경위도 좌표는 도엽코드로 계산할 수 있습니다.

맵매칭에서 가장 지루한  Grond Cotrol Point 작업은 최대한 적게 하기 위해서

qgis 지오레퍼런싱에서 사용할 .points 파일을 생성하는 파이썬 코드을 짰습니다.