fix ease index in ease pfaf tile file and remove pole tiles#1430
fix ease index in ease pfaf tile file and remove pole tiles#1430weiyuan-jiang wants to merge 5 commits into
Conversation
| do ii=1,np | ||
| if(type_tile(ii)==100)then | ||
| tile_id=map_tile(xi_tile(ii)+1,ny-yi_tile(ii)) ! note: EASE indices are 0-based | ||
| tile_id=map_tile(xi_tile(ii)+1, yi_tile(ii)+1) ! note: EASE indices are 0-based |
There was a problem hiding this comment.
@zyj8881357 : Is there a place elsewhere in the code (e.g., in preprocessing) where the orientation of the latitude indices (j_indg) from the EASE-Pfaf file needs to be flipped?
There was a problem hiding this comment.
Yes. for the file GEOSsurface_GridComp/Utils/Raster/preproc/routing_model/get_num_sub_catchment.f90.
The line 54: yi_tile = yi_tile + 1
should be changed to
yi_tile = nlatE - yi_tile
Also, add nlatE in line 6 as:
use routing_model_constants, only : nc, nlon, nlat, np=>np_tot, nlatE
That is all.
suggested by Yujin
| ! This is a special case for river-routing. The ocean grid is also EASE just for convenience | ||
| logical :: two_EASE | ||
| integer :: ip_keep, k | ||
| real(REAL64), parameter :: EASE_LAT_MAX = 85.04459_REAL64 |
There was a problem hiding this comment.
@weiyuan-jiang : The min/max latitudes for the EASE grid depend on the resolution. For consistency, we should probably get this from MAPL_ease_extent()
| xi_tile = xi_tile + 1 | ||
| call read_ncfile_int1d(trim(file_path1),"j_indg",yi_tile,np) | ||
| yi_tile = yi_tile + 1 | ||
| yi_tile = nlatE - yi_tile |
There was a problem hiding this comment.
@zyj8881357 : I tried to understand why the lat index is flipped here. It looks like yi_tile is needed to get the climatological runoff from the file SMAPL4_NRv11p4_M09_runoff_mean_2001_2023.nc, which is read in get_Qr_clmt.f90
It looks like the lat index orientation in the nc file above is consistent with the EASE convention.
But after reading the data, the lat index is flipped here:
It looks like if we don't flip the lat index at all, the processing should work out. But I didn't inspect every little bit of the code and didn't test it. There may be other processing involved that requires the lat index.
If my hunch is right, though, we should avoid flipping the lat index twice.
A tile file for the EASE grid that respects the Pfafstetter hydrological catchments for land tiles ("EASE-Pfaf" tile file) was introduced recently to support river routing. The original implementation increased the latitude (j) index of the tiles from south to north, which contradicts the EASE grid convention, which is North to South.
This PR fixes the orientation of the EASE grid lat indices in the to be consistent with the EASE grid convention.
The PR is trivially 0-diff for all current GCM tests and most GEOSldas tests.
The PR should be 0-diff for the new LDAS_[*]GLOBAL/model tests that include river routing.
The PR is not 0-diff for make_bcs when bcs on the EASE tile space are generated. In this case, the EASE-Pfaf file is expected to differ in the "j_indg" field.