corrected glacier gridpoints on flattened files (#1547)

* corrected glacier gridpoints on flattened files

* maximum distance <=0.25, not <0.25

* renamed tas_std to temp_std in flattened file+tests

* Small changes

---------

Co-authored-by: Fabien Maussion <fabien.maussion@uibk.ac.at>

oggm/shop/gcm_climate.py

@@ -224,7 +224,7 @@ def process_monthly_isimip_data(gdir, output_filesuffix='',
     else:
         gcm_server = 'https://cluster.klima.uni-bremen.de/~oggm/'
 
-    path = f'{gcm_server}/cmip6/isimip3b/flat/monthly/'
+    path = f'{gcm_server}/cmip6/isimip3b/flat/2023.2/monthly/'
     add = '_global_monthly_flat_glaciers.nc'
 
     fpath_spec = path + '{}_w5e5_'.format(member) + '{ssp}_{var}' + add

oggm/shop/w5e5.py

@@ -23,7 +23,7 @@
 
 GSWP3_W5E5_SERVER = 'https://cluster.klima.uni-bremen.de/~oggm/climate/'
 
-_base = 'gswp3-w5e5/flattened/monthly/'
+_base = 'gswp3-w5e5/flattened/2023.2/monthly/'
 
 BASENAMES = {
     'GSWP3_W5E5': {

oggm/tests/conftest.py

@@ -174,7 +174,6 @@ def secure_url_retrieve(url, *args, **kwargs):
                         'no_match/RGI62/b_040/{}/RGI60-15/RGI60-15.13.tar')
 
     assert ('github' in url or
-            'cluster.klima.uni-bremen.de/~oggm/ref_mb_params' in url or
             'cluster.klima.uni-bremen.de/~oggm/test_gdirs/' in url or
             'cluster.klima.uni-bremen.de/~oggm/demo_gdirs/' in url or
             'cluster.klima.uni-bremen.de/~oggm/test_climate/' in url or

oggm/tests/test_shop.py

@@ -253,8 +253,53 @@ def test_get_gswp3_w5e5_file(self):
         with pytest.raises(ValueError):
             w5e5.get_gswp3_w5e5_file(d, 'zoup')
 
-        # check if W5E5 and GSWP3_W5E5 are equal over common time period
-        # this is done in the flattening notebook
+    def test_glacier_gridpoint_selection(self):
+
+        from oggm.shop import w5e5
+        d = 'GSWP3_W5E5'
+        # test is only done for the `inv` file, as the other files are only
+        # downloaded for the HEF gridpoints as they would be too large otherwise.
+        # However, the same test and other tests are done for all files
+        # (also ISIMIP3b) and all glaciers in this notebook:
+        # https://nbviewer.org/urls/cluster.klima.uni-bremen.de/
+        # ~lschuster/example_ipynb/flatten_glacier_gridpoint_tests.ipynb
+        with xr.open_dataset(w5e5.get_gswp3_w5e5_file(d, 'inv')) as dinv:
+            dinv = dinv.load()
+
+        # select three glaciers where two failed in the
+        # previous gswp3_w5e5 version
+        for coord in [(10.7584, 46.8003),  # HEF
+                      (-70.8931 + 360, -72.4474),  # RGI60-19.00124
+                      (51.495, 30.9010),  # RGI60-12.01691
+                      (0, 0)  # random gridpoint not near to a glacier
+                      ]:
+            lon, lat = coord
+            # get the distances to the glacier coordinate
+            c = (dinv.longitude - lon) ** 2 + (dinv.latitude - lat) ** 2
+            c = c.to_dataframe('distance').sort_values('distance')
+            # select the nearest climate point from the flattened glacier gridpoint
+            lat_near, lon_near, dist = c.iloc[0]
+            # for a randomly chosen gridpoint, the next climate gridpoint is far away
+            if coord == (0, 0):
+                with pytest.raises(AssertionError):
+                    assert np.abs(lat_near - lat) <= 0.25
+                    assert np.abs(lon_near - lon) <= 0.25
+                    assert dist <= (0.25 ** 2 + 0.25 ** 2) ** 0.5
+            # for glaciers the next gridpoint should be the nearest
+            # (GSWP3-W5E5 resolution is 0.5°)
+            else:
+                assert np.abs(lat_near - lat) <= 0.25
+                assert np.abs(lon_near - lon) <= 0.25
+                assert dist <= (0.25 ** 2 + 0.25 ** 2) ** 0.5
+
+        # this only contains data for two glaciers, let's still check some basics
+        # both glaciers are not at latitude or longitude 0
+        with xr.open_dataset(w5e5.get_gswp3_w5e5_file(d, 'temp_std')) as dtemp_std:
+            assert np.all(dtemp_std.latitude != 0)
+            assert np.all(dtemp_std.longitude != 0)
+            assert dtemp_std.isel(time=0).temp_std.std() > 0
+            assert dtemp_std.longitude.std() > 0
+            assert dtemp_std.latitude.std() > 0
 
     def test_process_w5e5_data(self, class_case_dir):