[解析雨量を Pythonで読み、作図する (試作版)]で 図化した雨量の分布を、地形とくらべてみるために、 地面の標高のデータを、雨量の図と同じ図法・同じ大きさの地図に表示してみようと思った。
対象地域が日本だから、日本の国土地理院の標高データがよいと思ったのだが、 つぎのような状況で、つかいにくかった。
世界の陸の標高と海底の水深をまとめたデータとしては、 アメリカ合衆国の NOAA の 旧 National Geophysical Data Center (現在は National Centers for Environmental Information の一部) がつくったものが知られている。 いま代表的なデータセットは ETOPO1 とよばれるもので、 「ETOPO1 Global Relief Model」 https://www.ngdc.noaa.gov/mgg/global/ のページにある。 氷床のところで、氷の表面を採用したものと、氷の底の岩石の表面を採用したものとの2種類の データセットをつくっている。 測地系は「WGS84」で、日本が2002年に採用した「世界測地系」と厳密に同じではないがかなり近い。 データのダウンロードは https://www.ngdc.noaa.gov/mgg/global/relief/ETOPO1/data/ からできる。 格子点の位置に値がある grid-registered と、格子線でかこまれた升目に値がある cell-registered の データがつくられているが、grid-registered のほうがもとである。 データ形式は、netCDF、geoTIFF などもえらべるが、ここでは単純バイナリにした。 grid-registered の単純バイナリのデータは、4バイト浮動小数点数、下位バイト先行 (little endian)で、 格子点数は 東西 21601、南北 10801 である。これをダウンロードした。
このうち、関東甲信越をふくむ 北緯 34.5~38.5度、東経 136.5~141.5度 の範囲を切り出した。 (切り出しは GrADS でやった [注]。まず手作業で GrADS の ctl ファイルをつくり、 GrADS の中で緯度経度範囲を指定して「set gxout fwrite」でファイルに書きだした。) もとの ETOPO1 のデータの格子点のならびは北からだったが、切り出したものは (GrADSの default の) 南からになっている。 格子点数は 東西 301、南北 241 になった。
【[注] バイナリデータとして読んでその一部分をバイナリデータとして書きだすだけなので、 Python でもできるのだが、わたしが切り出すべき部分を正しくかぞえられるか自信がなく、 GrADS で 中間結果を図示しながらやるほうがまちがえにくいと判断した。】
単純バイナリデータを Python でよむ方法は、つぎのウェブページを参考にした。 (big endian のデータむけに「>」となっているところを、little endian なので「<」に変えた。)
https://qiita.com/mhangyo/items/76db7c6a6ebba6cf4330 @mhangyoが2017年08月07日に作成 Pythonでbig endianのbinaryを読み込んでndarrayに変換する
Cartopy を利用した作図は、解析雨量のばあいと基本的に同じで、手本はつぎのウェブページである。
https://qiita.com/wm-ytakano/items/1e6a5d250fe0d079b08e @wm-ytakanoが2021年02月14日に更新 PythonでGSMのGPVデータから気温の水平分布を描画する
(海岸線は Cartopy によって Natural Earth というサイトから取得されたデータがつかわれているが、 ここでは、解析雨量の図につかった「50 m」ではあらすぎると判断して「10 m」を採用した。)
標高 0 m をさかいに色をかえたが、ETOPO1 には陸と海の区別がないので、 海面より低い陸地が海と同様な色になり、湖が部分的に海と同様な色になって、 平野部の地図としてはちょっと気もちのよくないところがある。 山地の分布をみるにはさしつかえないだろう。
[2021-12-09 補足] 実行環境によっては、海が青でなく陸の低地と同じ緑でぬられてしまうことがある。 そこでうごいている Cartopy のバージョンか、matplotlib のバージョンでは、 mappable.cmap.set_under(), mappable.cmap.set_over() が無効になっているらしい。 そのばあいの対策として、0 m 未満、1600 m 以上の区間も levels にふくめ、 それに対応する色も colors にふくめてぬりわけるように修正したバージョンをつくった。