module set_max_speed_m
implicit none
contains
subroutine set_max_speed(speed_cont, max_speed, cont_list, cont_list_proj, &
n_cont, extr_coord, extr_coord_proj, ssh, u, v, corner, step)
! This procedure defines speed_cont and max_speed. On return,
! speed_cont may be a null ssh contour and max_speed may be set to
! missing_speed. max_speed is never NaN.
! The length of the interval of longitudes of the domain, step(1)
! * (size(ssh, 1) - 1), should be < 180 degrees, so that the
! geometrical processing done here is non-ambiguous.
use, intrinsic:: IEEE_ARITHMETIC, only: IEEE_IS_NAN, ieee_value, &
IEEE_QUIET_NAN
! Libraries:
use jumble, only: rad_to_deg
use shapelib, only: shpt_polygon
use shapelib_03, only: shp_append_object_03, dbf_write_attribute_03
use contour_531, only: convert_to_reg_coord, polyline
use ccw_orient_m, only: ccw_orient
use complete_ssh_m, only: complete_ssh
use derived_types, only: eddy, null_ssh_contour, ssh_contour, missing_speed
use good_contour_m, only: good_contour
use input_ssh_m, only: corner_whole => corner
use mean_speed_m, only: mean_speed
use spher_polyline_area_m, only: spher_polyline_area
use cont_list_m, only: write_cont_list, hshp_cont_list, &
hshp_cont_list_proj, ifield_ssh, ifield_speed
type(ssh_contour), intent(out):: speed_cont
real, intent(out):: max_speed
! Average of azimuthal speed on cont_list(i_outer) or speed_cont,
! in m s-1. Positive means counterclockwise circulation. If
! speed_cont is not a null ssh contour then max_speed is the
! speed on speed_cont.
type(polyline), intent(inout):: cont_list(:) ! (n_max_cont)
! Contour list, with points given by longitude and latitude, in
! radians. Defined only for subscripts 1:n_cont. On entry, the
! outermost contour is element with subscript n_cont. The order of
! points in each contour (clockwise or counter-clockwise) is not
! specified.
type(ssh_contour), intent(inout):: cont_list_proj(:) ! (n_max_cont)
! Contour list, with points given by projection
! coordinates. Defined only for subscripts 1:n_cont. On entry, the
! outermost contour is element with subscript n_cont and the
! contours are in monotonic order of ssh. The contour with a given
! subscript in cont_list_proj corresponds to the contour with the
! same subscript in cont_list. More precisely, the point with
! subscript j in contour i of cont_list_proj corresponds to the
! point with subscript j in contour i of cont_list (this is
! important for mean_speed). For i <= n_cont,
! cont_list_proj(i)%area has its default initialization value and
! the order of points in cont_list_proj(i) (clockwise or
! counter-clockwise) is not specified.
integer, intent(inout):: n_cont
! Number of good contours found and stored in cont_list. On entry,
! 1 <= n_cont <= n_max_cont - 1 (not 0 because set_max_speed is
! only called when there is an outermost contour). On return, 1 <=
! n_cont <= n_max_cont.
real, intent(in):: extr_coord(:)
! (2) longitude and latitude of extremum , in rad
real, intent(in):: extr_coord_proj(:) ! (2)
! coordinates of extremum in projection space
real, intent(in):: ssh(:, :), u(:, :), v(:, :)
! The domain should be the bounding box of out_cont, because we do
! not exclude other extrema.
real, intent(in):: corner(:) ! (2)
! Coordinates in projection plane of the lower left corner of the
! grid, corresponding to ssh(1, 1). A priori, this is not the
! corner of the grid of the NetCDF file.
real, intent(in):: step(:) ! (2) longitude and latitude steps, in rad
! Local:
real, allocatable:: speed(:) ! (n_cont) speed on the contour
integer i, i_outer, ishape
!---------------------------------------------------------------
i_outer = n_cont
if (n_cont >= 2) call complete_ssh(cont_list_proj%ssh, n_cont)
! Now find the contours associated to the new values of SSH:
do i = i_outer + 1, n_cont
cont_list_proj(i)%polyline = good_contour(corner, ssh, &
cont_list_proj(i)%ssh, extr_coord_proj)
cont_list(i) = convert_to_reg_coord(cont_list_proj(i)%polyline, &
corner_whole, step)
end do
allocate(speed(n_cont))
do i = 1, n_cont
if (cont_list(i)%closed) then
speed(i) = mean_speed(u, v, cont_list(i), cont_list_proj(i)%points, &
extr_coord, corner)
else
! It is possible that good_contour returned a null polyline
! if there is a missing value of SSH inside the outermost
! contour.
speed(i) = ieee_value(0., IEEE_QUIET_NAN)
end if
end do
if (write_cont_list) then
do i = 1, n_cont
call shp_append_object_03(ishape, hshp_cont_list_proj, &
shpt_polygon, cont_list_proj(i)%points)
call dbf_write_attribute_03(hshp_cont_list_proj, ishape, ifield_ssh, &
fieldvalue = cont_list_proj(i)%ssh)
call shp_append_object_03(ishape, hshp_cont_list, shpt_polygon, &
cont_list(i)%points * rad_to_deg)
if (ieee_is_nan(speed(i))) then
call dbf_write_attribute_03(hshp_cont_list, ishape, ifield_speed, &
missing_speed)
! (Cannot write NaN to dbf file.)
else
call dbf_write_attribute_03(hshp_cont_list, ishape, ifield_speed, &
speed(i))
end if
end do
end if
i = maxloc(abs(speed), dim = 1, mask = .not. IEEE_IS_NAN(speed))
! (The speed may also be NaN when we are near land.)
if (i == 0) then
! All speed values are NaN
speed_cont = null_ssh_contour()
max_speed = missing_speed
else
max_speed = speed(i)
if (i == i_outer) then
! Maximum speed on the outermost contour
speed_cont = null_ssh_contour()
else
speed_cont%polyline = cont_list(i)
speed_cont%ssh = cont_list_proj(i)%ssh
speed_cont%area = spher_polyline_area(speed_cont%polyline)
call ccw_orient(speed_cont)
end if
end if
end subroutine set_max_speed
end module set_max_speed_m