stackprobe7s_memo

何処にも披露する見込みの無いものを書き落とす場所

奥沢は目黒区だと勘違いしていた話 (C# + 国土地理院の地図データ)

つい最近(と言っても数年前)まで東急目黒線奥沢駅は目黒区か大田区だと思っていたのですが、ひょんなきっかけで世田谷区だと知りました。
まぁ勝手に思い込んでいただけだったので多少驚いて素直に納得したのですが、むしろ何十年もそのことに気づかなかったという方が不思議というか驚きだったわけです。
で、何でそんな勘違いをしていたかという小ネタをやるために国土地理院の地図データをちょっとだけ触ってみようかと思います。



国土地理院が全国の地図データを公開していて、利用規約や制限などはありますが割と詳細な地図データが無料で手に入ります。
 
fgd.gsi.go.jp
 
ダウンロード ⇒ 基盤地図情報 基本項目 ⇒ ファイル選択へ
と進んで行き、今回必要なメッシュ 533935 のみダウンロードします。

ダウンロードした FG-GML-533935-ALL-20221001.zip を展開すると

FG-GML-533935-AdmArea-20221001-0001.xml
FG-GML-533935-AdmBdry-20221001-0001.xml
FG-GML-533935-AdmPt-20221001-0001.xml
FG-GML-533935-BldA-20221001-0001.xml
FG-GML-533935-BldA-20221001-0002.xml
FG-GML-533935-BldA-20221001-0003.xml
FG-GML-533935-BldA-20221001-0004.xml
FG-GML-533935-BldL-20221001-0001.xml
FG-GML-533935-Cntr-20221001-0001.xml
FG-GML-533935-CommBdry-20221001-0001.xml
FG-GML-533935-CommPt-20221001-0001.xml
FG-GML-533935-Cstline-20221001-0001.xml
FG-GML-533935-ElevPt-20221001-0001.xml
FG-GML-533935-GCP-20221001-0001.xml
FG-GML-533935-RailCL-20221001-0001.xml
FG-GML-533935-RdCompt-20221001-0001.xml
FG-GML-533935-RdEdg-20221001-0001.xml
FG-GML-533935-SBAPt-20221001-0001.xml
FG-GML-533935-SBBdry-20221001-0001.xml
FG-GML-533935-WA-20221001-0001.xml
FG-GML-533935-WL-20221001-0001.xml
FG-GML-533935-WStrA-20221001-0001.xml
FG-GML-533935-WStrL-20221001-0001.xml
fmdid22-1001.xml

ファイル命名規則https://www.gsi.go.jp/common/000093391.pdf

といろいろ入っていますが、ダウンロードページにある FGDV というビューアで確認した上で今回必要なファイルを抜き出すと...

ファイル 内容
FG-GML-533935-AdmArea-20221001-0001.xml 行政区画
FG-GML-533935-BldL-20221001-0001.xml 建築物の外周線
FG-GML-533935-RailCL-20221001-0001.xml (鉄道の)軌道の中心線
FG-GML-533935-RdEdg-20221001-0001.xml 道路

...となるので、これらのファイルだけ見ていきます。
 
先ずそこそこファイルがでかいので私のショボいPCで読み込めるのか確認します。

public void Test01()
{
	LoadXMLTest(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-AdmArea-20221001-0001.xml");
	LoadXMLTest(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-BldL-20221001-0001.xml");
	LoadXMLTest(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RailCL-20221001-0001.xml");
	LoadXMLTest(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RdEdg-20221001-0001.xml");
}

private List<XMLNode> LoadedXMLs = new List<XMLNode>();

private void LoadXMLTest(string file)
{
	Console.WriteLine("Load " + file);
	DateTime stTm = DateTime.Now;
	LoadedXMLs.Add(XMLNode.LoadFromFile(file));
	DateTime edTm = DateTime.Now;
	Console.WriteLine("OK " + (edTm - stTm).TotalSeconds.ToString("F3") + "s");
}

Output

Load C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-AdmArea-20221001-0001.xml
OK 0.016s
Load C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-BldL-20221001-0001.xml
OK 9.249s
Load C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RailCL-20221001-0001.xml
OK 0.094s
Load C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RdEdg-20221001-0001.xml
OK 2.085s

読めました。
建築物の外周線はでかいだけあって重い。PCがショボいってのもあるけど・・・。
 
続いてフォーマットの確認。
昔この辺のファイルを触った記憶があるのだけど、もう詳細は忘れました。
どこかに書式に関する説明があるだろうしそれを読むべきだろうけど、ぱっと見それほど複雑ではなさそうなので強引にデータを読み解くことにします。

public void Test02()
{
	PrintXMLPaths(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-AdmArea-20221001-0001.xml");
	PrintXMLPaths(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-BldL-20221001-0001.xml");
	PrintXMLPaths(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RailCL-20221001-0001.xml");
	PrintXMLPaths(@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RdEdg-20221001-0001.xml");
}

private void PrintXMLPaths(string file)
{
	HashSet<string> xmlPaths = new HashSet<string>();
	XMLNode root = XMLNode.LoadFromFile(file);
	root.Search((xmlPath, node) => xmlPaths.Add(xmlPath));

	Console.WriteLine(file);
	Console.WriteLine();

	List<string> lines = new List<string>();

	lines.Add("XMLパス*最小の個数*最大の個数");
	lines.Add("----");

	foreach (string xmlPath in xmlPaths.OrderBy(SCommon.Comp))
		lines.Add(xmlPath
			+ "*" + GetNodeCount(root, xmlPath, int.MaxValue, Math.Min)
			+ "*" + GetNodeCount(root, xmlPath, int.MinValue, Math.Max));

	Common.ToConsoleTable(lines);

	foreach (string line in lines)
		Console.WriteLine(line);

	Console.WriteLine();
}

private int GetNodeCount(XMLNode root, string xmlPath, int count, Func<int, int, int> chooser)
{
	int p = xmlPath.LastIndexOf('/');

	if (p == -1) // ? ルートTag
		return 1;

	string parentXmlPath = xmlPath.Substring(0, p);
	string name = xmlPath.Substring(p + 1);

	root.Search((xp, node) =>
	{
		if (xp == parentXmlPath)
			count = chooser(count, node.Children.Where(n => n.Name == name).Count());

		return true;
	});

	return count;
}

Output

C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-AdmArea-20221001-0001.xml

XMLパス                                                                                                               最小の個数  最大の個数
--------------------------------------------------------------------------------------------------------------------------------------------
Dataset                                                                                                               1           1
Dataset/AdmArea                                                                                                       14          14
Dataset/AdmArea/admCode                                                                                               1           1
Dataset/AdmArea/area                                                                                                  1           1
Dataset/AdmArea/area/Surface                                                                                          1           1
Dataset/AdmArea/area/Surface/id                                                                                       1           1
Dataset/AdmArea/area/Surface/patches                                                                                  1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch                                                                     1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior                                                            1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring                                                       1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember                                           1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve                                     1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve/id                                  1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve/segments                            1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve/segments/LineStringSegment          1           1
Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve/segments/LineStringSegment/posList  1           1
Dataset/AdmArea/area/Surface/srsName                                                                                  1           1
Dataset/AdmArea/devDate                                                                                               1           1
Dataset/AdmArea/devDate/id                                                                                            1           1
Dataset/AdmArea/devDate/timePosition                                                                                  1           1
Dataset/AdmArea/fid                                                                                                   1           1
Dataset/AdmArea/id                                                                                                    1           1
Dataset/AdmArea/lfSpanFr                                                                                              1           1
Dataset/AdmArea/lfSpanFr/id                                                                                           1           1
Dataset/AdmArea/lfSpanFr/timePosition                                                                                 1           1
Dataset/AdmArea/name                                                                                                  1           1
Dataset/AdmArea/orgGILvl                                                                                              1           1
Dataset/AdmArea/type                                                                                                  1           1
Dataset/AdmArea/vis                                                                                                   1           1
Dataset/description                                                                                                   1           1
Dataset/gml                                                                                                           1           1
Dataset/id                                                                                                            1           1
Dataset/name                                                                                                          1           1
Dataset/schemaLocation                                                                                                1           1
Dataset/xlink                                                                                                         1           1
Dataset/xmlns                                                                                                         1           1
Dataset/xsi                                                                                                           1           1

C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-BldL-20221001-0001.xml

XMLパス                                                    最小の個数  最大の個数
---------------------------------------------------------------------------------
Dataset                                                    1           1
Dataset/BldL                                               384932      384932
Dataset/BldL/devDate                                       1           1
Dataset/BldL/devDate/id                                    1           1
Dataset/BldL/devDate/timePosition                          1           1
Dataset/BldL/fid                                           1           1
Dataset/BldL/id                                            1           1
Dataset/BldL/lfSpanFr                                      1           1
Dataset/BldL/lfSpanFr/id                                   1           1
Dataset/BldL/lfSpanFr/timePosition                         1           1
Dataset/BldL/loc                                           1           1
Dataset/BldL/loc/Curve                                     1           1
Dataset/BldL/loc/Curve/id                                  1           1
Dataset/BldL/loc/Curve/segments                            1           1
Dataset/BldL/loc/Curve/segments/LineStringSegment          1           1
Dataset/BldL/loc/Curve/segments/LineStringSegment/posList  1           1
Dataset/BldL/loc/Curve/srsName                             1           1
Dataset/BldL/orgGILvl                                      1           1
Dataset/BldL/type                                          1           1
Dataset/description                                        1           1
Dataset/gml                                                1           1
Dataset/id                                                 1           1
Dataset/name                                               1           1
Dataset/schemaLocation                                     1           1
Dataset/xlink                                              1           1
Dataset/xmlns                                              1           1
Dataset/xsi                                                1           1

C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RailCL-20221001-0001.xml

XMLパス                                                      最小の個数  最大の個数
-----------------------------------------------------------------------------------
Dataset                                                      1           1
Dataset/RailCL                                               3910        3910
Dataset/RailCL/devDate                                       1           1
Dataset/RailCL/devDate/id                                    1           1
Dataset/RailCL/devDate/timePosition                          1           1
Dataset/RailCL/fid                                           1           1
Dataset/RailCL/id                                            1           1
Dataset/RailCL/lfSpanFr                                      1           1
Dataset/RailCL/lfSpanFr/id                                   1           1
Dataset/RailCL/lfSpanFr/timePosition                         1           1
Dataset/RailCL/loc                                           1           1
Dataset/RailCL/loc/Curve                                     1           1
Dataset/RailCL/loc/Curve/id                                  1           1
Dataset/RailCL/loc/Curve/segments                            1           1
Dataset/RailCL/loc/Curve/segments/LineStringSegment          1           1
Dataset/RailCL/loc/Curve/segments/LineStringSegment/posList  1           1
Dataset/RailCL/loc/Curve/srsName                             1           1
Dataset/RailCL/orgGILvl                                      1           1
Dataset/RailCL/type                                          1           1
Dataset/RailCL/vis                                           1           1
Dataset/description                                          1           1
Dataset/gml                                                  1           1
Dataset/id                                                   1           1
Dataset/name                                                 1           1
Dataset/schemaLocation                                       1           1
Dataset/xlink                                                1           1
Dataset/xmlns                                                1           1
Dataset/xsi                                                  1           1

C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RdEdg-20221001-0001.xml

XMLパス                                                     最小の個数  最大の個数
----------------------------------------------------------------------------------
Dataset                                                     1           1
Dataset/RdEdg                                               46633       46633
Dataset/RdEdg/admOffice                                     1           1
Dataset/RdEdg/devDate                                       1           1
Dataset/RdEdg/devDate/id                                    1           1
Dataset/RdEdg/devDate/timePosition                          1           1
Dataset/RdEdg/fid                                           1           1
Dataset/RdEdg/id                                            1           1
Dataset/RdEdg/lfSpanFr                                      1           1
Dataset/RdEdg/lfSpanFr/id                                   1           1
Dataset/RdEdg/lfSpanFr/timePosition                         1           1
Dataset/RdEdg/loc                                           1           1
Dataset/RdEdg/loc/Curve                                     1           1
Dataset/RdEdg/loc/Curve/id                                  1           1
Dataset/RdEdg/loc/Curve/segments                            1           1
Dataset/RdEdg/loc/Curve/segments/LineStringSegment          1           1
Dataset/RdEdg/loc/Curve/segments/LineStringSegment/posList  1           1
Dataset/RdEdg/loc/Curve/srsName                             1           1
Dataset/RdEdg/orgGILvl                                      1           1
Dataset/RdEdg/type                                          1           1
Dataset/RdEdg/vis                                           1           1
Dataset/description                                         1           1
Dataset/gml                                                 1           1
Dataset/id                                                  1           1
Dataset/name                                                1           1
Dataset/schemaLocation                                      1           1
Dataset/xlink                                               1           1
Dataset/xmlns                                               1           1
Dataset/xsi                                                 1           1

ぱっと見の予想どおりシンプルな構造でした。
posList というタグにポリゴン (緯度・軽度の組み合わせのリスト) が入っているようなので、とりあえず抽出してみます。
行政区画には行政区名 (/Dataset/AdmArea/name) が入っているけど、建物名・路線名は入っていない模様。
なので線だけで何とかしてみます。

public void Test03()
{
	ExportPosList(
		@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-AdmArea-20221001-0001.xml",
		@"C:\temp\Area.txt",
		"Dataset/AdmArea/area/Surface/patches/PolygonPatch/exterior/Ring/curveMember/Curve/segments/LineStringSegment/posList"
		);
	ExportPosList(
		@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-BldL-20221001-0001.xml",
		@"C:\temp\Building.txt",
		"Dataset/BldL/loc/Curve/segments/LineStringSegment/posList"
		);
	ExportPosList(
		@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RailCL-20221001-0001.xml",
		@"C:\temp\Rail.txt",
		"Dataset/RailCL/loc/Curve/segments/LineStringSegment/posList"
		);
	ExportPosList(
		@"C:\temp\FG-GML-533935-ALL-20221001\FG-GML-533935-RdEdg-20221001-0001.xml",
		@"C:\temp\Road.txt",
		"Dataset/RdEdg/loc/Curve/segments/LineStringSegment/posList"
		);
}

private void ExportPosList(string mapFile, string polyFile, string polyXmlPath)
{
	XMLNode root = XMLNode.LoadFromFile(mapFile);

	using (StreamWriter writer = new StreamWriter(polyFile, false, Encoding.ASCII))
	{
		root.Search((xmlPath, node) =>
		{
			if (xmlPath == polyXmlPath)
			{
				writer.WriteLine(node.Value);
				writer.WriteLine(); // 空行 -- データの区切りとして
			}
			return true;
		});
	}
}

C:\temp\Area.txt , C:\temp\Building.txt , C:\temp\Rail.txt , C:\temp\Road.txt が作成されます。
中身は↓こんな感じです。

Building.txt
35.584445222 139.628468000
35.584447750 139.628575778
35.584368250 139.628578528
35.584365722 139.628470833
35.584445222 139.628468000

35.588461500 139.628454222
35.588464944 139.628520750
35.588304417 139.628533194
35.588301444 139.628475944
35.588307222 139.628475472
35.588306750 139.628466222
35.588461500 139.628454222

35.585365083 139.628463556
35.585365944 139.628491694
35.585374139 139.628491222
35.585375750 139.628545833
35.585312750 139.628548528
35.585310278 139.628466000
35.585365083 139.628463556
...

 
続いてこのデータの表示(画像ファイルへの変換)を試みます。
表示範囲としては以下の駅が入っていればよい想定です。(緯度経度は Wikipedia から)

場所 緯度 経度
九品仏駅 35.605417 139.661000
自由が丘駅 35.607500 139.668889
緑ヶ丘駅 35.606389 139.679361
大岡山駅 35.607500 139.685639
田園調布駅 35.596889 139.667306
奥沢駅 35.603833 139.672306

緯度経度を距離に変換することを考えます。

地球は楕円体なのですが、Wikipedia によると
赤道面での直径 12756.274 km (半径 6378.137 km)
極半径 6356.752 km
と、ほとんど変わらないので平均をとって半径 6367.4445 km の球として考えます。

表示範囲はだいたい北緯 35.6 度付近なので、経度については北緯 35.6 度における東西方向一周の長さから換算します。

double ER = 6367444.5; // 地球の半径(meter)
double LAT = 35.6; // 基準とする北緯(degree)

double latToMeter = ER * Math.PI / 180.0;
double lonToMeter = ER * Math.Cos(LAT * Math.PI / 180.0) * Math.PI / 180.0;

Console.WriteLine("緯度の 1 度を " + latToMeter.ToString("F9") + " メートルとする。");
Console.WriteLine("経度の 1 度を " + lonToMeter.ToString("F9") + " メートルとする。");

Output

緯度の 1 度を 111132.871463004 メートルとする。
経度の 1 度を 90362.222363910 メートルとする。

 
以上を踏まえて画像への変換を行います。

public void Test05()
{
	// 表示範囲(緯度経度)
	double NORTH_LAT = 35.607500; // 自由が丘駅
	double SOUTH_LAT = 35.596889; // 田園調布駅
	double WEST_LON = 139.661000; // 九品仏駅
	double EAST_LON = 139.685639; // 大岡山駅

	// 表示範囲のマージン(meter)
	double MARGIN = 300.0;

	// 出力画像の幅(pixel)
	int IMAGE_WIDTH = 1600;

	OutputImage(NORTH_LAT, SOUTH_LAT, WEST_LON, EAST_LON, MARGIN, IMAGE_WIDTH);
}

private void OutputImage(double nLat, double sLat, double wLon, double eLon, double margin, int imageWidth)
{
	double LAT_TO_METER = 111132.871463004;
	double LON_TO_METER = 90362.222363910;

	// マージン適用
	{
		nLat += margin / LAT_TO_METER;
		sLat -= margin / LAT_TO_METER;
		wLon -= margin / LON_TO_METER;
		eLon += margin / LON_TO_METER;
	}

	double xLon = eLon - wLon;
	double yLat = nLat - sLat;

	double xMeter = xLon * LON_TO_METER;
	double yMeter = yLat * LAT_TO_METER;

	double meterToPixel = imageWidth / xMeter;

	double latToPixel = LAT_TO_METER * meterToPixel;
	double lonToPixel = LON_TO_METER * meterToPixel;

	int w = imageWidth;
	int h = SCommon.ToInt(yMeter * meterToPixel);

	Bitmap image = new Bitmap(w, h);

	using (Graphics g = Graphics.FromImage(image))
	{
		g.FillRectangle(new SolidBrush(Color.White), 0, 0, w, h);

		Action<string, Pen> drawPolyFromFile = (polyFile, pen) =>
		{
			string[] polyLines = File.ReadAllLines(polyFile, Encoding.ASCII);
			int polyLineIndex = 0;

			while (polyLineIndex < polyLines.Length)
			{
				List<double[]> poly = new List<double[]>();

				while (polyLines[polyLineIndex] != "")
					poly.Add(polyLines[polyLineIndex++].Split(' ').Select(v => double.Parse(v)).ToArray());

				polyLineIndex++;

				foreach (double[] p in poly)
				{
					p[0] = (nLat - p[0]) * latToPixel; // Lat to pixel
					p[1] = (p[1] - wLon) * lonToPixel; // Lon to pixel
				}

				for (int polyIndex = 0; polyIndex + 1 < poly.Count; polyIndex++)
				{
					double[] p = poly[polyIndex + 0];
					double[] q = poly[polyIndex + 1];

					if (
						(0.0 <= p[0] && p[0] < h) ||
						(0.0 <= p[1] && p[1] < w) ||
						(0.0 <= q[0] && q[0] < h) ||
						(0.0 <= q[1] && q[1] < w)
						)
					{
						g.DrawLine(pen, (float)p[1], (float)p[0], (float)q[1], (float)q[0]);
					}
				}
			}
		};

		drawPolyFromFile(
			@"C:\temp\Building.txt",
			new Pen(new SolidBrush(Color.FromArgb(128, 255, 128)), 1.0F)
			);
		drawPolyFromFile(
			@"C:\temp\Road.txt",
			new Pen(new SolidBrush(Color.FromArgb(128, 128, 128)), 1.0F)
			);
		drawPolyFromFile(
			@"C:\temp\Rail.txt",
			new Pen(new SolidBrush(Color.FromArgb(0, 0, 0)), 1.0F)
			);
		drawPolyFromFile(
			@"C:\temp\Area.txt",
			new Pen(new SolidBrush(Color.FromArgb(64, 255, 0, 0)), 5.0F)
			);
	}
	image.Save(SCommon.NextOutputPath() + ".png");
}

Output

奥沢付近


で、何がやりたいかと云うと...
昔は Google map なんて便利なものがなく地図を眺める趣味もなかったので、古い人間の私は駅の住所から周辺の住所や行政区域をなんとなく把握していました。
東横線大井町線はよく使っていたので各駅がどの区にあるか把握していたのですが、目黒線(昔は目蒲線)は滅多に使わなかったので駅名すら曖昧な状態。
つまり奥沢はよくしらないけど、その周辺の駅はよくしっているという状態だったわけです。

src = https://github.com/soleil-taruto/StoreH/blob/main/Hatena/a20230108/Hatena20230108/Claes20200001/Claes20200001/Tests/Test0002.cs

奥沢付近

...のように駅があって、それぞれのどの区にあるかというと...

src = https://github.com/soleil-taruto/StoreH/blob/main/Hatena/a20230108/Hatena20230108/Claes20200001/Claes20200001/Tests/Test0003.cs

奥沢付近

...こうしてみると奥沢って目黒区か大田区かなと思っちゃうじゃないですか!!!
というだけの話。
念のため、奥沢は世田谷区です。

src = https://github.com/soleil-taruto/StoreH/blob/main/Hatena/a20230108/Hatena20230108/Claes20200001/Claes20200001/Tests/Test0004.cs

奥沢付近

区境を入れるとこんな感じ。


ここで使ってる src とか
https://github.com/soleil-taruto/StoreH/tree/main/Hatena/a20230108/Hatena20230108/Claes20200001/Claes20200001