ERA-5 remapping weights files flip longitude and latitude axes

[rmholmes]

My first test run with the https://github.com/COSIMA/1deg_era5_iaf repo has winds (and probably all the forcing fields) flipped in the longitude and latitude direction.

Jan 1981 zonal wind stress for ERA-5 (left) and JRA-55 (right);

This is likely a remapping weights issue, perhaps with the https://github.com/COSIMA/esmgrids/blob/master/esmgrids/era5_grid.py script?

[aekiss]

Actually, it appears from the phantom Australia in the North Atlantic and Africa in the Pacific that the data is flipped only in latitude, and offset (not flipped) in longitude.

ERA5 files, e.g. /g/data/rt52/era5/single-levels/reanalysis/10u/1979/10u_era5_oper_sfc_19790101-19790131.nc, have the latitude running in reverse from 90 to -90, which this np.flipud is apparently intended to deal with, but that doesn't seem the right thing to do if the data is not also being flipped.

esmgrids/esmgrids/era5_grid.py

Line 17 in 514b560

y_t = np.flipud(f.variables['latitude'][:])

ERA5 longitude starts at -180 (whereas JRA55-do starts at 0), which seems about right to explain the offset.

[rmholmes]

It looks like some reengineering of make_remap_weights.py and era5_grid.py is needed. Nic must have a version of make_remap_weights.py with the ERA5 options sitting around somewhere - but I can't find it. The changes look reasonably straightforward though - I can have a bash if no-one else is working on this?

For the latitude issue I could try just not flipping the latitude in era5_grid.py. For the longitude changes I could try modifying x_t in era5_grid.py. Otherwise understanding the ESMF_RegridWeightGen script will be neccessary.

[aekiss]

@nichannah's weights in /g/data/ik11/inputs/access-om2/input_20210915/common_1deg_era5 appear to have been generated by /g/data/v45/nah599/access-om2/tools/make_remap_weights.*.
I've pushed these into this PR: COSIMA/access-om2#258

[aekiss]

I'll be busy with other things today so feel free to have a go at fixing this issue. The first thing I'd try is omitting np.flipud and seeing what happens.

[aekiss]

CI tests failed on the PR. https://github.com/COSIMA/access-om2/runs/5964911440?check_suite_focus=true
I guess Nic didn't commit it because it had bugs?
I'll see if this is easily fixed.

[aekiss]

fixed bug that caused CI to fail.

[aekiss]

there are some tips on remapping here https://github.com/COSIMA/access-om2/wiki/Technical-documentation#creating-remapping-weights
and here https://github.com/COSIMA/access-om2/wiki/Tutorials#3-changing-the-oasis-remapping-files

[aekiss]

might need to be careful what version of ESMF_RegridWeightGen we use - see COSIMA/access-om2#216

[rmholmes]

Thanks @aekiss.

I've successfully generated remapping weights using module load esmf (which I think according to that issue should be ok for anything but the 1/10th-degree) without changing anything and am running a test (TBC).

If I try to remove the np.flipud on the latitude then a run into a check made in the make_corners routine in base_grid.py - it appears that the corner finding routines need latitude to be ordered south -> north.

We obviously don't want to derive new ERA-5 forcing files with latitude flipped. So it seems like the next step would be to rewrite the corner finding routines so that they work in both directions?

[aekiss]

thanks @rmholmes, that sounds sensible.

FYI I've put an executable using Russ' tripole bug fix here:
/g/data/ik11/inputs/access-om2/bin/ESMF_RegridWeightGen_f536c3e12d
and updated make_remap_weights.* in COSIMA/access-om2#258 to use it (as yet untested!)

[rmholmes]

Just a recap from my progress today:

If I remove the np.flipud in era5_grid.py and comment out the assert checks in make_corners (in base_grid.py) then the code succesfully creates remapping weight files, but nothing is changed. Similarly, if I add 90, 180 or 360 to x_t in era5_grid the files are created successfully but nothing is changed.

The code in esmgrids/base_grid.py seems apriori designed for dealing with latitude ordered south to north. Without this ordering, dy is negative which doesn't seem to be what is wanted. However, I don't have a good enough understanding of what the code is doing to figure out what should be done to fix things.

It is frustrating that the ERA-5 latitude data is ordered this way - doesn't seem to serve any practical purpose except to unneccesarily make things difficult.

[aekiss]

Thanks for all these experiments. I'll raise it at today's TWG meeting.

[aidanheerdegen]

Try this one
/scratch/v45/aph502/era5_mom10_wts.nc
I've tested it, and it gives the correct orientation.

[rmholmes]

Thanks @aidanheerdegen!!

That file has different variable names to the ones I've been working with (e.g. /g/data/e14/rmh561/access-om2/input/ERA-5/lat_flipud/ERA5_MOM1_patch.nc), so not sure it'll work out of the box?

[aidanheerdegen]

ncdump: /g/data/e14/rmh561/access-om2/input/ERA-5/lat_flipud/ERA5_MOM1_patch.nc: Permission denied

[rmholmes]

Fixed now.

[aidanheerdegen]

I've applied the same ncrename command to my weights file. Give it a burl.

[rmholmes]

Ok thanks. In setting this up I've found this issue, which I think invalidates all my checks yesterday. [EDIT: Nope - I was just using too many symbolic links]

But I can get around that now I know what it is. So I'll get back to you.

[rmholmes]

No luck, that gives me exactly the same result, again!

The run is at /home/561/rmh561/access-om2/1deg_era5_iaf if you want to check...

What did your testing consist of @aidanheerdegen? I feel like I'm missing something basic here...

[aidanheerdegen]

Hey Ryan

The remapping file we're altering is for the OASIS coupler, i.e. it is the one defined in namcouple

$ grep INPUT namcouple
../INPUT/ERA5_MOM1_conserve.nc dst
../INPUT/ERA5_MOM1_conserve.nc dst
../INPUT/ERA5_MOM1_conserve.nc dst
../INPUT/ERA5_MOM1_conserve.nc dst
../INPUT/ERA5_MOM1_patch.nc dst
../INPUT/ERA5_MOM1_conserve.nc dst
../INPUT/ERA5_MOM1_patch.nc dst
../INPUT/ERA5_MOM1_patch.nc dst
../INPUT/ERA5_MOM1_patch.nc dst
../INPUT/ERA5_MOM1_patch.nc dst

The INPUT directory referred to there is, rather confusingly, the one in work/INPUT, which is the input directory referred to in the access-om2 model, i.e. the top level

# Model configuration
name: common
model: access-om2
input: /g/data/ik11/inputs/access-om2/input_20210915/common_1deg_era5/

Looking at the git diff of your control directory compared to the previous run:

$ git diff HEAD^
diff --git a/manifests/input.yaml b/manifests/input.yaml
index f3653c6..4c7ed1c 100644
--- a/manifests/input.yaml
+++ b/manifests/input.yaml
@@ -30917,20 +30917,30 @@ work/atmosphere/INPUT/2021/sp_era5_oper_sfc_20211101-20211130.nc:                            
     binhash: 87a70249955c71776c9759a701b0f6f3
     md5: ba7376809074b024f3cb8f92177365ee
 work/atmosphere/INPUT/ERA5_MOM1_conserve.nc:
-  fullpath: /g/data/ik11/inputs/access-om2/input_20210915/common_1deg_era5/ERA5_MOM1_conserve.nc                      
+  fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/ERA5_MOM1_conserve.nc                                           
   hashes:
-    binhash: 9af5ca3052f6d5ed5d9139a23f121cf5
-    md5: 1bb0abf04369ea4ffe1d9da343d20f1f
+    binhash: b9d28101ccc3821fab8320579626812b
+    md5: e5c1c2ca736b4d394fdac397c46e7a21
 work/atmosphere/INPUT/ERA5_MOM1_patch.nc:
-  fullpath: /g/data/ik11/inputs/access-om2/input_20210915/common_1deg_era5/ERA5_MOM1_patch.nc                         
+  fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/ERA5_MOM1_patch.nc                                              
   hashes:
-    binhash: 2bab214c23d2eebc6d368e003fb641e6
-    md5: 8352ce46dfb3cf6f1cc6b118b8b5dde2
+    binhash: e450e9fdae1056617f2363522d41c38e
+    md5: e5c1c2ca736b4d394fdac397c46e7a21
 work/atmosphere/INPUT/PET0.RegridWeightGen.Log:
   fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/PET0.RegridWeightGen.Log                                        
   hashes:
     binhash: ffd3c47e6c1417cc8595e470ddae1fe3
     md5: 08c82ec3a08e47ccc8b2b5b29117d3bc
+work/atmosphere/INPUT/common_1deg_era5/rmp_jrar_to_cict_CONSERV.nc:                                                   
+  fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/common_1deg_era5/rmp_jrar_to_cict_CONSERV.nc                    
+  hashes:
+    binhash: 6563b5e196a4ccb9a189dc533893bbb2
+    md5: 10f0b5ea6b8102a03ad1140e5163a0f7
+work/atmosphere/INPUT/era5_mom10_wts.nc:
+  fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/era5_mom10_wts.nc                                               
+  hashes:
+    binhash: a881ab6afd66ab659457c62abd66a7b9
+    md5: e5c1c2ca736b4d394fdac397c46e7a21
 work/atmosphere/INPUT/lat_flipud/ERA5_MOM1_conserve.nc:
   fullpath: /g/data/e14/rmh561/access-om2/input/ERA-5/lat_flipud/ERA5_MOM1_conserve.nc                                
   hashes:

given the checksum of the updated weights file:

$ md5sum /scratch/v45/aph502/era5_mom10_wts.nc                                     
e5c1c2ca736b4d394fdac397c46e7a21  /scratch/v45/aph502/era5_mom10_wts.nc   

The only place the updated weights file shows up is in work/atmosphere/INPUT.

So you'll need to change

input: /g/data/ik11/inputs/access-om2/input_20210915/common_1deg_era5/

to another location that has the updated remapping weights files.

I guess none of the changed remapping weights were used if the same error was made in your other runs. Maybe those remapping weights files were ok? Possibly worth checking.

As an aside, the ERA5 documentation states:

When you download ERA5 data you have the option to interpolate the data to a custom grid and horizontal resolution (eg. 'grid':'0.5/0.5'). The default interpolation method is bilinear for continuous parameters (e.g. Temperature) and nearest neighbour for discrete parameters (eg. Vegetation).

If you request data in NetCDF format ('format':'netcdf'), interpolation to a regular grid is mandatory, because ECMWF's NetCDF implementation only supports regular grids. Specify the desired lat/lon grid with the 'grid' keyword, for example 'grid':'0.25/0.25'.

the ECMWF interpolation software does not conserve area integrals

So the data we're using isn't conservatively interpolated from the native ERA5 model output, so there seems little point in us doing so, in which case I would recommend changing all those entries in namcouple to ERA5_MOM1_patch.nc.

[rmholmes]

I need a face palm emogy...

Looks like the common_1deg_era5 was needed under both common (for the namcouple remaps) and atmosphere (for the rmp_jrar.. file). I just needed to look more carefully.

Thanks Aidan. I'll get on it.

[aekiss]

Thanks @aidanheerdegen - interesting point about the non-conservative interpolation. The patch interpolation gives much smoother (so more physically reasonable) results, but I wonder whether it's justifiable to add to the conservation error by using a second non-conservative interpolation. Not sure.

[aidanheerdegen]

It is justifiable to use conservation on fields that were never conserved? I think patch is more intellectually honest, and does not give a false impression.

[aekiss]

My thinking here is that ideally we'd want global integrals of fields like precip to be identical to the ERA5 native output to have the same water balance (likewise for radiation). The non-conservative interpolation to NetCDF introduces an error we can't avoid, but an additional patch interpolation could take it even further from native ERA5. I don't have a feel for the size of the non-conservation error in patch though. Maybe it's negligible, in which case the smoothness of patch is certainly appealing. I wince every time I look at the blocky conservative JRA55-do forcing we apply to the 0.1deg model.

[rmholmes]

Ok I now have success:

This is using the patch remapping weights file at /g/data/e14/rmh561/access-om2/input/ERA-5/ERA5_MOM1_patch.nc. I generated this with ../../tools/make_remap_weights.py --accessom2_input_dir ../ --atm_forcing_file /g/data/rt52/era5/single-levels/reanalysis/2t/1980/2t_era5_oper_sfc_19800101-19800131.nc --atm ERA5 --ocean MOM1 --npes 1 in that directory. I used the access-om2 tools code in this branch, that do not include Russ's updated ESMF code described at COSIMA/access-om2#258 - but I don't think that will be an issue until we go to the 1/10th.

@aidanheerdegen your file did not work for me, spitting out some CICE errors which I didn't investigate any further.

All I had to do was remove the np.flipud command in era5_grid.py and the assert statements in base_grid.py that checked the corners were correctly ordered as described in this PR. Is it worth altering these checks to do something sensible when latitude is flipped?

The offset in longitude appears to be an issue only with Nic's original files. The longitude appears correct for mine without any changes needed.

I am now proceeding with a 3 year run to compare with JRA55.