[GRASS-user] Grass, Python, mapcalc
Paul Meems
bontepaarden at gmail.com
Fri Dec 28 11:28:21 PST 2012
I understand I should use GRASS with Ubuntu to get the best results, but
I'm not only new with GRASS and Python but also with Ubuntu ;)
Ubuntu server I know a little and I have learned to appreciate it very
much, but Ubuntu desktop is completely new for me.
That is the reason I wanted to start on Window, because I feel more
comfortable with it.
It seems getting the annual solar radiation is working now. I need to
update my script for using r.series.
I'm going to run a full year now. I have the ArcMap results for the same
area as well and I first want to compare those.
After that I'm going to look at the rest of the ArcMap model, since the
solar radiation is only one of the steps I need to convert to GRASS/Python.
When I have reproduced all steps in GRASS I'm going to create a Python
script of it and then I can run the complete area on my Ubuntu server.
Thanks all for now and I will be back ;)
Paul
*Paul Meems *
Release manager, configuration manager
and forum moderator of MapWindow GIS.
www.mapwindow.org
Owner of MapWindow.nl - Support for
Dutch speaking users.
www.mapwindow.nl
*
*
2012/12/28 Daniel Lee <lee at isi-solutions.org>
> Hi Paul,
>
> I think one of your problems is... using windows ;) I have heard that the
> Windows version of GRASS works a LOT better than it used to, but GRASS was
> actually the reason I switched to Linux a couple of years ago. Back then it
> just ran better in Linux - a lot better. Since then, I've never gone back
> to Windows because I'm very satisfied with Linux, not just with GRASS but
> with everything else I do with the computer. Of course, Linux vs. Windows
> is a whole different story and doesn't belong in this list, but that was my
> experience.
>
> If you're wanting to use GRASS in a headless Ubuntu server anyway, why
> don't you set up a Linux machine or at least a Linux installation in a
> virtual machine and play around with it there? There everything should work
> well and the development environment would match your production
> environment. OSGeo has a great live installation including all of the
> latest Geo software, Python, etc. - everything the;heart could desire.
> since I'm out and about at the moment I can't send you the address, but
> just Google for an OSGeo live installation and you'll see it right away.
>
> Hope that helps,
> Daniel
> Am 28.12.2012 15:28 schrieb "Paul Meems" <bontepaarden at gmail.com>:
>
> Hi Sotiris,
>>
>> 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
>> parameters.
>>
>> 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)
>>
>> # Finalize:
>> 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*")
>> print "Done"
>>
>> Thanks
>>
>> Paul
>>
>> 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.
>>> best,
>>> sotiris
>>>
>>> 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
>>> with
>>> > 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
>>>
>>>
>>
>> _______________________________________________
>> 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...
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20121228/ca8bff69/attachment-0001.html>
More information about the grass-user
mailing list