Skip to content

quad_utils_mod regional grids that cross 0 longitude. #1124

@hkershaw-brown

Description

@hkershaw-brown

☠️ skeleton.

quad_utils_mod for regional grids that cross 0 longitude I do not think is correct.

! for a global grid, the initial values have already been set in
! the derived type. otherwise, find the min/max of lons and lats.
if (.not. h%opt%global_grid) then
h%ii%min_lon = minval(h%ii%lons_2d)
h%ii%max_lon = maxval(h%ii%lons_2d)
h%ii%lon_width = h%ii%max_lon - h%ii%min_lon ! FIXME: wrap?
if (h%ii%lon_width < 0) then
if(h%opt%spans_lon_zero) then
h%ii%lon_width = h%ii%lon_width + 360.0_r8
else
write(string1,*)'min_lon, max_lon, lon_width, spans_lon_zero: ', &
h%ii%min_lon, h%ii%max_lon, h%ii%lon_width, h%opt%spans_lon_zero
call error_handler(E_ERR,routine,'regional grid with bad longitudes', &
source, revision, revdate, text2=string1)
endif
endif
h%ii%min_lat = minval(h%ii%lats_2d)
h%ii%max_lat = maxval(h%ii%lats_2d)
h%ii%lat_width = h%ii%max_lat - h%ii%min_lat
endif

For lon -10 to 10 (350, 10]

max_lon = 350, min_lon = 10

h%ii%min_lon = minval(h%ii%lons_2d)
h%ii%max_lon = maxval(h%ii%lons_2d)
h%ii%lon_width = h%ii%max_lon - h%ii%min_lon 

lon_width = 350 - 10 = 340 rather than 20 ( -10 to 10 )

spans_lon_zero is being used as equivalent to cyclic which it is not.
spans_lon_zero is being used to find the pole (even if there is no pole), then fail boxes that go through the "pole" (not a pole).

! poles? span?
cyclic = h%opt%spans_lon_zero
pole = h%opt%pole_wrap

if (cyclic) then
! Begin by finding the quad that contains the pole for the dipole t_grid.
! To do this locate the u quad with the pole on its right boundary. This is on
! the row that is opposite the shifted pole and exactly follows a lon circle.
pole_x = nx / 2;
! Search for the row at which the longitude flips over the pole
pole_row_lon = h%ii%lons_2d(pole_x, 1);
do i = 1, ny
pindex = i
if(h%ii%lons_2d(pole_x, i) /= pole_row_lon) exit
enddo
! Pole boxes for u have indices pole_x or pole_x-1 and index - 1;
! (it's right before the flip).
u_pole_y = pindex - 1;
! Locate the T dipole quad that contains the pole.
! We know it is in either the same lat quad as the u pole edge or one higher.
! Figure out if the pole is more or less than halfway along
! the u quad edge to find the right one.
if(h%ii%lats_2d(pole_x, u_pole_y) > h%ii%lats_2d(pole_x, u_pole_y + 1)) then
t_pole_y = u_pole_y;
else
t_pole_y = u_pole_y + 1;
endif
endif

Which is then used to fail if "point is in one of the U boxes that go through the pole".

! Fail if point is in one of the U boxes that go through the
! pole (this could be fixed up if necessary)
if (lat_bot == u_pole_y .and. &
(lon_bot == pole_x -1 .or. lon_bot == pole_x)) then
istatus = 4
return
endif

Note I think there is some confusion in the code with u & t grids looks like used to deal with both. I think T (scalar) interpolation is being skipped at the pole in the above code. I will probably put this in as a separate issue since it affects global as well as regional grids.

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type
    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions