Create `loop_shp_tr_append.py`

Replace `loop_shp_append.py` by `loop_shp_tr_append.py`. Instead of
looping on years, we assume we are working on a single year. We could
call this script inside a loop on years in the shell. The advantage is
that we need to make fewer assumptions on the arborescence of
files. Also, we can make `loop_shp_tr_append.py` more general by
starting and ending at arbitrary days. We drop the timing that was in
`loop_shp_append.py`. We add the creation of file `ishape_last.txt`.

In script `shp_tr_append.sh`, assume that the directories contain
files `ishape_last.txt` and modify `ishape_last.txt` in the
destination directory.

In program `test_overlap`, read grid description from a file
`grid_nml.txt` in SHP triplet directory, instead of from standard
input. This is the same modification we made for program
`eddy_graph`. Modify input files accordingly.

README.md

 This code implements ideas and methods described in the following articles:
-[Chaigneau et al., 2009](https://doi.org/10.1016/j.pocean.2009.07.012)
-[Chaigneau et al., 2011](https://doi.org/10.1029/2011JC007134)
-[Pegliasco et al., 2015](https://doi.org/10.1002/2015jc010950)
-[Laxenaire et al., 2018](https://doi.org/10.1029/2018JC014270)
+[Chaigneau et al., 2009](https://doi.org/10.1016/j.pocean.2009.07.012);
+[Chaigneau et al., 2011](https://doi.org/10.1029/2011JC007134);
+[Pegliasco et al., 2015](https://doi.org/10.1002/2015jc010950);
+[Laxenaire et al., 2018](https://doi.org/10.1029/2018JC014270).

Tests/Input/Global/Extraction_eddies_periodic/SHP_triplet/grid_nml.txt

+&grid_nml CORNER_deg  = 0.000000       , -81.00000      ,
+ STEP_deg    = 20.00000       , 18.00000       ,
+ NLON    = 18,
+ NLAT    = 10,
+/

Tests/Input/Global/NSO_in.txt

-&MAIN_NML CORNER_DEG=0.125,-89.875, NLON    =1440 NLAT    =720,
-max_delta = 3
- k_test_1 = 20454, k_test_2 = 20455
-/
+&MAIN_NML max_delta = 3, k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/Global/successive_overlap_in.txt

-&MAIN_NML CORNER_DEG=0.125,-89.875, NLON    =1440 NLAT    =720
- k_test_1 = 20454, k_test_2 = 20455
-/
+&MAIN_NML k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/NSO_delta_in.txt

- &MAIN_NML
- CORNER_DEG      = 5.875000       , -32.12500      ,
- NLON    = 26,
- NLAT    = 22,
- MAX_DELTA       = 2
- k_test_1 = 20454, k_test_2 = 20455
- /
+ &MAIN_NML MAX_DELTA       = 2, k_test_1 = 20454, k_test_2 = 20455 /
 ,,,,,1
 ,,,,,,1

Tests/Input/Region_2/successive_overlap_in.txt

-&MAIN_NML CORNER_DEG = 1.625, -38.375, NLON = 57, NLAT = 33
- k_test_1 = 20454, k_test_2 = 20455
-/
+&MAIN_NML k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/Region_4/Eddies_region_4/grid_nml.txt

+&grid_nml CORNER_deg  = 16.12500       , -38.87500      ,
+ NLON    = 20,
+ NLAT    = 20
+/

Tests/Input/Region_4/NSODS_in.txt

-&MAIN_NML
- CORNER_deg  = 16.12500       , -38.87500      ,
- NLON    = 20,
- NLAT    = 20, max_delta = 3
- k_test_1 = 20454, k_test_2 = 20455
- /
+&MAIN_NML max_delta = 3, k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/Region_4/non_successive_overlap_in.txt

-&MAIN_NML
- CORNER_deg  = 16.12500       , -38.87500      ,
- NLON    = 20,
- NLAT    = 20, max_delta = 3
- k_test_1 = 20454, k_test_2 = 20454
- /
+&MAIN_NML max_delta = 3, k_test_1 = 20454, k_test_2 = 20454 /
 /
 /

Tests/Input/Region_5/NSO_in.txt

- &MAIN_NML
- CORNER_deg  = 0.12500       , -59.875,
- NLON    = 120,
- NLAT    = 180, max_delta = 3
- k_test_1 = 20454, k_test_2 = 20455
- /
+ &MAIN_NML max_delta = 3, k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/Region_5/SHP_triplet/grid_nml.txt

+&grid_nml CORNER_deg  = 0.12500       , -59.875,
+ NLON    = 120,
+ NLAT    = 180
+/

Tests/Input/Region_5/SOR5OD_in.txt

- &MAIN_NML
- CORNER_deg  = 0.12500       , -59.875,
- NLON    = 120,
- NLAT    = 180
- k_test_1 = 20454, k_test_2 = 20454
- /
+ &MAIN_NML k_test_1 = 20454, k_test_2 = 20454 /
 /
 /

Tests/Input/Region_5/SOR5_in.txt

- &MAIN_NML
- CORNER_deg  = 0.12500       , -59.875,
- NLON    = 120,
- NLAT    = 180
- k_test_1 = 20454, k_test_2 = 20455
- /
+ &MAIN_NML k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/Region_6_SHP_triplet/grid_nml.txt

+&grid_nml CORNER_DEG      = 5.875000       , -32.12500      ,
+ NLON    = 26,
+ NLAT    = 22 /

Tests/Input/SODS_in.txt

- &MAIN_NML
- CORNER_deg  = 16.12500       , -38.87500      ,
- NLON    = 20,
- NLAT    = 20
- k_test_1 = 20454, k_test_2 = 20455
- /
+ &MAIN_NML k_test_1 = 20454, k_test_2 = 20455 /
 /
 /

Tests/Input/successive_overlap_in.txt

- &MAIN_NML
- CORNER_deg  = 16.12500       , -38.87500      ,
- NLON    = 20,
- NLAT    = 20
- k_test_1 = 20454, k_test_2 = 20454
- /
+ &MAIN_NML k_test_1 = 20454, k_test_2 = 20454 /
 /
 /

Tests/Input/successive_overlap_periodic_in.txt

- &MAIN_NML
- CORNER_deg  = 0.000000       , -81.00000      ,
- STEP_deg    = 20.00000       , 18.00000       ,
- NLON    = 18,
- NLAT    = 10,
- DIST_LIM        = 4
- k_test_1 = 20454, k_test_2 = 20454
- /
+ &MAIN_NML DIST_LIM        = 4 k_test_1 = 20454, k_test_2 = 20454 /
 /
 /

Tests/test_overlap.f90

@@ -41,8 +41,8 @@ program test_overlap
   ! We look for an overlapping eddy at dist_lim (in grid points) of
   ! the first extremum.
 
-  namelist /main_nml/ corner_deg, step_deg, nlon, nlat, dist_lim, max_delta, &
-       k_test_1, k_test_2
+  namelist /grid_nml/ corner_deg, step_deg, nlon, nlat
+  namelist /main_nml/ dist_lim, max_delta, k_test_1, k_test_2
 
   !-------------------------------------------------------------------------
 
@@ -57,7 +57,11 @@ program test_overlap
 
   call get_command_arg_dyn(1, shp_tr_dir, &
        "Required argument: SHP-triplet-directory")
-
+  call new_unit(unit)
+  open(unit, file = shp_tr_dir // "/grid_nml.txt", status = "old", &
+       action = "read", position = "rewind")
+  read(unit, nml = grid_nml)
+  close(unit)
   write(unit = error_unit, nml = main_nml)
   write(unit = error_unit, fmt = *) "Enter namelist main_nml."
   read(unit = *, nml = main_nml)
@@ -83,16 +87,15 @@ program test_overlap
   call read_snapshot(flow(max_delta + 1), k_test_2, hshp, corner, step, nlon, &
        nlat, copy, k1, ishape_last)
   CALL shp_tr_close(hshp)
-  print *, "Enter flow(1)%list_vis%delta_out:"
+  print *, "Enter flow(1)%list_vis%delta_out (array):"
   read *, flow(1)%list_vis%delta_out
   print *, flow(1)%list_vis%delta_out
-  print *, "Enter flow(max_delta + 1)%list_vis%delta_in:"
+  print *, "Enter flow(max_delta + 1)%list_vis%delta_in (array):"
   read *, flow(max_delta + 1)%list_vis%delta_in
   print *, flow(max_delta + 1)%list_vis%delta_in
 
   flow(2:max_delta)%number_eddies = 20000
 
-  call new_unit(unit)
   open(unit, file = "edgelist.csv", status = "replace", &
        action  = "write")
 

loop_shp_append.py

-#!/usr/bin/env python3
-
-import datetime
-import subprocess
-from os import path
-import os
-import time
-import glob
-import shutil
-
-initial_dir = os.getcwd()
-
-for year in range(2011, 2017):
-    my_path = path.join(initial_dir, str(year))
-    os.chdir(my_path)
-
-    # The first day of the year is the base, so there is nothing to
-    # concatenate, we just copy it:
-    for src in glob.iglob(f"{year}_01/{year}_01_01/*"):
-        shutil.copy(src, os.curdir)
-
-    my_date = datetime.date(year, 1, 2)
-    final_date = datetime.date(year, 12, 31)
-    t0 = time.perf_counter()
-
-    with  open("perf_report.csv", "w") as perf_report:
-        perf_report.write("elapsed time, in s\n")
-
-        while my_date <= final_date:
-            my_dir = my_date.strftime("%Y_%m/%Y_%m_%d")
-
-            for filename in ["extremum", "max_speed_contour", 
-                             "outermost_contour"]:
-                from_file = path.join(my_dir, filename)
-
-                for command in ["dbfcat", "shpcat"]:
-                    subprocess.run([command, from_file, filename], check = True)
-
-            t1 = time.perf_counter()
-            perf_report.write(str(t1 - t0) + "\n")
-            t0 = t1
-            my_date += datetime.timedelta(1)

loop_shp_tr_append.py

+#!/usr/bin/env python3
+
+"""Concatenate a series of daily SHP triplet directories, in a given
+year, assuming their names are SHP_triplet_%Y_%m_%d."""
+
+import datetime
+import subprocess
+from os import path
+import os
+import shutil
+import sys
+import shapefile
+
+os.mkdir(sys.argv[1])
+year = 2006
+my_date = datetime.date(year, 1, 3)
+final_date = datetime.date(year, 1, 15)
+
+# The first day of the year is the base, so there is nothing to
+# concatenate, we just copy it:
+from_dir = my_date.strftime("SHP_triplet_%Y_%m_%d")
+for src in os.scandir(from_dir): shutil.copy(src, sys.argv[1])
+
+file = path.join(sys.argv[1], "ishape_last.txt")
+
+with open(file, "w") as ishape_last:
+    to_file = path.join(sys.argv[1], "extremum")
+    with shapefile.Reader(to_file) as shpr: n_shapes = len(shpr)
+    ishape_last.write(str(n_shapes - 1) + "\n")
+    my_date += datetime.timedelta(1)
+
+    while my_date <= final_date:
+        from_dir = my_date.strftime("SHP_triplet_%Y_%m_%d")
+
+        for filename in ["extremum", "max_speed_contour", "outermost_contour"]:
+            from_file = path.join(from_dir, filename)
+            to_file = path.join(sys.argv[1], filename)
+
+            for command in ["dbfcat", "shpcat"]:
+                subprocess.run([command, from_file, to_file], check = True)
+
+        with shapefile.Reader(to_file) as shpr: n_shapes = len(shpr)
+        ishape_last.write(str(n_shapes - 1) + "\n")
+        my_date += datetime.timedelta(1)