Move Sources/Tests to Tests.

test_good_contour has more general input.

test_inside tests only inside, no loop on radius.

GNUmakefile

@@ -3,7 +3,7 @@
 # 1. Source files and libraries
 
 makefile_dir = .
-VPATH = ${makefile_dir}/Sources ${makefile_dir}/Sources/Tests
+VPATH = ${makefile_dir}/Sources ${makefile_dir}/Tests
 
 src_test_local_extrema = test_local_extrema.f local_extrema.f
 

Sources/set_max_speed.f

@@ -24,7 +24,11 @@ contains
     
     integer, intent(in):: ind_targ_extr(:) ! (2) indices of the target extremum
     integer, intent(in):: extr_map(:, :)
+
     real, intent(in):: ssh(:, :), u(:, :), v(:, :)
+    ! The domain is allowed to be as small as the bounding box of
+    ! outermost_contour. If the domain is larger, this subroutine will
+    ! work too, but with wasted computation.
 
     real, intent(in):: corner(:) ! (2)
     ! longitude and latitude of the lower left corner of the grid, in rad

Sources/set_outermost_contour.f

@@ -7,6 +7,8 @@ contains
   subroutine set_outermost_contour(e, ind_targ_extr, innermost_level, &
        extr_map, ssh, corner, step, noise_around)
 
+    ! Defines e%outermost_contour. Method: dichotomy on level of ssh.
+
     use contour_531, only: polyline, convert_to_reg_coord
     use derived_types, only: eddy
     use good_contour_m, only: good_contour

Sources/Tests/test_get_snapshot.py → Tests/test_get_snapshot.py

@@ -9,7 +9,7 @@ corner = [5.125, -36.375]
 step = 0.25
 longitude = np.arange(29) * step + corner[0]
 latitude = np.arange(17) * step + corner[1]
-lon2d, lat_2d = np.meshgrid(longitude, latitude)
+lon_2d, lat_2d = np.meshgrid(longitude, latitude)
 
 reader_out = shapefile.Reader("outermost_contour_1")
 reader_m_s = shapefile.Reader("max_speed_contour_1")
@@ -30,7 +30,7 @@ for shape_out, shape_m_s in zip(reader_out.iterShapes(),
         plt.plot(points[:, 0], points[:, 1], color = lines[0].get_color(),
                  marker = marker)
 
-plt.plot(lon2d.reshape(-1), lat_2d.reshape(-1), "+", color="gray")
+plt.plot(lon_2d.reshape(-1), lat_2d.reshape(-1), "+", color="gray")
 plt.legend()
 
 reader = shapefile.Reader("extremum_1")

Sources/Tests/test_good_contour.f → Tests/test_good_contour.f

@@ -11,19 +11,29 @@ program test_good_contour
   type(polyline) c
   real x_min, y_min, step_x, step_y, level
   real, allocatable:: Z(:, :)
+  real:: inside_point(2) = [- 0.7, 0.]
   real outside_points(2, 2)
   TYPE(shpfileobject) shphandle
   integer field_number, shape_number
-  
+
+  namelist /main_nml/ inside_point, outside_points
+
   !-------------------------------------------------------------------
 
-  call define_example_reg(x_min, y_min, step_x, step_y, z)
-  level = 5.
+  ! Default values:
   outside_points(:, 1) = [- 0.5, - 1.]
   outside_points(:, 2) = [0.3, 0.]
-  c = good_contour([x_min, y_min], [step_x, step_y], z, level, [- 0.7, 0.], &
+
+  write(unit = *, nml = main_nml)
+  print *, "Enter namelist main_nml."
+  read(unit = *, nml = main_nml)
+  write(unit = *, nml = main_nml)
+
+  call define_example_reg(x_min, y_min, step_x, step_y, z)
+  level = 5.
+  c = good_contour([x_min, y_min], [step_x, step_y], z, level, inside_point, &
        outside_points)
-  
+
   if (c%n_points /= 0) then
      call shp_create_03("test_good_contour", shpt_polygon, shphandle)
      call dbf_add_field_03(field_number, shphandle, 'level', ftdouble, &
@@ -37,5 +47,5 @@ program test_good_contour
   else
      print *, "No good contour found."
   end if
-  
+
 end program test_good_contour

Sources/Tests/test_inside.f → Tests/test_inside.f

 program test_inside
 
   use inside_m, only: inside
+  use jumble, only: get_command_arg_dyn
   use netcdf, only: nf90_nowrite
   use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_gw_var
 
   implicit none
 
+  character(len = :), allocatable:: filename
   integer ncid, varid
   real, allocatable:: ssh(:, :, :) ! sea-surface height, in m
-  integer radius
+  integer:: radius = 1
+  integer:: ind_extr(2) = [37, 218]
+  real:: ssh_outermost = 0.640455
+  logical:: cyclone = .false.
+  logical i
+
+  namelist /main_nml/ radius, ind_extr, ssh_outermost, cyclone
 
   !-----------------------------------------------------------------------
 
-  call nf95_open("h_2011_01_12.nc", nf90_nowrite, ncid)
+  write(unit = *, nml = main_nml)
+  print *, "Enter namelist main_nml."
+  read(unit = *, nml = main_nml)
+  write(unit = *, nml = main_nml)
+
+  call get_command_arg_dyn(1, filename)
+  print *, "Reading from ", filename, "..."
+  call nf95_open(filename, nf90_nowrite, ncid)
 
   call nf95_inq_varid(ncid, "adt", varid)
   call nf95_gw_var(ncid, varid, ssh)
 
   call nf95_close(ncid)
-
-  radius = 1
-
-  do while(inside(ssh(:, :, 1), radius, ind_extr = [37, 218], &
-       ssh_outermost = 0.640455, cyclone = .false.))
-     radius = radius + 1
-  end do
-
-  print *, "radius = ", radius
+  i = inside(ssh(:, :, 1), radius, ind_extr, ssh_outermost, cyclone)
+  print *, "inside: ", i
 
 end program test_inside

Sources/Tests/test_inside.py → Tests/test_inside.py

@@ -3,6 +3,9 @@
 import netCDF4
 import matplotlib.pyplot as plt
 import numpy as np
+import sys
+
+if len(sys.argv) != 2: sys.exit("Required argument: ADT-file")
 
 # Lower left corner (llc):
 ilon_llc = 24
@@ -11,7 +14,7 @@ ilat_llc = 207
 ilon_urc = 50
 ilat_urc = 223
 
-f= netCDF4.Dataset("h_2011_01_12.nc")
+f= netCDF4.Dataset(sys.argv[1])
 longitude = f.variables["lon"][ilon_llc:ilon_urc + 1]
 latitude = f.variables["lat"][ilat_llc:ilat_urc + 1]
 ssh = f.variables["adt"][0, ilat_llc:ilat_urc + 1, ilon_llc:ilon_urc + 1]

Sources/Tests/test_local_extrema.f → Tests/test_local_extrema.f

@@ -47,6 +47,7 @@ program test_local_extrema
   allocate(longitude(nlon), latitude(nlat))
   allocate(ssh(nlon, nlat), extr_map(nlon, nlat))
   
+  print *, "Reading from ", filename, "..."
   call nf95_open(filename, nf90_nowrite, ncid)
 
   call nf95_inq_varid(ncid, "lon", varid)