はじめに

気象・海洋・気候変動などの分野で国際的に広く使われているデータフォーマットであるnetCDFを、GNU Octaveというソフトウェアを使って、簡単に可視化させることを目的としています。なおGNU Octaveはオープンソースであり、GPLの条件下で無料で使用することができます。

今回は、地形データをダウンロードして、世界地図を描画させていきます。プログラミングの学習に時間をかけず、簡単にデータを可視化させることに重点を置いています。機械学習や高度な解析が必要であれば、次のステップとして、Pythonなどを学習するとよいでしょう。

インストール

ここでは windows10を使います。

GNU Octave のインストールについては、こちらの記事を参照してください。

標高データの準備

ETOPO1地形データを使用して世界地図を描画させます。

ETOPO1はNOAA（アメリカ海洋大気庁）で公開されている陸地・海底の標高データです。緯度経度1分間隔の2次元格子データとなっています。

下記のページからダウンロードしてみましょう。

ETOPO1ダウンロードページ

今回は南極とグリーンランドの氷床下の地球表面のデータを選択

ETOPO1_Bed_g_gmt4.grd.gz をダウンロード

gzip圧縮してあるので解凍を行い、解凍後の「ETOPO1_Bed_g_gmt4.grd」ファイルを、作業フォルダに移動しておきましょう。ここでは作業フォルダを「 D:\octave\ 」とします（※作業フォルダは好みの場所に作成してください）。

世界地図を描画

Octaveを起動して「現在のディレクトリ」を、作業フォルダにして、ファイルブラウザに「ETOPO1_Bed_g_gmt4.grd」が表示されることを確認しましょう。

以下のコマンドをOctaveのコマンドウィンドウに入力してエンターを押しましょう。

# netcdf を取り扱うパッケージを読み込み
pkg load netcdf

# ETOPO1読み込み
topo_file = "ETOPO1_Bed_g_gmt4.grd";
x = ncread(topo_file,"x");
y = ncread(topo_file,"y");
z = ncread(topo_file,"z");
z=flipud(rot90(z));
[xm,ym] = meshgrid(x,y);

# 1分格子(21601x10801)を1度格子(360x180)に変換
lon=[-179.5:1:179.5];
lat=[-89.5:1:89.5];
[xi,yi]=meshgrid(lon,lat);
zi = interp2(xm,ym,z,xi,yi);

#西経180°から始まるデータを 東経0°から始まるデータに順番を入れ替える
zi=[zi(:,lon>=0) zi(:,lon<0)];
lon=[lon(lon>=0) lon(lon<0)];
lon(lon<0)=lon(lon<0)+360;
topo=zi;

# 描画 (標高0mの等高線)
contour(lon,lat,topo,[0 0]);
daspect([1 1]);

【補足】#記号からはじまる行は、コメントと呼ばれるメモで、プログラムとして認識されない行です。

Figure1ウィンドウが新たに開き、世界地図を描画したと思います。

描画された世界地図

【補足】1度格子の地形データを再利用するのであれば、このデータを保存するのが便利です。

# 作成した1度格子の地形データを etopo_1deg.mat というファイル名で保存
save("etopo_1deg.mat","lon","lat","topo");