Search icon CANCEL
Subscription
0
Cart icon
Your Cart (0 item)
Close icon
You have no products in your basket yet
Save more on your purchases! discount-offer-chevron-icon
Savings automatically calculated. No voucher code required.
Arrow left icon
Explore Products
Best Sellers
New Releases
Books
Videos
Audiobooks
Learning Hub
Free Learning
Arrow right icon
Geospatial Development By Example with Python
Geospatial Development By Example with Python

Geospatial Development By Example with Python: Build your first interactive map and build location-aware applications using cutting-edge examples in Python

eBook
$9.99 $43.99
Paperback
$54.99
Subscription
Free Trial
Renews at $19.99p/m

What do you get with Print?

Product feature icon Instant access to your digital eBook copy whilst your Print order is Shipped
Product feature icon Paperback book shipped to your preferred address
Product feature icon Download this book in EPUB and PDF formats
Product feature icon Access this title in our online reader with advanced features
Product feature icon DRM FREE - Read whenever, wherever and however you want
OR
Modal Close icon
Payment Processing...
tick Completed

Shipping Address

Billing Address

Shipping Methods
Table of content icon View table of contents Preview book icon Preview Book

Geospatial Development By Example with Python

Chapter 1. Preparing the Work Environment

Working with a programming language as a tool for geoprocessing provides the opportunity to construct a personalized application that can more optimally perform the task required by the user. This means that repetitive tasks can be automated, file inputs and outputs can be customized, and processes can be tuned to perform exactly what you want to be done.

Python is a powerful programming language that is gaining special attention as a tool for geoprocessing and scientific analysis. A number of factors may have contributed to its popularization, and three among them are worth mentioning: it's a scripting language, it's flexible and easy to learn, and it has a wide range of libraries available as open source software.

The number of available libraries and packages allow users to spend less time in programming basic functionalities and more in building processes and workflows to reach their goals.

In this first chapter, we will go through the process of installing all the libraries that you will need to go through the examples; it's likely that these same libraries will also satisfy most of your needs in real-world applications. Then, we will set up an Integrated Development Environment (IDE) that will help organize code and avoid mistakes. Finally, we will write a sample program with one of the libraries. Therefore, here are the topics that will be covered:

  • Installing Python and the packages that you need for the examples in this book
  • Learning the basics of the packages that you will use
  • Installing an IDE to write and organize your code
  • Creating a project for this book
  • Writing your first code

Installing Python

For this book, we suggest using Python 2.7; this version of Python is fully compatible with the libraries and packages that we will use in the examples and also has precompiled binary files available on the Internet for Windows users. We will keep all the examples as compatible as possible with Python 3.4 so that it would be easy to port future upgrades.

Windows users may find compatibility problems with the 64-bit packages, so we recommend the 32-bit version of Python for them.

For Linux users, we will show the installation procedures for Ubuntu Linux distribution and use package managers, so you don't need to worry about versions and requirements; the package managers will do this for you.

The libraries that you will install are written in Python and other languages, the most common being C and C++. These libraries can abstract classes, methods, and functions to Python objects or have an extra layer that makes the connection; when this happens, we say that the library has Python bindings.

Windows

Here are the steps to perform for the installation of Python on Windows:

  1. Go to https://www.python.org/downloads/windows/ and click on Download the latest Python 2.7 release for Windows.
  2. On the next page, roll down, and you will find a list of files; make sure that you download Windows x86 MSI installer.
  3. After the file is downloaded, open the installer by clicking on the file and following the instructions. Proceed with the default options by clicking on the Next button.

Ubuntu Linux

Ubuntu already comes with Python installed, so there is no need to install it. If for any reason, it's not available, you can install it with the following command:

sudo apt-get install python

Python 2.7.9 comes with Pip, but if you use an older version, you need to install Pip with the following command:

sudo apt-get install python-pip

Python packages and package manager

A Python package is a directory containing one or more Python files (that is, modules) plus one __init__.py file (this can be just an empty file). This file tells Python Interpreter that the directory is a package.

When writing Python code, we can import packages and modules and use them in our code. The Python community does this a lot; many packages use other packages and so on, forming an intricate network of requirements and dependencies.

In order to facilitate the installation of packages and all the requirements for it to run, Python has a package manager called pip.

Pip looks for packages in a central repository (or on a user-defined place), downloads it, then downloads its dependencies, and installs them. Some packages also use libraries in other languages, such as C. In these cases, these libraries need to be compiled during the installation. Ubuntu users don't have problem with this because many compilers are already installed on the system, but this won't work on Windows by default.

The repository of Python packages for Windows

Python makes it easy to install libraries and packages through pip. However, since Windows doesn't include any compiler by default, the installation of packages that needs the compilation of libraries fails. Instead of going through the process of installing a compiler, which is out of this book's scope, we can get the packages ready to use.

These packages come prebuilt for various types of system and don't need a compilation of its libraries. This type of package is called a wheel.

Note

Christoph Gohlke did a favor to all of us by building these packages and making them available for download at http://www.lfd.uci.edu/~gohlke/pythonlibs/.

Installing packages and required software

In this topic, we will go through the installation process of every package used in the book.

OpenCV

OpenCV is an optimized C/C++ library intended for video and image processing with hundreds of functions ranging from simple image resizing to object recognition, face detection, and so on. OpenCV is a big library, and we will use its capabilities of reading, transforming, and writing images. It's a good choice because its development is active, and it has a large user community and very good documentation.

Windows

Here is the installation procedure for Windows:

  1. Go to http://www.lfd.uci.edu/~gohlke/pythonlibs/.
  2. Press Ctrl + F to open the search dialog of your browser and then search for OpenCV.
  3. You will find a list of files; choose opencv_python‑2.4.11‑cp27‑none‑win32.whl or any OpenCV version that contains cp27 and win32. This means that this is the 32-bit version for Python 2.7.
  4. Save the downloaded file to a known location.
  5. Open Windows Command Prompt and run the following command:
    c:\Python27\scripts\pip install path_to_the_file_you_downloaded.whl
    
  6. You should see an output telling you that the installation was successful, as follows:
    Processing c:\downloads\opencv_python-2.4.12-cp27-none-win32.whl
    Installing collected packages: opencv-python
    Successfully installed opencv-python-2.4.12
    

    Tip

    You can drag and drop a file into the command prompt to enter its full path.

Ubuntu Linux

Here is the installation process for Ubuntu Linux:

  1. Open a new terminal with Ctrl + T.
  2. Then, enter the following command:
    sudo apt-get install python-opencv
    

Installing NumPy

NumPy is a package for scientific computing with Python. It handles multidimensional arrays of operations in a very efficient way. NumPy is required by OpenCV to run and will be used by many raster operations that we will perform in the examples. NumPy is also an efficient data container and will be our tool to calculate massive image data.

Windows

Repeat the same procedure as you did to install OpenCV; however, this time, search for NumPy and choose a file named numpy‑1.9.2+mkl‑cp27‑none‑win32.whl.

Ubuntu Linux

NumPy is automatically installed as a dependency of OpenCV on Ubuntu, but if you want to install it without OpenCV, follow these steps:

  1. Open a new terminal with Ctrl + T.
  2. Then, enter the following command:
    sudo pip install numpy
    

Installing GDAL and OGR

GDAL (Geospatial Data Abstraction Library) is composed of two packages that come together: OGR handles geospatial vector file formats, including coordinate system transformations and vector operations. GDAL is the raster part of the library, and in version 1.11, it comes packed with 139 drivers that can read, and some even create rasters. GDAL also comes packed with functions for raster transformations and calculations such as resizing, clipping, reprojecting, and so on.

In the following tables, there's an excerpt of the list of GDAL and OGR drivers with the most common formats that you may find:

Long format name

Code

Creation

Arc/Info ASCII Grid

AAIGrid

Yes

Arc/Info Export E00 GRID

E00GRID

No

ENVI .hdr Labelled Raster

ENVI

Yes

Generic Binary (.hdr Labelled)

GENBIN

No

Oracle Spatial GeoRaster

GEORASTER

Yes

GSat File Format

GFF

No

Graphics Interchange Format (.gif)

GIF

Yes

GMT Compatible netCDF

GMT

Yes

GRASS ASCII Grid

GRASSASCIIGrid

No

Golden Software ASCII Grid

GSAG

Yes

Golden Software Binary Grid

GSBG

Yes

Golden Software Surfer 7 Binary Grid

GS7BG

Yes

TIFF / BigTIFF / GeoTIFF (.tif)

GTiff

Yes

GXF (Grid eXchange File)

GXF

No

Erdas Imagine (.img)

HFA

Yes

JPEG JFIF (.jpg)

JPEG

Yes

NOAA Polar Orbiter Level 1b Data Set (AVHRR)

L1B

No

NOAA NGS Geoid Height Grids

NGSGEOID

No

NITF

NITF

Yes

NTv2 Datum Grid Shift

NTv2

Yes

PCI .aux Labelled

PAux

Yes

PCI Geomatics Database File

PCIDSK

Yes

PCRaster

PCRaster

Yes

Geospatial PDF

PDF

Yes

NASA Planetary Data System

PDS

No

Portable Network Graphics (.png)

PNG

Yes

R Object Data Store

R

Yes

Raster Matrix Format (*.rsw, .mtw)

RMF

Yes

RadarSat2 XML (product.xml)

RS2

No

Idrisi Raster

RST

Yes

SAGA GIS Binary format

SAGA

Yes

USGS SDTS DEM (*CATD.DDF)

SDTS

No

SGI Image Format

SGI

Yes

SRTM HGT Format

SRTMHGT

Yes

Terragen Heightfield (.ter)

TERRAGEN

Yes

USGS ASCII DEM / CDED (.dem)

USGSDEM

Yes

ASCII Gridded XYZ

XYZ

Yes

The following table describes the OGR drivers:

Format name

Code

Creation

Arc/Info Binary Coverage

AVCBin

No

Arc/Info .E00 (ASCII) Coverage

AVCE00

No

AutoCAD DXF

DXF

Yes

Comma Separated Value (.csv)

CSV

Yes

ESRI Shapefile

ESRI Shapefile

Yes

GeoJSON

GeoJSON

Yes

Géoconcept Export

Geoconcept

Yes

GeoRSS

GeoRSS

Yes

GML

GML

Yes

GMT

GMT

Yes

GPSBabel

GPSBabel

Yes

GPX

GPX

Yes

GPSTrackMaker (.gtm, .gtz)

GPSTrackMaker

Yes

Hydrographic Transfer Format

HTF

No

Idrisi Vector (.VCT)

Idrisi

No

KML

KML

Yes

Mapinfo File

MapInfo File

Yes

Microstation DGN

DGN

Yes

OpenAir

OpenAir

No

ESRI FileGDB

OpenFileGDB

No

PCI Geomatics Database File

PCIDSK

Yes

Geospatial PDF

PDF

Yes

PDS

PDS

No

PostgreSQL SQL dump

PGDump

Yes

U.S. Census TIGER/Line

TIGER

No

Note

You can find the full GDAL and OGR API documentation and the complete list of drivers at http://gdal.org/python/.

Windows

Again, we will use a wheel for the installation. Repeat the same procedure as before:

  1. Go to http://www.lfd.uci.edu/~gohlke/pythonlibs/.
  2. Now, search for GDAL and download the file named GDAL‑1.11.3‑cp27‑none‑win32.whl.
  3. Finally, install it with pip, as we did before.

Ubuntu Linux

Perform the following steps:

  1. Go to the terminal or open a new one.
  2. Then, enter the following command:
    sudo apt-get install python-gdal
    

Installing Mapnik

Mapnik is a map rendering package. It is a free toolkit to develop mapping applications. It produces high-quality maps and is used on many applications, including OpenStreetMaps.

Windows

Mapnik isn't available for installation as other libraries are. Instead, you need to go to http://mapnik.org/ and follow the download link:

  1. Download the Windows 32-bit package of Mapnik 2.2.
  2. Extract the mapnik-v2.2.0 to C:\ folder.
  3. Then, rename the extracted folder c:\mapnik.
  4. Now, add Mapnik to your PATH.
  5. Open Control Panel and go to System.
  6. Click on the Advanced System Settings link in the left-hand side column.
  7. In the System Properties window, click on the Advanced tab.
  8. Next, click on the Environment Variables button.
  9. In the System variables section, highlight the PATH variable and click on Edit. Add the following paths to the end of the list, each separated with a semicolon, as follows:
    c:\mapnik\bin;c:\mapnik\lib
    
  10. Now, click on the New button; then, set the variable name to PYTHONPATH and value to c:\mapnik\python\2.7\site-packages.

Ubuntu Linux

For this, perform the following:

  1. Go to the terminal or open a new one.
  2. Then, enter the following command:
    sudo apt-get install mapnik
    

Installing Shapely

Shapely is a package for the manipulation and analysis of two dimensional geometries. It can perform operations such as union and subtraction of geometries. It also can perform tests and comparisons, such as when a geometry intersects other geometries.

Windows

Here's what you need to do:

  1. As before, download the prebuilt wheel; this time, look for a file named Shapely‑1.5.13‑cp27‑none‑win32.whl.
  2. Then, install it with pip.

Ubuntu Linux

Here are the steps you need to perform:

  1. Go to the terminal or open a new one with Ctrl + T.
  2. Enter the following command:
    sudo pip install shapely
    

Installing other packages directly from pip

Some packages do not require compilation steps. For Windows users, these are easier to install because they can be obtained and installed directly with pip with a single command.

Windows

You need to simply type the following command in your Command Prompt:

c:\Python27\scripts\pip install django tabulate requests xmltodict psycopg2

Ubuntu Linux

In the terminal, type the following command:

sudo pip install django tabulate requests xmltodict psycopg2

For each package, you should see the progress of the installation, similar to the following:

Collecting django
  Downloading Django-1.9-py2.py3-none-any.whl (6.6MB)
    100% |################################| 6.6MB 43kB/s
Installing collected packages: django
Successfully installed django-1.9

Installing an IDE

IDEs are fancy text editors with tools and inspections regarding programming languages. You can surely use any text editor or IDE of your preference; none of the tasks in this book depends on the IDE, but an IDE will facilitate our work a lot because the suggested configuration will help you avoid mistakes and save time on typing, running, and debugging your code. The IDE checks the code for you and detects underlying errors; it even guesses what you are typing and completes the statements for you, runs the code with a simple command, and if there are exceptions, it provides links to the place where the exception occurred. For Windows or Linux, go to http://www.jetbrains.com/pycharm/ and click on the big orange button Get Pycharm Now. On the next page, select the free community edition.

Windows

Here are the steps you need to perform:

  1. After the download finishes, open the downloaded file; the Setup Wizard will pop up.
  2. Click on Next, and in the installation options, check both of the boxes: Create Desktop shortcut and Create associations.
  3. Click on Next and continue the installation.

Linux

Perform the following steps:

  1. Unpack the downloaded file in a directory.
  2. To open PyCharm, run pycharm.sh from the bin subdirectory. You can create a shortcut to it if you wish.

Tip

Downloading the example code

You can download the example code files for all Packt books you have purchased from your account at http://www.packtpub.com. If you purchased this book elsewhere, you can visit http://www.packtpub.com/support and register to have the files e-mailed directly to you.

Creating the book project

Perform the following steps:

  1. After installation, open Pycharm, and you will be prompted to create your first project:
    Creating the book project
  2. Click on create new project and then choose c:\geopy as your project location. In Linux, you can put the project inside your home folder—for example, /home/myname/geopy. Click on Create to create the project.
    Creating the book project
  3. In Windows, you will receive a security alert; this is Pycharm trying to access the Internet. It's recommended that you allow it so that you can later check for updates or download plugins:
    Creating the book project
  4. Finally, you should see the following window on your project workspace. Take some time to explore the menus and buttons, try right-clicking on your project folder to see the options:
    Creating the book project

Programming and running your first example

Now that we have all we need installed, we will go through our first example. In this example, we will test the installation and then see a glimpse of OGR's capabilities.

To do this, we will open a vector file containing the boundaries of all the countries in the world and make a list of country names. The objective of this simple example is to present the logic behind OGR objects and functions and give an understanding of how geospatial files are represented. Here's how:

  1. First, you need to copy the sample data provided with the book to your project folder. You can do this by dragging and dropping the data folder into the geopy folder. Make sure that the data folder is named data and that it's inside the geopy folder.
  2. Now, create a new directory for this chapter code, inside PyCharm. With your geopy project opened, right-click on the project folder and select New | Directory. Name it Chapter1.
  3. Create a new Python file. To do this, right-click on the Chapter1 folder and select New | Python File. Name it world_borders.py, and PyCharm will automatically open the file for editing.
  4. Type the following code in this file:
    import ogr
    # Open the shapefile and get the first layer.
    datasource = ogr.Open("../data/world_borders_simple.shp") 
    layer = datasource.GetLayerByIndex(0)
    print("Number of features: {}".format(layer.GetFeatureCount()))
  5. Now, run the code; in the menu bar, navigate to Run | Run, and in the dialog, choose world_borders. An output console will open at the bottom of the screen, and if everything goes fine, you should see this output:
    C:\Python27\python.exe C:/geopy/world_borders.py
    Number of features: 246
    
    Process finished with exit code 0
    

Congratulations! You successfully opened a Shapefile and counted the number of features inside it. Now, let's understand what this code does.

The first line imports the ogr package. From this point on, all the functions are available as ogr.FunctionName(). Note that ogr doesn't follow the Python naming conventions for functions.

The line after the comment opens the OGR datasource (this opens the shapefile containing the data) and assigns the object to the datasource variable. Note that the path, even on Windows, uses a forward slash (/) and not a backslash.

The next line gets the first layer of the data source by its index (0). Some data sources can have many layers, but this is not the case of a Shapefile, which has only one layer. So, when working with Shapefiles, we always know that the layer of interest is layer 0.

In the final line, the print statement prints the number of features returned by layer.GetFeatureCount(). Here, we will use Python's string formatting, where the curly braces are replaced by the argument passed to format().

Now, perform the following steps:

  1. In the same file, let's type the next part of our program:
    # Inspect the fields available in the layer.
    feature_definition = layer.GetLayerDefn()
    for field_index in range(feature_definition.GetFieldCount()):
        field_definition = feature_definition.GetFieldDefn(field_index)
        print("\t{}\t{}\t{}".format(field_index,
                                    field_definition.GetTypeName(),
                                    field_definition.GetName()))
  2. Rerun the code; you can use the Shift + F10 shortcut for this. Now, you should see the number of features as before plus a pretty table showing information on all the fields in the shapefile, as follows:
    Number of features: 246
     0 String  FIPS
     1 String  ISO2
     2 String  ISO3
     3 Integer UN
     4 String  NAME
     5 Integer POP2005
     6 Integer REGION
     7 Integer SUBREGION
    
    Process finished with exit code 0
    

    What happens in this piece of code is that feature_definition = layer.GetLayerDefn() gets the object that contains the definition of the features. This object contains the definition for each field and the type of geometry.

    In the for loop, we will get each field definition and print its index, name, and type. Note that the object returned by layer.GetLayerDefn() is not iterable, and we can't use for directly with it. So first, we will get the number of fields and use it in the range() function so that we can iterate through the indexes of the fields:

  3. Now, type the last part, as follows:
    # Print a list of country names.
    layer.ResetReading()
    for feature in layer:
        print(feature.GetFieldAsString(4))
  4. Run the code again and check the results in the output:
    Number of features: 246
     0 String FIPS
     1 String ISO2
     2 String ISO3
     3 Integer UN
     4 String NAME
     5 Integer POP2005
     6 Integer REGION
     7 Integer SUBREGION
    Antigua and Barbuda
    Algeria
    Azerbaijan
    Albania
    Armenia
    Angola
    ...
    Saint Barthelemy
    Guernsey
    Jersey
    South Georgia South Sandwich Islands
    Taiwan
    
    Process finished with exit code 0
    

The layer is iterable, but first, we need to ensure that we are at the beginning of the layer list with layer.ResetReading() (this is one of OGR's "gotcha" points).

The feature.GetFieldAsString(4) method returns the value of field 4 as a Python string. There are two ways of knowing whether the country names are in field 4:

  • Looking at the data's DBF file (by opening it with LibreOffice or Excel)
  • Looking at the table that we printed in the first part of the code

Your complete code should look similar to the following:

import ogr


# Open the shapefile and get the first layer.
datasource = ogr.Open("../data/world_borders_simple.shp")
layer = datasource.GetLayerByIndex(0)
print("Number of features: {}".format(layer.GetFeatureCount()))

# Inspect the fields available in the layer.
feature_definition = layer.GetLayerDefn()
for field_index in range(feature_definition.GetFieldCount()):
    field_definition = feature_definition.GetFieldDefn(field_index)
    print("\t{}\t{}\t{}".format(field_index,
                                field_definition.GetTypeName(),
                                field_definition.GetName()))

# Print a list of country names.
layer.ResetReading()
for feature in layer:
    print(feature.GetFieldAsString(4)) 

Transforming the coordinate system and calculating the area of all countries

Now, the objective is to know how much area is occupied by each country. However, the coordinates of country borders are expressed in latitude and longitude, and we can't calculate areas in this coordinate system. We want the area to be in the metric system, so first we need to convert the spatial reference system of the geometries.

Let's also take a step further in the programming techniques and start using functions to avoid the repetition of code. Perform the following steps:

  1. Create a new file in the Chapter1 directory, name this file world_areas.py, and program this first function:
    import ogr
    
    
    def open_shapefile(file_path):
        """Open the shapefile, get the first layer and returns
        the ogr datasource.
        """
        datasource = ogr.Open(file_path)
        layer = datasource.GetLayerByIndex(0)
        print("Opening {}".format(file_path))
        print("Number of features: {}".format(  
        layer.GetFeatureCount()))
        return datasource
  2. Run the code, go to Run | Run... in the menu, and select world_areas. If everything is correct, nothing should happen. This is because we are not calling our function. Add this line of code at the end and outside the function:
    datasource = open_shapefile("../data/world_borders_simple.shp") 
  3. Now, run the code again with Shift + F10 and check the output, as follows:
    Opening ../data/world_borders_simple.shp
    Number of features: 246
    
    Process finished with exit code 0
    

That's wonderful! You just created a piece of very useful and reusable code. You now have a function that can open any shapefile, print the number of features, and return the ogr datasource. From now on, you can reuse this function in any of your projects.

You are already familiar with how this code works, but there are a few novelties here that deserve an explanation. The def statement defines a function with the def function_name(arguments): syntax.

Remember when I told you that OGR doesn't follow Python's naming convention? Well, the convention is that function names should all be in lowercase with an underscore between words. A good hint for names is to follow the verb_noun rule.

Tip

These conventions are described in a document called PEP-8, where PEP stands for Python Enhancement Program. You can find this document at https://www.python.org/dev/peps/pep-0008/.

Right after the function's definition, you can see a description between triple quotes; this is a docstring, and it is used to document the code. It's optional but very useful for you to know what the function does.

Now, let's get back to our code. The second important thing to point out is the return statement. This makes the function return the values of the variables listed after the statement—in this case, the datasource.

Note

It's very important that all pieces of the OGR objects flow together through the program. In this case, if we return only the layer, for example, we will get a runtime error later in our program. This happens because in OGR internals, the layer has a reference to the data source, and when you exit a Python function, all objects that don't exit the function are trashed, and this breaks the reference.

Now, the next step is to create a function that performs the transformation. In OGR, the transformation is made in the feature's geometry, so we need to iterate over the features, get the geometry, and transform its coordinates. We will do this using the following steps:

  1. Add the following function to your world_areas.py file just after the open_shapefile function:
    def transform_geometries(datasource, src_epsg, dst_epsg):
        """Transform the coordinates of all geometries in the
     first layer.
     """
        # Part 1
        src_srs = osr.SpatialReference()
        src_srs.ImportFromEPSG(src_epsg)
        dst_srs = osr.SpatialReference()
        dst_srs.ImportFromEPSG(dst_epsg)
        transformation = osr.CoordinateTransformation(src_srs, dst_srs)
        layer = datasource.GetLayerByIndex(0)
    
        # Part 2
        geoms = []
        layer.ResetReading()
        for feature in layer:
            geom = feature.GetGeometryRef().Clone()
            geom.Transform(transformation)
            geoms.append(geom)
         return geoms

    The function takes three arguments: the ogr layer, the EPSG code of the coordinate system of the file, and the EPSG code for the transformation output.

    Here, it created an osr.CoordinateTransformation object; this object contains the instructions to perform the transformation.

    Probably by now, Pycharm should be complaining that osr is an unresolved reference; osr is the part of GDAL that deals with coordinate systems.

  2. Now, import the module by adding this line at the top of your code:
    import osr
    

    Here, the code iterates over all features, gets a reference to the geometry, and performs the transformation. As we don't want to change the original data, the geometry is cloned, and the transformation is made on the clone.

    Python lists are ordered; this means that the elements are in the same order in which they are appended to the list, and this order is always kept. This allows us to create a list of geometries in the same order of the features that are in the data source. This means that the geometries in the list and the features have the same index and can be related in the future by the index.

  3. Now, let's test the code; add the following lines at the end of the file (the first line is the one that you already added before):
    datasource = open_shapefile("../data/world_borders_simple.shp")
    layer = datasource.GetLayerByIndex(0)
    feature = layer.GetFeature(0)
    print("Before transformation:")
    print(feature.GetGeometryRef())
    transformed_geoms = transform_geometries(datasource, 4326, 3395)
    print("After transformation:")
    print(transformed_geoms[0])
  4. Finally, before you run the code, add one more import at the beginning of the program. It should be the first statement of your code, as follows:
    from __future__ import print_function

    This import allows us to use the print() function from Python 3 with the desired behavior, thus maintaining the compatibility.

  5. The complete code should look similar to this:
    from __future__ import print_function
    import ogr
    import osr
    
    
    def open_shapefile(file_path):
        ...
    
    def transform_geometries(datasource, src_epsg, dst_epsg):
        ...
        
    datasource = open_shapefile("../data/world_borders_simple.shp")
    layer = datasource.GetLayerByIndex(0)
    feature = layer.GetFeature(0)
    print("Before transformation:")
    print(feature.GetGeometryRef())
    transformed_geoms = transform_geometries(datasource, 4326, 3395)
    print("After transformation:")
    print(transformed_geoms[0])
  6. Run your program again by pressing Shift + F10. In the output, note the difference in the coordinates before and after the transformation:
    Opening ../data/world_borders_simple.shp
    Number of features: 246
    Before transformation:
    MULTIPOLYGON (((-61.686668 17.024441000000138 ... )))
    After transformation:
    MULTIPOLYGON (((-6866928.4704937246 ... )))
    
    Process finished with exit code 0
    
  7. Now, add another function. This function will calculate the area in square meters (because we will use the geometries that have coordinates in meters), convert the value (or not) to square kilometers or square miles, and store the values in another list with the same order as before. Execute the following code:
    def calculate_areas(geometries, unity='km2'):
        """Calculate the area for a list of ogr geometries."""
        # Part 1
        conversion_factor = {
            'sqmi': 2589988.11,
            'km2': 1000000,
            'm': 1}
        # Part2
        if unity not in conversion_factor:
            raise ValueError(
                "This unity is not defined: {}".format(unity))
        # Part 3
        areas = []
        for geom in geometries:
            area = geom.Area()
            areas.append(area / conversion_factor[unity])
        return areas

    Firstly, note that in the function definition, we use unity='km2'; this is a keyword argument, and when you call the functions, this argument is optional.

    In Part 1, a dictionary is used to define a few conversion factors for the area unit. Feel free to add more units if you wish. By the way, Python doesn't care if you use single or double quotes.

    In Part 2, a verification is made to check whether the passed unity exists and whether it is defined in conversion_factor. Another way of doing this is catching the exception later; however, for now, let's opt for readability.

    In Part 3, the code iterates the ogr geometries, calculates the area, converts the values, and puts it on a list.

  8. Now, to test the code, edit your first line, including division to the future imports. This will ensure that all divisions return floating point numbers and not integers. Here's how it should look:
    from __future__ import print_function, division
  9. Then, update the testing part of your code to the following:
    datasource = open_shapefile("../data/world_borders_simple.shp")
    transformed_geoms = transform_geometries(datasource, 4326, 3395)
    calculated_areas = calculate_areas(transformed_geoms, unity='sqmi')
    print(calculated_areas)
  10. Run it, change the unity, then run again, and note how the results change.

Very well, unity conversion is another very important procedure in geoprocessing, and you just implemented it in your calculate_areas function.

However, having a list of numbers as the output is not very useful to us. So, it's time to combine everything that we did so far in order to extract valuable information.

Sort the countries by area size

You programmed three functions so far; now, let's add another one to our list by converting the code that generated a list of country names to a function and add this function to world_areas.py, as follows:

def get_country_names(datasource):
    """Returns a list of country names."""
    layer = datasource.GetLayerByIndex(0)
    country_names = []
    layer.ResetReading()
    for feature in layer:
        country_names.append(feature.GetFieldAsString(4))
    return country_names

Now, we have four functions, which are:

  • open_shapefile
  • transform_geometries
  • calculate_areas
  • get_country_names

All these functions return iterables, with each item sharing the same index on all of them, thus making it easy to combine the information.

So, let's take advantage of this feature to sort the countries by area size and return a list of the five biggest countries and their areas. For this, add another function, as follows:

def get_biggest_countries(countries, areas, elements=5):
    """Returns a list of n countries sorted by area size."""    
    countries_list = [list(country) 
                      for country in zip(areas, countries)]

    sorted_countries = sorted(countries_list, 
                              key=itemgetter(0), reverse=True) 
    return sorted_countries[:5]

In the first line, the two lists are zipped together, producing a list of country-area pairs. Then, we used the Python list's sorted method, but as we don't want the lists to be sorted by both values, we will define the key for sorting. Finally, the list is sliced, returning only the desired number of values.

  1. In order to run this code, you need to import the itemgetter function and put it at the beginning of the code but after from __future__ imports, as follows:
    from operator import itemgetter
  2. Now, edit the testing part of your code to look similar to the following:
    datasource = open_shapefile("../data/world_borders_simple.shp")
    transformed_geoms = transform_geometries(datasource, 4326, 3395)
    country_names = get_country_names(datasource)
    country_areas = calculate_areas(transformed_geoms)
    biggest_countries = get_biggest_countries(country_names, 
                                              country_areas)
    for item in biggest_countries:
        print("{}\t{}".format(item[0], item[1])) 
  3. Now, run the code and take a look at the results, as follows:
    Opening ../data/world_borders_simple.shp
    Number of features: 246
    82820725.1423 Russia
    51163710.3726 Canada
    35224817.514 Greenland
    21674429.8403 United States
    14851905.8596 China
    
    Process finished with exit code 0

Summary

In this chapter, we had a brief introduction to the libraries and packages that we will use in this book. By installing these libraries, you also learned the general procedure of how to search and install Python packages. You can use this procedure in other cases whenever you feel the need for other libraries in your applications.

Then, we wrote code that made use of the OGR library to open a shapefile and perform area calculation and sorting. These simple procedures showed a little bit of the internal organization of OGR, how it handles geographic data, and how it is possible to extract information from them. In the next chapter, we will use some of the techniques learned here to read data and process vector points.

Left arrow icon Right arrow icon

Key benefits

  • • Learn the full geo-processing workflow using Python with open source packages
  • • Create press-quality styled maps and data visualization with high-level and reusable code
  • • Process massive datasets efficiently using parallel processing

Description

From Python programming good practices to the advanced use of analysis packages, this book teaches you how to write applications that will perform complex geoprocessing tasks that can be replicated and reused. Much more than simple scripts, you will write functions to import data, create Python classes that represent your features, and learn how to combine and filter them. With pluggable mechanisms, you will learn how to visualize data and the results of analysis in beautiful maps that can be batch-generated and embedded into documents or web pages. Finally, you will learn how to consume and process an enormous amount of data very efficiently by using advanced tools and modern computers’ parallel processing capabilities.

Who is this book for?

Geospatial Development By Example with Python is intended for beginners or advanced developers in Python who want to work with geographic data. The book is suitable for professional developers who are new to geospatial development, for hobbyists, or for data scientists who want to move into some simple development.

What you will learn

  • • Prepare a development environment with all the tools needed for geo-processing with Python
  • • Import point data and structure an application using Python's resources
  • • Combine point data from multiple sources, creating intuitive and functional representations of geographic objects
  • • Filter data by coordinates or attributes easily using pure Python
  • • Make press-quality and replicable maps from any data
  • • Download, transform, and use remote sensing data in your maps
  • • Make calculations to extract information from raster data and show the results on beautiful maps
  • • Handle massive amounts of data with advanced processing techniques
  • • Process huge satellite images in an efficient way
  • • Optimize geo-processing times with parallel processing
Estimated delivery fee Deliver to Ecuador

Standard delivery 10 - 13 business days

$19.95

Premium delivery 3 - 6 business days

$40.95
(Includes tracking information)

Product Details

Country selected
Publication date, Length, Edition, Language, ISBN-13
Publication date : Jan 30, 2016
Length: 340 pages
Edition : 1st
Language : English
ISBN-13 : 9781785282355
Category :
Languages :
Tools :

What do you get with Print?

Product feature icon Instant access to your digital eBook copy whilst your Print order is Shipped
Product feature icon Paperback book shipped to your preferred address
Product feature icon Download this book in EPUB and PDF formats
Product feature icon Access this title in our online reader with advanced features
Product feature icon DRM FREE - Read whenever, wherever and however you want
OR
Modal Close icon
Payment Processing...
tick Completed

Shipping Address

Billing Address

Shipping Methods
Estimated delivery fee Deliver to Ecuador

Standard delivery 10 - 13 business days

$19.95

Premium delivery 3 - 6 business days

$40.95
(Includes tracking information)

Product Details

Publication date : Jan 30, 2016
Length: 340 pages
Edition : 1st
Language : English
ISBN-13 : 9781785282355
Category :
Languages :
Tools :

Packt Subscriptions

See our plans and pricing
Modal Close icon
$19.99 billed monthly
Feature tick icon Unlimited access to Packt's library of 7,000+ practical books and videos
Feature tick icon Constantly refreshed with 50+ new titles a month
Feature tick icon Exclusive Early access to books as they're written
Feature tick icon Solve problems while you work with advanced search and reference features
Feature tick icon Offline reading on the mobile app
Feature tick icon Simple pricing, no contract
$199.99 billed annually
Feature tick icon Unlimited access to Packt's library of 7,000+ practical books and videos
Feature tick icon Constantly refreshed with 50+ new titles a month
Feature tick icon Exclusive Early access to books as they're written
Feature tick icon Solve problems while you work with advanced search and reference features
Feature tick icon Offline reading on the mobile app
Feature tick icon Choose a DRM-free eBook or Video every month to keep
Feature tick icon PLUS own as many other DRM-free eBooks or Videos as you like for just $5 each
Feature tick icon Exclusive print discounts
$279.99 billed in 18 months
Feature tick icon Unlimited access to Packt's library of 7,000+ practical books and videos
Feature tick icon Constantly refreshed with 50+ new titles a month
Feature tick icon Exclusive Early access to books as they're written
Feature tick icon Solve problems while you work with advanced search and reference features
Feature tick icon Offline reading on the mobile app
Feature tick icon Choose a DRM-free eBook or Video every month to keep
Feature tick icon PLUS own as many other DRM-free eBooks or Videos as you like for just $5 each
Feature tick icon Exclusive print discounts

Frequently bought together


Stars icon
Total $ 164.97
Python Geospatial Development
$54.99
Python Geospatial Analysis Cookbook
$54.99
Geospatial Development By Example with Python
$54.99
Total $ 164.97 Stars icon
Banner background image

Table of Contents

11 Chapters
1. Preparing the Work Environment Chevron down icon Chevron up icon
2. The Geocaching App Chevron down icon Chevron up icon
3. Combining Multiple Data Sources Chevron down icon Chevron up icon
4. Improving the App Search Capabilities Chevron down icon Chevron up icon
5. Making Maps Chevron down icon Chevron up icon
6. Working with Remote Sensing Images Chevron down icon Chevron up icon
7. Extract Information from Raster Data Chevron down icon Chevron up icon
8. Data Miner App Chevron down icon Chevron up icon
9. Processing Big Images Chevron down icon Chevron up icon
10. Parallel Processing Chevron down icon Chevron up icon
Index Chevron down icon Chevron up icon

Customer reviews

Rating distribution
Full star icon Full star icon Full star icon Full star icon Full star icon 5
(4 Ratings)
5 star 100%
4 star 0%
3 star 0%
2 star 0%
1 star 0%
Winston Apr 04, 2016
Full star icon Full star icon Full star icon Full star icon Full star icon 5
Great book!!! This book takes the reader down a well structured path using Python to develop geospatial app. First the author has the reader download and install various Python packages like NumPy and Mapnik combined with multiple datasources to create fully-functional apps. Must read.
Amazon Verified review Amazon
ruben Apr 06, 2016
Full star icon Full star icon Full star icon Full star icon Full star icon 5
This is an excellent book for programming with Python very interesting book for users who want to get involve in the languange programming.Much more than simple scripts, you will write functions to import data, create Python classes that represent your features, and learn how to combine and filter them.
Amazon Verified review Amazon
Bill Jones Apr 08, 2016
Full star icon Full star icon Full star icon Full star icon Full star icon 5
Awesome book, well written, and structured nicely. I'm new to Python but in that regard I didn't care if it was 3 or 2.7 until I started to want to use other modules or play around with some libraries. So one of the reviews did mention that but I didn't have any issues with the content in this book. I would highly recommend this book if you want to do some Geospatial Development with Python.
Amazon Verified review Amazon
Tim Crothers Apr 07, 2016
Full star icon Full star icon Full star icon Full star icon Full star icon 5
I really enjoyed this book. Prior I had a lot of experience with Python but only limited experience with various geospatial packages. By the end of the book I had a solid foundation for how to use several libraries I was unfamiliar with and how to apply them to various uses related to geospatial applications. Be aware the book is written in Python 2.7 so if you are a Python 3 developer you'll need to translate and many of the libraries aren't available yet in Python 3 so the value might be a bit less
Amazon Verified review Amazon
Get free access to Packt library with over 7500+ books and video courses for 7 days!
Start Free Trial

FAQs

What is the delivery time and cost of print book? Chevron down icon Chevron up icon

Shipping Details

USA:

'

Economy: Delivery to most addresses in the US within 10-15 business days

Premium: Trackable Delivery to most addresses in the US within 3-8 business days

UK:

Economy: Delivery to most addresses in the U.K. within 7-9 business days.
Shipments are not trackable

Premium: Trackable delivery to most addresses in the U.K. within 3-4 business days!
Add one extra business day for deliveries to Northern Ireland and Scottish Highlands and islands

EU:

Premium: Trackable delivery to most EU destinations within 4-9 business days.

Australia:

Economy: Can deliver to P. O. Boxes and private residences.
Trackable service with delivery to addresses in Australia only.
Delivery time ranges from 7-9 business days for VIC and 8-10 business days for Interstate metro
Delivery time is up to 15 business days for remote areas of WA, NT & QLD.

Premium: Delivery to addresses in Australia only
Trackable delivery to most P. O. Boxes and private residences in Australia within 4-5 days based on the distance to a destination following dispatch.

India:

Premium: Delivery to most Indian addresses within 5-6 business days

Rest of the World:

Premium: Countries in the American continent: Trackable delivery to most countries within 4-7 business days

Asia:

Premium: Delivery to most Asian addresses within 5-9 business days

Disclaimer:
All orders received before 5 PM U.K time would start printing from the next business day. So the estimated delivery times start from the next day as well. Orders received after 5 PM U.K time (in our internal systems) on a business day or anytime on the weekend will begin printing the second to next business day. For example, an order placed at 11 AM today will begin printing tomorrow, whereas an order placed at 9 PM tonight will begin printing the day after tomorrow.


Unfortunately, due to several restrictions, we are unable to ship to the following countries:

  1. Afghanistan
  2. American Samoa
  3. Belarus
  4. Brunei Darussalam
  5. Central African Republic
  6. The Democratic Republic of Congo
  7. Eritrea
  8. Guinea-bissau
  9. Iran
  10. Lebanon
  11. Libiya Arab Jamahriya
  12. Somalia
  13. Sudan
  14. Russian Federation
  15. Syrian Arab Republic
  16. Ukraine
  17. Venezuela
What is custom duty/charge? Chevron down icon Chevron up icon

Customs duty are charges levied on goods when they cross international borders. It is a tax that is imposed on imported goods. These duties are charged by special authorities and bodies created by local governments and are meant to protect local industries, economies, and businesses.

Do I have to pay customs charges for the print book order? Chevron down icon Chevron up icon

The orders shipped to the countries that are listed under EU27 will not bear custom charges. They are paid by Packt as part of the order.

List of EU27 countries: www.gov.uk/eu-eea:

A custom duty or localized taxes may be applicable on the shipment and would be charged by the recipient country outside of the EU27 which should be paid by the customer and these duties are not included in the shipping charges been charged on the order.

How do I know my custom duty charges? Chevron down icon Chevron up icon

The amount of duty payable varies greatly depending on the imported goods, the country of origin and several other factors like the total invoice amount or dimensions like weight, and other such criteria applicable in your country.

For example:

  • If you live in Mexico, and the declared value of your ordered items is over $ 50, for you to receive a package, you will have to pay additional import tax of 19% which will be $ 9.50 to the courier service.
  • Whereas if you live in Turkey, and the declared value of your ordered items is over € 22, for you to receive a package, you will have to pay additional import tax of 18% which will be € 3.96 to the courier service.
How can I cancel my order? Chevron down icon Chevron up icon

Cancellation Policy for Published Printed Books:

You can cancel any order within 1 hour of placing the order. Simply contact [email protected] with your order details or payment transaction id. If your order has already started the shipment process, we will do our best to stop it. However, if it is already on the way to you then when you receive it, you can contact us at [email protected] using the returns and refund process.

Please understand that Packt Publishing cannot provide refunds or cancel any order except for the cases described in our Return Policy (i.e. Packt Publishing agrees to replace your printed book because it arrives damaged or material defect in book), Packt Publishing will not accept returns.

What is your returns and refunds policy? Chevron down icon Chevron up icon

Return Policy:

We want you to be happy with your purchase from Packtpub.com. We will not hassle you with returning print books to us. If the print book you receive from us is incorrect, damaged, doesn't work or is unacceptably late, please contact Customer Relations Team on [email protected] with the order number and issue details as explained below:

  1. If you ordered (eBook, Video or Print Book) incorrectly or accidentally, please contact Customer Relations Team on [email protected] within one hour of placing the order and we will replace/refund you the item cost.
  2. Sadly, if your eBook or Video file is faulty or a fault occurs during the eBook or Video being made available to you, i.e. during download then you should contact Customer Relations Team within 14 days of purchase on [email protected] who will be able to resolve this issue for you.
  3. You will have a choice of replacement or refund of the problem items.(damaged, defective or incorrect)
  4. Once Customer Care Team confirms that you will be refunded, you should receive the refund within 10 to 12 working days.
  5. If you are only requesting a refund of one book from a multiple order, then we will refund you the appropriate single item.
  6. Where the items were shipped under a free shipping offer, there will be no shipping costs to refund.

On the off chance your printed book arrives damaged, with book material defect, contact our Customer Relation Team on [email protected] within 14 days of receipt of the book with appropriate evidence of damage and we will work with you to secure a replacement copy, if necessary. Please note that each printed book you order from us is individually made by Packt's professional book-printing partner which is on a print-on-demand basis.

What tax is charged? Chevron down icon Chevron up icon

Currently, no tax is charged on the purchase of any print book (subject to change based on the laws and regulations). A localized VAT fee is charged only to our European and UK customers on eBooks, Video and subscriptions that they buy. GST is charged to Indian customers for eBooks and video purchases.

What payment methods can I use? Chevron down icon Chevron up icon

You can pay with the following card types:

  1. Visa Debit
  2. Visa Credit
  3. MasterCard
  4. PayPal
What is the delivery time and cost of print books? Chevron down icon Chevron up icon

Shipping Details

USA:

'

Economy: Delivery to most addresses in the US within 10-15 business days

Premium: Trackable Delivery to most addresses in the US within 3-8 business days

UK:

Economy: Delivery to most addresses in the U.K. within 7-9 business days.
Shipments are not trackable

Premium: Trackable delivery to most addresses in the U.K. within 3-4 business days!
Add one extra business day for deliveries to Northern Ireland and Scottish Highlands and islands

EU:

Premium: Trackable delivery to most EU destinations within 4-9 business days.

Australia:

Economy: Can deliver to P. O. Boxes and private residences.
Trackable service with delivery to addresses in Australia only.
Delivery time ranges from 7-9 business days for VIC and 8-10 business days for Interstate metro
Delivery time is up to 15 business days for remote areas of WA, NT & QLD.

Premium: Delivery to addresses in Australia only
Trackable delivery to most P. O. Boxes and private residences in Australia within 4-5 days based on the distance to a destination following dispatch.

India:

Premium: Delivery to most Indian addresses within 5-6 business days

Rest of the World:

Premium: Countries in the American continent: Trackable delivery to most countries within 4-7 business days

Asia:

Premium: Delivery to most Asian addresses within 5-9 business days

Disclaimer:
All orders received before 5 PM U.K time would start printing from the next business day. So the estimated delivery times start from the next day as well. Orders received after 5 PM U.K time (in our internal systems) on a business day or anytime on the weekend will begin printing the second to next business day. For example, an order placed at 11 AM today will begin printing tomorrow, whereas an order placed at 9 PM tonight will begin printing the day after tomorrow.


Unfortunately, due to several restrictions, we are unable to ship to the following countries:

  1. Afghanistan
  2. American Samoa
  3. Belarus
  4. Brunei Darussalam
  5. Central African Republic
  6. The Democratic Republic of Congo
  7. Eritrea
  8. Guinea-bissau
  9. Iran
  10. Lebanon
  11. Libiya Arab Jamahriya
  12. Somalia
  13. Sudan
  14. Russian Federation
  15. Syrian Arab Republic
  16. Ukraine
  17. Venezuela