J-EGG500の水深データをGeoTiffに変換しました(補間がいまいち)
秋田港周辺の海底地形を作成したかったので、まず"海底地形データ フリー"で検索したところ
のリンクがありJ-EGG500というデータが日本周辺では定番のようです。早速、日本海洋データセンターのサイトに行ってデータをダウンロードしてきました。
データはこんな感じで説明もあります。
0 39.10205 139.00005 672 0 39.09754 139.00021 736 0 39.09303 139.00038 704 0 39.08852 139.00055 676 0 39.08401 139.00070 669 0 39.07949 139.00087 662 0 39.07498 139.00104 661
点のデータとしてQGISで直接読み込んでみると
それらしい位置に表示されますね。
これを扱いやすいようにGeoTiffに変換しようと検索したところ、結局GDALという地理関連の処理ライブラリで変換するのが主な方法のようです。
このYoutubeビデオ
https://www.youtube.com/watch?v=zLNLG0j13Cw&ab_channel=MakingSenseRemotely
#%% from osgeo import gdal import pandas as pd import numpy as np import os filepath="jodc-depth500mesh-20210728053445.txt" df = pd.read_csv(filepath, sep="\s+", header=None) # Keep only xyz df.drop(columns=[0], inplace=True) df.columns=['x','y','value'] df = df.sort_values(by=['y','x'], ascending=[False, True]) df.to_csv("dem.xyz", index=False, header=None, sep=" ",float_format='%.8f') #demn = gdal.Translate('dem.tif', 'dem.xyz', outputSRS='EPSG:32719') df.to_csv("dem.csv", index=False) if os.path.exists('dem.vrt'): os.remove('dem.vrt') with open("dem.vrt", 'w') as f: f.write("<OGRVRTDataSource>\n \ <OGRVRTLayer name=\"dem\">\n \ <SrcDataSource>dem.csv</SrcDataSource>\n \ <GeometryType>wkbPoint</GeometryType>\n \ <GeometryField encoding=\"PointFromColumns\" x=\"y\" y=\"x\" z=\"value\"/>\n \ </OGRVRTLayer>\n \ </OGRVRTDataSource>") g = gdal.Grid("dem.tif", "dem.vrt", outputSRS="EPSG:4326") g = None
それらしい位置に表示されました。補間の影響と思われる人工的なパターンが見られますがおそらく補間の方法とパラメータの設定で消えそうな感じです。
スクリプトを作成したというより、Youtubeのビデオにあるスクリプトをコピーしただけですね。作者に感謝です。また、こちらの記事も参考になりました。 500mメッシュ水深データをGISで利用する方法 - 自然環境保全のための周辺技術