[GRASS5] r.mapcalc: x^y and exp(x,y)
Glynn Clements
glynn.clements at virgin.net
Wed Feb 6 13:26:11 EST 2002
Markus Neteler wrote:
> as you are currently on r.mapcalc:
> Do you know the difference between x^y and exp(x,y)?
> It appears to be doing the same thing. Do you know what is the difference
> between those two?
Both call pow(), but the "^" operator (implemented by power_x() in
math.c) explicitly checks for the case of a negative number raised to
a non-integral power, and returns NULL (using G_set_d_null_value())
However, in this case, pow() will return NaN (and set errno to EDOM),
so there shouldn't be any difference. However, there is; the two cases
return different NaN values: exp(x,y) returns 0x8000000000000 while
x^y returns 0xfffffffffffff. The former isn't NULL according to
G_is_d_null_value(), while the latter is.
BTW, the current behaviour of r.mapcalc3 matches that of exp() in for
both exp() and the "^" operator when x and y are FCEL/DCELL values.
However, another difference is that the "^" operator in r.mapcalc3
returns CELL if both arguments are CELL, while exp() always returns
DCELL.
I'll fix both r.mapcalc and r.mapcalc3 to ensure that the case of a
negative number raised to a non-integral power is always NULL.
However, I suspect that it might be better to use:
errno = 0;
result = pow(x, y);
if (errno != 0)
SETNULL_D(&result);
rather than the existing test:
if (x < 0.0 && y != ceil(y))
The latter can misbehave due to differences in precision.
More generally, this raises the question of whether the existing
implementation of G_is_[cd]_null_value() is correct, i.e. whether it
should treat all NaN values as NULL rather than one specific bit
pattern.
--
Glynn Clements <glynn.clements at virgin.net>
More information about the grass-dev
mailing list