[GRASS-dev] r.profile tests

Hi,

I was trying to fix failing test of r.profile:
http://fatra.cnr.ncsu.edu/grassgistests/reports_for_date-2014-10-19-07-00/report_for_nc_spm_08_grass7_nc/raster/r.profile/test_profile_ncspm/index.html

which is working on 32-bit machine. I committed this change:
http://trac.osgeo.org/grass/changeset/62299

which fixed it but I was wondering, if r.profile shouldn’t use atan2 instead of atan? It seems that the division in atan could possibly lead to overflow? But nobody complained about it yet… Tests work with atan2.

Thanks,

Anna

Anna Petrášová wrote:

I was wondering, if r.profile shouldn't use atan2
instead of atan? It seems that the division in atan could possibly lead to
overflow? But nobody complained about it yet... Tests work with atan2.

Using atan2() shouldn't be necessary, as each quadrant is handled
separately. However, the due south case is problematic:

      if (rows > 0 && cols >= 0) {
          /* SW Quad or due south */
          AZI = atan((rows / cols));

as "cols" could equal zero.

Floating-point division by zero won't result in overflow, it will
result in infinity (which is representable). The problem is that it's
not clear whether the argument to atan() will be positive infinity
(resulting in AZI=pi/2) or negative infinity (resulting in AZI=-pi/2).

If cols is -0.0, the argument will (typically) be negative infinity
(IIRC, C implementations aren't required to preserve the sign of a
zero).

Given that AZI is only used for calculating X and Y, I would suggest
using neither atan() nor atan2() but simply normalising the vector,
e.g.

  double k = res / hypot(rows, cols);
  double X = k * rows;
  double Y = k * cols;

Note that this doesn't need a separate case for each quadrant.

Technically, we shouldn't be using hypot() as it isn't part of C89 (it
was added in C99). But it's already being used in:

  lib/vector/diglib/prune.c
  lib/vector/Vlib/line.c
  lib/vector/Vlib/poly.c
  lib/gmath/la.c
  lib/gis/distance.c

--
Glynn Clements <glynn@gclements.plus.com>

On Mon, Oct 20, 2014 at 9:45 PM, Glynn Clements <glynn@gclements.plus.com>
wrote:

Anna Petrášová wrote:

> I was wondering, if r.profile shouldn't use atan2
> instead of atan? It seems that the division in atan could possibly lead
to
> overflow? But nobody complained about it yet... Tests work with atan2.

Using atan2() shouldn't be necessary, as each quadrant is handled
separately. However, the due south case is problematic:

            if (rows > 0 && cols >= 0) {
                /* SW Quad or due south */
                AZI = atan((rows / cols));

as "cols" could equal zero.

Floating-point division by zero won't result in overflow, it will
result in infinity (which is representable). The problem is that it's
not clear whether the argument to atan() will be positive infinity
(resulting in AZI=pi/2) or negative infinity (resulting in AZI=-pi/2).

If cols is -0.0, the argument will (typically) be negative infinity
(IIRC, C implementations aren't required to preserve the sign of a
zero).

Given that AZI is only used for calculating X and Y, I would suggest
using neither atan() nor atan2() but simply normalising the vector,
e.g.

        double k = res / hypot(rows, cols);
        double X = k * rows;
        double Y = k * cols;

Thanks for the suggestion, I committed it in 62327. It passed my tests.

Anna

Note that this doesn't need a separate case for each quadrant.

Technically, we shouldn't be using hypot() as it isn't part of C89 (it
was added in C99). But it's already being used in:

        lib/vector/diglib/prune.c
        lib/vector/Vlib/line.c
        lib/vector/Vlib/poly.c
        lib/gmath/la.c
        lib/gis/distance.c

--
Glynn Clements <glynn@gclements.plus.com>