Move the computation of nearby_extr out of procedure get_1_outerm, to

procedure set_all_outerm. The idea is that, in this branch,
nearby_extr will eventually include extrema depending not only on them
having sufficient amplitude but also on their cyclonicity. Both pieces
of information cannot be simply encoded in extr_map.

GNUmakefile

@@ -48,8 +48,7 @@ execut = test_good_contour test_inside_4 test_get_1_outerm test_local_extrema te
 
 # 3. Compiler-dependent part
 
-mode = debug
-include ${general_compiler_options_dir}/${FC}_${mode}.mk
+include ${general_compiler_options_dir}/settings.mk
 
 ifeq (${FC},gfortran)
 # gfortran bug:
@@ -59,7 +58,6 @@ endif
 # 4. Rules
 
 all: ${execut} log
-include ${general_compiler_options_dir}/settings.mk
 test_get_1_outerm: ${obj_test_get_1_outerm}
 test_set_max_speed: ${obj_test_set_max_speed}
 test_max_speed_contour_ssh: max_speed_contour_ssh.o

Tests/test_get_1_outerm.f

@@ -2,16 +2,20 @@ program test_get_1_outerm
 
   use, intrinsic:: ISO_FORTRAN_ENV
 
-  use derived_types, only: ssh_contour
+  ! Libraries:
+  use contour_531, only: convert_to_reg_coord
+  use jumble,  only: argwhere
   use netcdf, only: nf90_nowrite
   use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, &
        nf95_gw_var
   use nr_util, only: pi
-  use get_1_outerm_m, only: get_1_outerm
   use shapelib, only: shpt_polygon, shpfileobject, ftdouble, shpclose
   use shapelib_03, only: shp_create_03, dbf_add_field_03, &
        shp_append_object_03, dbf_write_attribute_03
 
+  use derived_types, only: ssh_contour
+  use get_1_outerm_m, only: get_1_outerm
+
   implicit none
 
   integer ncid, varid
@@ -30,6 +34,10 @@ program test_get_1_outerm
   TYPE(shpfileobject) shphandle
   integer field_number, shape_number
 
+  real, allocatable:: nearby_extr(:, :) ! (2, :) longitude and
+  ! latitude, in rad, of all the significant extrema, except the
+  ! target extremum
+
   namelist /main_nml/ ind_targ_extr, innermost_level, cyclone
 
   !----------------------------------------------------------------
@@ -69,11 +77,13 @@ program test_get_1_outerm
   call nf95_close(ncid)
 
   step = [longitude(2) - longitude(1), latitude(2) - latitude(1)]
+  nearby_extr = convert_to_reg_coord(argwhere(extr_map > 0 .and. extr_map &
+       /= extr_map(ind_targ_extr(1), ind_targ_extr(2))), corner &
+       = [longitude(1), latitude(1)], step = step)
   outermost_contour = get_1_outerm(cyclone, &
        coord_extr = [longitude(1), latitude(1)] + (ind_targ_extr - 1) * step, &
-       i = extr_map(ind_targ_extr(1), ind_targ_extr(2)), &
-       innermost_level = innermost_level, extr_map = extr_map, ssh = ssh, &
-       corner = [longitude(1), latitude(1)], step = step)
+       innermost_level = innermost_level, nearby_extr = nearby_extr, &
+       ssh = ssh, corner = [longitude(1), latitude(1)], step = step)
 
   if (outermost_contour%closed) then
      print *, "Radius of disk of equal area: ", &

Tests/test_get_snapshot.f

@@ -11,7 +11,7 @@ program test_get_snapshot
   use netcdf, only: nf90_nowrite
   use netcdf95, only: nf95_open, find_coord, nf95_inquire_dimension, &
        nf95_get_var, nf95_close
-  use nr_util, only: pi, assert, deg_to_rad
+  use nr_util, only: assert, deg_to_rad
   use shapelib, only: shpfileobject, shpclose
 
   implicit none

depend.mk

 dispatch_snapshot.o : write_eddy.o send_snapshot.o derived_types.o 
-get_1_outerm.o : spherical_polyline_area.o outermost_possible_level.o good_contour.o derived_types.o 
+get_1_outerm.o : spherical_polyline_area.o good_contour.o derived_types.o 
 get_snapshot.o : set_all_outerm.o set_max_speed.o receive_snapshot.o derived_types.o 
 read_eddy.o : derived_types.o 
 read_snapshot.o : read_eddy.o derived_types.o 

get_1_outerm.f

@@ -4,8 +4,8 @@ module get_1_outerm_m
 
 contains
 
-  type(ssh_contour) function get_1_outerm(cyclone, coord_extr, i, &
-       innermost_level, extr_map, ssh, corner, step)
+  type(ssh_contour) function get_1_outerm(cyclone, coord_extr, &
+       innermost_level, nearby_extr, ssh, corner, step)
 
     ! Gets one outermost contour, if it can find one. Method:
     ! dichotomy on level of ssh. If the procedure cannot find an
@@ -13,27 +13,24 @@ contains
     ! contour is not found if, and only if, there is no good contour
     ! at innermost level.
 
-    use contour_531, only: polyline, convert_to_reg_coord
+    use contour_531, only: polyline
     use derived_types, only: ssh_contour, missing_ssh
     use good_contour_m, only: good_contour
-    use jumble, only: argwhere, pack_indices
     use nr_util, only: assert
-    use outermost_possible_level_m, only: outermost_possible_level
     use spherical_polyline_area_m, only: spherical_polyline_area
 
     logical, intent(in):: cyclone
     real, intent(in):: coord_extr(:) ! (2)
 
-    integer, intent(in):: i
-    ! identifying number of the target extremum, that is, the extremum
-    ! at coord_extr
-
     real, intent(in):: innermost_level
     ! ssh level of the innermost contour around the target extremum,
     ! in m. By construction, innermost_level < extremum for a maximum,
     ! > extremum for a minimum.
 
-    integer, intent(in):: extr_map(:, :)
+    real, intent(in):: nearby_extr(:, :) ! (2, :) longitude and
+    ! latitude, in rad, of all the significant extrema, except the
+    ! target extremum
+
     real, intent(in):: ssh(:, :) ! in m
 
     real, intent(in):: corner(:) ! (2)
@@ -44,10 +41,6 @@ contains
 
     ! Local:
 
-    real, allocatable:: nearby_extr(:, :) ! (2, :) longitude and
-    ! latitude, in rad, of all the significant extrema, except the
-    ! target extremum
-
     real level_try, level_good, level_bad ! in m
     type(polyline) tentative_contour
     real, parameter:: precision = 5e-4 ! in m
@@ -55,8 +48,6 @@ contains
 
     !-----------------------------------------------------------------
 
-    nearby_extr = convert_to_reg_coord(argwhere(extr_map > 0 &
-         .and. extr_map /= i), corner, step)
     get_1_outerm%polyline = good_contour(corner, step, ssh, innermost_level, &
          coord_extr, nearby_extr)
 

set_all_outerm.f

@@ -10,8 +10,10 @@ contains
     ! contours in snapshot. Not a function because snapshot is not
     ! completely defined.
 
+    use contour_531, only: convert_to_reg_coord
     use derived_types, only: snapshot
     use get_1_outerm_m, only: get_1_outerm
+    use jumble, only: argwhere
     use local_extrema_m, only: local_extrema
 
     type(snapshot), intent(out):: s
@@ -55,6 +57,10 @@ contains
     real, allocatable:: corner_window(:, :) ! (2, s%number_vis_eddies)
     ! longitude and latitude, in rad
 
+    real, allocatable:: nearby_extr(:, :) ! (2, :) longitude and
+    ! latitude, in rad, of all the significant extrema, except the
+    ! target extremum
+
     !--------------------------------------------------------------
 
     nlon = size(ssh, 1)
@@ -86,9 +92,12 @@ contains
 
     do i = 1, s%number_vis_eddies
        if (flat_extr(i)) then
+          nearby_extr &
+               = convert_to_reg_coord(argwhere(s%extr_map(llc(1, i):urc(1, i), &
+               llc(2, i):urc(2, i)) > 0 .and. s%extr_map(llc(1, i):urc(1, i), &
+               llc(2, i):urc(2, i)) /= i), corner_window(:, i), step)
           s%list_vis(i)%out_cont = get_1_outerm(s%list_vis(i)%cyclone, &
-               s%list_vis(i)%coord_extr, i, innermost_level(i), &
-               s%extr_map(llc(1, i):urc(1, i), llc(2, i):urc(2, i)), &
+               s%list_vis(i)%coord_extr, innermost_level(i), nearby_extr, &
                ssh(llc(1, i):urc(1, i), llc(2, i):urc(2, i)), &
                corner_window(:, i), step)
           if (s%list_vis(i)%out_cont%n_points == 0) then
@@ -123,9 +132,12 @@ contains
     do i = 1, s%number_vis_eddies
        if (s%list_vis(i)%suff_amp .and. noise_around(i) &
             .or. .not. flat_extr(i)) then
+          nearby_extr &
+               = convert_to_reg_coord(argwhere(s%extr_map(llc(1, i):urc(1, i), &
+               llc(2, i):urc(2, i)) > 0 .and. s%extr_map(llc(1, i):urc(1, i), &
+               llc(2, i):urc(2, i)) /= i), corner_window(:, i), step)
           s%list_vis(i)%out_cont = get_1_outerm(s%list_vis(i)%cyclone, &
-               s%list_vis(i)%coord_extr, i, innermost_level(i), &
-               s%extr_map(llc(1, i):urc(1, i), llc(2, i):urc(2, i)), &
+               s%list_vis(i)%coord_extr, innermost_level(i), nearby_extr, &
                ssh(llc(1, i):urc(1, i), llc(2, i):urc(2, i)), &
                corner_window(:, i), step)
           s%list_vis(i)%suff_amp = s%list_vis(i)%out_cont%n_points /= 0