ジュリアでHTTPプロトコルを通じてHYCOMから海洋データをダウンロードする方法 📂データ取得

ジュリアでHTTPプロトコルを通じてHYCOMから海洋データをダウンロードする方法

概要

ハイコンでデータを一つずつ手でクリックして取得・整理するのは容易ではない。条件を指定するURLを与えて取得する方法を紹介する。

コード

`HTTP.jl` ¹

for y = 2019:2024
    url = replace("https://ncss.hycom.org/thredds/ncss/GLBy0.08/expt_93.0/ts3z/$y?var=water_temp
    &north=36.88
    &west=129.52
    &east=130.24
    &south=36.32
    &disableProjSubset=on
    &horizStride=1
    &time_start=$y-01-01T00%3A00%3A00Z
    &time_end=$y-12-31T21%3A00%3A00Z
    &timeStride=1
    &vertCoord=0
    &accept=netcdf4", '\n' => "", " " => "")
    @time response = HTTP.get(url)

    open("data_$y.nc4", "w") do f
        write(f, response.body)
    end
end

例えば上のコードは GLBy0.08 の93番実験で、2019年から2024年までの海面水温(SST, Sea Surface Temperature)を緯度36.32〜36.88、経度129.52〜130.24の範囲で3時間間隔のデータを取得するコードだ。

`NCDatasets.jl` ²

data_ = []
for y = 2019:2024
    ds  = NCDataset("data_$y.nc4")
    push!(data_, DataFrame([ds["time"][:] reshape(ds["water_temp"][:, :, 1, :], 150, :)'], ["t"; repeat("i" .* string.(1:10), outer = 15) .* repeat("j" .* string.(1:15), inner = 10)]))
    close(ds)
end

ダウンロードされるファイルは NetCDF-4 (Network Common Data Form version 4) を意味する *.nc4 ファイルだ。高次元データを扱うために広く使われているが、このガイドでは DataFrames.jl でデータフレームを作成している。

保存

_data_ = vcat(data_...)
_data_.t = Date.(_data_.t)
_data_ = combine(groupby(_data_, :t), names(_data_, Not(:t)) .=> mean .=> names(_data_, Not(:t)))
CSV.write("data_GLBy0.08_expt_93.0.csv", _data_)
tnsr = reshape(Matrix(_data_[:, 2:end])', 10, 15, :)

@save "data_tnsr.jld2" tnsr
# @load "data_tnsr.jld2"; tnsr

JLD2.jl を通じてテンソルそのままを保存することもできるし、書き出しが簡単な CSV.jl で保存する方法もある。

全体のコード

using HTTP, NCDatasets, DataFrames, Dates, CSV, JLD2, StatsBase
for y = 2019:2024
    url = replace("https://ncss.hycom.org/thredds/ncss/GLBy0.08/expt_93.0/ts3z/$y?var=water_temp
    &north=36.88
    &west=129.52
    &east=130.24
    &south=36.32
    &disableProjSubset=on
    &horizStride=1
    &time_start=$y-01-01T00%3A00%3A00Z
    &time_end=$y-12-31T21%3A00%3A00Z
    &timeStride=1
    &vertCoord=0
    &accept=netcdf4", '\n' => "", " " => "")
    @time response = HTTP.get(url)

    open("data_$y.nc4", "w") do f
        write(f, response.body)
    end
end

# ds  = NCDataset("data_$y.nc4")
# # ds.attrib
# # keys(ds)
# # ds["water_temp"][:, :, 1, :]
# CSV.write("lon.csv", DataFrame(i = 1:10, x = round.(ds["lon"][:], digits = 2)))
# CSV.write("lat.csv", DataFrame(j = 1:15, y = round.(ds["lat"][:], digits = 2)))
# close(ds)

data_ = []
for y = 2019:2024
    ds  = NCDataset("data_$y.nc4")
    push!(data_, DataFrame([ds["time"][:] reshape(ds["water_temp"][:, :, 1, :], 150, :)'], ["t"; repeat("i" .* string.(1:10), outer = 15) .* repeat("j" .* string.(1:15), inner = 10)]))
    close(ds)
end
_data_ = vcat(data_...)
_data_.t = Date.(_data_.t)
_data_ = combine(groupby(_data_, :t), names(_data_, Not(:t)) .=> mean .=> names(_data_, Not(:t)))
CSV.write("data_GLBy0.08_expt_93.0.csv", _data_)
tnsr = reshape(Matrix(_data_[:, 2:end])', 10, 15, :)

@save "data_tnsr.jld2" tnsr
# @load "data_tnsr.jld2"; tnsr

ニーニョ 3.4

次は ニーニョ 3.4^{Niño 3.4} と呼ばれる³、緯度5°S〜5°N、経度170°W〜120°Wの海面水温の平均値を取得できるコードだ。

using HTTP, NCDatasets, DataFrames, Dates, CSV, JLD2, StatsBase, ProgressMeter
using Base.Threads

@showprogress @threads for day in Date(2019, 1, 1):Date(2024, 9, 5)
try
    str_day = string(day)
    ("data_$str_day.nc4" |> isfile) && continue
    url = replace("https://ncss.hycom.org/thredds/ncss/GLBy0.08/expt_93.0/ts3z/$(Year(day).value)?var=water_temp
    &north=5
    &west=120
    &east=170
    &south=-5
    &disableProjSubset=on
    &horizStride=1
    &time_start=$(day)T00%3A00%3A00Z
    &time_end=$(day)T21%3A00%3A00Z
    &timeStride=1
    &vertCoord=0
    &accept=netcdf4", '\n' => "", " " => "")
    # @time response = HTTP.get(url)
    response = HTTP.get(url)

    open("data_$str_day.nc4", "w") do f
        write(f, response.body)
    end
catch e
    @warn e
end
end

nc4_ = filter(endswith(".nc4"), readdir())
@time data_ = CSV.read("data_form.csv", DataFrame)
@showprogress for k in eachindex(nc4_)
    NCDataset(nc4_[k]) do ds
        df = DataFrame([ds["time"][:] reshape(ds["water_temp"][:, :, 1, :], 251*626, :)'], ["t"; repeat("i" .* string.(1:626), outer = 251) .* repeat("j" .* string.(1:251), inner = 626)])
        df.t = Date.(df.t)
        data_[k, :] = combine(groupby(df, :t), names(df, Not(:t)) .=> mean .=> names(df, Not(:t)))[1, :]
    end
end

CSV.write("data_GLBy0.08_expt_93.0.csv", data_)
tnsr = reshape(Matrix(data_[:, 2:end])', 626, 251, :)

@save "data_tnsr.jld2" tnsr