gutmann
commented
6 years ago Thanks, I thought this was fixed in df28968b3e34e4d906cd12405c732015ccbbebbd but it seems to still be stuck. This seems to come from a relatively subtle bug related to the precision with which one tests for a point to be on a line vs inside a polygon more generally in different parts of https://github.com/NCAR/icar/blob/develop/src/utilities/geo_reader.f90
Essentially it is testing that a point falls right on the line between two polygons and at one point it decides that the point is close enough to the line that is should be considered "on" it, but in another part of the code it decides it is too far away. In this case it selects the wrong subset of points because it doesn't find any valid polygons, then it later recognizes that the subset it is using is not valid when it attempts to compute the interpolation weights.
This kind of problem seems to be more common on highly skewed / rotated grids, though I think it is a more general problem.
When using WRF output (27 km as low res, 9 km as high res) as forcing data for ICAR, I get the following error:
Initializing Boundaries Interpreting geopotential height as residing between model layers Setting up domain geographic Look Up Tables ERROR finding triangle 39 39 38 38 61 60 61 60 33.12846 238.2691
33.35105 238.3338
33.14337 238.4846
33.22448 238.0857
33.01720 238.2369
ERROR: Point not located in bounding triangle 238.2691 33.12846
Triangle vertices 238.3338 33.35105
238.4846 33.14337
238.2852 33.18402
-0.3326202 2.0896307E-04 1.332411
Point in poly returns: F Triangulation is broken