位置情報データの距離測定と座標変換

import pandas as pd
import geopandas as gpd
# 鉄道駅の位置情報データ（国土数値情報）
# https://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-N05-2023.html
file_path = '../../00_common_data/01_国土数値情報/02_鉄道ライン/2023/N05-23_GML/N05-23_Station2.shp'
# シェープファイルの読み込み
gdf_station = gpd.read_file(file_path, encoding='shift-jis')

# 座標参照系の確認
print(gdf_station.crs)
# epsg:6668

# 位置情報データの読み込み
file_path = '../00_data/GPS_241109_SJIS_LF_新橋渋谷.CSV'
df = pd.read_csv(file_path, encoding='cp932')
df['geometry'] = gpd.points_from_xy(df['経度'], df['緯度'])
gdf = gpd.GeoDataFrame(df, geometry='geometry')
# WGS84として座標参照系を設定（EPSG:4326）
gdf = gdf.set_crs(epsg=4326) 
print(gdf.crs)
# epsg:4326

# 座標変換の実施
crs_after = 'epsg:6677' # JGD2011 Ⅸ(9)系
gdf_crs = gdf.to_crs(crs=crs_after)
print(gdf_crs.crs)
# epsg:6677

# 渋谷駅との距離を計算
# 渋谷駅の位置情報データを抽出
gdf_station = gdf_station.rename(columns={
    'N05_001':'事業者種別', 'N05_002':'路線名', 'N05_003':'運営会社', 'N05_004':'供用開始年', 'N05_005b':'設置期間（設置開始）', 'N05_005e':'設置期間（設置終了）',
    'N05_006':'関係ID', 'N05_007':'変遷ID', 'N05_008':'変遷備考', 'N05_009':'備考', 'N05_011':'駅名',
    })
gdf_station = gdf_station[(gdf_station['路線名']=='山手線')&(gdf_station['駅名']=='渋谷')]
display(gdf_station)
# POINT (139.70165 35.65809)

# 座標変換の実施（地理座標系から平面直角座標系へ）
crs_after = 'epsg:6677' # JGD2011 Ⅸ(9)系 ※人流データと座標系を揃える
gdf_station_crs = gdf_station.to_crs(crs=crs_after)
print(gdf_station_crs.crs)
display(gdf_station_crs)
# POINT (-11923.357 -37924.796)

# 渋谷駅の位置座標（座標変換後）
loc_shibuya = gdf_station_crs['geometry'].iloc[0]
print(loc_shibuya)

# 渋谷駅との距離を測定
gdf_dist = gdf_crs.copy()
gdf_dist['渋谷駅との距離'] = gdf_dist.distance(loc_shibuya)
display(gdf_dist.head())
print(gdf_dist['渋谷駅との距離'].head(6))
# epsg:6677 # JGD2011 Ⅸ(9)系
# 0    78.943404
# 1    78.943404
# 2    78.943404
# 3    78.943404
# 4    88.596504
# 5    88.596504

print(gdf_dist['渋谷駅との距離'].head(6))
# epsg:3100 UTM54N系
# 0    78.933013
# 1    78.933013
# 2    78.933013
# 3    78.933013
# 4    88.584843
# 5    88.584843

from shapely.geometry import Point, LineString, Polygon
from geopy.distance import geodesic

gdf_dist = gdf.copy()
# 距離測定用前処理
gdf_dist['経度緯度'] = gdf_dist.apply(lambda x: Point(x['経度'], x['緯度']), axis=1)
# 渋谷駅の位置座標
print(loc_shibuya.y, loc_shibuya.x)
# 35.658093 139.701645

# 渋谷駅との距離を測定
gdf_dist['渋谷駅との距離'] = [geodesic((loc.y, loc.x), (loc_shibuya.y, loc_shibuya.x)).meters for loc in gdf_dist['geometry']]
print(gdf_dist['渋谷駅との距離'].head(10))
# 0    78.951161
# 1    78.951161
# 2    78.951161
# 3    78.951161
# 4    88.605210
# 5    88.605210