From e8581d244840246a759890483a61f140c478d0c0 Mon Sep 17 00:00:00 2001
From: Lionel GUEZ <guez@lmd.ens.fr>
Date: Fri, 9 Feb 2018 14:10:24 +0100
Subject: [PATCH] Circumventing gfortran bug: gfortran makes no difference
 between quiet and signalling NaN so we cannot use the -ffpe-trap=invalid
 option.

Accept that we might not be able to compute the mean speed on the
outermost contour: there might be undefined values of velocity
around. In that case just use the contour coming from
max_speed_contour_ssh.

Add quiver_key in plot_snapshot.py.

New test for region 5.
---
 GNUmakefile            |  4 ++++
 Tests/plot_snapshot.py |  7 +++++--
 Tests/tests.json       | 12 ++++++++++--
 mean_speed.f           |  4 ----
 set_max_speed.f        |  9 ++++++++-
 5 files changed, 27 insertions(+), 9 deletions(-)

diff --git a/GNUmakefile b/GNUmakefile
index 4b07fa03..bb36ba89 100644
--- a/GNUmakefile
+++ b/GNUmakefile
@@ -35,7 +35,11 @@ execut = test_good_contour test_inside_4 test_get_1_outerm test_local_extrema te
 mode = debug
 include ${general_compiler_options_dir}/${FC}_${mode}.mk
 comma=,
+
+ifeq (${FC},gfortran)
+# gfortran bug:
 FFLAGS := $(subst invalid${comma},,${FFLAGS})
+endif
 
 # 4. Rules
 
diff --git a/Tests/plot_snapshot.py b/Tests/plot_snapshot.py
index 13045832..1b9814fd 100755
--- a/Tests/plot_snapshot.py
+++ b/Tests/plot_snapshot.py
@@ -92,8 +92,11 @@ if args.velocity:
         u_rotated, v_rotated = m.rotate_vector(f["u"][0], f["v"][0], longitude,
                                                latitude)
         
-    m.quiver(lon_2d, lat_2d, u_rotated, v_rotated, latlon = True)
-    
+    quiver_return = m.quiver(lon_2d, lat_2d, u_rotated, v_rotated,
+                             latlon = True, scale = 5, scale_units = "width")
+    quiverkey(quiver_return, 0.75, 0.9, 0.1, r"0.1 m s$^{-1}$",
+              coordinates = "figure")
+
 m.drawmeridians(np.linspace(llcrnrlon, urcrnrlon, 6), labels = [0, 0, 0, 1])
 m.drawparallels(np.linspace(llcrnrlat, urcrnrlat, 6), labels = [1, 0, 0, 0])
 m.drawcoastlines()
diff --git a/Tests/tests.json b/Tests/tests.json
index e7ad8dd0..a8958fb1 100644
--- a/Tests/tests.json
+++ b/Tests/tests.json
@@ -118,14 +118,14 @@
     },
     {
 	"args": "$compil_prod_dir/test_get_snapshot",
-	"title": "Get_snapshot",
+	"title": "Get_snapshot_region_1",
 	"required": [["$input_dir/2006_01/h_region_1.nc", "h.nc"],
 		     ["$input_dir/2006_01/uv_region_1.nc", "uv.nc"]],
 	"input": "0.\n"
     },
     {
 	"args": "$compil_prod_dir/test_get_snapshot",
-	"title": "Get_snapshot_noise",
+	"title": "Get_snapshot_region_1_noise",
 	"required": [["$input_dir/2006_01/h_region_1.nc", "h.nc"],
 		     ["$input_dir/2006_01/uv_region_1.nc", "uv.nc"]],
 	"input": "1e-3\n"
@@ -177,5 +177,13 @@
 		     ["$input_dir/2006_01/uv_region_4.nc", "uv.nc"]],
 	"input": "1e-3\n",
 	"description": "Part of the domain has missing values."
+    },
+    {
+	"args": "$compil_prod_dir/test_get_snapshot",
+	"title": "Get_snapshot_region_5",
+	"required": [["$input_dir/2006_01/h_region_5.nc", "h.nc"],
+		     ["$input_dir/2006_01/uv_region_5.nc", "uv.nc"]],
+	"input": "1e-3\n",
+	"description": "Same as Get_snapshot_region_4 with larger domain."
     }
 ]
diff --git a/mean_speed.f b/mean_speed.f
index c56168d5..d3737069 100644
--- a/mean_speed.f
+++ b/mean_speed.f
@@ -9,10 +9,7 @@ contains
     ! Interpolates the wind at each point of the input polygon and
     ! computes the mean azimuthal speed at interpolation points.
 
-    use, intrinsic:: IEEE_ARITHMETIC, only: IEEE_IS_NAN
-
     use contour_531, only: polyline
-    use nr_util, only: assert
     use numer_rec_95, only: bilinear_interp2_reg
 
     real, intent(in):: u(:, :), v(:, :) ! velocity
@@ -49,7 +46,6 @@ contains
     end do
 
     mean_speed = sum(v_azim) / ni
-    call assert(.not. IEEE_IS_NAN(mean_speed), "mean_speed: nan")
 
   end function mean_speed
 
diff --git a/set_max_speed.f b/set_max_speed.f
index 11e588a9..147f3e65 100644
--- a/set_max_speed.f
+++ b/set_max_speed.f
@@ -8,6 +8,8 @@ contains
 
     ! Defines the components max_speed_contour and max_speed of argument e.
 
+    use, intrinsic:: IEEE_ARITHMETIC, only: IEEE_IS_NAN
+
     use contour_531, only: convert_to_reg_coord, convert_to_ind, null_polyline
     use derived_types, only: eddy
     use geometry, only: polygon_point_dist_2d
@@ -87,15 +89,20 @@ contains
             e%max_speed_contour%ssh, e%coord_extr, nearby_extr)
        e%max_speed = mean_speed(u, v, e%max_speed_contour%polyline, &
             e%coord_extr, corner, step)
+
        ! Might the outermost contour be a better choice?
        speed_outerm = mean_speed(u, v, e%outermost_contour%polyline, &
             e%coord_extr, corner, step)
-       
+
+       if (IEEE_IS_NAN(speed_outerm)) speed_outerm = 0.
+       ! (This may happen when the outermost contour is near land.)
+
        if (abs(speed_outerm) > abs(e%max_speed)) then
           ! Abandon the contour coming from max_speed_contour_ssh:
           e%max_speed_contour = e%outermost_contour
           e%max_speed = speed_outerm
        else
+          ! Stick to the contour coming from max_speed_contour_ssh.
           e%max_speed_contour%area &
                = spherical_polygon_area(e%max_speed_contour%polyline)
        end if
-- 
GitLab