qgis3でジオリファレンスする方法

#ジオリファレンスとは
　おそらく「Geo + Reference」で地理情報を扱う際に、位置情報が無い地図に位置情報を設定する作業のこと。
　朝倉書店刊『地理情報科学事典』には、ジオリファレンスという項目は無いが58～61ページの「5－2位置合わせ」が該当すると思います。
　パスコの[用語集とGISの使い方](https://www.pasco.co.jp/recommend/word/word089/)によると、イメージファイルにGISで使用するために位置情報等を設定する作業のことと解説されていますね。
　
　こういった文章ではイメージしにくいのでたとえ話。
　例えば、田舎のおじいちゃんの家に遊びにいったとき、土蔵から古い絵地図を見つけたとする。内容を推測するに、どうもご先祖様が宝を埋めた地図のようだが、そのままでは、どこに埋まっているかわからない。
　宝探しをするためには、まず地図がどこの範囲を示しているのか特定する必要がある。
　次に宝のマークがある位置を特定するには、絵地図のランドマークが、実際の世界でどこの何を示しているのかの対応をとり、そのランドマークを頼りにめどをつける必要があるわけですね。
　で、この作業をする際、まず地図上で絵地図の内容を確認していくわけですが、現在の地図と絵地図を重ね合わせて、ぴったり重なることはまずありません。重ね合わせるためには、絵地図を描かれている内容を確認しつつ、変形して重ねることになりますが、「現在の地図のように、一定の規則で描かれていない」「測量技術が未熟で包囲や距離などが不正確」といったことが問題になるわけです。
　しかし、絵地図などをスキャナなどでデータすると、幾何変換などの作業はコンピュータが行ってくれるわけです。
　そとて、ここで幾何変換をどのように行うのかを指示するための指示や実際の変換作業がジオリファレンスということになります。
　qgisでは、ジオリファレンスを行うための機能があり、基準になる地図と変換したい地図を比較し、地図上で一致するポイントを指示し、変換方法を指定して変換するという流れになります。
　
#題材「戦前の地図」
　宝探しの地図といった魅力的な題材が無いため、スタンフォード大学が公開している戦前の日本の測量図を現在の地図と重ね合わせてみます。
　スタンフォード大学では、アメリカが日本を占領していた際にアメリカ軍が接収した日本陸軍の測量図を[ネット上で公開](https://stanford.maps.arcgis.com/apps/SimpleViewer/index.html?appid=733446cc5a314ddf85c59ecc10321b41)しています。　
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/776ac9bb-74f9-0078-07ba-2fab1ff0dac9.png)

　ここから、秋田市の地図をダウロンードします。
　地図が表示されているウィンドウの右上に、矢印が飛び出そうとしているアイコンがあるので、これをクリックすると解像度を選んでダウンロードできます。
　![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/5fb8ada1-0733-cf23-30e3-a442b26fb530.png)
　オリジナル解像度を選ぶと、JPEG2000形式で約30MBytes位のファイルがダウンロードできました。
　ちなみにライセンスは、下記の通りでパブリック・ドメインです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/545245f6-0564-7e8b-611e-2e22142f8b40.png)

#ジオリファレンスの手順
　ジオリファレンスを行う前に、qgis3.10ではジオリファレンサーというプラグインを登録する必要があります。
　その後、変換したい地図をジオリファレンサーで読み込み、変換の基準となるポイントを数カ所以上登録することとなります。
##基準となる地図の準備
　もし、基準となる地図が無い場合、先に登録しておきます。
　今回は[地理院地図](https://maps.gsi.go.jp/development/ichiran.html#std)を登録したプロジェクトを作っておいた。

##プラグインの登録
　qgisを起動し、【ラスタ(R)】を開いて「ジオリファレンサ」というメニューがあれば、ジオリファレンスを行える。
　無い場合、以下の手順でプラグインを登録する必要がある。
　まず、【プラグイン(P)】というメニューを開く。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/c976352d-c705-0efa-8110-b5c3d42c1343.png)
「プラグインの管理とインストール」というメニューを実行する。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/d1dfcd32-a39f-23de-b6fa-4cf51e5e1bef.png)

　検索のテキストボックスに「ジオ」と打ち込むと「GDALジオリファレンサ」が選択されてくるのでチェックを入れる。必要なファイルがセットアップされるとジオリファレンスができるようになる。

##基準点の登録
　【ラスタ(R)】から「ジオリファレンサ」を実行する。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/00790217-fa85-6988-9f75-2d7dd5e28d7a.png)
　ジオリファレンサのメニューから【ファイル(F)】を開き、「ラスタを開く】を実行すると、ファイル選択のダイアログが開くので、読み込ませるファイルを選択する。ジオリファレンサはjpeg2000も読み込めるので、今回の地図ファイルも処理できる。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/fdc74e74-1cf4-facd-fc26-98e090835e00.png)
　
　qgisのジオリファレンサでは、地図を変換するための基準点をGCPと読んでいる(Geo Control Point？)。
　GCPは二枚の地図で同一となる地点をしているわけで、作業効率や正確性を考えると、年月がたっても同じ位置であることが明確であり、かつ地図上で特定しやすいところを選ぶべきですね。
　年月がたっても位置が変わらないポイントは、以外とないものです。
　自然の地形の場合、年月が経ると災害などで変わってしまうことがあります。特に川の流域などは災害が無くても、水の流れで流路が変わることがよくりあますし、河川改修や水路建設などで変わることもあります。
　道路など人工物も位置が変わることがありますが、こちらはある程度、変わった時期などが特定できると思います。
　今回は変換対象が地図であるため、重ね合わせやすいことから、GCPの数はそんなに要らないだろうと見込み、鉄道の線路の分岐点などからGCPを選ぶことにしました。　
　
　まず、北西のGCPとして男鹿線と羽越本線の分岐点を選びます。
　地理院地図とジオリファレンサーの地図で同じあたりを表示しておくと作業がしやすいです。
　ジオリファレンサーの地図を動かすには、ツールバーから「手」のアイコンを選択します。次にGCPを設定するには、次のアイコンを選択します。![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/78c02602-6241-6ade-c2dd-71626edb3ab1.png)
　この状態でジオリファレンサー上の地図で、任意の位置をポイントすると、次のようなダイアログがポップアップするので、「地図キャンパスから選ぶ」のボタンをクリックします(座標がわかっている場合、座標を入力しても良いと思いますが・・・)。
　![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/d634784f-1e30-15c7-aca6-449d15b1d196.png)
　すると、qgisの地図が表示されるので、先ほど選択したポイントと対応するポイントを選択しGCPを登録します。
　これを数カ所分行います。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/e8d5d691-f1f1-9660-c001-f95de1201933.png)
　今回はGCPとして、男鹿線との分岐点と合わせ4カ所(残りは補陀寺、羽州街道と奥羽本線の交点、羽越本線の鉄橋)を選びました(地図上の赤点)。

##地図の変換
　ジオリファレンサーのメニューバーにある「実行」![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/ed604752-afff-d018-bbb7-122fe2a3302e.png)を押すと次のようなダイアログが開きます。
　![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/feb81d9e-64e0-e93a-663e-3fa1403048b6.png)
　今回は変換する地図自体が、きちんとした地図なのでシンプルな変換で良いと思われるため「線形」を選んでいます。
　GISのお約束として「変換先SRS」を正しくしておきます。秋田市の地図なので、EPSG:6678を指定しています。
　出力ラスタには、変換した画像ファイルの名前を入力します。
　完了時にqgisにロードするようにしておくと、自動的に重ねてくれます。
　今回気がついたのですが、「PDFマップの生成」という項目があるので、ジオPDFも作れるようです(これは後で検証予定)。
　ここまで設定したら変換実行です。
　再度、実行を押します。
　自動的に変換した地図が表示されるので、レイヤの透明度を調整して地理院地図と重ねてみます。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/347462/3025bf07-32bc-8ba6-06bd-b4359405b38f.png)
　微妙にずれてる気もしますが、実用上は問題なさそうです。
#今後の課題
##変換方法の研究
　変換方法が何種類か選べるので、その意味を調べる予定です。
##生成されるPDFの確認
　ジオPDFか確認します。