[GRASS-dev] question about any July changes in hydrology functions

Beginning sometime in the 1st half of July, we’ve started getting some weird artifacts in surface process models that we’ve been developing using r.flow, r.slop.aspect, and r.mapcalc. The results are linear zones of deposits that seem to follow cell grid lines. This wasn’t happening previous to early July. We’ve scouring our scripts and so far can’t find anything to account for this changed behavior. It acts like some kind of dynamic resolution change artifact, but we are not doing anything to cause this AFAICT.

It is probably something squirrelly in our scripts that we just haven’t seen. But I thought I go for a long shot and ask if there have been any changes that might affect the behavior of r.flow or r.slope.aspect recently. I didn’t see anything obvious in the commit list archives, but may have missed a library change of some kind.

Thanks
Michael


Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics and Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

Beginning sometime in the 1st half of July, we¹ve started getting some weird
artifacts in surface process models that we¹ve been developing using r.flow,
r.slop.aspect, and r.mapcalc. The results are linear zones of deposits that
seem to follow cell grid lines. This wasn¹t happening previous to early
July. We¹ve scouring our scripts and so far can¹t find anything to account
for this changed behavior. It acts like some kind of dynamic resolution
change artifact, but we are not doing anything to cause this AFAICT.

It is probably something squirrelly in our scripts that we just haven¹t
seen. But I thought I go for a long shot and ask if there have been any
changes that might affect the behavior of r.flow or r.slope.aspect recently.
I didn¹t see anything obvious in the commit list archives, but may have
missed a library change of some kind.

$ cvs diff -D 2007-07-01 -D 2007-07-15 raster/r.mapcalc raster/r.slope.aspect raster/r.flow
cvs server: Diffing raster/r.mapcalc
cvs server: Diffing raster/r.slope.aspect
cvs server: Diffing raster/r.flow

Assuming that your dates are correct, it isn't the programs. Extending
the comparison to lib/{gis,vector,bitmap,btree,rowio,segment} (the
only libraries which those programs use directly) doesn't show
anything which looks remotely relevant.

I can only suggest checking out and compiling "before" and "after"
versions, and examining any intermediate data produced by your
scripts.

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

Thanks Glynn,

We've been diffing our scripts to the extent possible--i.e., once we started
committing them to the SVN and haven't found anything obvious yet.

Helena also mentioned r.neighbors as a possible culprit, as I'd forgotten to
mention that we also use this to smooth results.

Michael

On 8/1/07 12:31 AM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

Beginning sometime in the 1st half of July, we¹ve started getting some weird
artifacts in surface process models that we¹ve been developing using r.flow,
r.slop.aspect, and r.mapcalc. The results are linear zones of deposits that
seem to follow cell grid lines. This wasn¹t happening previous to early
July. We¹ve scouring our scripts and so far can¹t find anything to account
for this changed behavior. It acts like some kind of dynamic resolution
change artifact, but we are not doing anything to cause this AFAICT.

It is probably something squirrelly in our scripts that we just haven¹t
seen. But I thought I go for a long shot and ask if there have been any
changes that might affect the behavior of r.flow or r.slope.aspect recently.
I didn¹t see anything obvious in the commit list archives, but may have
missed a library change of some kind.

$ cvs diff -D 2007-07-01 -D 2007-07-15 raster/r.mapcalc raster/r.slope.aspect
raster/r.flow
cvs server: Diffing raster/r.mapcalc
cvs server: Diffing raster/r.slope.aspect
cvs server: Diffing raster/r.flow

Assuming that your dates are correct, it isn't the programs. Extending
the comparison to lib/{gis,vector,bitmap,btree,rowio,segment} (the
only libraries which those programs use directly) doesn't show
anything which looks remotely relevant.

I can only suggest checking out and compiling "before" and "after"
versions, and examining any intermediate data produced by your
scripts.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

Helena also mentioned r.neighbors as a possible culprit, as I'd forgotten to
mention that we also use this to smooth results.

Neither r.neighbors or lib/stats changed during the first half of July.

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

Tests over the last couple days suggest that r.neighbors may be the, or one
of, the causes. We lose most of the artifacts if we turn off smoothing using
r.neighbors, and the artifacts are much worse with a neighborhood of 7x7
than 3x3.

We're probably wrong about the date, however. This seems to only show up
clearly in very long runs (a simulation of 50 recursive models) and is most
pronounced with larger smoothing neighborhoods. Previously we'd done a small
neighborhood of 3x3 and done most of our tests for no more than 10
iterations. We only did a couple of long ones and were looking more at stats
from the output than the maps themselves. Now we are doing a number of 50+
iteration runs (the most recent one ran for nearly 600 years simulated
time).

using a median smoother gives much worse results than a mean smoother,
though a median ought to be better the larger the neighborhood is, since it
should not be affected as much by extreme values.

Our next test it to find out if there is some kind of issue with region
setting that is interacting with the smoothing. Probably not, but we need to
make sure.

I'll report more later after. If anyone is interested in seeing what I've
tried to describe, I've posted images of the effects of different smoothing
parameters at:

<http://www.public.asu.edu/~cmbarton/files/LandDyn&gt;

Michael

On 8/1/07 7:34 PM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

Helena also mentioned r.neighbors as a possible culprit, as I'd forgotten to
mention that we also use this to smooth results.

Neither r.neighbors or lib/stats changed during the first half of July.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

That is what I call a gridding effect -
check it out here - and this was done only with simwe and r.mapcalc - no r.flow, r.neighbors, etc.
(some of the problem actually may be hidden in r.mapcalc computations)

http://skagit.meas.ncsu.edu/~helena/gmslab/oxfordms/oxf.html

you can see how it evolves here (you start seeing it after several iterations - it is there from the start
but it is small at the beginning and just grows over time until you start seeing it)

http://skagit.meas.ncsu.edu/~helena/gmslab/oxfordms/eledhedge.gif

I thought that smoothing will take care of it but apparently just the opposite
We need to find out what to do about it - I think that
  it is more related to the modeling approach than GRASS itself,

BTW as I understand it median is computed by dividing the values into a finite number of classes
so the median from a continuous field would not change continuously from cell to cell and
the effect you are getting could be expected especially for larger number of cells.
Mean is computed directly from the floating point values - is that right?

Helena

Helena Mitasova
Dept. of Marine, Earth and Atm. Sciences
1125 Jordan Hall, NCSU Box 8208,
Raleigh NC 27695
http://skagit.meas.ncsu.edu/~helena/

On Aug 1, 2007, at 11:39 PM, Michael Barton wrote:

Tests over the last couple days suggest that r.neighbors may be the, or one
of, the causes. We lose most of the artifacts if we turn off smoothing using
r.neighbors, and the artifacts are much worse with a neighborhood of 7x7
than 3x3.

We're probably wrong about the date, however. This seems to only show up
clearly in very long runs (a simulation of 50 recursive models) and is most
pronounced with larger smoothing neighborhoods. Previously we'd done a small
neighborhood of 3x3 and done most of our tests for no more than 10
iterations. We only did a couple of long ones and were looking more at stats
from the output than the maps themselves. Now we are doing a number of 50+
iteration runs (the most recent one ran for nearly 600 years simulated
time).

using a median smoother gives much worse results than a mean smoother,
though a median ought to be better the larger the neighborhood is, since it
should not be affected as much by extreme values.

Our next test it to find out if there is some kind of issue with region
setting that is interacting with the smoothing. Probably not, but we need to
make sure.

I'll report more later after. If anyone is interested in seeing what I've
tried to describe, I've posted images of the effects of different smoothing
parameters at:

<http://www.public.asu.edu/~cmbarton/files/LandDyn&gt;

Michael

On 8/1/07 7:34 PM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

Helena also mentioned r.neighbors as a possible culprit, as I'd forgotten to
mention that we also use this to smooth results.

Neither r.neighbors or lib/stats changed during the first half of July.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

On 8/1/07 10:24 PM, "Helena Mitasova" <hmitaso@unity.ncsu.edu> wrote:

That is what I call a gridding effect -
check it out here - and this was done only with simwe and r.mapcalc
- no r.flow, r.neighbors, etc.
(some of the problem actually may be hidden in r.mapcalc computations)

http://skagit.meas.ncsu.edu/~helena/gmslab/oxfordms/oxf.html

you can see how it evolves here (you start seeing it after several
iterations - it is there from the start
but it is small at the beginning and just grows over time until you
start seeing it)

http://skagit.meas.ncsu.edu/~helena/gmslab/oxfordms/eledhedge.gif

These displays match what we are getting. Apparently, this only begins to
show up with long, recursive simulations.

One characteristic of recursive cellular automata simulations, of course, is
emergent patterns that are not easily predicted by the original
processes--including cascades of deviation amplifying effects. I didn't
expect this to pop up in this part of our simulation, but in retrospect, it
is not all that surprising. But is is annoying and we'll have to correct it
somehow.

I thought that smoothing will take care of it but apparently just the
opposite

This does seem strange. Even odder is the fact that the anomalies get bigger
the larger the smoothing neighborhood. Usually, increasing the neighborhood
creates greater smoothing.

It might be the way we are using smoothing. We check to see if an
erosion/deposition value exceeds the smoothed value by some user-determined
amount. It is does, the smoothed value is used instead of the original
value. The goal was to get rid of extreme values but leave the rest alone.

Perhaps we should just use the smoothed value for all cells instead of just
extreme ones.

We need to find out what to do about it - I think that
  it is more related to the modeling approach than GRASS itself,

Is this perhaps a function of dx and dy calculations?

BTW as I understand it median is computed by dividing the values into
a finite number of classes
so the median from a continuous field would not change continuously
from cell to cell and
the effect you are getting could be expected especially for larger
number of cells.

I don't know how r.neighbors computes the median. However, it is supposed to
be a value such that half the neighbor cells have larger values and half
have smaller values. I would have thought that the way to compute would be
to order the neighborhood values from smallest to largest and find the
middle. Since the neighborhood is always an odd number, the middle is always
a single cell.

Mean is computed directly from the floating point values - is that
right?

This ought to be the sum of the neighbor cell values divided by the number
of cells.

The median is less affected by a single extreme cell than the mean, which is
why we preferred the median. However, if the median calculation algorithm is
problematic, then we'll have to switch to the mean.

Thanks for looking at this. I'll be very interested to hear any further
ideas you have.

Michael

Helena

Helena Mitasova
Dept. of Marine, Earth and Atm. Sciences
1125 Jordan Hall, NCSU Box 8208,
Raleigh NC 27695
http://skagit.meas.ncsu.edu/~helena/

On Aug 1, 2007, at 11:39 PM, Michael Barton wrote:

Tests over the last couple days suggest that r.neighbors may be
the, or one
of, the causes. We lose most of the artifacts if we turn off
smoothing using
r.neighbors, and the artifacts are much worse with a neighborhood
of 7x7
than 3x3.

We're probably wrong about the date, however. This seems to only
show up
clearly in very long runs (a simulation of 50 recursive models) and
is most
pronounced with larger smoothing neighborhoods. Previously we'd
done a small
neighborhood of 3x3 and done most of our tests for no more than 10
iterations. We only did a couple of long ones and were looking more
at stats
from the output than the maps themselves. Now we are doing a number
of 50+
iteration runs (the most recent one ran for nearly 600 years simulated
time).

using a median smoother gives much worse results than a mean smoother,
though a median ought to be better the larger the neighborhood is,
since it
should not be affected as much by extreme values.

Our next test it to find out if there is some kind of issue with
region
setting that is interacting with the smoothing. Probably not, but
we need to
make sure.

I'll report more later after. If anyone is interested in seeing
what I've
tried to describe, I've posted images of the effects of different
smoothing
parameters at:

<http://www.public.asu.edu/~cmbarton/files/LandDyn&gt;

Michael

On 8/1/07 7:34 PM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

Helena also mentioned r.neighbors as a possible culprit, as I'd
forgotten to
mention that we also use this to smooth results.

Neither r.neighbors or lib/stats changed during the first half of
July.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

Tests over the last couple days suggest that r.neighbors may be the, or one
of, the causes. We lose most of the artifacts if we turn off smoothing using
r.neighbors, and the artifacts are much worse with a neighborhood of 7x7
than 3x3.

We're probably wrong about the date, however. This seems to only show up
clearly in very long runs (a simulation of 50 recursive models) and is most
pronounced with larger smoothing neighborhoods. Previously we'd done a small
neighborhood of 3x3 and done most of our tests for no more than 10
iterations. We only did a couple of long ones and were looking more at stats
from the output than the maps themselves. Now we are doing a number of 50+
iteration runs (the most recent one ran for nearly 600 years simulated
time).

using a median smoother gives much worse results than a mean smoother,
though a median ought to be better the larger the neighborhood is, since it
should not be affected as much by extreme values.

The nonlinear nature of the median can exacerbate discontinuities. The
median chooses a single value from the set of inputs; it doesn't
average them.

This can make it a poor choice for a smoothing filter, particularly if
you are going to calculate derivatives (e.g. r.slope.aspect) on the
result.

If you want to eliminate outliers, a better option might be an average
over some quantile range (e.g. 20%-80%), or a weighted average where
values near the median get a larger weight.

If you're trying to remove spatial noise, you normally want a filter
which gives greater weight to the values nearer the centre. To avoid
"contamination" by noise, one method is to use a variety of linear
filters and take the median of the results.

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

Michael Barton wrote:

> I thought that smoothing will take care of it but apparently just the
> opposite

This does seem strange. Even odder is the fact that the anomalies get bigger
the larger the smoothing neighborhood. Usually, increasing the neighborhood
creates greater smoothing.

That's the case if the neighbours are being averaged, but it doesn't
hold for a quantile.

It might be the way we are using smoothing. We check to see if an
erosion/deposition value exceeds the smoothed value by some user-determined
amount. It is does, the smoothed value is used instead of the original
value. The goal was to get rid of extreme values but leave the rest alone.

That will typically cause gradient discontinuities as you switch
between filtered and unfiltered data. That wouldn't necessarily be an
issue in itself, but it will be if you start calculating derivatives
on the result.

Actually, ignore my earlier suggestion of using a
median-of-filtered-values approach. That will have similar problems
regarding derivatives.

Moreover, any mechanism that has a discrete "cut" (i.e. if <criterion>
then <value1> else <value2>) will have the same issue. You need to
find some way to introduce a continuous transition. E.g. rather than
selecting filtered vs unfiltered based upon a threshold, interpolate
between filtered and unfiltered based upon the distance from the
threshold, turning the discrete step into an ogive (S-curve).

Even then, any kind of nonlinearity risks introducing artifacts when
it occurs as part of a complex process, particularly if it's
iterative.

In general, iterative processes won't end up being stable just because
you want them to be stable. If there's any kind of "amplification"
involved, instability is likely to be the "default" behaviour;
stability has to be forced.

> BTW as I understand it median is computed by dividing the values into
> a finite number of classes
> so the median from a continuous field would not change continuously
> from cell to cell and
> the effect you are getting could be expected especially for larger
> number of cells.

I don't know how r.neighbors computes the median. However, it is supposed to
be a value such that half the neighbor cells have larger values and half
have smaller values. I would have thought that the way to compute would be
to order the neighborhood values from smallest to largest and find the
middle. Since the neighborhood is always an odd number, the middle is always
a single cell.

This is essentially what happens, although if there are any null
cells, the median may be computed over an even number of cells
resulting in the mean of the two middle cells. The actual
implementation (from lib/stats/c_median.c) is:

  void c_median(DCELL *result, DCELL *values, int n)
  {
    n = sort_cell(values, n);
  
    if (n < 1)
      G_set_d_null_value(result, 1);
    else
      *result = (values[(n-1)/2] + values[n/2]) / 2;
  }

where sort_cell() sorts the cell values into ascending order.

[If the number of cells is odd, (n-1)/2 and n/2 will be equal, so the
middle value is averaged with itself.]

> Mean is computed directly from the floating point values - is that
> right?

This ought to be the sum of the neighbor cell values divided by the number
of cells.

The median is less affected by a single extreme cell than the mean, which is
why we preferred the median. However, if the median calculation algorithm is
problematic, then we'll have to switch to the mean.

It's problematic if you're going to differentiate the result, because
it's a discontinuous function; it's essentially a tree of if/then/else
choices.

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

Glynn and Helena,

Thanks very much for the information and advice.

So the upshot is, for this kind of simulation modeling, we avoid using a
cut-off function and avoid using a median filter. Both create slight but
abrupt discontinuities in the landscape which can amplify in recursive runs.

I think we'll try just using the r.neighbors smoothed map (using a mean
smoother) rather than some function of the original and smoothed map.

I have a couple further questions before we launch into revising these
scripts.

1. Is there any speed advantage/disadvantage to using r.mfilter over
r.neighbors? At the moment, because we are trying to get rid of occasional
spikes, I don't think that a weighted mean will give better results than an
unweighted mean. But the only way to get functions like a weighted mean is
to use r.mfilter, so it would be good to know if there are any known
performance hits or if it is faster.

2. We are starting with a depressionless DEM, but now I'm wondering if it
would be a good idea to fill any depressions between each iteration of the
simulation. If so, it looks like r.terraflow builds this into its routine
and could save considerable time over running r.fill.dir each time. Any
thoughts along those lines?

Thanks again
Michael

On 8/2/07 5:14 AM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

I thought that smoothing will take care of it but apparently just the
opposite

This does seem strange. Even odder is the fact that the anomalies get bigger
the larger the smoothing neighborhood. Usually, increasing the neighborhood
creates greater smoothing.

That's the case if the neighbours are being averaged, but it doesn't
hold for a quantile.

It might be the way we are using smoothing. We check to see if an
erosion/deposition value exceeds the smoothed value by some user-determined
amount. It is does, the smoothed value is used instead of the original
value. The goal was to get rid of extreme values but leave the rest alone.

That will typically cause gradient discontinuities as you switch
between filtered and unfiltered data. That wouldn't necessarily be an
issue in itself, but it will be if you start calculating derivatives
on the result.

Actually, ignore my earlier suggestion of using a
median-of-filtered-values approach. That will have similar problems
regarding derivatives.

Moreover, any mechanism that has a discrete "cut" (i.e. if <criterion>
then <value1> else <value2>) will have the same issue. You need to
find some way to introduce a continuous transition. E.g. rather than
selecting filtered vs unfiltered based upon a threshold, interpolate
between filtered and unfiltered based upon the distance from the
threshold, turning the discrete step into an ogive (S-curve).

Even then, any kind of nonlinearity risks introducing artifacts when
it occurs as part of a complex process, particularly if it's
iterative.

In general, iterative processes won't end up being stable just because
you want them to be stable. If there's any kind of "amplification"
involved, instability is likely to be the "default" behaviour;
stability has to be forced.

BTW as I understand it median is computed by dividing the values into
a finite number of classes
so the median from a continuous field would not change continuously
from cell to cell and
the effect you are getting could be expected especially for larger
number of cells.

I don't know how r.neighbors computes the median. However, it is supposed to
be a value such that half the neighbor cells have larger values and half
have smaller values. I would have thought that the way to compute would be
to order the neighborhood values from smallest to largest and find the
middle. Since the neighborhood is always an odd number, the middle is always
a single cell.

This is essentially what happens, although if there are any null
cells, the median may be computed over an even number of cells
resulting in the mean of the two middle cells. The actual
implementation (from lib/stats/c_median.c) is:

void c_median(DCELL *result, DCELL *values, int n)
{
n = sort_cell(values, n);

if (n < 1)
G_set_d_null_value(result, 1);
else
*result = (values[(n-1)/2] + values[n/2]) / 2;
}

where sort_cell() sorts the cell values into ascending order.

[If the number of cells is odd, (n-1)/2 and n/2 will be equal, so the
middle value is averaged with itself.]

Mean is computed directly from the floating point values - is that
right?

This ought to be the sum of the neighbor cell values divided by the number
of cells.

The median is less affected by a single extreme cell than the mean, which is
why we preferred the median. However, if the median calculation algorithm is
problematic, then we'll have to switch to the mean.

It's problematic if you're going to differentiate the result, because
it's a discontinuous function; it's essentially a tree of if/then/else
choices.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics and Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

1. Is there any speed advantage/disadvantage to using r.mfilter over
r.neighbors? At the moment, because we are trying to get rid of occasional
spikes, I don't think that a weighted mean will give better results than an
unweighted mean. But the only way to get functions like a weighted mean is
to use r.mfilter, so it would be good to know if there are any known
performance hits or if it is faster.

r.mfilter currently works with integers. AFAICT, there's no
fundamental reason why it can't be upgraded to FP.

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

Then, as is it won't work for us. If it can be upgraded to integers we can
try it.

Thanks for the info. Saved us some time in fruitless testing.

Now that we understand what is going on, we're getting better results. We
will still need to conduct more tests, but will report back to the list with
the results.

Michael

On 8/2/07 5:31 PM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

1. Is there any speed advantage/disadvantage to using r.mfilter over
r.neighbors? At the moment, because we are trying to get rid of occasional
spikes, I don't think that a weighted mean will give better results than an
unweighted mean. But the only way to get functions like a weighted mean is
to use r.mfilter, so it would be good to know if there are any known
performance hits or if it is faster.

r.mfilter currently works with integers. AFAICT, there's no
fundamental reason why it can't be upgraded to FP.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

Michael Barton wrote:

>> 1. Is there any speed advantage/disadvantage to using r.mfilter over
>> r.neighbors? At the moment, because we are trying to get rid of occasional
>> spikes, I don't think that a weighted mean will give better results than an
>> unweighted mean. But the only way to get functions like a weighted mean is
>> to use r.mfilter, so it would be good to know if there are any known
>> performance hits or if it is faster.
>
> r.mfilter currently works with integers. AFAICT, there's no
> fundamental reason why it can't be upgraded to FP.

Then, as is it won't work for us. If it can be upgraded to integers we can
try it.

I've added an FP version, r.mfilter.fp.

I didn't modify r.mfilter as it turns out that there is some behaviour
which is specific to integer maps, and to the zero-is-null behaviour
of 4.x. It wasn't practical to add support for FP and null values
without breaking compatibility.

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

Thanks Glynn. We'll certainly give it a try.

Michael

On 8/3/07 3:18 AM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

1. Is there any speed advantage/disadvantage to using r.mfilter over
r.neighbors? At the moment, because we are trying to get rid of occasional
spikes, I don't think that a weighted mean will give better results than an
unweighted mean. But the only way to get functions like a weighted mean is
to use r.mfilter, so it would be good to know if there are any known
performance hits or if it is faster.

r.mfilter currently works with integers. AFAICT, there's no
fundamental reason why it can't be upgraded to FP.

Then, as is it won't work for us. If it can be upgraded to integers we can
try it.

I've added an FP version, r.mfilter.fp.

I didn't modify r.mfilter as it turns out that there is some behaviour
which is specific to integer maps, and to the zero-is-null behaviour
of 4.x. It wasn't practical to add support for FP and null values
without breaking compatibility.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton

I wanted to report that we managed to get rid of the artifacts and produce a
very nice looking map, with a 50-year recursive landscape evolution
simulation.

Surprisingly, we ended up using a median filter after all, but without the
cutoff function. I realized that we needed to 'despeckle' the map in order
to remove the spikes that were the original problem we were trying to fix
when we generated the gridded peaks artifacts.

We're running more tests to see what neighborhood size works best. 7x7
looked really nice for most things, but was a little rough in areas where
grazing animals had randomly denuded cells of vegetation and it had
regenerated to a varying amount (though still far better than any other
result we'd had so far). Trying 9x9 now.

Also running it with the coupled agent based model of village farming over
the weekend to see how it performs.

I'll try to post some pictures of the results next week.

This updated version of r.landscape.evol is now on the GRASS Addons SVN
server.

Thanks to both of you for helping us to sort out this mysterious simulation
effect.

Michael

On 8/3/07 3:18 AM, "Glynn Clements" <glynn@gclements.plus.com> wrote:

Michael Barton wrote:

1. Is there any speed advantage/disadvantage to using r.mfilter over
r.neighbors? At the moment, because we are trying to get rid of occasional
spikes, I don't think that a weighted mean will give better results than an
unweighted mean. But the only way to get functions like a weighted mean is
to use r.mfilter, so it would be good to know if there are any known
performance hits or if it is faster.

r.mfilter currently works with integers. AFAICT, there's no
fundamental reason why it can't be upgraded to FP.

Then, as is it won't work for us. If it can be upgraded to integers we can
try it.

I've added an FP version, r.mfilter.fp.

I didn't modify r.mfilter as it turns out that there is some behaviour
which is specific to integer maps, and to the zero-is-null behaviour
of 4.x. It wasn't practical to add support for FP and null values
without breaking compatibility.

__________________________________________
Michael Barton, Professor of Anthropology
Director of Graduate Studies
School of Human Evolution & Social Change
Center for Social Dynamics & Complexity
Arizona State University

phone: 480-965-6213
fax: 480-965-7671
www: http://www.public.asu.edu/~cmbarton