GISデータの基本

# この記事について

仕事で GIS データを扱うことになったので、基本から学ぶことにしました。
以下の講座を受けながら、最初に出てきたジオメトリオブジェクトの基本をメモ。

https://www.udemy.com/course/python_gis/

# GISデータの基本

## ジオメトリオブジェクト

Point(点)、LineString(線)、Polygon(面)の3つのオブジェクトタイプが基本となる。
いずれも項目名として `type`, `coordinates` を持つ。

`type` にはそれぞれのオブジェクトタイプの文字列が入り、それに応じて `coordinates` の持ち方が変わる。

### Point

平面上の点をあらわす。
ある施設の場所など。


```json: Pointタイプのジオメトリオブジェクト
{"type": "Point",
  "coordinates": [30, 10]
}
```

`coordinates` は、緯度と経度をとる配列。

#### shapely を使った　Point オブジェクトの作成と描画

```python: Point オブジェクトの作成
from shapely.geometry import Point
point = Point(0.0, 1.0)
```

上記の `point` を表示すると、以下のような点が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/ad1cc012-7478-38bb-a053-e4c5946573e8.png)

#### shapely を使った　Point オブジェクトの型やデータ内容の表示

```python: Point オブジェクトの表示
print(point)     # POINT (0 1) のように表示される
point.geom_type  # 'Point' と表示される
```

### LineString

平面上の線をあらわす。
道路など。

```json: LineStringタイプのジオメトリオブジェクト
{"type": "LineString",
  "coordinates": [[30, 10], [10, 30], [40,40]]
}
```

`coordinates` の中身は、緯度と経度をとる配列（Pointと同じ）が複数入った配列。

#### shapely を使った　LineString オブジェクトの作成と描画

```python: LineString オブジェクトの作成
from shapely.geometry import LineString
line = LineString([(0,0),(1,1),(3,0),(5,2)])
```

上記の `line` を表示すると、以下のような線が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/1c67f762-f5d4-2a68-394d-5094a51c8aad.png)


#### shapely を使った　LineString オブジェクトの型やデータ内容の表示

```python: LineString オブジェクトの表示
print(line)     # LINESTRING (0 0, 1 1, 3 0, 5 2) のように表示される
line.geom_type  # 'LineString' と表示される
```

#### LineString オブジェクト(線)の長さを求める

`length` を使います。

```python: LineString オブジェクトの長さを求める
line.length # 6.4787086646190755
```

### Polygon

一定の面を表す。
丁目、街区など。

```json: Polygonタイプのジオメトリオブジェクト
{"type": "Polygon",
  "coordinates": [
                  [[35, 10], [45, 45], [15, 40], [10, 20], [35, 10]],
                  [[20, 30], [35, 35], [30, 20], [20, 30]]
                  ]
}
```

`coordinates` の中身は、最初と最後が同じ座標を持つ LineString のような配列を複数持つ配列となっている。（上記の図は、四角形の中に三角形が入っているような図形）

#### shapely を使った　Polygon オブジェクトの作成と描画

```python: Polygon オブジェクトの作成
from shapely.geometry import Polygon
polygon = Polygon([(0,0),(0,1),(1,1),(1,0),(0,0)])
```

上記の `polygon` を表示すると、以下のような図形が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/68849689-7182-98bf-bff7-9de62de8f44f.png)


#### shapely を使った　Polygon オブジェクトの型やデータ内容の表示

```python: Polygon オブジェクトの表示
print(polygon)     # POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0)) のように表示される
polygon.geom_type  # 'Polygon' と表示される
```

#### Polygon オブジェクト(面)の面積を求める

`area` を使います。

```python: Polygon オブジェクトの面積を求める
polygon.area # 1.0
```

### Multi○○

複数の `Point`, `LineString`, `Polygon` を持つ `MultiPoint`, `MultiLineString`, `MultiPolygon`, `MultiCollection` といったタイプもある。

ただこれらは基本の3つを集めたものなので、基本は `Point`, `LineString`, `Polygon` の3つ。

#### shapely を使った　MultiPoint オブジェクトの作成と描画

```python: MultiPoint オブジェクトの作成
from shapely.geometry import MultiPoint
points = MultiPoint([(0,0),(0,1),(1,1),(1,0),(0,0)])
```

上記の `points` を表示すると、以下のような4つの点が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/6d718e21-3326-87f2-0cc6-424f8c5dcb8f.png)


#### shapely を使った　MultiPoint オブジェクトの型やデータ内容の表示

```python: MultiPoint オブジェクトの表示
print(points)     # MULTIPOINT (0 0, 0 1, 1 1, 1 0, 0 0) のように表示される
points.geom_type  # 'MultiPoint' と表示される
```

#### shapely を使った　MultiLineString オブジェクトの作成と描画

```python: MultiLineString オブジェクトの作成
from shapely.geometry import MultiLineString

# lines に突っ込む線の座標をあらかじめ定義
line_coords = [((0,0),(1,1)),((-1,0), (-2,0))]

# MultiLineString オブジェクトの生成
lines = MultiLineString(line_coords)
```

上記の `lines` を表示すると、以下のような2つの線が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/c92621c6-5213-44a9-9842-ee5c5acca0a4.png)


#### shapely を使った　MultiLineString オブジェクトの型やデータ内容の表示

```python: MultiPoint オブジェクトの表示
print(lines)     # MULTILINESTRING ((0 0, 1 1), (-1 0, -2 0)) のように表示される
lines.geom_type  # 'MultiLineString' と表示される
```

#### shapely を使った　MultiPolygon オブジェクトの作成と描画

```python: MultiPolygon オブジェクトの作成
from shapely.geometry import MultiPolygon

# polygons に突っ込む図形をあらかじめ定義
# 四角形
polygon_a = Polygon([(0,0),(0,1),(1,1),(1,0),(0,0)])
# 三角形
polygon_b = Polygon([(-1,-1),(-2,-2),(1,-2),(-1,-1)])

# MultiLinePolygon オブジェクトの生成
polygons = MultiPolygon([polygon_a, polygon_b])
```

上記の `polygons` を表示すると、以下のような2つの線が表示される。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/587845/2f2719ff-8117-e647-c47c-35d2909242a3.png)


#### shapely を使った　MultiLineString オブジェクトの型やデータ内容の表示

```python: MultiPoint オブジェクトの表示
print(polygons)     # MULTIPOLYGON (((0 0, 0 1, 1 1, 1 0, 0 0)), ((-1 -1, -2 -2, 1 -2, -1 -1))) のように表示される
polygons.geom_type  # 'MultiPolygon' と表示される
```

## GISデータのファイル形式

| ファイル名称 | 拡張子 | 備考 |
|:-|:-|:-|
|[シェープファイル](https://www.esrij.com/gis-guide/esri-dataformat/shapefile/)| - |以下の拡張子を持つファイルの集合体<br>（他にも任意のファイル形式あり）<br>各ファイルサイズは2GB以内、フィールド名は10バイトまで|
|〃| shp|図形の座標|
|〃| dbf|属性情報|
|〃| shx|図形と属性の対応関係|
|[GeoJSON](https://ja.wikipedia.org/wiki/GeoJSON)|geojson|JSON形式|
|XML|xml|マークアップ形式|
|[KML](https://developers.google.com/kml/documentation/kml_tut?hl=ja)|kml|Google Earth に使われている|
|[GML](https://ja.wikipedia.org/wiki/Geography_Markup_Language)|gml|XMLに基づく、Geoに特化したもの|

## CRS - 座標参照系

位置情報の表現方法の規格というか体系として、「座標参照系」というものがある。

- 測地系 - 地球を球体としてとらえた位置を定義
- 座標参照系 - 地球（の一部）を平面としてとらえた位置を定義
    - 地理座標系 - 緯度と経度の組み合わせによる表現
    - 投影座標系 - 位置情報を二次元座標上に投影

### 座標参照系の例

geopandas では、以下のように標準のサンプルデータを利用することができる。

```python
import geopandas as gpd
import folium
import matplotlib.pyplot as plt
from geodatasets import get_path # 標準データセットをダウンロードし、パスを取得する

path = get_path('nybb')

gdf = gpd.read_file(path)
```

こうして取得したデータ `gdf` の CRSを表示するには、以下のようにする。

```python
# 座標参照系(CRS)の情報を表示
gdf.crs
```

すると、以下のような情報が返ってくる。

```
<Projected CRS: EPSG:2263>
Name: NAD83 / New York Long Island (ftUS)
Axis Info [cartesian]:
- X[east]: Easting (US survey foot)
- Y[north]: Northing (US survey foot)
Area of Use:
- name: United States (USA) - New York - counties of Bronx; Kings; Nassau; New York; Queens; Richmond; Suffolk.
- bounds: (-74.26, 40.47, -71.8, 41.3)
Coordinate Operation:
- name: SPCS83 New York Long Island zone (US survey foot)
- method: Lambert Conic Conformal (2SP)
Datum: North American Datum 1983
- Ellipsoid: GRS 1980
- Prime Meridian: Greenwich
```

いろいろ細かいことが書いてあるが、最上段の `<Projected CRS: EPSG:2263>` を見ることで、 このデータはCRSの `EPSG2263` という体系に則って表現されているんだな、ということがわかる。

https://epsg.io/2263

⬆️によれば、EPSG2263はニューヨークのロングアイランド周辺の位置情報をけっこう細かめに表現する体系らしい。

以下のようなコードで、CRSの体系を変換することもできる。

```python
gdf_converted = gdf.to_crs(epsg=4326)
```

:::note
`inplace=True` オプションを指定してインプレースでの変換もできる
:::

https://epsg.io/4326

⬆️によれば、EPSG4326は世界地図を平面としてとらえ、GPSなどでも使われている体系のようです。

### よく使う座標系

日本の領土・領海に特化した座標系としては、EPSG6668が有名とのことです。

https://epsg.io/6668

⬆️は単位が度ですが、メートル単位の座標系として EPSG6691 もまぁまぁ使うらしいです。

https://epsg.io/6691