ネットにある地図では緯度経度を正確に作図するという用途には向いていません。国土地理院のデータから高精度に描画する地図を作成する方法を紹介します。他の緯度経度データと合わせてデータビジュアライゼーションで活用することを目的としています。また専用のGISソフトウェアではなくPHPとAdobe Illustratorでの処理を前提としています。
1 データの入手先
データは国土数値地理情報から入手します。国土地理情報のダウンロードサイト
ここにある行政区域のデータを使います。行政区域データは各都道府県毎にダウンロードできます。市区町村の名称と境界線のデータがXML形式で入っています。このXMLを読み取り境界線を画面に描画するためのプログラムをPHPで作成します。描画はSVG形式で書き出します。SVGは座標を記述したテキストデータの様な物でブラウザでの表示やIllustratorなどのグラフィックソフトでも開く事ができます。
2 データのダウンロードと保存
47都道府県のデータをダウンロードします。北から番号順になっています。ダウンロードしたデータを解凍するとフォルダが出来上がります。そのフォルダを開くと下図のようなファイルが入っています。
この中で拡張子がxmlのファイルが目的のデータが入ったファイルです。このファイルから座標を読み取り描画を行います。
3.ダウンロードしたXMLの仕組み
PHPではxmlを読み取るための仕組みがあります。ネットにも色々とサンプルがあります。そこらへんを参考にコードを書いてあげます。それほど難しくはありませんがダウンロードしたデータは名前空間付きのデータになっています。名前空間が付いているとちょっと手間が増えます。
まずはデータが入っているXMLをみてみます。つぎのようなタグの中にデータが入っています。
<?xml version="1.0" encoding="UTF-8"?> <ksj:Dataset gml:id= "N03Dataset" xmlns:ksj="http://nlftp.mlit.go.jp/ksj/schemas/ksj-app" xmlns:gml="http://www.opengis.net/gml/3.2" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://nlftp.mlit.go.jp/ksj/schemas/ksj-app KsjAppSchema-N03-v2_2.xsd"> <gml:description>国土数値情報 行政区域 インスタンス文書</gml:description> <!-- データ提供範囲 --> <gml:boundedBy> <gml:EnvelopeWithTimePeriod srsName="JGD2011 / (B, L)" frame="GC / JST"> <gml:lowerCorner>20.0 120.0</gml:lowerCorner> <gml:upperCorner>46.0 154.0</gml:upperCorner> <gml:beginPosition calendarEraName="西暦">1900</gml:beginPosition> <gml:endPosition indeterminatePosition="unknown"/> </gml:EnvelopeWithTimePeriod> </gml:boundedBy> <!-- 図形 --> <gml:Curve gml:id="cv0_0"> <gml:segments> <gml:LineStringSegment> <gml:posList> 39.73879062 140.05369904 39.73887229 140.05365484 39.73890323 140.05363629 39.73885796 140.05349189 39.73874256 140.05355737 39.73879062 140.05369904 </gml:posList> </gml:LineStringSegment> </gml:segments> </gml:Curve> <gml:Surface gml:id="sf0"> <gml:patches> <gml:PolygonPatch> <gml:exterior> <gml:Ring> <gml:curveMember xlink:href="#cv0_0"/> </gml:Ring> </gml:exterior> </gml:PolygonPatch> </gml:patches> </gml:Surface> ...中略 <!-- 属性 --> <ksj:AdministrativeBoundary gml:id="gy0"> <ksj:bounds xlink:href="#sf0"/> <ksj:prefectureName>秋田県</ksj:prefectureName> <ksj:subPrefectureName></ksj:subPrefectureName> <ksj:countyName></ksj:countyName> <ksj:cityName>秋田市</ksj:cityName> <ksj:administrativeAreaCode codeSpace="AdministrativeAreaCode.xml">05201</ksj:administrativeAreaCode> </ksj:AdministrativeBoundary> ...中略 </ksj:Dataset>
ざっくりと説明すると<ksj:Dataset>〜</ksj:Dataset>にデータが入っています。データは3つの種類に分かれています。
<gml:Curve gml:id=”cv0_0″>が座標データ
<gml:Surface gml:id=”sf0″>がサーフェスの説明
<ksj:AdministrativeBoundary gml:id=”gy0″>が行政区域名が入った部分。
例えば上記の例では<ksj:AdministrativeBoundary gml:id=”gy0″>を見ると秋田県、秋田市という事がわかります。そして<ksj:bounds xlink:href=”#sf0″/>に#sf0とあるので、サーフェスはsf0を見れば良いとわかります。そこで<gml:Surface gml:id=”sf0″>を見ると<gml:curveMember xlink:href=”#cv0_0″/>でcv0_0とあるので、座標（緯度、東経）はcv0_0を見ます。そこに<gml:posList>〜 </gml:posList>として
<gml:posList> 39.73879062 140.05369904 39.73887229 140.05365484 39.73890323 140.05363629 39.73885796 140.05349189 39.73874256 140.05355737 39.73879062 140.05369904 </gml:posList>
が入っているとなります。緯度と東経が文字列で並んでいます。この部分は改行コードで区切ってあり緯度と東経の区切りは半角スペース区切りです。ちなみにこの場所はどこなのかGoogle Mapで確認すると
ここになってます。google map上では表示されない小さな島（岩？）のようです。数値情報地図はここまで細かいところまで載っています。
4.座標データの取り出し
さてここからPHPでデータを取り出す方法に入ります。まずはダウンロードしたXMLファイルから座標データを読み出すためにファイルを開きます。
$dom = new DOMDocument('1.0', 'UTF-8');
$dom->preserveWhiteSpace = false;
$dom->formatOutput = true;
$dom->load($loadXML);
でXMLファイルを読み込みます。$loadXMLはxmlファイルのファイル名（ファイルパス）です。これで変数$domにXMLの全体が入ります。
つぎにksjの部分を取り出します。ここらへんの方法はサンプルが色々とサイトにも紹介されていますが名前空間付きだとなんか処理が違う様で混乱しました。とりあえず＊で指定して中身を取り出します。
$root = $dom->getElementsByTagName("*")->item(0)
次に上記の１行で変数$rootにksjの中身を代入します。ここから名前空間付きのXMLデータを処理します。次のforeachループで名前空間Curveのデータをとりだします。座標データが欲しいだけなので他のデータは使いません。
foreach( $root->getElementsByTagNameNS("*","Curve") as $itemGML ){
}
このループでCurveごとのデータが取り出せます。ループの中では
$dataArray = explode("\n",$itemGML->nodeValue);//改行コードで配列にわける
として改行コードで配列にして緯度と東経のデータを取り出します。配列$dataArrayに文字列で緯度東経が入ります。そこから緯度東経を納めるための配列$coordinateに値を入れます。その時に前後にスペースやタブがあると邪魔なのでtrimで綺麗な文字列にしておきます。
$j=0;
foreach ($dataArray as $item) {//各要素ごとに繰り返す
$item = trim($item);
if($item != ""){//空白文字列の時には処理をしない。
$coordinate[$i][$j]=$item;
$j++;
}
}
$i++;
配列$coordinateは２次元の配列ですが、最初の添え字で島もしくは市区町村の境界線、２番目の添え字が座標データの並びになります。
5.緯度経度の取り出し
$coordinateには文字列で緯度経度が入っています。これをスペース区切りで取り出すと緯度と経度（東経）になります。
$latilongArray = explode(" ",$coordinate[$i][$j]);//空白で区切る
$latitude1String = $latilongArray[0];//各要素を取り出す
$longitude1String = $latilongArray[1];
そしてその文字列を数値（float）に変換します。数値への変換は投影の計算のところでやっています。
6.地図の投影
平面の地図を作るためにメルカトル図法で投影します。投影するまえに作図の範囲を40cm x 40cm。仮のスクリーンサイズを4000 x 4000pxとしています。ブラウザで表示するためにはpxが重要ですが今回はIllustratorなどで扱うので適当にしています。縮尺を求めるために日本全体を描画するために日本の最大サイズ（緯度経度）を求めておきます。
$screenWidth = 4000;
$screenHeight = 4000;
$north_max = 46.0;//最北端：北緯45度33分28秒 東経148度45分14秒 46.0
$south_max = 20.0;//最南端：北緯20度25分30.6585秒 東経136度4分11.1766秒 20.0
$east_max = 154.0;//最東端：北緯24度17分12秒 東経153度58分50秒 154.0
$west_max = 122.0;// 最西端：北緯24度26分58秒 東経122度56分01秒 122.0
$y_max = log(tan(pi()/4.0 + deg2rad($north_max)/2.0));
$y_min = log(tan(pi()/4.0 + deg2rad($south_max)/2.0));
$x_max = deg2rad($east_max);//横方向
$x_min= deg2rad($west_max);
$cLo = ($x_max - $x_min)/2 + $x_min;//描画中央になる東経
//描画するための縮尺を求める
$dx = $x_max - $x_min;
$dy = $y_max - $y_min;
$mX = $screenWidth/$dx;
$mY = $screenHeight/$dy;
//描画するための縮尺で縦方向と横方向で大きい方を使う。
if ($mX < $mY){
$multiple = $mX;
} else {
$multiple = $mY;
}
縦方向、横方向の縮尺を求めてその大きい方を描画の時の倍率$multipleとして使います。
7. メルカトル図法の計算
上記の地図の投影でもすでに出てきていますが緯度経度から座標に計算する部分は次の様な計算になります。
$radian = deg2rad(floatval($latitude1String));
$x = deg2rad(floatval($longitude1String)) - $cLo;//横方向
$y = log(tan(pi()/4.0 + $radian/2.0));//縦方向
$x = $multiple * $x + $screenWidth/2;
$y = $multiple * (1.0 - $y) - ($multiple * (1.0-$yOffset));
ここで$xと$yが描画する座標になります。
8.SVGファイルの作成
SVGファイルは中身はテキストファイルなので座標をsvgのルールで書き出して作ります。polylineで座標を書き出してます。
$contents='<?xml version="1.0" standalone="no"?>';
$contents .="\n";
$contents .='<!DOCTYPE svg PUBLIC "-//W3C//DTD SVG 1.1//EN" "http://www.w3.org/Graphics/SVG/1.1/DTD/svg11.dtd">';
$contents .="\n";
$contents .='<svg width="40cm" height="40cm" viewBox="0 0 4000 4000" xmlns="http://www.w3.org/2000/svg" version="1.1">';
$contents .="\n";
$contents .='<desc>都道府県の境界線データ</desc>';
$contents .="\n";
$contents .='<polyline fill="none" stroke="black" stroke-width="1" points="0,0 4000,0 4000,4000,0,4000" />';//図形の書き出し
//境界線データの書き出し polylineでデータで書き出す。
$imax = count($coordinate);
//print("imax=$imax\n");
for ($i=0 ;$i<$imax ; $i++){
$jmax = count($coordinate[$i]);
//print("jmax=$jmax\n");
$svgData = "";
for ($j= 0;$j<$jmax ; $j++){
//メルカトル図法による計算。
$latilongArray = explode(" ",$coordinate[$i][$j]);//空白で区切る
//print("strring=" . $coordinate[$i][$j] . "\n");
$latitude1String = $latilongArray[0];//各要素を取り出す
$longitude1String = $latilongArray[1];
$radian = deg2rad(floatval($latitude1String));
$x = deg2rad(floatval($longitude1String)) - $cLo;//横方向
$y = log(tan(pi()/4.0 + $radian/2.0));//縦方向
$x = $multiple * $x + $screenWidth/2;
$y = $multiple * (1.0 - $y) - ($multiple * (1.0-$yOffset));
$svgData .= $x . "," . $y . " ";
}
$contents .='<polyline fill="none" stroke="blue" stroke-width="1" points="';//図形の書き出し
$contents .= $svgData . '" />' . "\n";
}
$contents .="</svg>";
file_put_contents($saveSVG, $contents);
9.Finder上でのファイルの配置
phpのプログラムとXMLの配置は次の様にします。
フォルダにダウンロードしたフォルダが入っています。
下の階層はこんな感じ
10.コードとスクリプトの起動
プログラム全文は次の通りです。
plotPrefData.php
<?php
ini_set('memory_limit', '256M');
mb_language("Japanese");//文字コードの設定
mb_internal_encoding("UTF-8");
$loadXML = "./国土数値情報行政区域/N03-170101_01_GML/N03-17_01_170101.xml";
$saveSVG= "hokkaido.svg";//SVGファイルのファイル名
loadXMLsaveSVGXML($loadXML,$saveSVG);
function loadXMLsaveSVGXML($loadXML,$saveSVG){
print("XML file loading...");
$dom = new DOMDocument('1.0', 'UTF-8');
$dom->preserveWhiteSpace = false;
$dom->formatOutput = true;
$dom->load($loadXML);
$root = $dom->getElementsByTagName("*")->item(0);//最初の要素の中身を取り出す
$i = 0;
foreach( $root->getElementsByTagNameNS("*","Curve") as $itemGML ){//これでgml:Curveのところが取れるらしい。
$dataArray = explode("\n",$itemGML->nodeValue);//改行コードで配列にわける
$j=0;
foreach ($dataArray as $item) {//各要素ごとに繰り返す
$item = trim($item);
if($item != ""){//空白文字列の時には処理をしない。
$coordinate[$i][$j]=$item;
$j++;
}
}
$i++;
}
//=========================SVGで書き出す
//スクリーンサイズ
$screenWidth = 4000;
$screenHeight = 4000;
$north_max = 46.0;//最北端：北緯45度33分28秒 東経148度45分14秒 46.0
$south_max = 20.0;//最南端：北緯20度25分30.6585秒 東経136度4分11.1766秒 20.0
$east_max = 154.0;//最東端：北緯24度17分12秒 東経153度58分50秒 154.0
$west_max = 122.0;// 最西端：北緯24度26分58秒 東経122度56分01秒 122.0
$y_max = log(tan(pi()/4.0 + deg2rad($north_max)/2.0));
$y_min = log(tan(pi()/4.0 + deg2rad($south_max)/2.0));
$x_max = deg2rad($east_max);//横方向
$x_min= deg2rad($west_max);
$cLo = ($x_max - $x_min)/2 + $x_min;//描画中央になる東経
//描画するための縮尺を求める
$dx = $x_max - $x_min;
$dy = $y_max - $y_min;
$mX = $screenWidth/$dx;
$mY = $screenHeight/$dy;
//描画するための縮尺で縦方向と横方向で大きい方を使う。
if ($mX < $mY){
$multiple = $mX;
} else {
$multiple = $mY;
}
//y軸の描画用の定数
$yOffset = $y_max;
//var_dump($coordinate);
//==========================SVGファイルを書き出す
$contents='<?xml version="1.0" standalone="no"?>';
$contents .="\n";
$contents .='<!DOCTYPE svg PUBLIC "-//W3C//DTD SVG 1.1//EN" "http://www.w3.org/Graphics/SVG/1.1/DTD/svg11.dtd">';
$contents .="\n";
$contents .='<svg width="40cm" height="40cm" viewBox="0 0 4000 4000" xmlns="http://www.w3.org/2000/svg" version="1.1">';
$contents .="\n";
$contents .='<desc>都道府県の境界線データ</desc>';
$contents .="\n";
$contents .='<polyline fill="none" stroke="black" stroke-width="1" points="0,0 4000,0 4000,4000,0,4000" />';//図形の書き出し
//境界線データの書き出し polylineでデータで書き出す。
$imax = count($coordinate);
//print("imax=$imax\n");
for ($i=0 ;$i<$imax ; $i++){
$jmax = count($coordinate[$i]);
//print("jmax=$jmax\n");
$svgData = "";
for ($j= 0;$j<$jmax ; $j++){
//メルカトル図法による計算。
$latilongArray = explode(" ",$coordinate[$i][$j]);//空白で区切る
//print("strring=" . $coordinate[$i][$j] . "\n");
$latitude1String = $latilongArray[0];//各要素を取り出す
$longitude1String = $latilongArray[1];
$radian = deg2rad(floatval($latitude1String));
$x = deg2rad(floatval($longitude1String)) - $cLo;//横方向
$y = log(tan(pi()/4.0 + $radian/2.0));//縦方向
$x = $multiple * $x + $screenWidth/2;
$y = $multiple * (1.0 - $y) - ($multiple * (1.0-$yOffset));
$svgData .= $x . "," . $y . " ";
}
$contents .='<polyline fill="none" stroke="blue" stroke-width="1" points="';//図形の書き出し
$contents .= $svgData . '" />' . "\n";
}
$contents .="</svg>";
file_put_contents($saveSVG, $contents);
}
このコードをplotPrefData.phpという名前で保存しておきます。
phpの起動にはターミナルが必要になります。
下図では省略していますがcdコマンドで作業中のフォルダ（ディレクトリ）に移動します。
cd でスペースを入れて目的のディレクトリのパスを入力してreturnキー。ディレクトリのパスを入力するのが面倒な時にはFinderからフォルダをターミナルにドラッグ＆ドロップでファイルパスに入力が出来ます。
php ./plotPrefData.php とターミナルに入力してリターンキーで作成したphpのスクリプトが動きます。
11.作成したSVGファイル
無事に処理が終了すれば同じフォルダにhokkaido.svgという画像ファイルがつくれれているはずです。この画像ファイルは41.6MBあり開くにはかなり時間がかかります。ブラウザで開いてしばらく反応がありませんが、じっとまっててください。
Illustratorではデータ数がありすぎて開くまで20分くらいかかりました。また、すべての行政区域の境界線が入っているので用途によっては使いにくいものになります。また下図では択捉島の一部が欠けています。これはIllustratorのパスの制限によるものです。ベジェ曲線の点が30,000点以上になると図形が表示されないようです。下図はIllustratorのアウトライン表示。
非常に細かい部分まで入っています。下図はえりも岬の先端部分。
とりあえずダウンロードしたデータを描画してSVGファイルにすることはできました。使いやすくするために以下の事を行います。
・県境以外の行政区域の区切り線を削除する。
・描画時にデータサイズを調整できるようにする。
・3万点以上のデータへの対応。
続きは次回で。