[GRASS-dev] [GRASS-user] Gauss filter in r.neighbors

Can someone determine whether there's a problem with
"r.neighbors gauss=..." in 6.x?

Martin Album Ytre-Eide wrote:

I followed you examples and it did not work in grass 6.4.4 so I
tried them in grass7.1 and they worked for me as well.

Tried my map - and it worked. I result looks good :slight_smile:

Thanks!

(maybe there is a problem in r.neighbors for grass 6.4.4?)

Martin

________________________________________
Fra: Glynn Clements [glynn@gclements.plus.com]
Sendt: 2. oktober 2014 14:27
Til: Martin Album Ytre-Eide
Kopi: grass-user@lists.osgeo.org
Emne: Re: [GRASS-user] Gauss filter in r.neighbors

Martin Album Ytre-Eide wrote:

> I am using the the r.neighbors function and I have an idea about
> what the 'Gauss' parameter will do for me, but I can't find any
> documentation for it. Experimenting with it, left me clueless.

Using gauss= is equivalent to using weight= with a weights file
corresponding to a Gaussian filter.

The exact weights which are used can be found in the
gaussian_weights() function in raster/r.neighbors/readweights.c. They
are effectively:

        ncb.weights[i][j] = exp(-(x*x+y*y)/(2*sigma2))/(2*M_PI*sigma2);

where x, y is the offset (in cells) from the centre of the
neighbourhood and sigma2 is the square of the value passed to the
gauss= option.

> What I won't to do:

> Apply a Gaussian weighting of cell neighbors and sum the values of
> the neighbors. If all cells=1, I would get only 1's back.

You need to use method=average to make the effective weights sum to
one (when weights are used, method=average divides the weighted sum of
the cell values by the sum of the weights).

IIRC, the weights are such that they would sum to one if you had an
infinitely-large neighbourhood, but the fact that you're cropping the
curve to a finite neighbourhood means that the sum will always be less
than one.

The weights are determined by sigma alone, and are unaffected by the
neighbourhood size.

> I have a map of zeros and ones. I tried this:
>
> r.neighbors -c input='input_map' output='output_map' method=sum gauss=10 size=33
>
> The result is a map with only zeros.

It works for me; e.g.

        $ r.mapcalc --o 'foo = rand(0,2)' seed=1
        $ r.neighbors --o -c in=foo out=bar method=average size=33 gauss=10
        $ r.info -r bar
        min=0.442173333781042
        max=0.568925311919939
        $ r.neighbors --o -c in=foo out=bar method=sum size=33 gauss=10
        $ r.info -r bar
        min=0.0985901293262625
        max=0.449425795562072

        $ r.mapcalc --o 'foo = 1'
        $ r.neighbors --o -c in=foo out=bar method=average size=33 gauss=10
        $ r.info -r bar
        min=1
        max=1
        $ r.neighbors --o -c in=foo out=bar method=sum size=33 gauss=10
        $ r.info -r bar
        min=0.221413499330273
        max=0.812157275619236

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

Glynn Clements wrote:

Can someone determine whether there's a problem with
"r.neighbors gauss=..." in 6.x?

There is no problem with the Gaussian filter, but a problem with the
output map type of r.neighbours in G64: the output type is in G64 the
same like the input type, whereas in G7 the output type depends on the
chosen function. That means, if you use r.neighbours method=average in
G64 with a CELL map, the output is also CELL.

Markus M

Martin Album Ytre-Eide wrote:

I followed you examples and it did not work in grass 6.4.4 so I
tried them in grass7.1 and they worked for me as well.

Tried my map - and it worked. I result looks good :slight_smile:

Thanks!

(maybe there is a problem in r.neighbors for grass 6.4.4?)

Martin

________________________________________
Fra: Glynn Clements [glynn@gclements.plus.com]
Sendt: 2. oktober 2014 14:27
Til: Martin Album Ytre-Eide
Kopi: grass-user@lists.osgeo.org
Emne: Re: [GRASS-user] Gauss filter in r.neighbors

Martin Album Ytre-Eide wrote:

> I am using the the r.neighbors function and I have an idea about
> what the 'Gauss' parameter will do for me, but I can't find any
> documentation for it. Experimenting with it, left me clueless.

Using gauss= is equivalent to using weight= with a weights file
corresponding to a Gaussian filter.

The exact weights which are used can be found in the
gaussian_weights() function in raster/r.neighbors/readweights.c. They
are effectively:

        ncb.weights[i][j] = exp(-(x*x+y*y)/(2*sigma2))/(2*M_PI*sigma2);

where x, y is the offset (in cells) from the centre of the
neighbourhood and sigma2 is the square of the value passed to the
gauss= option.

> What I won't to do:

> Apply a Gaussian weighting of cell neighbors and sum the values of
> the neighbors. If all cells=1, I would get only 1's back.

You need to use method=average to make the effective weights sum to
one (when weights are used, method=average divides the weighted sum of
the cell values by the sum of the weights).

IIRC, the weights are such that they would sum to one if you had an
infinitely-large neighbourhood, but the fact that you're cropping the
curve to a finite neighbourhood means that the sum will always be less
than one.

The weights are determined by sigma alone, and are unaffected by the
neighbourhood size.

> I have a map of zeros and ones. I tried this:
>
> r.neighbors -c input='input_map' output='output_map' method=sum gauss=10 size=33
>
> The result is a map with only zeros.

It works for me; e.g.

        $ r.mapcalc --o 'foo = rand(0,2)' seed=1
        $ r.neighbors --o -c in=foo out=bar method=average size=33 gauss=10
        $ r.info -r bar
        min=0.442173333781042
        max=0.568925311919939
        $ r.neighbors --o -c in=foo out=bar method=sum size=33 gauss=10
        $ r.info -r bar
        min=0.0985901293262625
        max=0.449425795562072

        $ r.mapcalc --o 'foo = 1'
        $ r.neighbors --o -c in=foo out=bar method=average size=33 gauss=10
        $ r.info -r bar
        min=1
        max=1
        $ r.neighbors --o -c in=foo out=bar method=sum size=33 gauss=10
        $ r.info -r bar
        min=0.221413499330273
        max=0.812157275619236

--
Glynn Clements <glynn@gclements.plus.com>
_______________________________________________
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

On Mon, Oct 6, 2014 at 9:34 AM, Markus Metz
<markus.metz.giswork@gmail.com> wrote:

Glynn Clements wrote:

Can someone determine whether there's a problem with
"r.neighbors gauss=..." in 6.x?

There is no problem with the Gaussian filter, but a problem with the
output map type of r.neighbours in G64: the output type is in G64 the
same like the input type, whereas in G7 the output type depends on the
chosen function. That means, if you use r.neighbours method=average in
G64 with a CELL map, the output is also CELL.

FYI - the related (bulk) backports from trunk to relbranch7 were these:
r59668, r59678, r60456

(in case of further backporting to GRASS 6)

markusN