Skip to content

Add flat noise weighted by lat to calibration #1154

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
wants to merge 1 commit into
base: main
Choose a base branch
from
Draft

Add flat noise weighted by lat to calibration #1154

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .buildkite/Project.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,7 @@ Thermodynamics = "b60c26fb-14c3-4610-9d3e-2d17fe7ff00c"

[compat]
ClimaAnalysis = "0.5.17"
ClimaCalibrate = "0.0.16"
ClimaCalibrate = "0.1.0"
ClimaDiagnostics = "0.2.13"
ClimaTimeSteppers = "0.7, 0.8"
EnsembleKalmanProcesses = "2.4.1"
Expand Down
7 changes: 5 additions & 2 deletions experiments/calibration/calibrate_land.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,9 +6,10 @@ const n_iterations = 10 # 1 iterations takes ~ 1.5 hour with current settings ((
const spinup_period = Year(1)
const start_date = DateTime(2008, 12, 01) # this is the start of the forward model spinup
@assert month(start_date + spinup_period) == 12 "The start of your calibration period should be December."
const nelements = (50, 15) # resolution - (horizontal elements (lon, lat), vertical elements (soil discretization))
const dirname = "land_swu_zenith_only" # ideally, describe your calibration in a few words
const nelements = (75, 15) # resolution - (horizontal elements (lon, lat), vertical elements (soil discretization))
const dirname = "land_swu_zenith_only_antarctica" # ideally, describe your calibration in a few words
const caldir = joinpath("/glade/campaign/univ/ucit0011/alexis/", dirname) # you might want to save somewhere else than login
const flat_noise = 5 # in unit of your variable_list. TODO: make this a vector, same dim as variable_list
# Don't forget to adjust your priors and forward_model files.


Expand Down Expand Up @@ -73,6 +74,8 @@ end
# Locations used for calibration (currently all coordinates on land):
include(joinpath(@__DIR__, "make_training_locations.jl"))
training_locations = make_training_locations(nelements)
training_locations = filter(loc -> loc[2] <= -60.0, training_locations) # Antartica only

# potentially we can add regional runs or specific lon lat bands or filter (e.g., regions with snow)

# NOTE: The noise is set in observationseries_era5.jl - adjust if needed.
Expand Down
33 changes: 24 additions & 9 deletions experiments/calibration/observationseries_era5.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ end

results = Dict(
k => Dict(v => [] for v in variable_list) for
k in (:all_seasons, :seasonal_vars, :seasonal_means)
k in (:all_seasons, :seasonal_vars, :seasonal_means, :flat_weigh)
)

"""
Expand All @@ -92,6 +92,11 @@ years_to_delete =
[y for y in unique(years) if y < year(start_date + spinup_period)]
valid_years = [y for y in unique(years) if !(y in years_to_delete)]

function flat_weigh_season(flat_noise, lat)
val = flat_noise^2 ./ max(cosd(lat), 0.1)
return fill(val, 4)
end

for (lon, lat) in training_locations
seasonal = Dict(
var_name => group_by_years(
Expand Down Expand Up @@ -123,13 +128,23 @@ for (lon, lat) in training_locations
var_name in variable_list
)

flat_weigh_dict = Dict(
var_name => flat_weigh_season(flat_noise, lat) for
var_name in variable_list
)

for var in variable_list
push!(results[:all_seasons][var], seasonal[var])
push!(results[:seasonal_vars][var], variances[var])
push!(results[:seasonal_means][var], means[var])
push!(results[:flat_weigh][var], flat_weigh_dict[var])
end
end

# results[:flat_weigh] contains a flat noise weighted by lats.
# for example, with flat = 5 W m^-2, variance from 250 (15.8 W m-2 at poles)
# to 25 (5 W m-2 at equator)

# Build noise
# Currently, replaces values in results[:seasonal_vars] with 25.
# For example, results[:seasonal_vars]["lhf"][1] is a 4 element vector of 25.0 values. (no annual variation)
Expand All @@ -143,13 +158,13 @@ end
#end

# latitudes = map(x -> x[2], training_locations)
for var in values(results[:seasonal_vars])
for i in 1:length(var)
var[i] .= 5^2 # flat noise
# IF scale by lat, # var[i] .= var[i] ./ max(cosd(lat), 0.1) # need to get lat
# IF add model error, add results[:seasonal_means][var] .* 0.05 # or some factor
end
end
#for var in values(results[:seasonal_vars])
# for i in 1:length(var)
# var[i] .= 5^2 ./ max(cosd(lat), 0.1)
# # IF scale by lat, # var[i] .= var[i] ./ max(cosd(lat), 0.1) # need to get lat
# # IF add model error, add results[:seasonal_means][var] .* 0.05 # or some factor
# end
#end

# Generate observation series
obs_y = [
Expand All @@ -159,7 +174,7 @@ obs_y = [
Dict(
"samples" => results[:all_seasons][var_name][i][y],
"covariances" => LinearAlgebra.Diagonal(
results[:seasonal_vars][var_name][i],
results[:flat_weigh][var_name][i],
),
"names" => "$(var_name)_$(lon)_$(lat)_$(y)",
),
Expand Down