use residuals and not χ2 for Gauss-Newton

This commit is contained in:
Gaspard Jankowiak 2024-05-16 09:39:28 +02:00
parent 339cdf749d
commit 1aa2170ada
2 changed files with 13 additions and 2 deletions

View file

@ -77,6 +77,10 @@ function lowpass_filter(i; σ=1)
return r return r
end end
function residuals(I_data, I_model, err)
return (I_data .- I_model ./ err) .^ 2
end
function chi2(I_data, I_model, err) function chi2(I_data, I_model, err)
return sum((I_data .- I_model ./ err) .^ 2) return sum((I_data .- I_model ./ err) .^ 2)
end end

11
test.jl
View file

@ -68,14 +68,21 @@ padding_factor[18] = 10.0
bounds = MH.boxconstraints(lb=lb_reduced, ub=ub_reduced) bounds = MH.boxconstraints(lb=lb_reduced, ub=ub_reduced)
function obj(P) function obj_χ2(P)
I_model, neg_H20 = intensity_reduced(P) I_model, neg_H20 = intensity_reduced(P)
χ2 = Utils.chi2(I_data, I_model, err) χ2 = Utils.chi2(I_data, I_model, err)
factor = (neg_H20 == 0 ? 1.0 : 5.0 * (neg_H20 + 1)) factor = (neg_H20 == 0 ? 1.0 : 5.0 * (neg_H20 + 1))
return factor * χ2 return factor * χ2
end end
barriered_obj = Utils.add_log_barriers(obj, simple_bounds; padding_factor=padding_factor) function obj_residuals(P)
I_model, neg_H20 = intensity_reduced(P)
residuals = Utils.residuals(I_data, I_model, err)
factor = (neg_H20 == 0 ? 1.0 : 5.0 * (neg_H20 + 1))
return factor * residuals
end
barriered_obj = Utils.add_log_barriers(obj_residuals, simple_bounds; padding_factor=padding_factor)
information = MH.Information(f_optimum=0.0) information = MH.Information(f_optimum=0.0)
information = MH.Information() information = MH.Information()