set_contours.f90

module set_contours_m

  implicit none

contains

  subroutine set_contours(out_cont, speed_cont, max_speed, cyclone, extr, &
       innermost_level, step, ssh, u, v, corner, min_area, outside_points)

    ! This procedure sets contours in an eddy.

    ! Libraries:
    use contour_531, only: convert_to_ind

    use config_m, only: min_amp
    use derived_types, only: extremum, ssh_contour
    use get_1_outerm_m, only: get_1_outerm
    use set_max_speed_m, only: set_max_speed

    type(ssh_contour), intent(out):: out_cont ! outermost contour

    type(ssh_contour), intent(out):: speed_cont
    ! contour with maximum average azimuthal speed

    real, intent(out):: max_speed
    ! Average of azimuthal speed on out_cont 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.

    logical, intent(in):: cyclone
    type(extremum), intent(in):: extr

    real, intent(in):: innermost_level
    ! SSH level of the innermost contour around each extremum, in
    ! m. By construction, innermost_level < extremum for a maximum, >
    ! extremum for a minimum.

    real, intent(in):: step(:) ! (2) longitude and latitude steps, in rad

    real, intent(in), dimension(:, :):: ssh, u, v
    ! sea-surface height, in m, and wind, in m s-1

    real, intent(in):: corner(:) ! (2) longitude and latitude of the
    ! corner of the grid, in rad

    real, intent(in):: min_area ! minimum area of an outermost contour, in m2

    real, intent(in):: outside_points(:, :) ! (2, :) longitude and
    ! latitude, in rad, of all the significant extrema, except the
    ! target extremum

    ! Local:

    real innermost_level_2

    ! Window around each extremum:

    integer llc(2) ! indices in grid of lower left corner
    integer urc(2) ! indices in grid of upper right corner

    real corner_window(2) ! longitude and latitude of the window
    ! around each extremum, in rad

    integer, parameter:: n_max_cont = 31 ! must be >= 3

    type(ssh_contour) cont_list(n_max_cont)
    ! Contour list. Defined only for subscripts 1:n_cont. If n_cont >=
    ! 1, the outermost contour is element with subscript n_cont. The
    ! contours are in monotonic order of ssh. If i /= n_cont - 1,
    ! cont_list(i)%area has its default initialization value and the
    ! order of points in cont_list(i)%polyline (clockwise or
    ! counter-clockwise) is not specified.

    integer n_cont
    ! number of good contours found and stored in cont_list, 0 <=
    ! n_cont <= n_max_cont

    !--------------------------------------------------------------

    ! Only look for good contours with amplitudes >= min_amp:
    innermost_level_2 = merge(extr%ssh + merge(min_amp, - min_amp, cyclone), &
         innermost_level, abs(extr%ssh - innermost_level) < min_amp)

    call get_1_outerm(out_cont, cont_list, n_cont, cyclone, extr%coord, &
         innermost_level_2, outside_points, ssh, corner, step, min_area)

    ! Done with outermost contour, now let us take care of the
    ! max-speed contour.
    if (out_cont%closed) then
       ! Restrict the field to the outermost contour:

       llc = floor(convert_to_ind(minval(out_cont%points, dim = 2), corner, &
            step))
       urc = ceiling(convert_to_ind(maxval(out_cont%points, dim = 2), corner, &
            step))

       ! Should have no effect except because of roundup error:
       urc(2) = min(urc(2), size(ssh, 2))
       urc(1) = min(urc(1), size(ssh, 1))

       corner_window = corner + (llc - 1) * step

       ! Done restricting field

       call set_max_speed(speed_cont, max_speed, cont_list, n_cont, &
            extr%coord, ssh(llc(1):urc(1), llc(2):urc(2)), &
            u(llc(1):urc(1), llc(2):urc(2)), &
            v(llc(1):urc(1), llc(2):urc(2)), corner_window, step)
    end if

  end subroutine set_contours

end module set_contours_m