[GRASS-user] Grass, Python, mapcalc
bontepaarden at gmail.com
Fri Dec 28 06:27:42 PST 2012
Thanks for your help.
I will have a look at the r.series.
Meanwhile I've discovered the graphical modeler. It is not working 100% in
all cases but it does generates a Python script with all the necessary
This is my script so far, which seems to be working:
for x in range(1, 10, 1):
print "Working on day %d" % (x)
# Create radiation for each day:
glob_rad = 'global' + str(x)
grass.run_command('r.sun', flags = 's', overwrite = True, elevin =
'w001001 at temp', aspin = 'aspect at temp', slopein = 'slope at temp', glob_rad =
glob_rad, day = x)
# Combine the radiation file with the previous one:
if x == 2:
grass.run_command("r.mapcalculator", overwrite = True, amap =
"global1", bmap = "global2", formula = "A+B",
outfile = "suntemp2")
elif x > 2:
# Use temp name or else it won't work:
amap = "suntemp" + str(x-1)
bmap = "global" + str(x) + "@temp"
outfile = "suntemp" + str(x)
grass.run_command("r.mapcalculator", overwrite = True, amap = amap,
bmap = bmap, formula = "A+B",
outfile = outfile)
grass.mapcalc('sun=' + outfile, quiet=False, verbose=True, overwrite=True)
# Clean up:
grass.run_command('g.mremove', flags = 'f', rast = "suntemp*")
grass.run_command('g.mremove', flags = 'f', rast = "global*")
2012/12/28 S. Koukoulas (lists) <sotkouklistes at gmail.com>
> Dear Paul,
> having said all that for r.mapcalc, I wonder whether you really need
> mapcalc. In your case,
> couldn't you just use the * r.series * command with Sum (and perhaps
> g.mlist can help you make a list of all the daily maps).
> I think that it might be faster as well.
> On 12/28/2012 03:08 PM, Paul Meems wrote:
> > Thanks Sotiris,
> > I've tried your suggestion but converted it to Python because I'm on
> > Windows.
> > This is my script:
> > grass.run_command('g.mremove', rast = ('out*'))
> > grass.mapcalc('sun=global1', quiet=False, verbose=False, overwrite=True)
> > for x in range(2, 10, 1):
> > exp = 'out=sun+global' + str(x)
> > print exp
> > grass.mapcalc(exp, quiet=False, verbose=False, overwrite=True)
> > grass.run_command('g.rename', flag='o', rast = ('out','sun'))
> > This doesn't result in the expected outcome.
> > sun is equal to global1 and no values are added
> > out is global1+global2 but is no more values are added. But mapcalc is
> > called several times.
> > I've tested this with:
> > r.what --v -f -n input=global1 east_north=156287.225281,383024.737026
> > r.what --v -f -n input=sun east_north=156287.225281,383024.737026
> > r.what --v -f -n input=out east_north=156287.225281,383024.737026
> > Not only am I a new user to GRASS but I'm also still a novice in
> > Python. Most likely it is something in my script.
> > If somebody has a different approach to get the solar radiation of a
> > whole year I'm happy to try that.
> > Thanks,
> > Paul
> > Date: Fri, 28 Dec 2012 14:36:15 +0200
> > From: "S. Koukoulas (lists)" <sotkouklistes at gmail.com
> > <mailto:sotkouklistes at gmail.com>>
> > To: grass-user at lists.osgeo.org <mailto:grass-user at lists.osgeo.org>
> > Subject: Re: [GRASS-user] Grass, Python, mapcalc
> > Message-ID: <50DD923F.806 at gmail.com <mailto:50DD923F.806 at gmail.com>>
> > Content-Type: text/plain; charset=windows-1252
> > To my experience with r.mapcalc you cannot do what yo do with
> > programming languages... e.g. i=i+1 will not work directly.
> > Instead you
> > could try a more indirect way, something like (here is an example
> > bash shell script):
> > #g.mremove rast=out* (to remove previous stuff, if needed)
> > r.mapcalc sun=global1
> > for a in ?seq 2 10?; do
> > r.mapcalc out = sun + global$a
> > g.rename --o rast=out,sun
> > done
> > basically the main idea is to use an intermediate map and then before
> > the loop ends, pass it (with g.rename) to the map variable you want.
> > There might be a more elegant way with python, but hopefully this
> > shell
> > script example might help you.
> > HTH,
> > sotiris
> > _______________________________________________
> > grass-user mailing list
> > grass-user at lists.osgeo.org
> > http://lists.osgeo.org/mailman/listinfo/grass-user
-------------- next part --------------
An HTML attachment was scrubbed...
More information about the grass-user