[GRASS-user] calculate mode value on moving window (with mapalgebra) with null value cells

I need to assign values to a the cells on the “border” of a raster. The inside and the outside are distinguished by having or not having null values assigned.
I also need to keep the other cell values (internals) untouched.
I thought to use a mapcalc expression like the following (knowing the the external cells have value 9999)

if(A==9999,max(A[1,1],A[1,0],A[1,-1],A[0,-1],A[-1,-1],A[-1,0],A[-1,1],A[0,1]),A)

In this case I’m using max, and I alwasy obtain null values as result, as if nulls (from the cells surrounding the external side of the border) were considered maximum values respect to numeric values (of the internal cells).
I could substiture null values for the external cells with a very low negative value, but I would prefer being able to exlude (or manage somehow) null values in my filter.
Any suggestion?

giovanni

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster. The inside and the outside are distinguished by having or not having null values assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in the manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you want with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 && isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


Hi Marcello,
thanks, but my problem is not finding the border. I already have it (with its own category) but assigning it the value from the surrounding (not null) value, e.g. max/min/etc.
The problem with kernel filters/moving windows is that they do not filter out null values…

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster. The inside and the outside are distinguished by having or not having null values assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in the manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you want with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 && isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


OK, sorry for that. So your problem is the famous edge effect that affects us all.

Well, I think you have a real problem there.

If you find a good solution, please share.

Best,
Marcello.

On Tue, Jun 19, 2012 at 11:39 AM, G. Allegri <giohappy@gmail.com> wrote:

Hi Marcello,
thanks, but my problem is not finding the border. I already have it (with its own category) but assigning it the value from the surrounding (not null) value, e.g. max/min/etc.
The problem with kernel filters/moving windows is that they do not filter out null values…

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster. The inside and the outside are distinguished by having or not having null values assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in the manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you want with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 && isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


On 19/06/12 10:15, G. Allegri wrote:

I need to assign values to a the cells on the "border" of a raster. The
inside and the outside are distinguished by having or not having null
values assigned.
I also need to keep the other cell values (internals) untouched.
I thought to use a mapcalc expression like the following (knowing the
the external cells have value 9999)

if(A==9999,max(A[1,1],A[1,0],A[1,-1],A[0,-1],A[-1,-1],A[-1,0],A[-1,1],A[0,1]),A)

How about using an isnull() function on each cell, i.e.:

if(A==9999,max(if(isnull(A[1,1]),0,A[1,1]),if(isnull(A[1,0]),0,A[1,0]), etc, ),A)

?

Moritz

On Tue, Jun 19, 2012 at 12:39 PM, G. Allegri <giohappy@gmail.com> wrote:

Hi Marcello,
thanks, but my problem is not finding the border. I already have it (with
its own category) but assigning it the value from the surrounding (not null)
value, e.g. max/min/etc.
The problem with kernel filters/moving windows is that they do not filter
out null values...

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

HTH,

Markus M

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the "border" of a raster. The
inside and the outside are distinguished by having or not having null values
assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in the
manual to find the borders of your raster. Something like:

# creates an inverted raster from your raster
> r.mapcalc "raster_inverted=if(isnull(raster,1,null())"
# grow this inverted raster by one cell
> r.grow in=raster_inverted out=raster_inverted_grown
# now both rasters overlap at the border, so you can do whatever you want
with it, for instance, extract it
> r.mapcalc "border=if(raster_inverted_grown==1 &&
> isnull(raster)==0,raster,null())"

Hope it helps.

Cheers,
Marcello.

_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Hi Moritz,
this is a partial solution (as assigning unlikely high/low values to the null cells), but it doesn’t solve, for example, the mode() function.
Anyway, this was the way I already followed :wink:

giovanni

2012/6/19 Moritz Lennert <mlennert@club.worldonline.be>

On 19/06/12 10:15, G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster. The
inside and the outside are distinguished by having or not having null
values assigned.
I also need to keep the other cell values (internals) untouched.
I thought to use a mapcalc expression like the following (knowing the
the external cells have value 9999)

if(A==9999,max(A[1,1],A[1,0],A[1,-1],A[0,-1],A[-1,-1],A[-1,0],A[-1,1],A[0,1]),A)

How about using an isnull() function on each cell, i.e.:

if(A==9999,max(if(isnull(A[1,1]),0,A[1,1]),if(isnull(A[1,0]),0,A[1,0]), etc, ),A)

?

Moritz

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

I have considered r.neighbors but I need to apply the filter only to some specific categories (border cells in my case). A MASK wouldn’t solve my problem, because it would mask the neighbors…

giovanni

HTH,

Markus M

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster. The
inside and the outside are distinguished by having or not having null values
assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in the
manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you want

with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 &&
isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

On Tue, Jun 19, 2012 at 6:15 PM, G. Allegri <giohappy@gmail.com> wrote:

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

I have considered r.neighbors but I need to apply the filter only to some
specific categories (border cells in my case). A MASK wouldn't solve my
problem, because it would mask the neighbors...

You could replace 9999 (according to your first post the category of
border cells) with NULL, then run r.neighbors, without a MASK?

Markus M

giovanni

HTH,

Markus M

>
> giovanni
>
>
> 2012/6/19 Marcello Gorini <gorini@gmail.com>
>>
>>
>>
>> G. Allegri wrote:
>>
>>> I need to assign values to a the cells on the "border" of a raster.
>>> The
>>> inside and the outside are distinguished by having or not having null
>>> values
>>> assigned.
>>> I also need to keep the other cell values (internals) untouched.
>>>
>>
>> Hey,
>>
>> Check out r.grow. I am pretty sure you can modify the example given in
>> the
>> manual to find the borders of your raster. Something like:
>>
>> # creates an inverted raster from your raster
>> > r.mapcalc "raster_inverted=if(isnull(raster,1,null())"
>> # grow this inverted raster by one cell
>> > r.grow in=raster_inverted out=raster_inverted_grown
>> # now both rasters overlap at the border, so you can do whatever you
>> want
>> with it, for instance, extract it
>> > r.mapcalc "border=if(raster_inverted_grown==1 &&
>> > isnull(raster)==0,raster,null())"
>>
>> Hope it helps.
>>
>> Cheers,
>> Marcello.
>>
>>
>>
>
>
> _______________________________________________
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
>

Yes, with r.neighbors I woul chenge the values to NULL, but I don’t know how to solve the problem of applying the filter only to the cells with a specific value (involving only the values of their neighbors).
My border cells have assigned a 10000 value. I’m working on a DEM., so this fictious value is a way to code my borders. I want to apply the filter only on these cells, while leaving the others (not on the border) unaltered.
With r.mapcalc I use an if() statement, but with r.neighbors I would need a MASK.

giovanni

2012/6/19 Markus Metz <markus.metz.giswork@googlemail.com>

On Tue, Jun 19, 2012 at 6:15 PM, G. Allegri <giohappy@gmail.com> wrote:

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

I have considered r.neighbors but I need to apply the filter only to some
specific categories (border cells in my case). A MASK wouldn’t solve my
problem, because it would mask the neighbors…

You could replace 9999 (according to your first post the category of
border cells) with NULL, then run r.neighbors, without a MASK?

Markus M

giovanni

HTH,

Markus M

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster.
The
inside and the outside are distinguished by having or not having null
values
assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given in
the
manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you

want
with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 &&
isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

On Tue, Jun 19, 2012 at 6:36 PM, G. Allegri <giohappy@gmail.com> wrote:

Yes, with r.neighbors I woul chenge the values to NULL, but I don't know how
to solve the problem of applying the filter only to the cells with a
specific value (involving only the values of their neighbors).

r.neighbors followed by r.patch, no MASK? This will keep the valid
cells of the original map and only replace cells that are NULL in the
original map with values from the filtered map. E.g.

g.region -p rast=DEM
r.neighbors in=DEM out=DEM_filt method=mode
r.patch in=DEM,DEM_filt out=DEM_patched

My border cells have assigned a 10000 value. I'm working on a DEM., so this
fictious value is a way to code my borders. I want to apply the filter only
on these cells, while leaving the others (not on the border) unaltered.
With r.mapcalc I use an if() statement, but with r.neighbors I would need a
MASK.

BTW, if the border problem arouse from reprojecting with a method
other than nearest, you can use in trunk as resampling method one of
bilinear_f,cubic_f,lanczos_f with r.proj to avoid the border effect.

Markus M

giovanni

2012/6/19 Markus Metz <markus.metz.giswork@googlemail.com>

On Tue, Jun 19, 2012 at 6:15 PM, G. Allegri <giohappy@gmail.com> wrote:
>> AFAICT, r.neighbors ignores NULL values and assigns the new value from
>> the surrounding non-NULL values. You could then patch your original
>> map with the output of r.neighbors.
>
>
> I have considered r.neighbors but I need to apply the filter only to
> some
> specific categories (border cells in my case). A MASK wouldn't solve my
> problem, because it would mask the neighbors...

You could replace 9999 (according to your first post the category of
border cells) with NULL, then run r.neighbors, without a MASK?

Markus M

>
> giovanni
>
>
>>
>>
>> HTH,
>>
>> Markus M
>>
>> >
>> > giovanni
>> >
>> >
>> > 2012/6/19 Marcello Gorini <gorini@gmail.com>
>> >>
>> >>
>> >>
>> >> G. Allegri wrote:
>> >>
>> >>> I need to assign values to a the cells on the "border" of a raster.
>> >>> The
>> >>> inside and the outside are distinguished by having or not having
>> >>> null
>> >>> values
>> >>> assigned.
>> >>> I also need to keep the other cell values (internals) untouched.
>> >>>
>> >>
>> >> Hey,
>> >>
>> >> Check out r.grow. I am pretty sure you can modify the example given
>> >> in
>> >> the
>> >> manual to find the borders of your raster. Something like:
>> >>
>> >> # creates an inverted raster from your raster
>> >> > r.mapcalc "raster_inverted=if(isnull(raster,1,null())"
>> >> # grow this inverted raster by one cell
>> >> > r.grow in=raster_inverted out=raster_inverted_grown
>> >> # now both rasters overlap at the border, so you can do whatever you
>> >> want
>> >> with it, for instance, extract it
>> >> > r.mapcalc "border=if(raster_inverted_grown==1 &&
>> >> > isnull(raster)==0,raster,null())"
>> >>
>> >> Hope it helps.
>> >>
>> >> Cheers,
>> >> Marcello.
>> >>
>> >>
>> >>
>> >
>> >
>> > _______________________________________________
>> > grass-user mailing list
>> > grass-user@lists.osgeo.org
>> > http://lists.osgeo.org/mailman/listinfo/grass-user
>> >
>
>

g.region -p rast=DEM
r.neighbors in=DEM out=DEM_filt method=mode
r.patch in=DEM,DEM_filt out=DEM_patched

Right Markus :wink:
This is the simplest approach.

BTW, if the border problem arouse from reprojecting with a method
other than nearest, you can use in trunk as resampling method one of
bilinear_f,cubic_f,lanczos_f with r.proj to avoid the border effect.

I didn’t know about these “fallbacks” methods in trunk. Very good.
I was thinking about that since a long time, because it’s a featue I find very useful in ArcGIS. GRASS gives a plus: the option to use fallbacks or not :wink:

It would be useful to have them in r.resamp.interp too. This is the command that “eats” my border.

giovanni

Markus M

giovanni

2012/6/19 Markus Metz <markus.metz.giswork@googlemail.com>

On Tue, Jun 19, 2012 at 6:15 PM, G. Allegri <giohappy@gmail.com> wrote:

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

I have considered r.neighbors but I need to apply the filter only to
some
specific categories (border cells in my case). A MASK wouldn’t solve my
problem, because it would mask the neighbors…

You could replace 9999 (according to your first post the category of
border cells) with NULL, then run r.neighbors, without a MASK?

Markus M

giovanni

HTH,

Markus M

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster.
The
inside and the outside are distinguished by having or not having
null
values
assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given
in
the
manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you

want
with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 &&
isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

FYI, Glynn has answered me right now on grass-dev saying that in trunk nmax, nmin, nmode, etc. are available, which discard NULL values.

giovanni

2012/6/20 G. Allegri <giohappy@gmail.com>

g.region -p rast=DEM
r.neighbors in=DEM out=DEM_filt method=mode
r.patch in=DEM,DEM_filt out=DEM_patched

Right Markus :wink:
This is the simplest approach.

BTW, if the border problem arouse from reprojecting with a method
other than nearest, you can use in trunk as resampling method one of
bilinear_f,cubic_f,lanczos_f with r.proj to avoid the border effect.

I didn’t know about these “fallbacks” methods in trunk. Very good.
I was thinking about that since a long time, because it’s a featue I find very useful in ArcGIS. GRASS gives a plus: the option to use fallbacks or not :wink:

It would be useful to have them in r.resamp.interp too. This is the command that “eats” my border.

giovanni

Markus M

giovanni

2012/6/19 Markus Metz <markus.metz.giswork@googlemail.com>

On Tue, Jun 19, 2012 at 6:15 PM, G. Allegri <giohappy@gmail.com> wrote:

AFAICT, r.neighbors ignores NULL values and assigns the new value from
the surrounding non-NULL values. You could then patch your original
map with the output of r.neighbors.

I have considered r.neighbors but I need to apply the filter only to
some
specific categories (border cells in my case). A MASK wouldn’t solve my
problem, because it would mask the neighbors…

You could replace 9999 (according to your first post the category of
border cells) with NULL, then run r.neighbors, without a MASK?

Markus M

giovanni

HTH,

Markus M

giovanni

2012/6/19 Marcello Gorini <gorini@gmail.com>

G. Allegri wrote:

I need to assign values to a the cells on the “border” of a raster.
The
inside and the outside are distinguished by having or not having
null
values
assigned.
I also need to keep the other cell values (internals) untouched.

Hey,

Check out r.grow. I am pretty sure you can modify the example given
in
the
manual to find the borders of your raster. Something like:

creates an inverted raster from your raster

r.mapcalc “raster_inverted=if(isnull(raster,1,null())”

grow this inverted raster by one cell

r.grow in=raster_inverted out=raster_inverted_grown

now both rasters overlap at the border, so you can do whatever you

want
with it, for instance, extract it

r.mapcalc “border=if(raster_inverted_grown==1 &&
isnull(raster)==0,raster,null())”

Hope it helps.

Cheers,
Marcello.


grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user