PyWPS / GRASS analýza viditelnosti / Verze 1

Z FreeGIS portál
Přejít na: navigace, hledání

Tato stránka popisuje první verzi WPS služby pro výpočet viditelnosti. Implementováno pomocí PyWPS a GRASS GIS. GRASS lokace 'eu-dem-merc' je ke stažení zde.

Tato verze nemá definovány žádné parametry, výstupem je výměra (v mapových jednotkách) území viditelného z daného bodu.

Skript pro poskytování služby

Skript v jazyku Python je úmístěn v adresáři processes. Více o psaní skriptů pro systém GRASS zde.

Definice třídy Process

Začneme deficí třídy Process odvozené z třídy WPSProcess, která je definována aplikací PyWPS. Soubor se skriptem pojmenujeme viewshed-1.py, výsledek je ke stažení zde.

# -*- coding: utf-8 -*-

from pywps.Process import WPSProcess

class Process(WPSProcess):
     def __init__(self, location='eu-dem-merc'):
          WPSProcess.__init__(self,
                              identifier="viewshed-1", # odpovídá názvu souboru
                              version = "0.1",       # verze služby
                              title="GRASS GIS: Výpočet viditelnosti",
                              abstract="Ukázka WPS služby pro výpočet viditelnosti GRASS nástrojem r.viewshed. " \
                                   "Služba nepřijímá vstupní parametry a vrací výměru viditelného území v hektarech.",
                              grassLocation=location)

          # proměnná do které se uloží výsledek analýzy
          self.answer=self.addLiteralOutput(identifier = "visible_area",
                                            title = "Výměra viditelného území v hektarech")

Test definice služby

Definici služby můžeme otestovat buď z výstupu GetCapabilities anebo přímo pomocí DescribeProcess.

Příklad z příkazové řádky:

 /usr/local/bin/pywps.py "service=wps&version=1.0.0&request=GetCapabilities"

Příklad z webového prohlížeče:

http://geo102.fsv.cvut.cz/cgi-bin/pywps.cgi?service=wps&version=1.0.0&request=GetCapabilities
...
<wps:Process wps:processVersion="0.1">
 <ows:Identifier>viewshed-1</ows:Identifier>
 <ows:Title>GRASS GIS: Výpočet viditelnosti</ows:Title>
 <ows:Abstract>Ukázka WPS služby pro výpočet viditelnosti GRASS nástrojem r.viewshed. Služba nepřijímá vstupní parametry a vrací výměru viditelného území v hektarech</ows:Abstract>
</wps:Process>
...
 http://geo102.fsv.cvut.cz/cgi-bin/pywps.cgi?service=wps&version=1.0.0&request=DescribeProcess&identifier=viewshed-1
...
<ProcessDescription wps:processVersion="0.1" storeSupported="false" statusSupported="true">
 <ows:Identifier>viewshed-1</ows:Identifier>
 <ows:Title>GRASS GIS: Výpočet viditelnosti</ows:Title>
 <ows:Abstract>Ukázka WPS služby pro výpočet viditelnosti GRASS nástrojem r.viewshed. Služba nepřijímá vstupní parametry a vrací výměru viditelného území v hektarech</ows:Abstract>
 <ProcessOutputs>
  <Output>
   <ows:Identifier>visible_area</ows:Identifier>
   <ows:Title>Výměra viditelného území v hektarech</ows:Title>
   <LiteralOutput>
    <ows:DataType ows:reference="http://www.w3.org/TR/xmlschema-2/#integer">integer</ows:DataType>
   </LiteralOutput>
  </Output>
 </ProcessOutputs>
</ProcessDescription>
...

Rozšíření skriptu o výpočet viditelnosti

# -*- coding: utf-8 -*-

import os
import sys

### gisbase = os.environ['GISBASE'] = "/home/martin/src/grass_trunk/dist.x86_64-unknown-linux-gnu/"
### gisbase = os.environ['GISBASE'] = "/usr/local/grass-7.0.svn"
gisbase = os.environ['GISBASE'] = "/home/student/grass-geo102/grass_trunk/dist.i686-pc-linux-gnu/"
sys.path.append(os.path.join(os.environ["GISBASE"], "etc", "python"))

from grass.script.core import run_command, read_command

from pywps.Process import WPSProcess

class Process(WPSProcess):
     def __init__(self, location='eu-dem-merc'):
          WPSProcess.__init__(self,
                              identifier="viewshed-1", # odpovídá názvu souboru
                              version = "0.1",       # verze služby
                              title="GRASS GIS: Výpočet viditelnosti",
                              abstract="Ukázka WPS služby pro výpočet viditelnosti GRASS nástrojem r.viewshed. " \
                                   "Služba nepřijímá vstupní parametry a rací výměru viditelného území v hektarech.",
                              grassLocation=location)

          # proměnná do které se uloží výsledek analýzy
          self.answer=self.addLiteralOutput(identifier = "visible_area",
                                            title = "Výměra viditelného území v hektarech")
          
          os.environ['GRASS_VERBOSE'] = '0'
          os.environ['GRASS_SKIP_MAPSET_OWNER_CHECK'] = '1'

          # název výstupní rastrové mapy
          self.view_map = 'view_%d' % os.getpid()
          
     # výpočet viditelnosti v daném výpočetním regionu (offset 5km) a observačního bodu (viewpoint)
     def execute(self, raster='dmt', region = { 'north' : 6547933.83, 'south' : 6537908.83,
                                                'west' : 1545808.36, 'east' : 1555833.36 },
                 viewpoint=(1550820.86,6542921.33)):
          # nastavení regionu
          self.set_region(raster, region)
          # výpočet viditelnosti
          self.viewshed(raster, viewpoint)
          # nastavení výstupní hodnoty
          self.answer.setValue(self.get_visible_area())

     # nastavení výpočetního regionu (rozlišení převzato ze vstupního rastru)
     def set_region(self, raster, region):
          run_command('g.region', rast=raster, n=region['north'], s=region['south'],
                      w=region['west'], e=region['east'])

     # výpočet viditelnosti
     def viewshed(self, elev, coord):
          run_command('r.viewshed', input=elev, output=self.view_map,
                      coord=coord, flags='b')
          
     # určení výměry viditelného území (v mapových jednotkách)
     def get_visible_area(self):
          try:
               value = "%d" % round(float(read_command('r.stats', input=self.view_map, flags='a').splitlines()[0].split(' ', 1)[1]) / 1e4)
          except Exception, e:
               print >> sys.stderr, e
               value = "-1"
          
          return value

# pro účely testování skriptu z příkazové řádky
if __name__ == "__main__":
     process = Process()
     process.execute()
     print "RESULT: %s" % process.get_visible_area()

Ladění služby

TODO

Test skriptu

Příklad spuštění skriptu z adresáře, kde je umístěn.

python viewshed-1.py

TODO: přidat výstup skriptu

Test služby

Z příkazové řádky

export PYWPS_CFG=/usr/local/pywps/pywps.cfg
export PYWPS_PROCESSES=/usr/local/pywps/processes
/usr/local/bin/pywps.py "service=wps&version=1.0.0&request=Execute&identifier=viewshed-1"

Z webové prohlížeče

http://geo102.fsv.cvut.cz/services/yfsgwps?service=wps&version=1.0.0&request=Execute&identifier=viewshed-1

Výsledek

...
    <wps:Status creationTime="2014-01-23T12:51:42Z">
        <wps:ProcessSucceeded>PyWPS Process viewshed successfully calculated</wps:ProcessSucceeded>
    </wps:Status>
    <wps:ProcessOutputs>
        <wps:Output>
            <ows:Identifier>visible_area</ows:Identifier>
            <ows:Title>Výměra viditelného území v hektarech</ows:Title>
            <wps:Data>
                <wps:LiteralData dataType="integer">1028</wps:LiteralData>
            </wps:Data>
        </wps:Output>
    </wps:ProcessOutputs>
...