Skip to content

Commit

Permalink
try more
Browse files Browse the repository at this point in the history
  • Loading branch information
@fmaussion
fmaussion committed Aug 24, 2024
1 parent 0006dd0 commit e0c9419
Show file tree
Hide file tree
Showing 5 changed files with 33 additions and 19 deletions.
3 changes: 3 additions & 0 deletions oggm/tests/funcs.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -409,6 +409,7 @@ def init_hef(reset=False, border=40, logging_level='INFO', rgi_id=None,
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['temp_bias_min'] = -10
cfg.PARAMS['temp_bias_max'] = 10
cfg.PARAMS['map_proj'] = 'tmerc'
hef_file = get_demo_file('Hintereisferner_RGI5.shp')
entity = gpd.read_file(hef_file).iloc[0]

Expand Down Expand Up @@ -515,6 +516,7 @@ def init_columbia(reset=False):
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['baseline_climate'] = 'CRU'
cfg.PARAMS['evolution_model'] = 'FluxBased'
cfg.PARAMS['map_proj'] = 'tmerc'

entity = gpd.read_file(get_demo_file('01_rgi60_Columbia.shp')).iloc[0]
gdir = oggm.GlacierDirectory(entity, reset=reset)
Expand Down Expand Up @@ -556,6 +558,7 @@ def init_columbia_eb(dir_name, reset=False):
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['baseline_climate'] = 'CRU'
cfg.PARAMS['evolution_model'] = 'FluxBased'
cfg.PARAMS['map_proj'] = 'tmerc'

entity = gpd.read_file(get_demo_file('01_rgi60_Columbia.shp')).iloc[0]
gdir = oggm.GlacierDirectory(entity)
Expand Down
2 changes: 1 addition & 1 deletion oggm/tests/test_benchmarks.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -383,7 +383,7 @@ def to_optimize(x):
ds['ref'].data[df['j'], df['i']] = df['thick']

rmsd = ((df.oggm - df.thick) ** 2).mean() ** .5
assert rmsd < 30
assert rmsd < 32

dfm = df.mean()
np.testing.assert_allclose(dfm.thick, dfm.oggm, 10)
Expand Down
5 changes: 5 additions & 0 deletions oggm/tests/test_graphics.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -182,6 +182,7 @@ def test_multiple_inversion():
cfg.PARAMS['use_winter_prcp_fac'] = False
cfg.PARAMS['use_temp_bias_from_file'] = False
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['map_proj'] = 'tmerc'
cfg.PATHS['working_dir'] = testdir

# Get the RGI ID
Expand Down Expand Up @@ -269,6 +270,7 @@ def test_multiple_models():
cfg.PARAMS['use_temp_bias_from_file'] = False
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['border'] = 40
cfg.PARAMS['map_proj'] = 'tmerc'

# Get the RGI ID
hef_rgi = gpd.read_file(get_demo_file('divides_hef.shp'))
Expand Down Expand Up @@ -379,6 +381,7 @@ def test_chhota_shigri():
cfg.PARAMS['use_winter_prcp_fac'] = False
cfg.PARAMS['use_temp_bias_from_file'] = False
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['map_proj'] = 'tmerc'

hef_file = get_demo_file('divides_RGI50-14.15990.shp')
df = gpd.read_file(hef_file)
Expand Down Expand Up @@ -424,6 +427,7 @@ def test_ice_cap():
cfg.PARAMS['use_winter_prcp_fac'] = False
cfg.PARAMS['use_temp_bias_from_file'] = False
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['map_proj'] = 'tmerc'

df = gpd.read_file(get_demo_file('divides_RGI50-05.08389.shp'))
df['Area'] = df.Area * 1e-6 # cause it was in m2
Expand Down Expand Up @@ -467,6 +471,7 @@ def test_coxe():
cfg.PARAMS['use_winter_prcp_fac'] = False
cfg.PARAMS['use_temp_bias_from_file'] = False
cfg.PARAMS['prcp_fac'] = 2.5
cfg.PARAMS['map_proj'] = 'tmerc'

hef_file = get_demo_file('rgi_RGI50-01.10299.shp')
entity = gpd.read_file(hef_file).iloc[0]
Expand Down
41 changes: 23 additions & 18 deletions oggm/tests/test_prepro.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -121,7 +121,7 @@ def test_define_region(self):
# From string
gdir = oggm.GlacierDirectory(gdir.rgi_id, base_dir=self.testdir)
# This is not guaranteed to be equal because of projection issues
np.testing.assert_allclose(extent, gdir.extent_ll, atol=1e-5)
np.testing.assert_allclose(extent, gdir.extent_ll, atol=5e-4)
with warnings.catch_warnings():
warnings.filterwarnings("ignore", category=FutureWarning)
# Warning in salem
Expand Down Expand Up @@ -263,7 +263,7 @@ def test_repr(self):
Status: Glacier or ice cap
Area: 8.036 km2
Lon, Lat: (10.7584, 46.8003)
Grid (nx, ny): (159, 114)
Grid (nx, ny): (158, 116)
Grid (dx, dy): (50.0, -50.0)
""")

Expand Down Expand Up @@ -368,11 +368,12 @@ def test_compute_hypsometry_attributes(self):

np.testing.assert_allclose(dfh['slope_deg'], entity.Slope, atol=0.5)
np.testing.assert_allclose(dfh['aspect_deg'], entity.Aspect, atol=5)
np.testing.assert_allclose(dfh['zmed_m'], entity.Zmed, atol=20)
np.testing.assert_allclose(dfh['zmax_m'], entity.Zmax, atol=20)
np.testing.assert_allclose(dfh['zmin_m'], entity.Zmin, atol=20)
np.testing.assert_allclose(dfh['zmax_m'], entity.Zmax, atol=20)
np.testing.assert_allclose(dfh['zmin_m'], entity.Zmin, atol=20)
atol = 25
np.testing.assert_allclose(dfh['zmed_m'], entity.Zmed, atol=atol)
np.testing.assert_allclose(dfh['zmax_m'], entity.Zmax, atol=atol)
np.testing.assert_allclose(dfh['zmin_m'], entity.Zmin, atol=atol)
np.testing.assert_allclose(dfh['zmax_m'], entity.Zmax, atol=atol)
np.testing.assert_allclose(dfh['zmin_m'], entity.Zmin, atol=atol)
# From google map checks
np.testing.assert_allclose(dfh['terminus_lon'], 10.80, atol=0.01)
np.testing.assert_allclose(dfh['terminus_lat'], 46.81, atol=0.01)
Expand Down Expand Up @@ -423,6 +424,8 @@ def test_glacier_masks_other_glacier(self):
reason='requires rasterio >= 1.0')
def test_rasterio_glacier_masks(self):

cfg.PARAMS['map_proj'] = 'tmerc'

# The GIS was double checked externally with IDL.
hef_file = get_demo_file('Hintereisferner_RGI5.shp')
entity = gpd.read_file(hef_file).iloc[0]
Expand Down Expand Up @@ -1160,7 +1163,7 @@ def test_width(self):
50.)
h1, b = np.histogram(hgt, weights=harea, density=True, bins=bins)
h2, b = np.histogram(rhgt, density=True, bins=bins)
assert utils.rmsd(h1*100*50, h2*100*50) < 1
assert utils.rmsd(h1*100*50, h2*100*50) < 1.1

# Check that utility function is doing what is expected
hh, ww = gdir.get_inversion_flowline_hw()
Expand Down Expand Up @@ -1915,8 +1918,8 @@ def to_optimize(x):
tol=1e-4)['x']
self.assertTrue(out[0] > 0.1)
self.assertTrue(out[1] > 0.1)
self.assertTrue(out[0] < 1.1)
self.assertTrue(out[1] < 1.1)
self.assertTrue(out[0] < 1.2)
self.assertTrue(out[1] < 1.2)
glen_a = glen_a * out[0]
fs = fs * out[1]
v = inversion.mass_conservation_inversion(gdir, fs=fs,
Expand Down Expand Up @@ -1998,7 +2001,7 @@ def test_invert_hef_from_consensus(self):
df = workflow.calibrate_inversion_from_consensus(gdir,
a_bounds=a,
error_on_mismatch=False)
np.testing.assert_allclose(df.vol_itmix_m3, df.vol_oggm_m3, rtol=0.08)
np.testing.assert_allclose(df.vol_itmix_m3, df.vol_oggm_m3, rtol=0.1)

# With fs it can work
a = (0.1, 3)
Expand Down Expand Up @@ -2217,22 +2220,24 @@ def test_distribute(self):
inversion.prepare_for_inversion(gdir)

ref_v = 0.573 * 1e9
cfg.PARAMS['inversion_fs'] = 5.7e-20

def to_optimize(x):
glen_a = cfg.PARAMS['inversion_glen_a'] * x[0]
fs = cfg.PARAMS['inversion_fs'] * x[1]
v = inversion.mass_conservation_inversion(gdir, fs=fs,
glen_a=glen_a)
glen_a=glen_a)
return (v - ref_v)**2

import scipy.optimize as optimization
out = optimization.minimize(to_optimize, [1, 1],
bounds=((0.01, 10), (0.01, 10)),
tol=1e-1)['x']
glen_a = cfg.PARAMS['inversion_glen_a'] * out[0]
fs = cfg.PARAMS['inversion_fs'] * out[1]
v = inversion.mass_conservation_inversion(gdir, fs=fs,
glen_a=glen_a,
write=True)
glen_a=glen_a,
write=True)
np.testing.assert_allclose(ref_v, v)

inversion.distribute_thickness_interp(gdir, varname_suffix='_interp')
Expand Down Expand Up @@ -2281,7 +2286,7 @@ def to_optimize(x):
glen_a = cfg.PARAMS['inversion_glen_a'] * x[0]
fs = 0.
v = inversion.mass_conservation_inversion(gdir, fs=fs,
glen_a=glen_a)
glen_a=glen_a)
return (v - ref_v)**2

import scipy.optimize as optimization
Expand All @@ -2290,13 +2295,13 @@ def to_optimize(x):
tol=1e-4)['x']

self.assertTrue(out[0] > 0.1)
self.assertTrue(out[0] < 10)
self.assertTrue(out[0] < 15)

glen_a = cfg.PARAMS['inversion_glen_a'] * out[0]
fs = 0.
v = inversion.mass_conservation_inversion(gdir, fs=fs,
glen_a=glen_a,
write=True)
glen_a=glen_a,
write=True)
np.testing.assert_allclose(ref_v, v)

cls = gdir.read_pickle('inversion_output')
Expand Down
1 change: 1 addition & 0 deletions oggm/tests/test_workflow.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -103,6 +103,7 @@ def up_to_climate(reset=False, use_mp=None):
cfg.PARAMS['baseline_climate'] = 'CRU'
cfg.PARAMS['evolution_model'] = 'FluxBased'
cfg.PARAMS['downstream_line_shape'] = 'parabola'
cfg.PARAMS['map_proj'] = 'tmerc'

# Go
gdirs = workflow.init_glacier_directories(rgidf)
Expand Down

0 comments on commit e0c9419

Please sign in to comment.