[GRASS5] r.random broken

I just tried to use r.random to generate a set of random points across the Spearfish dem file—to illustrate interpolation.

When I display the points, they are all clustered along the top (i.e., north) edge. This happens with both raster and vector points, and when I run it against both elevation.dem and elevation.10m. It looks something is wrong with the algorithm to generate the y values.

Michael


Michael Barton, Professor of Anthropology
School of Human Evolution and Social Change
Arizona State University
Tempe, AZ 85287-2402
USA

voice: 480-965-6262; fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

I just tried to use r.random to generate a set of random points across the
Spearfish dem file‹to illustrate interpolation.

When I display the points, they are all clustered along the top (i.e.,
north) edge. This happens with both raster and vector points, and when I run
it against both elevation.dem and elevation.10m. It looks something is wrong
with the algorithm to generate the y values.

r.random totally won't work if either __CYGWIN__ or __APPLE__ are
defined, due to the following in raster/r.random/creat_rand.c:

  #if defined(__CYGWIN__) || defined(__APPLE__)
  #define lrand48() rand()/32767.0
  #define srand48(sv) (srand((unsigned)(sv)))
  #else

lrand48() is supposed to return random numbers between 0 and 2^31. The
above lrand48() macro will end up returning values which are far too
small, so, the random test will return true too often, resulting in
the desired number of random cells being reached far too early.

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

Glynn,

Thanks for identifying the cause of the behavior I observed.

This is weird. Why is this module programmed to work incorrectly on Apple
and Cygwin systems??? Can it be corrected?

Michael

On 6/4/05 12:10 AM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

I just tried to use r.random to generate a set of random points across the
Spearfish dem file‹to illustrate interpolation.

When I display the points, they are all clustered along the top (i.e.,
north) edge. This happens with both raster and vector points, and when I run
it against both elevation.dem and elevation.10m. It looks something is wrong
with the algorithm to generate the y values.

r.random totally won't work if either __CYGWIN__ or __APPLE__ are
defined, due to the following in raster/r.random/creat_rand.c:

#if defined(__CYGWIN__) || defined(__APPLE__)
#define lrand48() rand()/32767.0
#define srand48(sv) (srand((unsigned)(sv)))
#else

lrand48() is supposed to return random numbers between 0 and 2^31. The
above lrand48() macro will end up returning values which are far too
small, so, the random test will return true too often, resulting in
the desired number of random cells being reached far too early.

____________________
C. Michael Barton, Professor of Anthropology
School of Human Evolution and Social Change
PO Box 872402
Arizona State University
Tempe, AZ 85287-2402
USA

Phone: 480-965-6262
Fax: 480-965-7671
www: <www.public.asu.edu/~cmbarton>

>> I just tried to use r.random to generate a set of random points
>> across the Spearfish dem file_to illustrate interpolation.
>>
>> When I display the points, they are all clustered along the top
>> (i.e., north) edge. This happens with both raster and vector
>> points, and when I run it against both elevation.dem and
>> elevation.10m. It looks something is wrong with the algorithm to
>> generate the y values.
>
> r.random totally won't work if either __CYGWIN__ or __APPLE__ are
> defined, due to the following in raster/r.random/creat_rand.c:
>
> #if defined(__CYGWIN__) || defined(__APPLE__)
> #define lrand48() rand()/32767.0
> #define srand48(sv) (srand((unsigned)(sv)))
> #else
>
> lrand48() is supposed to return random numbers between 0 and 2^31.
> The above lrand48() macro will end up returning values which are far
> too small, so, the random test will return true too often, resulting
> in the desired number of random cells being reached far too early.

..

Thanks for identifying the cause of the behavior I observed.

This is weird. Why is this module programmed to work incorrectly on
Apple and Cygwin systems??? Can it be corrected?

just replace 32767.0 with 2^31. Maybe try-

#define lrand48() rand()/2147483648.;

or
rand()/2.1475e+09;

or
#include <math.h> // & link with -lm
..
rand()/exp2(31);

?

Hamish

> >> I just tried to use r.random to generate a set of random points
> >> across the Spearfish dem file_to illustrate interpolation.
> >>
> >> When I display the points, they are all clustered along the top
> >> (i.e., north) edge. This happens with both raster and vector
> >> points, and when I run it against both elevation.dem and
> >> elevation.10m. It looks something is wrong with the algorithm to
> >> generate the y values.
> >
> > r.random totally won't work if either __CYGWIN__ or __APPLE__ are
> > defined, due to the following in raster/r.random/creat_rand.c:
> >
> > #if defined(__CYGWIN__) || defined(__APPLE__)
> > #define lrand48() rand()/32767.0
> > #define srand48(sv) (srand((unsigned)(sv)))
> > #else
> >
> > lrand48() is supposed to return random numbers between 0 and 2^31.
> > The above lrand48() macro will end up returning values which are far
> > too small, so, the random test will return true too often, resulting
> > in the desired number of random cells being reached far too early.
..

just replace 32767.0 with 2^31. Maybe try-

sorry, ignore that.

better:

#define lrand48() rand()*65536.

Hamish

Hamish wrote:

> > > lrand48() is supposed to return random numbers between 0 and 2^31.
> > > The above lrand48() macro will end up returning values which are far
> > > too small, so, the random test will return true too often, resulting
> > > in the desired number of random cells being reached far too early.
> ..
>
> just replace 32767.0 with 2^31. Maybe try-

sorry, ignore that.

better:

#define lrand48() rand()*65536.

Still wrong.

lrand48() returns a value between 0 and 2^31.

rand() returns a value between 0 and RAND_MAX, where RAND_MAX is an
implementation-defined constant.

To emulate lrand48() using rand(), you would need something like:

  #define lrand48() (rand() * (1<<31) / RAND_MAX)

Except that the above will overflow. To prevent that, convert to
double then back to long, e.g.:

  #define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))

Or, if your compiler supports it, use "long long":

  #define lrand48() ((long)(rand() * (1LL<<31) / RAND_MAX))

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

A month ago, I noted that r.random produces points all clustered along
the top of a region, instead of randomly distributed (see attached
output from r.random for 1000 points). You quickly identified the
problem. But it looks like it hasn¹t been fixed yet.

Following fix applied to 6.1-cvs.

In order to keep things generic, I applied this version:
  (I have no idea what compilers people use)

#define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))

Could this be fixed for 6.0.1?

Could someone please test it on Mac & Cygwin, if it works on both we can
apply it to the 6.0 release branch too.

Hamish

======================================================================
http://grass.itc.it/pipermail/grass5/2005-June/018528.html
[here slighly cropped -HB]

From: Glynn Clements
Date: Sun, 5 Jun 2005 01:12:05 +0100

> > > > lrand48() is supposed to return random numbers between 0 and
> > > > 2^31. The above lrand48() macro will end up returning values
> > > > which are far too small, so, the random test will return true
> > > > too often, resulting in the desired number of random cells
> > > > being reached far too early.

..

lrand48() returns a value between 0 and 2^31.

rand() returns a value between 0 and RAND_MAX, where RAND_MAX is an
implementation-defined constant.

To emulate lrand48() using rand(), you would need something like:

  #define lrand48() (rand() * (1<<31) / RAND_MAX)

Except that the above will overflow. To prevent that, convert to
double then back to long, e.g.:

  #define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))

Or, if your compiler supports it, use "long long":

  #define lrand48() ((long)(rand() * (1LL<<31) / RAND_MAX))

Hamish,

Thanks very much. I don't have the setup to compile to test. Though if a new
binary comes out in the next few days I will give it a go.

Michael

On 7/4/05 10:59 PM, "Hamish" <hamish_nospam@yahoo.com> wrote:

A month ago, I noted that r.random produces points all clustered along
the top of a region, instead of randomly distributed (see attached
output from r.random for 1000 points). You quickly identified the
problem. But it looks like it hasn¹t been fixed yet.

Following fix applied to 6.1-cvs.

In order to keep things generic, I applied this version:
  (I have no idea what compilers people use)

#define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))

Could this be fixed for 6.0.1?

Could someone please test it on Mac & Cygwin, if it works on both we can
apply it to the 6.0 release branch too.

Hamish

======================================================================
http://grass.itc.it/pipermail/grass5/2005-June/018528.html
[here slighly cropped -HB]

From: Glynn Clements
Date: Sun, 5 Jun 2005 01:12:05 +0100

lrand48() is supposed to return random numbers between 0 and
2^31. The above lrand48() macro will end up returning values
which are far too small, so, the random test will return true
too often, resulting in the desired number of random cells
being reached far too early.

..

lrand48() returns a value between 0 and 2^31.

rand() returns a value between 0 and RAND_MAX, where RAND_MAX is an
implementation-defined constant.

To emulate lrand48() using rand(), you would need something like:

#define lrand48() (rand() * (1<<31) / RAND_MAX)

Except that the above will overflow. To prevent that, convert to
double then back to long, e.g.:

#define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))

Or, if your compiler supports it, use "long long":

#define lrand48() ((long)(rand() * (1LL<<31) / RAND_MAX))

____________________
C. Michael Barton, Professor of Anthropology
School of Human Evolution and Social Change
PO Box 872402
Arizona State University
Tempe, AZ 85287-2402
USA

Phone: 480-965-6262
Fax: 480-965-7671
www: <www.public.asu.edu/~cmbarton>

Hi,

there is probably another canditate:

~/grass61/raster/r.mapcalc > grep lrand *
xrand.c:#define lrand48() ((long)rand())
xrand.c: long x = lrand48();

Markus

On 7/4/05 10:59 PM, "Hamish" <hamish_nospam@yahoo.com> wrote:

>> A month ago, I noted that r.random produces points all clustered along
>> the top of a region, instead of randomly distributed (see attached
>> output from r.random for 1000 points). You quickly identified the
>> problem. But it looks like it hasn¹t been fixed yet.
>
> Following fix applied to 6.1-cvs.
>
> In order to keep things generic, I applied this version:
> (I have no idea what compilers people use)
>
> #define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))
>
>
>> Could this be fixed for 6.0.1?
>
> Could someone please test it on Mac & Cygwin, if it works on both we can
> apply it to the 6.0 release branch too.
>
>
>
> Hamish
>
> ======================================================================
> http://grass.itc.it/pipermail/grass5/2005-June/018528.html
> [here slighly cropped -HB]
>
>
> From: Glynn Clements
> Date: Sun, 5 Jun 2005 01:12:05 +0100
>
>>>>>> lrand48() is supposed to return random numbers between 0 and
>>>>>> 2^31. The above lrand48() macro will end up returning values
>>>>>> which are far too small, so, the random test will return true
>>>>>> too often, resulting in the desired number of random cells
>>>>>> being reached far too early.
> ..
>> lrand48() returns a value between 0 and 2^31.
>>
>> rand() returns a value between 0 and RAND_MAX, where RAND_MAX is an
>> implementation-defined constant.
>>
>> To emulate lrand48() using rand(), you would need something like:
>>
>> #define lrand48() (rand() * (1<<31) / RAND_MAX)
>>
>> Except that the above will overflow. To prevent that, convert to
>> double then back to long, e.g.:
>>
>> #define lrand48() ((long)((double) rand() * (1<<31) / RAND_MAX))
>>
>> Or, if your compiler supports it, use "long long":
>>
>> #define lrand48() ((long)(rand() * (1LL<<31) / RAND_MAX))

____________________
C. Michael Barton, Professor of Anthropology
School of Human Evolution and Social Change
PO Box 872402
Arizona State University
Tempe, AZ 85287-2402
USA

Phone: 480-965-6262
Fax: 480-965-7671
www: <www.public.asu.edu/~cmbarton>

Markus Neteler wrote:

there is probably another canditate:

~/grass61/raster/r.mapcalc > grep lrand *
xrand.c:#define lrand48() ((long)rand())
xrand.c: long x = lrand48();

That one is less of a problem, as it doesn't matter what the range of
lrand48() is, so long as it's larger than the interval between the
arguments to r.mapcalc's rand() function. The macro could still be a
problem if rand() only returns e.g. 16-bit values.

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