[GRASS-user] r.mapcalc median filters entire map instead of only null values

Patton, Eric epatton at nrcan.gc.ca
Tue Jan 16 16:46:37 EST 2007


I'm trying to use r.neighbors and r.mapcalc in a script to fill holes in a raster where there are null values. A temporary raster is created with r.neighbors that contains a 3x3 median-interpolated surface of the primary DEM. r.mapcalc should only use the value in this temporary raster where a null value occurs in the primary DEM.

Here's the relevant section from my script:

# Procedure for filling nulls with multiple median filter passes

if [ "$PASSES" -gt 1 ] ; then

 COUNT=1
 while [ "$COUNT" -le "$PASSES" ] ; do
					
  r.neighbors --o input=${MAP} output=fill.tmp method=median size=3 
  r.mapcalc "${MAP}_fill = if(isnull(${MAP}), fill.tmp, ${MAP})"
				
  MAP=${MAP}_fill
  COUNT=$(($COUNT + 1))
  g.remove rast=fill.tmp 

 done

else
 echo ""
 r.neighbors --o input=${MAP} output=fill.tmp method=median size=3
 r.mapcalc "${MAP}_fill = if(isnull("$MAP"), fill.tmp, "$MAP")"
 MAP=${MAP}_fill
 g.remove rast=fill.tmp 

fi

My problem is that when my script is run, the output $MAP_fill raster has it's entire surface smoothed with the median-filtered one, not just the null values like I intended. But if I run the r.neighbors and r.mapcalc lines in isolation outside of the script, I get the intended null-filled DEM with no unwanted smoothing. Am I not quoting the r.mapcalc line properly? 

FWIW, I've also tried different quoting methods for the r.mapcalc line, with the same problem occurring each time:

r.mapcalc ${MAP}_fill = "if(isnull($"MAP"), fill.tmp, $"MAP")"
r.mapcalc "${MAP}_fill = if(isnull(${MAP}), fill.tmp, ${MAP})"
r.mapcalc ${MAP}_fill = "if(isnull($MAP), fill.tmp, $MAP)"

Any ideas?

~ Eric.




More information about the grass-user mailing list