diff --git a/Documentation_texfol/Graphiques/periodicity.odg b/Documentation_texfol/Graphiques/periodicity.odg
index 36b9ecb48d5d7003dde801837ef641061a348dfa..0aff6989bddf09d78b7c73784e04df28c3d41ba9 100644
Binary files a/Documentation_texfol/Graphiques/periodicity.odg and b/Documentation_texfol/Graphiques/periodicity.odg differ
diff --git a/Documentation_texfol/documentation.tex b/Documentation_texfol/documentation.tex
index 0a29d868e5f2026b308fd13e07ae2ebe94aafdff..9c912b3cf0de3a909004a880656214f7f0acf567 100644
--- a/Documentation_texfol/documentation.tex
+++ b/Documentation_texfol/documentation.tex
@@ -1100,14 +1100,15 @@ Cf. figure (\ref{fig:periodicity}).
\centering
\includegraphics{periodicity}
\caption{Cas de périodicité dans local\_extrema. Le tableau
- représenté est extr\_map. Pour définir la case (1, j), il faut
- avoir copié la case (m, j - 1) (qui peut être non nulle) dans
- local\_map. Pour définir la case (m, j), il faut avoir copié les
- cases (1, j - 1) et (1, j) (dont l'une peut être non nulle) dans
- local\_map.}
+ représenté est extr\_map. Pour définir la case (2, j), il faut
+ avoir copié la case (m - 1, j - 1) (qui peut être non nulle). Pour
+ définir la case (m - 1, j), il faut avoir copié la case (2, j)
+ (qui peut être non nulle). (La case (2, j - 1) doit aussi avoir
+ été copiée mais cela a été fait au passage sur la colonne
+ précédente.)}
\label{fig:periodicity}
\end{figure}
-Cas de périodicité : field(i + m, :) = field(i, :).
+Cas de périodicité : field(i + m - 2, :) = field(i, :).
\subsection{get\_1\_outerm}
diff --git a/Tests/test_get_snapshot.f b/Tests/test_get_snapshot.f
index 15b547571ae50f687d9952bf5b7af1cf5f635d44..cdf6a7fd4fab42a3df062075d470ca1c9ac6ea6d 100644
--- a/Tests/test_get_snapshot.f
+++ b/Tests/test_get_snapshot.f
@@ -3,16 +3,18 @@ program test_get_snapshot
use, intrinsic:: ieee_arithmetic, only: IEEE_SUPPORT_DATATYPE, &
ieee_support_nan
- use derived_types, only: snapshot
- use dispatch_snapshot_m, only: dispatch_snapshot
- use get_snapshot_m, only: get_snapshot
- use init_shapefiles_m, only: init_shapefiles
+ ! Libraries:
use jumble, only: new_unit
use netcdf, only: nf90_nowrite
use netcdf95, only: nf95_open, find_coord, nf95_inquire_dimension, &
nf95_get_var, nf95_close
- use nr_util, only: assert, deg_to_rad
+ use nr_util, only: assert, deg_to_rad, twopi
use shapelib, only: shpfileobject, shpclose
+
+ use derived_types, only: snapshot
+ use dispatch_snapshot_m, only: dispatch_snapshot
+ use get_snapshot_m, only: get_snapshot
+ use init_shapefiles_m, only: init_shapefiles
implicit none
@@ -32,6 +34,9 @@ program test_get_snapshot
! Just the first two values to get the corner and step:
real longitude(2), latitude(2) ! in degrees
+ real step(2) ! longitude and latitude steps, in rad
+ logical periodic ! grid is periodic in longitude
+
!--------------------------------------------------------------
call assert(IEEE_SUPPORT_DATATYPE(0.), ieee_support_nan(0.), &
@@ -52,12 +57,18 @@ program test_get_snapshot
call nf95_get_var(ncid, varid, latitude)
call nf95_close(ncid)
+
+ step = [longitude(2) - longitude(1), latitude(2) - latitude(1)] * deg_to_rad
+
+ ! As we are requiring the grid spacing to be uniform, the value of
+ ! "periodic" must be consistent with the values of step(1) and nlon:
+ periodic = nint(twopi / step(1)) == nlon
+ print *, "periodic = ", periodic
call get_snapshot(s, min_amp, m = 1, n_proc = 1, k_end = 1, max_delta = 4, &
nlon = nlon, nlat = nlat, k = 1, max_radius = [20, 20], &
- corner = [longitude(1), latitude(1)] * deg_to_rad, &
- step = &
- [longitude(2) - longitude(1), latitude(2) - latitude(1)] * deg_to_rad)
+ corner = [longitude(1), latitude(1)] * deg_to_rad, step = step, &
+ periodic = periodic)
call init_shapefiles(hshp_extremum, hshp_outermost, hshp_max_speed)
diff --git a/Tests/test_local_extrema.f b/Tests/test_local_extrema.f
index b5792c0a515f1e7a176c3ec7807c487590c44bac..fe893ffe073894957a1b7726c7ab0ec5571bd495 100644
--- a/Tests/test_local_extrema.f
+++ b/Tests/test_local_extrema.f
@@ -3,6 +3,8 @@ program test_local_extrema
! Reads ADT from a NetCDF file and writes the list of extrema to a
! CSV file. Also creates a NetCDF file containing the map of extrema.
+ ! Note: the grid spacing does not need to be uniform.
+
! Libraries:
use jumble, only: new_unit, get_command_arg_dyn
use netcdf, only: nf90_nowrite, nf90_clobber, NF90_FLOAT, NF90_INT
@@ -21,8 +23,12 @@ program test_local_extrema
character(len = :), allocatable:: filename
real, allocatable:: longitude(:) ! (nlon) in degrees
real, allocatable:: latitude(:) ! (nlat) in degrees
- real, allocatable:: ssh(:, :) ! (nlon, nlat) sea-surface height, in m
- integer, allocatable:: extr_map(:, :) ! (nlon, nlat) map of extrema
+
+ real, allocatable:: ssh(:, :) ! (1 - copy:nlon + copy, nlat)
+ ! sea-surface height, in m
+
+ integer, allocatable:: extr_map(:, :) ! (1 - copy:nlon + copy, nlat)
+ ! map of extrema
integer, allocatable:: ind_extr(:, :) ! (2, n_extr)
! indices in the two dimensions of extrema
diff --git a/Tests/test_max_speed_contour_ssh.f b/Tests/test_max_speed_contour_ssh.f
index c4a40f5158cbeda4bf35b505d98f9ff7e7df7add..10a18029ea8e1f102d8770cf2b5f3305971a9b79 100644
--- a/Tests/test_max_speed_contour_ssh.f
+++ b/Tests/test_max_speed_contour_ssh.f
@@ -2,18 +2,20 @@ program test_max_speed_contour_ssh
use, intrinsic:: ISO_FORTRAN_ENV
- use max_speed_contour_ssh_m, only: max_speed_contour_ssh
+ ! Libraries:
+ use jumble, only: get_command_arg_dyn
use netcdf, only: nf90_nowrite
use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var
use nr_util, only: assert
+ use max_speed_contour_ssh_m, only: max_speed_contour_ssh
+
implicit none
integer:: ilon_llc = 21, ilat_llc = 215 ! lower left corner
integer:: ilon_urc = 49, ilat_urc = 231 ! upper right corner
integer n_lon, n_lat
- integer length, status
character(len = :), allocatable:: adt_file, velocity_file
integer ncid, varid
real, allocatable:: ssh(:, :) ! (n_lon, n_lat) sea-surface height, in m
@@ -28,14 +30,8 @@ program test_max_speed_contour_ssh
call assert(COMMAND_ARGUMENT_COUNT() == 2, &
"Required arguments: ADT-file velocity-file")
-
- call get_command_argument(number = 1, length = length, status = status)
- allocate(character(len = length):: adt_file)
- call get_command_argument(1, adt_file)
-
- call get_command_argument(number = 2, length = length, status = status)
- allocate(character(len = length):: velocity_file)
- call get_command_argument(2, velocity_file)
+ call get_command_arg_dyn(1, adt_file)
+ call get_command_arg_dyn(2, velocity_file)
write(unit = error_unit, nml = main_nml)
write(unit = error_unit, fmt = *) "Enter namelist main_nml."
diff --git a/Tests/test_set_all_outerm.f b/Tests/test_set_all_outerm.f
index dfc52246f429f96a8f7486801a7b110912613552..53f22a905eac82aa475b5f2029997029434652d7 100644
--- a/Tests/test_set_all_outerm.f
+++ b/Tests/test_set_all_outerm.f
@@ -1,17 +1,19 @@
program test_set_all_outerm
+ ! Libraries:
use contour_531, only: null_polyline
- use derived_types, only: snapshot
use jumble, only: new_unit, get_command_arg_dyn
use netcdf, only: nf90_nowrite
use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, &
find_coord, nf95_inquire_dimension, nf95_get_att
use nr_util, only: pi
- use set_all_outerm_m, only: set_all_outerm
use shapelib, only: shpt_point, shpt_polygon, shpfileobject, ftdouble, &
shpclose, ftinteger
use shapelib_03, only: shp_create_03, dbf_add_field_03, shp_append_point_03, &
dbf_write_attribute_03, shp_append_null_03, shp_append_object_03
+
+ use derived_types, only: snapshot
+ use set_all_outerm_m, only: set_all_outerm
implicit none
@@ -82,7 +84,7 @@ program test_set_all_outerm
call set_all_outerm(s, min_amp, max_radius, &
corner = [longitude(1), latitude(1)] * deg_over_rad, step &
= [longitude(2) - longitude(1), latitude(2) - latitude(1)] &
- * deg_over_rad, ssh = ssh)
+ * deg_over_rad, periodic = .false., ssh = ssh)
call shp_create_03("extremum_1", shpt_point, hshp_extremum)
call dbf_add_field_03(ifield, hshp_extremum, 'ssh', ftdouble, nwidth = 13, &
diff --git a/Tests/tests.json b/Tests/tests.json
index 8d6c9f5e012c7c8ecb8c29f3286e48603d3540fd..670164bc1bd95fc72bc1da5597116d7a92ff0fa6 100644
--- a/Tests/tests.json
+++ b/Tests/tests.json
@@ -5,7 +5,8 @@
"$input_dir/example.nc"],
"title" : "Good_contour",
"required": [["$input_dir/outside_points_1.csv", "outside_points.csv"]],
- "description": "3 contours at that level. 2 contain the inside point, one of which contains the 2 outside points."
+ "description":
+ "3 contours at that level. 2 contain the inside point, one of which contains the 2 outside points."
},
{
"stdin_filename" : "$input_dir/good_contour_2.txt",
@@ -13,7 +14,8 @@
"$input_dir/example.nc"],
"title" : "Good_contour_2",
"required": [["$input_dir/outside_points_2.csv", "outside_points.csv"]],
- "description": "Select another good contour. Case where one of the contours tested does not contain inside point."
+ "description":
+ "Select another good contour. Case where one of the contours tested does not contain inside point."
},
{
"stdin_filename" : "$input_dir/no_good_contour.txt",
@@ -70,13 +72,15 @@
},
{
"args" : ["$compil_prod_dir/test_max_speed_contour_ssh",
- "$large_input_dir/h_2006_01_01.nc", "$large_input_dir/uv_2006_01_01.nc"],
+ "$large_input_dir/h_2006_01_01.nc",
+ "$large_input_dir/uv_2006_01_01.nc"],
"title" : "Max_speed_contour_ssh",
"input" : "&main_nml /\n"
},
{
"args" : ["$compil_prod_dir/test_max_speed_contour_ssh",
- "$large_input_dir/h_2006_01_01.nc", "$large_input_dir/uv_2006_01_01.nc"],
+ "$large_input_dir/h_2006_01_01.nc",
+ "$large_input_dir/uv_2006_01_01.nc"],
"title" : "Max_speed_contour_ssh_north",
"stdin_filename" : "$input_dir/max_speed_contour_ssh_nml.txt"
},
diff --git a/derived_types.f b/derived_types.f
index 7790e0561a291429af6f4eca51b014a9de73ffbf..8e1cc019f8aee57aa48e0cac74883597326a857a 100644
--- a/derived_types.f
+++ b/derived_types.f
@@ -42,9 +42,11 @@ module derived_types
integer number_vis_eddies ! number of visible eddies
- integer, allocatable:: extr_map(:, :) ! (nlon, nlat)
- ! At a point of extremum SSH: identification number or this
- ! extremum. 0 at other points.
+ integer, allocatable:: extr_map(:, :)
+ ! (1 - max_radius(1):nlon + max_radius(1), nlat) if the grid is
+ ! periodic in longitude, else (nlon, nlat). At a point of
+ ! extremum SSH: identification number or this extremum. 0 at
+ ! other points.
integer, allocatable:: ind_extr(:, :) ! (2, number_vis_eddies)
! List of coordinates of ssh extrema in the global grid, in index
diff --git a/get_1_outerm.f b/get_1_outerm.f
index 0689f934b7e0f87ad6a223d419076d4848bbb08f..a4d51c503001d85e79570e22012d9c390e6c675b 100644
--- a/get_1_outerm.f
+++ b/get_1_outerm.f
@@ -14,10 +14,20 @@ contains
! at innermost level. If a contour is found then its level is
! farther from the extremum than innermost level.
+ ! 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. The
+ ! longitudes of outside points and the target extremum should be
+ ! in the interval [corner(1), corner(1) + step(1) * (size(ssh, 1)
+ ! - 1)] as there will be no automatic shifting by a mulitple of
+ ! 360 degrees.
+
+ ! Libraries:
use contour_531, only: polyline
+ use nr_util, only: assert
+
use derived_types, only: ssh_contour, missing_ssh
use good_contour_m, only: good_contour
- use nr_util, only: assert
use spherical_polyline_area_m, only: spherical_polyline_area
logical, intent(in):: cyclone
diff --git a/get_snapshot.f b/get_snapshot.f
index fbf9c299c456039d92de79f022a836753b6ba527..8f693f3b0109585cdeda17e49f9f51ad32283b6d 100644
--- a/get_snapshot.f
+++ b/get_snapshot.f
@@ -5,7 +5,7 @@ module get_snapshot_m
contains
subroutine get_snapshot(s, min_amp, m, n_proc, k_end, max_delta, nlon, &
- nlat, k, max_radius, corner, step)
+ nlat, k, max_radius, corner, step, periodic)
use, intrinsic:: ieee_arithmetic, only: ieee_value, IEEE_QUIET_NAN
@@ -44,12 +44,15 @@ contains
! corner of the whole grid, in rad
real, intent(in):: step(:) ! (2) longitude and latitude steps, in rad
+ logical, intent(in):: periodic ! grid is periodic in longitude
! Local:
integer ncid
- real ssh(nlon, nlat) ! sea-surface height, in m
- real u(nlon, nlat), v(nlon, nlat) ! wind, in m s-1
+ real ssh(1 - merge(max_radius(1), 0, periodic):nlon + merge(max_radius(1), &
+ 0, periodic), nlat) ! sea-surface height, in m
+ real, dimension(1 - merge(max_radius(1), 0, periodic):nlon &
+ + merge(max_radius(1), 0, periodic), nlat):: u, v ! wind, in m s-1
integer i
! Window around each extremum:
@@ -76,7 +79,7 @@ contains
! and v when we compute the max-speed contours.)
call nf95_close(ncid)
- call set_all_outerm(s, min_amp, max_radius, corner, step, ssh)
+ call set_all_outerm(s, min_amp, max_radius, corner, step, periodic, ssh)
! Done with outermost contours, now let us take care of the
! max-speed contours.
@@ -88,11 +91,12 @@ contains
llc = floor(convert_to_ind(minval(s%list_vis(i)%out_cont%points, &
dim = 2), corner, step))
- urc = min(ceiling(convert_to_ind(maxval( &
- s%list_vis(i)%out_cont%points, dim = 2), corner, step)), &
- [nlon, nlat])
- ! (min should have no effect except because of roundup error)
-
+ urc = ceiling(convert_to_ind(maxval( &
+ s%list_vis(i)%out_cont%points, dim = 2), corner, step))
+
+ if (.not. periodic) urc = min(urc, [nlon, nlat])
+ ! (should have no effect except because of roundup error)
+
call set_max_speed(s%list_vis(i), s%ind_extr(:, i) - llc + 1, &
nearby_extr(s%extr_map(llc(1):urc(1), llc(2):urc(2)), &
s%list_vis, i), ssh(llc(1):urc(1), llc(2):urc(2)), &
@@ -119,7 +123,7 @@ contains
use netcdf95, only: nf95_inq_varid, nf95_get_var, nf95_get_att
real, intent(out):: values(:, :)
- ! (1 - merge(max_radius, 0, periodic):nlon + merge(max_radius, 0,
+ ! (1 - merge(max_radius(1), 0, periodic):nlon + merge(max_radius(1), 0,
! periodic), nlat) ssh, u or v
character(len = *), intent(in):: name ! of NetCDF variable
@@ -138,6 +142,13 @@ contains
! Change the missing value:
where (values(1:nlon, :) == Fill_Value) values(1:nlon, :) = new_fill_value
+ if (periodic) then
+ ! Extend in longitude:
+ values(1 - max_radius(1):0, :) &
+ = values(nlon + 1 - max_radius(1):nlon, :)
+ values(nlon + 1:nlon + max_radius(1), :) = values(1:max_radius(1), :)
+ end if
+
end subroutine get_var
end subroutine get_snapshot
diff --git a/good_contour.f b/good_contour.f
index 69f41f8f92c4d12d68a2eb0d79af15db2ae3b82c..491220f7e43fe9fe669633385af7333b7e4f3447 100644
--- a/good_contour.f
+++ b/good_contour.f
@@ -13,6 +13,14 @@ contains
! first satisfying contour it finds. If such a contour does not
! exist, the procedure returns a null polyline.
+ ! The length of the interval of longitudes of the domain, step(1)
+ ! * (size(z, 1) - 1), should be < 180 degrees, so that the
+ ! geometrical processing done here is non-ambiguous. The
+ ! longitudes of outside points and inside point should be in the
+ ! interval [corner(1), corner(1) + step(1) * (size(z, 1) - 1)] as
+ ! there will be no automatic shifting by a mulitple of 360
+ ! degrees.
+
use contour_531, only: polyline, find_contours_reg_grid, null_polyline
use geometry, only: polygon_contains_point
diff --git a/local_extrema.f b/local_extrema.f
index 9906317260791251ace1b31337b4cefc4289a405..aba4ecd60f9637de71b3325c886854390ff937ab 100644
--- a/local_extrema.f
+++ b/local_extrema.f
@@ -17,6 +17,9 @@ contains
! in extr_map but there is no duplication in ind_extr,
! innermost_level and local_min.
+ ! Note that no coordinate grid is used here so there is no
+ ! assumption on the grid underlying "field".
+
use nr_util, only: assert_eq
integer, intent(out):: extr_map(:, :) ! (m, n) Map of
diff --git a/mean_speed.f b/mean_speed.f
index d3737069d6dc0c22bf6a3663ff6814775601b564..5823a291426af39266dd76d27e4e8de9cb9b27ae 100644
--- a/mean_speed.f
+++ b/mean_speed.f
@@ -9,6 +9,7 @@ contains
! Interpolates the wind at each point of the input polygon and
! computes the mean azimuthal speed at interpolation points.
+ ! Libraries:
use contour_531, only: polyline
use numer_rec_95, only: bilinear_interp2_reg
diff --git a/nearby_extr.f b/nearby_extr.f
index 45d66bacc554f84f76efdd0d3ff021293db765a5..6f69677f9faceae2d2f03430c88457937f74a35d 100644
--- a/nearby_extr.f
+++ b/nearby_extr.f
@@ -6,6 +6,9 @@ contains
pure function nearby_extr(extr_map, list_vis, i)
+ ! Returns a list of extrema that cannot be engulfed in a good
+ ! contour around the target extremum.
+
use derived_types, only: eddy
real, allocatable:: nearby_extr(:, :) ! (2, :) longitude and
diff --git a/set_all_outerm.f b/set_all_outerm.f
index 5d6f85737f21bd380fc4194399855256278986c8..893d14d23bb45f7616b8a94c2b8d7fe080a4857c 100644
--- a/set_all_outerm.f
+++ b/set_all_outerm.f
@@ -4,7 +4,7 @@ module set_all_outerm_m
contains
- subroutine set_all_outerm(s, min_amp, max_radius, corner, step, ssh)
+ subroutine set_all_outerm(s, min_amp, max_radius, corner, step, periodic, ssh)
! Extraction of eddies: find extrema and set all outermost
! contours in the snapshot "s". Not a function because "s" is not
@@ -13,6 +13,7 @@ contains
! Libraries:
use jumble, only: argwhere
use numer_rec_95, only: indexx
+ use nr_util, only: twopi
use derived_types, only: snapshot
use get_1_outerm_m, only: get_1_outerm
@@ -33,43 +34,67 @@ contains
! corner of the global grid, in rad
real, intent(in):: step(:) ! (2) longitude and latitude steps, in rad
+ logical, intent(in):: periodic ! grid is periodic in longitude
- real, intent(in):: ssh(:, :) ! (nlon, nlat) sea-surface height, in m
+ real, intent(in):: ssh(:, :)
+ ! (1 - max_radius(1):nlon + max_radius(1), nlat) if the grid is periodic
+ ! in longitude, else (nlon, nlat). Sea-surface height, in m.
! Local:
integer nlon, nlat
+
real, allocatable:: innermost_level(:) ! (s%number_vis_eddies)
- ! ssh level of the innermost contour around each extremum, in
+ ! SSH level of the innermost contour around each extremum, in
! m. By construction, innermost_level < extremum for a maximum, >
! extremum for a minimum.
logical, allocatable:: cyclone(:) ! (s%number_vis_eddies)
- integer i, l
+ integer i, l, copy
integer n_cycl ! number of cyclones
real, parameter:: min_area = 8e8
! minimum area of an outermost contour, in m2
integer, allocatable:: sorted_extr(:) ! (s%number_vis_eddies)
- ! identifying number of extrema, first minima sorted in order of
- ! decreasing SSH, followed by maxima sorted in order of increasing
- ! SSH
+ ! Sorted identifying number of extrema: first, minima sorted in
+ ! order of decreasing SSH, and second, maxima sorted in order of
+ ! increasing SSH.
integer, allocatable:: selection(:)
! identifying numbers of a selection of extrema
! Window around each extremum:
+
integer llc(2) ! indices in global grid of lower left corner
integer urc(2) ! indices in global grid of upper right corner
+ real corner_window(2) ! longitude and latitude of the window around each
+ ! extremum, in rad
+
+ real, allocatable:: outside_points(:, :) ! (2, :) longitude and
+ ! latitude, in rad, of all the significant extrema, except the
+ ! target extremum
+
!--------------------------------------------------------------
- nlon = size(ssh, 1)
+ copy = merge(max_radius(1), 0, periodic)
+ nlon = size(ssh, 1) - 2 * copy
nlat = size(ssh, 2)
- allocate(s%extr_map(nlon, nlat))
- call local_extrema(s%extr_map, s%ind_extr, innermost_level, cyclone, ssh, &
- periodic = .false., my_lbound = [1, 1])
+ allocate(s%extr_map(1 - copy:nlon + copy, nlat))
+ copy = merge(1, 0, periodic)
+ call local_extrema(s%extr_map(1 - copy:nlon + copy, :), s%ind_extr, &
+ innermost_level, cyclone, ssh(1 - copy:nlon + copy, :), periodic, &
+ my_lbound = [1 - copy, 1])
+
+ if (periodic) then
+ ! Extend in longitude:
+ s%extr_map(1 - max_radius(1):- 1, :) &
+ = s%extr_map(nlon + 1 - max_radius(1):nlon - 1, :)
+ s%extr_map(nlon + 2:nlon + max_radius(1), :) &
+ = s%extr_map(2:max_radius(1), :)
+ end if
+
s%number_vis_eddies = size(s%ind_extr, 2)
allocate(s%list_vis(s%number_vis_eddies), sorted_extr(s%number_vis_eddies))
@@ -96,19 +121,33 @@ contains
i = sorted_extr(l)
! Define the geographical window around each eddy extremum:
- llc = max(s%ind_extr(:, i) - max_radius, 1)
- urc = min(s%ind_extr(:, i) + max_radius, [nlon, nlat])
+
+ llc = s%ind_extr(:, i) - max_radius
+ urc = s%ind_extr(:, i) + max_radius
+
+ if (.not. periodic) then
+ llc = max(llc, 1)
+ urc = min(urc, [nlon, nlat])
+ end if
+
+ corner_window = corner + (llc - 1) * step
! No need to consider contours with amplitudes < min_amp:
if (abs(s%list_vis(i)%ssh_extr - innermost_level(i)) < min_amp) &
innermost_level(i) = s%list_vis(i)%ssh_extr &
+ merge(min_amp, - min_amp, s%list_vis(i)%cyclone)
+ outside_points = nearby_extr(s%extr_map(llc(1):urc(1), llc(2):urc(2)), &
+ s%list_vis, i)
+
+ if (periodic) outside_points(1, :) = outside_points(1, :) &
+ + ceiling((corner_window(1) - outside_points(1, :)) / twopi) * twopi
+ ! (Shift the longitude of each point to a value close to the
+ ! longitude of the target extremum.)
+
s%list_vis(i)%out_cont = get_1_outerm(s%list_vis(i)%cyclone, &
- s%list_vis(i)%coord_extr, innermost_level(i), &
- nearby_extr(s%extr_map(llc(1):urc(1), llc(2):urc(2)), s%list_vis, &
- i), ssh(llc(1):urc(1), llc(2):urc(2)), &
- corner = corner + (llc - 1) * step, step = step)
+ s%list_vis(i)%coord_extr, innermost_level(i), outside_points, &
+ ssh(llc(1):urc(1), llc(2):urc(2)), corner_window, step)
s%list_vis(i)%valid = s%list_vis(i)%out_cont%area >= min_area
end do
diff --git a/set_max_speed.f b/set_max_speed.f
index 5a458aff05bdbbb5bb16fab361f25fef00659807..42ede51daf38dd3cfb959dbb3760996569e2d138 100644
--- a/set_max_speed.f
+++ b/set_max_speed.f
@@ -9,6 +9,13 @@ contains
! Defines the components speed_cont, max_speed and radius4 of argument e.
+ ! 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. The
+ ! longitudes of outside points should be in the interval
+ ! [corner(1), corner(1) + step(1) * (size(ssh, 1) - 1)] as there
+ ! will be no automatic shifting by a mulitple of 360 degrees.
+
use, intrinsic:: IEEE_ARITHMETIC, only: IEEE_IS_NAN
! Lbraries: