[GRASS-user] Get non-NULL cells over multiple rasters

[Johannes Radinger]

Hi again...

-------- Original-Nachricht --------
> Datum: Thu, 01 Mar 2012 15:18:24 +0100
> Von: Moritz Lennert <mlennert at club.worldonline.be>
> An: Johannes Radinger <JRadinger at gmx.at>
> CC: grass-user at lists.osgeo.org
> Betreff: Re: [GRASS-user] Get non-NULL cells over multiple rasters

> Johannes Radinger wrote:
> > Hi,
> >
> > -------- Original-Nachricht --------
> >    
> >> Datum: Thu, 01 Mar 2012 13:18:30 +0100
> >> Von: Moritz Lennert<mlennert at club.worldonline.be>
> >> An: Johannes Radinger<JRadinger at gmx.at>
> >> CC: GRASS user list<grass-user at lists.osgeo.org>
> >> Betreff: Re: [GRASS-user] Get non-NULL cells over multiple rasters
> >>      
> >    
> >> Johannes Radinger wrote:
> >>      
> >>> Hi,
> >>>
> >>> I have a bunch of raster maps (~50 are named crast_* and 50 are named
> >>> rast_*). Is there a simple way (eg. using regular expressions) to
> >>> process them with the mapcalculator to get one resulting map which
> >>> shows only those cells (value 1) which are non-NULL cell in ALL maps?
> >>> Like I'd like the cells which all rasters have commonly populated.
> >>> Maybe there is a simple way to do this in the mapcalculator (or any
> >>> other modul) resp. in a python script?? Suggestions are mostly
> >>> welcome. Thanks,
> >>>        
> >> r.series with the -n flag ? From the man page:
> >>
> >> "With -n flag, any cell for which any of the corresponding input cells
> >> are NULL is automatically set to NULL (NULL propagation). The aggregate
> >> function is not called, so all methods behave this way with respect to
> >> the -n flag."
> >>
> >> See the man page for examples of how to use r.series in conjunction
> with
> >> g.mlist to work with a series of maps.
> >>      
> >
> > Thank you very  much, r.series is the right tool for my task. I tried to
> > use it in combination with g.mlist but somehow I don't manage it. I
> tried different methods e.g.:
> >
> > 1) in the GRASS command console I typed:
> > r.series -n input="'g.mlist pattern='rast_*|crast_*' sep=,'"
> output=Treene_tmp_mask method=count
> >    
> 
> Are you talking about the wxGUI command console, or the terminal ?

I was trying it with the wxGUI command console, but I now I'll focus
on the python script...

> 
> In the terminal the whole g.mlist part has to be in backticks (`) which 
> means "replace everything which is within the backticks with the results 
> of that command". I imagine that this does not work in the GUI command 
> console.
> 
> In the terminal your command should be:
> 
> r.series -n input=`g.mlist pattern="rast_*|crast_*" sep=,`
> output=Treene_tmp_mask method=count
> 
> 
> > 2) and I tried it in a python script with:
> > 	grass.run_command("r.series",
> > 			input = "'g.mlist pattern='rast_*|crast_*|Treene_raster' sep=,'",
> > 			overwrite=True,
> > 			flags = "n",
> > 			output = "Treene_tmp_mask",
> > 			method = "count")
> >    
> 
> 
> If you want to work in a python script, you have to first run 
> read_command(g.mlist), put the result in a comma-seperated string which 
> you then feed into r.series. There are possibly other ways that I don't 
> know about...

I followed your suggestions and tried to split the command up into a part
with read_command(g.mlist) and a part with grass.run_command("r.series"). This nearly works, there is just one problem with the read_command(g.mlist). It seems that the output after the last listed map includes a linebreak or something similar. Thus using the output of read_command causes a problem when r.series tries to load the last map.
I made a reproducable script for the spearfish dataset:

#!/usr/bin/env python

import sys

import grass.script as grass
import grass.script.setup as gsetup

def main():	

	gmlist = grass.read_command("g.mlist",
			pattern = "soils.*",
			sep=",")
	print gmlist

	grass.run_command("r.series",
					input = gmlist,
					overwrite=True,
					flags = "n",
					output = "test_output",
					method = "count")
	return 0

if __name__ == "__main__":
	sys.exit(main())


cheers,

/johannes

> 
> Moritz

-- 
NEU: FreePhone 3-fach-Flat mit kostenlosem Smartphone!                                  
Jetzt informieren: http://mobile.1und1.de/?ac=OM.PW.PW003K20328T7073a