-
Notifications
You must be signed in to change notification settings - Fork 0
/
DAMM_LsqFit.jl
42 lines (33 loc) · 1.58 KB
/
DAMM_LsqFit.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
# Least square error fit of the DAMM model to data
using LsqFit, CSV
# Load fixed param
#EaSx = 53
R = 8.314472*10^-3 # Universal gas constant
Dl = 3.17 # Diffusion coeff of substrate in liquid phase
Sxt = 0.0125 # Soil C content
p = 2.4*10^-2 # Fraction of soil C that is considered soluble
Dgas = 1.67 # Diffusion coefficient of oxygen in air
# Db = 1.53 # Soil bulk density !! DAMM complex if SWC > 0.365. old 1.5396
Db = 1.43
Dp = 2.52 # Soil particle density
# DAMM model, !for the algorithm to work, Param are scaled to be of similar magnitude
@. multimodel(Ind_var, Param) = (1e8*Param[1]*exp(-Param[2]/(R*(273.15+Ind_var[:, 1]))))*
(((p*Sxt)*Dl*Ind_var[:, 2]^3)/(1e-8*Param[3]+((p*Sxt)*Dl*Ind_var[:, 2]^3)))*
((Dgas*0.209*(1-Db/Dp-Ind_var[:, 2])^(4/3))/
(1e-3*Param[4]+(Dgas*0.209*(1-Db/Dp-Ind_var[:, 2])^(4/3))))*
10000*10/1000/12*1e6/60/60
# Example
Param = Param_ini = [1.0, 62.0, 3.46, 2.0] # Scaled Param, as explained above. AlphaSx, EaSx, kMsx, kMo2
Ind_var = [15.0 0.3; 8 0.3]
multimodel(Ind_var, Param)
# Load data
df = CSV.read("Input\\RSMmean.csv",dateformat="yyyy-mm-dd")
SWC = df.SWC; Tsoil = df.Tsoil; Rsoil = Dep_var = df.RSM
Ind_var = hcat(Tsoil, SWC)
# Fit DAMM to data
# Need to load data first... TO DO
fit = curve_fit(multimodel, Ind_var, Dep_var, Param_ini) # For DAMM, Ind_var is Tsoil and SWC, Dep_var is Rsoil
Param_fit = coef(fit)
Modeled_data = multimodel(Ind_var,Param_fit)
# NB. fit on manual Rsoil survey, daily average (64 locations)
# Need to get an actual estimate of 484 site parameters! (using EucFACE values here), e.g. soil particle density, soil C content, etc.