example_01.qmd

# BEPS

```{julia}
# import RTableTools: fread
# d = fread("examples/input/p1_meteo")
```

```{julia}
using BEPS
d = deserialize("data/p1_meteo")
lai = readdlm("examples/input/p1_lai.txt")[:]

par = (lon=120.5, lat=30.5, landcover=25, clumping=0.85,
  soil_type=8, Tsoil=2.2,
  soilwater=0.4115, snowdepth=0.0)
```

```{julia}
@time df_jl, df_ET_jl, Tg = besp_main(d, lai, par; version="julia");
@time df_c, df_ET_c, Tg = besp_main(d, lai, par; version="c");

sum(df_jl)
sum(df_c)

df_diff = abs.(df_jl .- df_c)
df_diff_perc = abs.(df_jl .- df_c) ./ df_c .* 100

max(df_diff)
max(df_diff_perc) # SH, 1.48%的误差
```

# 绘图展示结果

```{julia}
using Plots
using Ipaper
gr(framestyle=:box)

function plot_var(var=:SH)
  n = size(df_jl, 1)
  # n = 24 * 7
  inds = 1:n

  y_jl = df_jl[inds, var]
  y_c = df_c[inds, var]
  diff_perc = (y_jl .- y_c) ./ y_c .* 100
  plot(diff_perc; title="$var", label=nothing) # 
end

vars = names(df_jl)[[1:4; 8; 12:16]]
ps = map(plot_var, vars)

plot(ps..., layout=(3, 4), size=(1400, 800))
# write_fig("images/Figure1_bias_of_julia-version.png", 15, 8; show=false)
```


```{julia}
## 土壤温度的变化
begin
  using Plots
  gr(framestyle=:box)

  depths = [0.05, 0.10, 0.20, 0.40, 1.25] |> cumsum |> x -> round.(x, digits=4)

  p = plot(; size=(1400, 700), title="Soil temperature")
  plot!(d[:, :tem]; label="Tair")
  for i = 1:5
    label = "depth: $(depths[i])"
    plot!(Tg[:, i]; label)
  end
  p
end

# write_fig("images/Figure2_variation_of_Tg.png", 15, 8; show=false)
```