Work with Spatial Data¶
Last Updated: August 1, 2014
Up to this point in the tutorial, you have not been using the spatial functionality in GsshaPy. This is where most of the progress occurred in version 2.0. In this tutorial you will learn how to read a project using spatial database objects. This means that instead of storing the points, lines, polygons, and rasters as plain text, they will be stored using the spatial field types provided by PostGIS (raster and geometry). Once stored in PostGIS spatial fields, the data is exposed to over 1000 spatial database functions that can be used in queries to convert the data to different formats (e.g.: KML, WKT, GeoJSON, PNG), transform the coordinate reference system, and perform geoprocessing tasks (e.g.: buffer, intersect, union).
Spatial Project Read¶
To read the GSSHA files into spatial fields in the database we setup as we did before, creating a session and an
instance of the
>>> from gsshapy.lib import db_tools as dbt >>> from gsshapy.orm import ProjectFile >>> sqlalchemy_url = dbt.sqlalchemy_url = dbt.init_postgresql_db(username='gsshapy', host='localhost', database='gsshapy_tutorial', port='5432', password='pass') >>> spatialSession = dbt.create_session(sqlalchemy_url) >>> spatialProjectFile = ProjectFile()
Then we call the
readProject() method enabling spatial objects by setting the
spatial argument to
readProject() method looks like this:
>>> readDirectory = '/path_to/tutorial-data' >>> filename = 'parkcity.prj' >>> spatialProjectFile.readProject(directory=readDirectory, projectFileName=filename, session=spatialSession, spatial=True)
You will probably notice that reading the project into the database using spatial objects takes a little more time than when reading it without spatial objects. This delay is caused primarily by the conversion of the native GSSHA spatial formats to the PostGIS spatial formats. For example, GSSHA stores raster files in GRASS ASCII format while PostGIS stores rasters in the binary version of the Well Known Text format (Well Known Binary).
readOutput() methods attempt to automatically lookup the spatial reference
ID using the projection file that can be included with GSSHA models. This uses a web service, so an internet connection
is required. If the projection file is not included or no id can be found, the default spatial reference ID will be used:
4236 WGS 84. Warnings will be thrown in this case, but GsshaPy will not throw an error. If you want to force the spatial
reference ID, pass it in manually using the
Spatial Read for Individual Files¶
You can also apply the spatial read methodology to individual files. Instantiate the file object for the file you wish
to read into the database with spatial objects and call the
read() method with the same spatial arguments as
illustrated above. For example, let’s read in the mask map file.
The file object that is used to read in the mask map is the
gsshapy.orm.RasterMapFile object (use the
Supported Files page to determine what objects are used to read each file). Automatic spatial reference
id lookup is not a feature for the read methods of individual files. Use the
lookupSpatialReferenceID, too look up the srid or enter it manually. Create a new instance
of this object and invoke the
read() method pointing it at the mask map file (
>>> from gsshapy.orm import RasterMapFile, ProjectionFile >>> filename = 'parkcity.msk' >>> srid = ProjectionFile.lookupSpatialReferenceID('readDirectory, 'parkcity_prj.pro') >>> maskMap = RasterMapFile() >>> maskMap.read(directory=readDirectory, filename=filename, session=spatialSession, spatial=True, spatialReferenceID=srid)
Not all files have spatial data, so passing in the spatial arguments to the read methods of these objects has no effect on the reading process. Refer to the API Documentation for a file object to determine if it supports spatial objects.
Spatial Project Write¶
There is no change needed to write a project that has been read in spatially. Use the same write methods as illustrated in the previous tutorial. This will not be demonstrated here.
After a project or file object has been read into the database with spatial objects, it is exposed to a number of spatial methods that can be used to generate visualizations of the data in various formats.
To demonstrate how these methods can be used to generate spatial visualizations, we will use the
method of the
gsshapy.orm.ProjectFile. This method uses the mask map and the stream network to generate a
summary visualization of the GSSHA model. Define the path where you want to write the kml file to. Then, query the
database for the project file that was read in as part of the spatial reading (id = 3 if you have been following the
tutorial) and invoke the
getModelSummaryAsKml() method on it:
>>> from gsshapy.orm import ProjectFile >>> import os >>> kml_path = os.path.join(writeDirectory, 'model_summary.kml') >>> newSpatialProjectFile = spatialSession.query(ProjectFile).filter(ProjectFile.id == 3).one() >>> newSpatialProjectFile.getModelSummaryAsKml(session=spatialSession, path=kml_path)
You will find the
model_summary.kml file in your write directory. If you have the
Google Earth Desktop application, you can view the
visualization. KML can also be loaded into the Google Maps and Google Earth web viewers to embed it in a website.
You can experiment with the other spatial methods to understand how they work. Refer to the API Documentation
for details in how to use each method.
Spatial Methods Available¶
File objects that include spatial methods include:
gsshapy.base.GeometricObjectBase offers several general purpose methods for objects that inherit from
gsshapy.base.RasterObjectBase offers several general purpose methods for objects that inherit from it: