Verifygrid: check lon/lat bounds.

src/Verifygrid.cc

@@ -281,7 +281,7 @@ are_polygon_vertices_arranged_in_clockwise_order(double cell_area)
 }
 
 static void
-create_sorted_point_array(size_t gridsize, const double *grid_center_lon, const double *grid_center_lat, Varray<double> &center_point_array)
+create_sorted_point_array(size_t gridsize, const Varray<double> &grid_center_lon, const Varray<double> &grid_center_lat, Varray<double> &center_point_array)
 {
   for (size_t cell_no = 0; cell_no < gridsize; cell_no++)
     {
@@ -475,8 +475,9 @@ void copy_unique_corners(int actual_number_of_corners, const Varray<Point3D> &ce
 }
 
 static void
-verify_grid(int gridtype, size_t gridsize, size_t nx, int gridno, int ngrids, int ncorner, double *grid_center_lon, double *grid_center_lat,
-            double *grid_corner_lon, double *grid_corner_lat)
+verify_grid(size_t gridsize, size_t nx, int ncorner,
+            const Varray<double> &grid_center_lon, const Varray<double> &grid_center_lat,
+            const Varray<double> &grid_corner_lon, const Varray<double> &grid_corner_lat)
 {
   /*
      First, this function performs the following test:
@@ -508,14 +509,12 @@ verify_grid(int gridtype, size_t gridsize, size_t nx, int gridno, int ngrids, in
 
   std::vector<int> no_cells_with_a_specific_no_of_corners(ncorner, 0);
 
-  print_header(gridtype, gridsize, nx, gridno, ngrids);
-
-  /* For performing the first test, an array of all center point coordinates is built. */
+  // For performing the first test, an array of all center point coordinates is built.
 
   Varray<double> center_point_array(gridsize * 2);
   create_sorted_point_array(gridsize, grid_center_lon, grid_center_lat, center_point_array);
 
-  /* Now checking for the number of unique center point coordinates. */
+  // Now checking for the number of unique center point coordinates.
 
   const auto no_unique_center_points = get_no_unique_center_points(gridsize, nx, center_point_array);
 
@@ -734,11 +733,13 @@ print_lonlat(int dig, int gridID, size_t gridsize, const Varray<double> &lon, co
   cdiInqKeyString(gridID, CDI_XAXIS, CDI_KEY_NAME, name, &length);
   const auto xmm = varrayMinMax(gridsize, lon);
   cdoPrint("%10s : %.*g to %.*g degrees", name, dig, xmm.min, dig, xmm.max);
+  if (xmm.min < -720.0 || xmm.max > 720.0) cdoWarning("Longitudes out of bounds!");
 
   length = CDI_MAX_NAME;
   cdiInqKeyString(gridID, CDI_YAXIS, CDI_KEY_NAME, name, &length);
   const auto ymm = varrayMinMax(gridsize, lat);
   cdoPrint("%10s : %.*g to %.*g degrees", name, dig, ymm.min, dig, ymm.max);
+  if (ymm.min < -90.0 || ymm.max > 90.0) cdoWarning("Latitudes out of bounds!");
 }
 
 void *
@@ -757,7 +758,6 @@ Verifygrid(void *argument)
   for (int gridno = 0; gridno < ngrids; ++gridno)
     {
       bool lgeo = true;
-      bool luse_grid_corner = true;
 
       auto gridID = vlistGrid(vlistID, gridno);
       const auto gridtype = gridInqType(gridID);
@@ -790,7 +790,6 @@ Verifygrid(void *argument)
 
           if  (!gridHasCoordinates(gridID)) cdoAbort("Cell center coordinates missing!");
 
-          Varray<double> grid_corner_lat, grid_corner_lon;
           Varray<double> grid_center_lat(gridsize), grid_center_lon(gridsize);
 
           gridInqYvals(gridID, grid_center_lat.data());
@@ -800,24 +799,27 @@ Verifygrid(void *argument)
           cdo_grid_to_degree(gridID, CDI_XAXIS, gridsize, grid_center_lon.data(), "grid center lon");
           cdo_grid_to_degree(gridID, CDI_YAXIS, gridsize, grid_center_lat.data(), "grid center lat");
 
-          if (luse_grid_corner)
+          print_header(gridtype, gridsize, nx, gridno, ngrids);
+
+          if (gridHasBounds(gridID))
             {
               if (ncorner == 0) cdoAbort("Number of cell corners undefined!");
               const auto nalloc = ncorner * gridsize;
-              grid_corner_lat.resize(nalloc);
-              grid_corner_lon.resize(nalloc);
-
-              if (!gridHasBounds(gridID)) cdoAbort("Cell corner coordinates missing!");
+              Varray<double> grid_corner_lat(nalloc), grid_corner_lon(nalloc);
 
               gridInqYbounds(gridID, grid_corner_lat.data());
               gridInqXbounds(gridID, grid_corner_lon.data());
 
               cdo_grid_to_degree(gridID, CDI_XAXIS, ncorner * gridsize, grid_corner_lon.data(), "grid corner lon");
               cdo_grid_to_degree(gridID, CDI_YAXIS, ncorner * gridsize, grid_corner_lat.data(), "grid corner lat");
-            }
 
-          verify_grid(gridtype, gridsize, nx, gridno, ngrids, ncorner, grid_center_lon.data(), grid_center_lat.data(),
-                      grid_corner_lon.data(), grid_corner_lat.data());
+              verify_grid(gridsize, nx, ncorner, grid_center_lon, grid_center_lat,
+                          grid_corner_lon, grid_corner_lat);
+            }
+          else
+            {
+              cdoWarning("Grid cell corner coordinates missing!");
+            }
 
           const int dig = Options::CDO_flt_digits;
           print_lonlat(dig, gridID, gridsize, grid_center_lon, grid_center_lat);