Tuesday, July 10, 2018

QGIS Python Tutorial: Animate GPS Track

In this QGIS python tutorial series, I will explain about how to animate a GPS track on QGIS map. Animation is a really cool way to visualize spatial data, because our eyes and brain are more connected with a series of moving object in the screen. It made us to watch it until finish or at a point when our brain says enough.

For this tutorial, I will use a GPX file which is called ride.gpx. I will parse  the file using xml dom to extract the coordinate of track points. Then animate the points with markers with a defined time interval. The result of track animation can be seen in figure 1.

Tracking animation Python QGIS
Figure 1. Tracking animation
Let's start this tutorial with importing some python libraries such as xml dom and time.

import time
from xml.dom import minidom

Using the minidom library, then we read the gpx file and get the track element. It will return a track list for each point.

#READ GPX FILE
data=open('F:/ride.gpx')
xmldoc = minidom.parse(data)
track = xmldoc.getElementsByTagName('trkpt')
n_track=len(track)

Next, the list is parsed to get each point coordinate (longitude and latitude). Based on the coordinates, then a marker is created for each point and stored in a list(marker_list). The marker symbol also customized in this step by setting up the color, icon type, fill color, size, and line width. Finally the markers are displayed on QGIS map canvas with interval 1 second.  

#PARSING GPX ELEMENT AND CREATING MARKER
marker_list=[]
canvas=iface.mapCanvas()
for s in range(n_track):
    lon,lat=track[s].attributes['lon'].value,track[s].attributes['lat'].value
    x=float(lon)
    y=float(lat)
    m = QgsVertexMarker(canvas)
    m.setCenter(QgsPointXY(x,y))
    m.setColor(QColor(255,0,0))
    m.setFillColor(QColor(255,255,0))
    m.setIconSize(10)
    m.setIconType(QgsVertexMarker.ICON_CIRCLE)
    m.setPenWidth(3)
    marker_list.append(m)
    time.sleep(1)
    print(m)

Next, a function which is called hide_track is created to hide the GPS track. And also a play_track function is created to play the track animation. These functions will take each marker in marker list and hide or show it with 1 second time interval. Both functions can be used only after the code execution. To hide the track, just type in the python console hide_track() and play_track() to start the GPS track animation again.

#FUNCTION TO HIDE AND PLAY TRACK ANIMATION    
def hide_track():
    for i in range(n_track):
        marker_list[i].hide()

def play_track():
    hide_track()
    for j in range(n_track):
        marker_list[j].show()
        time.sleep(1) #Time interval. You can change it here
        print(marker_list[j])

The complete code  can be seen below. In figure 2 shows the code action in QGIS.

import time
from xml.dom import minidom

#READ GPX FILE
data=open('F:/ride.gpx')
xmldoc = minidom.parse(data)
track = xmldoc.getElementsByTagName('trkpt')
n_track=len(track)

#PARSING GPX ELEMENT AND CREATING MARKER
marker_list=[]
canvas=iface.mapCanvas()
for s in range(n_track):
    lon,lat=track[s].attributes['lon'].value,track[s].attributes['lat'].value
    x=float(lon)
    y=float(lat)
    m = QgsVertexMarker(canvas)
    m.setCenter(QgsPointXY(x,y))
    m.setColor(QColor(255,0,0))
    m.setFillColor(QColor(255,255,0))
    m.setIconSize(10)
    m.setIconType(QgsVertexMarker.ICON_CIRCLE)
    m.setPenWidth(3)
    marker_list.append(m)
    time.sleep(1)
    print(m)
   
#FUCNTION TO HIDE AND PLAY TRACK ANIMATION    
def hide_track():
    for i in range(n_track):
        marker_list[i].hide()

def play_track():
    hide_track()
    for j in range(n_track):
        marker_list[j].show()
        time.sleep(1)
        print(marker_list[j])

GPS track animation in QGIS
Figure 2. GPS track animation in QGIS
That's all the tutorial how to anime GPS track in QGIS using python. If you want to learn more about QGIS Python Programming, please see the QGIS Python Tutorial series which covers some other topics. 

Saturday, July 7, 2018

QGIS Python Tutorial: Drawing Polyline and Polygon

In another post has been discussed how to draw marker in QGIS map canvas. Maybe you have a question what about polyline and polygon? How to draw those in QGIS using Python? This tutorial will give the answer.

Drawing Polyline

To draw a polyline we are using QgsRubberBand class. But Firstly we need to initiate the canvas first using QGIS interface.

canvas=iface.mapCanvas()

Using QgsPoints class then a list of points coordinate for a polyline are defined. Then we set the geometry with QgsGeometry class and set the defined points as the argument. Finally customize  how the polyline looks like. Set it's color and the width of the polyline. The complete code as follow:

#Drawing Polyline
polyline = QgsRubberBand(canvas, False)  # False = not a polygon
points =[QgsPoint(-124,48), QgsPoint(-124,51 ), QgsPoint(-120,51),QgsPoint(-120,48)]
polyline.setToGeometry(QgsGeometry.fromPolyline(points), None)
polyline.setColor(QColor(255, 0, 0))
polyline.setWidth(3)

Drawing Polygon

Similar with polyline, we also use QgsRubberBand class to draw polygon. But in defining the point list we are using QgsPointXY class instead of QgsPoint. Then we customize the outline color of the polygon, it's widht and fill color. The complete code to draw polygon can be seen below.

#Drawing Polygon
polygon = QgsRubberBand(canvas)
p_points=[QgsPointXY(-123, 49), QgsPointXY(-123, 50),QgsPointXY(-121, 50),QgsPointXY(-121,49)]
polygon.setToGeometry(QgsGeometry.fromPolygonXY([p_points]),None)
polygon.setColor(QColor(0, 0, 255))
polygon.setFillColor(QColor(255,255,0))
polygon.setWidth(3)

Figure 1 shows the polyline and polygon from the code above. We can see a red polyline and a yellow polygon with blue outline in QGIS map canvas. To show and hide the drawing object use hide and show method. Type it in the console.

polyline.hide() # to hide polygon
Polyline.show() # to show polygon

Polyline and Polygon drawing QGIS
Figure 1. Polyline and Polygon drawing
That's all the tutorial how to add polyline and polygon with python in QGIS. Please see other tutorials which discuss other topics in the QGIS python programming tutorial series. To learn further about this topic, take a look the PygGIS Developer Cookbook.

Friday, July 6, 2018

Python QGIS Tutorial: Getting Started

Welcome to the first tutorial about QGIS Python Programming. This tutorial will discuss about the python console in QGIS. You can skip this tutorial if you already familiar to work with python console in QGIS.

Getting Started with QGIS Python Programming

To get started using python code in QGIS, we are using Python Console. Click the python console in the plugins menu as in figure 1.
Open Python Console
Figure 1. Open Python Console


The python console will be docked at the bottom of QGIS window as in figure 2.

Qgis Window with python console
Figure 2. Qgis Window with python console

In the console you can type the python syntax. for example to check the python version we can use the following code.

import sys
print (sys.version)

Type those code in the console and hit enter you will get the python version and some info of the machine like figure 3.

Python console
Figure 3. Python console
In the console we can just type a line of syntax or code and then execute it. Although we can make a function or a loop process in it, but it would be difficult to handle. That's why we use the console just for short purpose, such as to get a variable, object, list, etc. For a complex or long coding we are using the editor. To access the editor just click the paper icon with a pencil as in figure 4.

Python editor
Figure 4. Python editor
In the figure 4 can be seen a blank editor. In this editor we can write our coding and can be saved for the next purpose.

That's all the tutorial how to get started QGIS with python programming. In the other QGIS Python Programming post we will make a journey to explore more about this topic. Hopefully you can have fun with it.

Python QGIS Tutorial: Working with Vector Layer

In this tutorial we will explore more about how to work with vector layer like getting fields name, layer extent, getting feature, and many more. To work with vector layer with python we are using QgsVectorLayer class. In this tutorial we will use some methods that available in this class.

Adding Vector Layer 

Let's start this tutorial with adding a vector layer. (See the tutorial how to add vector layer with python to get more explanation about adding vector layer in QGIS with python)

park_layer=iface.addVectorLayer("F:\park_polygons.shp","park","ogr")

Checking Coordinate System

When getting a vector layer, usually we want to know the coordinate system of the layer. We can do this with crs method.

#Check Layer's Coordinate System
crs=park_layer.crs()
print(crs.description())

Output: NAD83 / UTM zone 10N

Get The Layer Extent

Next let's check the the extent of the layer. For this we use the extent method. It will return a QgsRectangle class. Then to get the extent coordinate we use minimum and maximum coordinate method.

#Check Layer's Extent
extent=park_layer.extent()
min_x=extent.xMinimum()
max_x=extent.xMaximum()
min_y=extent.yMinimum()
max_y=extent.yMaximum()
print (min_x,min_y,max_x,max_y)

Output:
483719.3637926 5449712.914727 498270.49187981 5462386.1670929

Counting Feature Number

To count number featureCount method is used. It will return the number of feature in the layer.

#Counting Feature
n_feature=park_layer.featureCount()
print("number feature: ", n_feature)

Output: number feature:  264

Getting Attributes Column

Next, we want to know what information that can be found in the layer by checking its attribute's column name. For this we use fields method in a loop and print the column name and field type.

#Get Fields Name and Type
for field in park_layer.fields():
    print (field.name(),field.typeName())

Output:
PARK_NAME String
PARK_ID Integer64
PARK_URL String
AREA_HEC String

Layer Metadata

Metadata is information about the data. We can see the metadata using htmlMetada method. The result of metada in html format will be displayed in the console.

#Get HTML Metadata
metadata=park_layer.htmlMetadata()
print (metadata) 

Output: (Something like this)
<html>
<body>
<h1>Information from provider</h1>
<hr>
<table class="list-view">
<tr><td class="highlight">Original</td><td>park</td></tr>
<tr><td class="highlight">Name</td><td>park</td></tr>
<tr><td class="highlight">Source</td><td>F:\park_polygons.shp</td></tr>
<tr><td class="highlight">Storage</td><td>ESRI Shapefile</td></tr>
<tr><td class="highlight">Comment</td><td></td></tr>
<tr><td class="highlight">Encoding</td><td>UTF-8</td></tr>
<tr><td class="highlight">Geometry</td><td>Polygon (MultiPolygon)</td></tr>
<tr><td class="highlight">CRS</td><td>EPSG:26910 - NAD83 / UTM zone 10N - Projected</td></tr>
<tr><td class="highlight">Extent</td><td>483719.3637926569208503,5449712.9147277921438217 : 498270.4918798131402582,5462386.1670929444953799</td></tr>
<tr><td class="highlight">Unit</td><td>meters</td></tr>
<tr><td class="highlight">Feature count</td><td>264</td></tr>
</table>
<br><br><h1>Identification</h1>
<hr>
................
................
</body>
</html>

Getting Attribute Value

To get attribute values of the layer firstly we have have to get the feature. Then in a loop we use attributes method to get each feature's attribute. The output will return as a list for each feature.

#Getting Attribute Value
features=park_layer.getFeatures()
for f in features:
    attr=f.attributes()
    print (attr)

Sample Outoput:
['Yaletown Park', 237, 'http://covapp.vancouver.ca/parkfinder/parkdetail.aspx?inparkid=237', '0.17']

There are many methods still available for QgsVectorLayer  class. If you want to explore more, take a look at QGIS Python API website. Visit other topics in QGIS Python Programming

Python QGIS Tutorial: Adding Vector Layer

Commonly adding or displaying geospatial data is the first step when starting working with GIS. Mainly there are two types of geospatial data format, vector and raster. Vector is stored in a coordinate or a list of coordinates in point, polyline or polygon geometry. Meanwhile the raster format storing the data in picture element (pixel).

This tutorial will discuss how to display vector data with python in QGIS. It will cover two methods to display vector data using QGIS interface and QGIS project instance. This tutorial is using QGIS 3.0.

Adding Vector Layer with QGIS Interface API

We are starting this tutorial with adding a vector layer in shp format using QGIS interface API. To add a vector layer with QGIS Interface API can be done with addVectorLayer method as in the following expression:
layer=iface.addVectorLayer(path:string,layer_name:string,library:string)

For example to add a street layer with the name public_streets.shp with ogr simple feature library will be:
layer=iface.addVectorLayer("F:\public_streets.shp","street","ogr")

Figure 1 shows the result when the expression was executed. The street data will be added into QGIS layer, and the data displayed in QGIS map canvas.

Adding a vector layer with qgis interface
figure 1. Adding a vector layer with qgis interface






  



Adding Vector Layer with Project Instance

The second method we will use project instance to add vector data. For this we use the following syntax expression:

QgsProject.instance().addMapLayer(layer:QgsVectorLayer) #to add a vector layer
QgsProject.instance().addMapLayers(layer list:QgsVectorLayer]) #to add a list of vector layer

From the syntax above can be seen that we can use addMapLayer method to add a single layer or addMapLayers to add more than one layer in a list.  Moreover we can also observe that the layer to be displayed must be defined as vector layer with QgsVectorLayer class. So firstly we have to define the layer using the class with the following syntax:

layer=QgsVectorLayer(path:string,layer_name:string,library:string)

The listing code below shows the example of adding layer with this method.

#Defining Layer with QgsVectorLayer
street_layer=QgsVectorLayer("F:\public_streets.shp","street","ogr")
park_layer=QgsVectorLayer("F:\park_polygons.shp","park","ogr")
crime_layer=QgsVectorLayer("F:\crime_shp_2016.shp","crime","ogr")
#Add a list of layer 
QgsProject.instance().addMapLayers([street_layer,park_layer])
#Add a single layer 
QgsProject.instance().addMapLayer(crime_layer)

Figure 2 shows the result of the code execution. Can be seen that those three layers have been added to the layer tree and displayed in QGIS map canvas.
Adding vector layer with QGIS Project instance
Figure 2. Adding vector layer with QGIS Project instance
That's all the tutorial how to add vector layer with python in QGIS. Please see other tutorials which discuss other topics in this QGIS python programming tutorial series.

Python QGIS Tutorial: Adding CSV Data

Comma Separated Value (CSV) data usually contains coordinates information that can be displayed on a map. We are frequently dealing with CSV data when working in GIS. Therefore it's very important we can work with this data type.

In QGIS we can display CSV data using Add Delimited Text Layer tool. If you haven't tried it, it can be found in Add Layer menu. But we will not discuss about it, because this is a series tutorial of QGIS Python programming. So we will discuss how to add a CSV data using python code in QGIS.

Checking CSV Data

Before adding a CSV data with python code in QGIS. I suggest to check the CSV data firstly.  The reason is, a CSV data may be structured in different way for any CSV data. It could be found the value is separated with comma(,), semicolon (;), tab, etc. So knowing how it is formatted will be a guidance for us how to handle the data. To check the data can be done with opening it in a text application like Notepad, Notepad++, etc. Any text application can be used, just choose your favorite one. 

In this tutorial I will use a CSV data from USGS Earthquake Catalog that contains earthquake location around the world. I opened the data with Notepad++ and it looks like figure 1.

Screenshot of earthquake csv data
Figure 1. Screenshot of earthquake CSV data
From the screenshot in figure 1, can be seen that the values are separated with comma(,) and there is no parentheses. So we will use this format in coding later.

Adding CSV Data

Now, let's adding the CSV data into QGIS with python. To add a CSV data we're using QGIS interface or QGIS project instance. But firstly we have to make the CSV file as a vector layer using QgsVectorLayer class with the following expression:
csv_layer=QgsVectorLayer(uri:string,layer name: string, library:string)

In the expression above can be seen that the first argument is we have to define an Uniform Resources Identifier (URI) for the CSV file. An URI could be a file address on a server or over internet. If the file is at local machine the URI can be written as file:///path. The second argument is layer name, give it what ever you like. The third one is library, for CSV we are using delimitedtext

Next we have to define some parameters in the URI such as:
  • encoding : data encoding (optional)
  • xField : Column name for longitude value
  • yField : Column name for latitude value
  • crs : Coordinate system in EPSG number
The complete URI for the earthquake CSV data will be like the expression below:

uri = "file:///F:/eq-data.csv?encoding=%s&delimiter=%s&xField=%s&yField=%s&crs=%s" % ("UTF-8",",", "longitude", "latitude","epsg:4326")

So making sure the URI is correct is the key to add a CSV data. After defining the URI, then we make a vector layer and add it to QGIS. The complete code can be seen as below.

uri = "file:///F:/eq-data.csv?encoding=%s&delimiter=%s&xField=%s&yField=%s&crs=%s" % ("UTF-8",",", "longitude", "latitude","epsg:4326")

#Make a vector layer
eq_layer=QgsVectorLayer(uri,"eq-data","delimitedtext")

#Check if layer is valid
if not eq_layer.isValid():
    print ("Layer not loaded")

#Add CSV data    
QgsProject.instance().addMapLayer(eq_layer)

Figure 2 shows the result when the code was executed.
Earthquake csv data added with python
Figure 2. Earthquake csv data added with python


That's all the tutorial how to add vector layer with python in QGIS. Please see other tutorials which discuss other topics in this QGIS python programming tutorial series.

Python QGIS Tutorial: Drawing Marker

Marker is a point on the map which is not a point feature. It just a graphic on the map that can be used to mark a location. In this tutorial we will see how to draw a marker on the QGIS map canvas.

Before drawing a marker we have to define the canvas using QGIS interface with the following expression. Otherwise the marker will not show up.

canvas=iface.mapCanvas()

To draw a marker, we are using QgsVertexMarker class. Therefore we need to create a marker using the class into the canvas.

m1 = QgsVertexMarker(canvas)

After creating the marker, then we can set the marker position with setCenter method. Moreover the color, icon and the width line of the marker can be customized with  other methods as in the following code.

m1.setCenter(QgsPointXY(260,104))
m1.setColor(QColor(255,0, 0)) #(R,G,B)
m1.setIconSize(10)
m1.setIconType(QgsVertexMarker.ICON_X)
m1.setPenWidth(3)

The icon that used in the code was X in red color. There are several type of icons avalailabe:
  • ICON_BOX
  • ICON_CIRCLE
  • ICON_CROSS
  • ICON_DOUBLE_TRIANGLE
  • ICON_NONE
  • ICON_X
If you want to change the icon, just change the icon type in the fourth line with a type above. The icon looks like as in figure 1. Must be remember, every time you want to create a maker, a new marker must be defined using the QgsVertexMarker class.

Marker icon type QGIS
Figure 1. Marker icon type
Lastly the marker can be hidden and shown using hide and show method. So the expression will be like this.

m1.hide() # to hide a marker
m1.show() # to show a marker

Take a look the QGIS Python API documentation to learn more about QgsVertexMarker class. 

Friday, May 11, 2018

How to Create Contour Lines in QGIS

This tutorial will discuss how to create contour lines in QGIS. But before we are discussing about drawing contour lines in QGIS, let's start talking a little bit about the concept of contour line, what contour line definition is, contour lines characteristics and and how the method in drawing contour line . In the last section we will see how to create contour lines from points and also from DEM in QGIS. For this tutorial I am using QGIS 3.0.

What is Contour Line?

Contour is a height isoline. Then what is an isoline? Isoline is a line on a map with a constant value. Therefore the contour definition is a line which has a constant elevation along the line. There are some characteristics of contour line such as contour interval, contour elevation and distance between contour lines. If we see a contour map we found many lines on it. Those lines are never crossing to each other and each line has a unique elevation value that define its elevation on the map. Contour interval defines height interval between the line. For example if a contour map has 5 m interval contour, it means the height different between the adjacent contour line is 5 m. If we know the distance between contour lines then we can calculate the slope. The rule of thumb in defining contour interval is 1/2000 x map scale.       



How to Draw Contour Lines

The basic method in drawing contour line is interpolation from a point dataset. The simplest interpolation method is linear interpolation which is frequently used in drawing contour lines manually.  In drawing a contour line with linear interpolation, at least three steps need to be done:
  1. Connecting each point to form a TIN (Triangular Irregular Network) structure.
  2. Do interpolation for each line in the TIN.
  3. Connecting the interpolation point with the same elevation.
Figure 1 shows step by step approach in drawing contour line from connecting each point to form TIN, do interpolation along the TIN line and drawing contour line by connecting the same elevation of interpolation points.
Contour line drawing
Figure 1. Contour line drawing
We already got some basic concepts about contour line. Now let's move to technical part how to create contour lines with QGIS. Firstly we will see how to create contour line from a set of height points.

Create Contours from Points

To create contour lines from points, we are using Contour plugin. As usual we can get the plugin from Plugins menu >> Manage and install plugins. In the search option type "Contour". You will see the Contour plugin appears as in figure 2. Check it and Install the plugin if you don't have it in your QGIS 3.

Installing Contour Plugin QGIS 3
Figure 2. Installing Contour Plugin
Now, let's make contour lines from points dataset with the following steps:

1. Add the points into QGIS map canvas as in figure3.
Points dataset
Figure 3. Points dataset
2. Open the Contour plugin by clicking the icon. The contour plugin window will appear as in figure 4.
3. Select the point's layer name and the elevation data field. In the contouring menu, we can choose  to create contour line, filled contour or both. Next, we have to define the number of contour lines. We can set any number, but I suggest the number of contour line follows the rule of thumb as mentioned above. For example if we want to make a contour map in 100.000 scale, then the contour interval will be 1/2000 x 100.000 = 50 m. Then the number of contour line will be: maximum elevation/contour interval. 3554/50=71 lines. When we set the calculated number, the contour interval can be seen in the right list. We can tweak the number if the interval is not correct. That's why I change the number to 72. Moreover the interval can be edited manually by double click the interval number.
4. Next we can specify the output name. The output will be written to memory, so we need to export to it to a GIS file to keep it permanently.
5. In the next setting we can specify the precision and unit for the contour label. Furthermore the color graduation for the contour lines can be changed by selecting a color ramp for the output contour lines.
6. When finish click the Add button. The result can be seen in the QGIS map canvas without closing the contour plugin window, therefore we can switch back to change a parameter if the result not looks as we want. Figure 5 shows the result of contour lines.     

Contour plugin window
Figure 4. Contour plugin window

Contour lines result from points
Figure 5. Contour lines result from points

Create Contours from DEM

In this part we will create contour lines from DEM. To create contour lines from DEM we can use the Contour tool in the Raster menu >> Extraction >> Contour. See figure 6.

QGIS Contour tool
Figure 6. Contour tool
To create the contour lines from DEM, can be done with several steps:
1. Add DEM data into QGIS map canvas, as in figure 7. If you don't have one. You can download DEM data directly from QGIS.

DEM data in QGIS map canvas
Figure 7. DEM data in QGIS map canvas
2. Open the Contour tool. The contour tool window will appear as in figure 8.

QGIS Raster contour tool
Figure 8. Raster contour tool
3. Select the DEM data. Then set the contour interval, in this case  I used 50 m. In the Contours option, we can specify the contour output path or save it to temporary file.
4. If all the settings already set. Click the Run in Background button. The result will shown in QGIS  map canvas as in figure 9.

Contour lines from DEM
Figure 9. Contour lines from DEM

Contour Smoothing

If we compare the result, the contour lines from DEM looks smoother than points. That's because the DEM contains more elevation information compare to points dataset. The DEM has information in each pixel, so it could consist of thousands pixel and elevation information. But don't worry. We can make the contour line from points smoother using Smooth tool as follow:

1.  In the processing toolbox type "smooth". The smooth tool will appear under  Vector Geometry, as in figure 10.
Vector geometry smooth tool
Figure 10. Vector geometry smooth tool
2. Open the smooth tool. Select the input contour layer. Set the number of iteration, offset and maximum angle to smooth. The explanation of each parameter can be seen in the help window on the right. See figure 11.

Contour smoothing
Figure 11. Contour smoothing
3. If all the parameters have been set. Click Run in Background button to start the smoothing process. Figure 12 shows the comparison of contour lines before and after smoothing.
Contour lines before and after smoothing
Figure 12. Contour lines before and after smoothing
That's all the tutorial how to create contour lines in QGIS both from points and DEM. I hope it's useful. Please give comment if you find another tool in QGIS or any improvement in creating contour lines using QGIS software.


Saturday, April 14, 2018

QGIS Tutorial: How to Create an Action

The previous tutorial discussed about how to perform spatial query in QGIS 3. In this tutorial, I will discuss about Actions in QGIS. What is Action? Action is a command or task to be executed based on a click trigger. For example if a user click a feature then an image will be opened. Here the click is the trigger and opening an image is the action.

We can define some types of actions in QGIS like opening an application or file in an operating system (Windows, Mac or UNIX), open a web page,  or create a specific action with a python script. In this tutorial I will cover how to create an open action and a python action. The action example is to open an URL or web link of traffic camera in Vancouver city. For that, we will use traffic camera location data that can be downloaded from Vancouver City Data Catalogue. For this tutorial I'm using QGIS 3.0.

Opening an URL Action

Let's start with creating a simple action to view a camera traffic when a camera point is clicked on QGIS map canvas. The action can be done with the following steps:

1. Add a basemap. For this case I'm using Open Street Map.
2. Add traffic camera data into QGIS map canvas. The camera location over the basemap will look like as in figure 1.

QGIS Action Vancouver City Traffic Camera Location
Figure 1. Vancouver City Traffic Camera Location
3. Open camera location layer properties and choose Actions icon. In this step, by clicking Create default actions, some actions will be listed in the Action List. To create a new one click the green + button as in figure 2.

Add a new action in QGIS
Figure 2. Add a new action
4. An Add action window will appear as in figure 3. In the Type option select Open. Give a description and a short name, for example I gave the name open traffic camera. In the Action Scope check Canvas and Feature Scope, cause we want the action will be executed for the feature in QGIS map canvas.
5. Next in the Action text, type [%URL%]. URL is a column name in the attribute table which contains a HTML link to each  traffic camera. You can see how the attribute looks like in figure 4.

Edit an action in QGIS
Figure 3. Add/Edit an action
Traffic camera attribute table QGIS Action
Figure 4. Traffic camera attribute table
6. Last step after clicking OK button both in the Add new action window and layer properties, let's run the action by selecting the action from Action Icons in the toolbar as in figure 5.
Select an Action QGIS
Figure 5. Select an Action
The mouse pointer will change into cross hair. Click a camera point. Soon the URL of the camera will be opened if the default web browser as in figure 6.
Traffic camera view in web browser
Figure 6. Traffic camera view in web browser

Python URL Action

In this part we will make an action with Python script for the same example. The action will open a QT web browser instead of the default web browser. Most steps are similar like above, except:
1. Select Python in the type option.
2. Give a new name for example: open traffic camera python.
3. In the Action text, type the listing code below. Figure 7 shows the new action with python.
  
from PyQt5 import QtCore
from PyQt5.QtCore import *
from PyQt5.QtWebKit import *
from PyQt5.QtWebKitWidgets import *
web = QWebView()
web.load(QUrl("[%URL%]"))
web.show()

QGIS Action with Python
Figure 7. QGIS Action with Python
When the action is executed by clicking a traffic camera point. A QT web browser will be opened as in figure 8.
Traffic camera view in QT web browser
Figure 8. Traffic camera view in QT web browser
Congratulation! Now you can observe the traffic condition around Vancouver City.
That's all the tutorial how to create an Action in QGIS. I think this is a great tool, cause gives more functionality to our map application. If you want to explore more, visit QGIS actions documentation.         

Monday, April 2, 2018

Creating GPX File Viewer in Python

GPX file is a GPS data that stored in XML format. There are some types of data that stored in GPS namely waypoint, route and track. Track is a type of data which is recorded regularly by  GPS in an interval. Therefore, with GPS tracker data we can visualize a trip that we did, which road we passed by, the length of track, time taken, etc. There are some tools or apps that available to view or visualize a GPS track file such as ArcGIS, QGIS, Google Earth and also some online GPX viewer app also available to visualize a GPS tracker data. Just do online searching and you will find many of them.

In this tutorial, I will discuss how to create a GPX track file viewer in Python. I'm using Jupyter Notebook with Python 3. So if you want to run the code make sure you have Jupyter Notebook installed in your machine. We will visualize the track itself, speed and elevation profile along the track. Moreover the visualization will be dynamic in an animation that plot each GPS tracking point as in the figure below.
To create the GPX file viewer we need some modules like matplotlib, time, IPython, xml dom and math. Import all the required modules as in the code below.

import matplotlib.pyplot as plt
import time
from IPython import display
from xml.dom import minidom
import math  

Then open a GPX file and parse the data using minidom. Here we are taking some important elements: 'trkpt', 'ele' and 'time' for track point, elevation and time.

#READ GPX FILE
data=open('F:/ride.gpx')
xmldoc = minidom.parse(data)
track = xmldoc.getElementsByTagName('trkpt')
elevation=xmldoc.getElementsByTagName('ele')
datetime=xmldoc.getElementsByTagName('time')
n_track=len(track)

After getting those elements, we need to get the value from each elements. Then we need to parse the elements to get latitude, longitude, elevation and time. At the end, the time is converted to second. It will be used later in speed calculation.

#PARSING GPX ELEMENT
lon_list=[]
lat_list=[]
h_list=[]
time_list=[]
for s in range(n_track):
    lon,lat=track[s].attributes['lon'].value,track[s].attributes['lat'].value
    elev=elevation[s].firstChild.nodeValue
    lon_list.append(float(lon))
    lat_list.append(float(lat))
    h_list.append(float(elev))
    # PARSING TIME ELEMENT
    dt=datetime[s].firstChild.nodeValue
    time_split=dt.split('T')
    hms_split=time_split[1].split(':')
    time_hour=int(hms_split[0])
    time_minute=int(hms_split[1])
    time_second=int(hms_split[2].split('Z')[0])
    total_second=time_hour*3600+time_minute*60+time_second
    time_list.append(total_second)

Now we have all parameters that are needed to create the GPX file viewer. To do further calculation like distance and speed, we create some functions namely geo2cart, distance and speed. geo2cart is a function to convert Geodetic coordinate into Cartesian coordinate with WGS-84 ellipsoid reference datum. We need the Cartesian coordinate because it will be used to calculate a distance in meter between two GPS tracking points. After getting the distance, then the speed in m/s can be calculated by divided with time difference.
 
#GEODETIC TO CARTERSIAN FUNCTION
def geo2cart(lon,lat,h):
    a=6378137 #WGS 84 Major axis
    b=6356752.3142 #WGS 84 Minor axis
    e2=1-(b**2/a**2)
    N=float(a/math.sqrt(1-e2*(math.sin(math.radians(abs(lat)))**2)))
    X=(N+h)*math.cos(math.radians(lat))*math.cos(math.radians(lon))
    Y=(N+h)*math.cos(math.radians(lat))*math.sin(math.radians(lon))
    return X,Y

#DISTANCE FUNCTION
def distance(x1,y1,x2,y2):
    d=math.sqrt((x1-x2)**2+(y1-y2)**2)
    return d

#SPEED FUNCTION
def speed(x0,y0,x1,y1,t0,t1):
    d=math.sqrt((x0-x1)**2+(y0-y1)**2)
    delta_t=t1-t0
    s=float(d/delta_t)
    return s

Next step, using the function we do calculation of distance and speed. The result are stored in distance and speed list.

#POPULATE DISTANCE AND SPEED LIST
d_list=[0.0]
speed_list=[0.0]
l=0
for k in range(n_track-1):
    if k<(n_track-1):
        l=k+1
    else:
        l=k
    XY0=geo2cart(lon_list[k],lat_list[k],h_list[k])
    XY1=geo2cart(lon_list[l],lat_list[l],h_list[l])
    
    #DISTANCE
    d=distance(XY0[0],XY0[1],XY1[0],XY1[1])
    sum_d=d+d_list[-1]
    d_list.append(sum_d)
    
    #SPEED
    s=speed(XY0[0],XY0[1],XY1[0],XY1[1],time_list[k],time_list[l])
    speed_list.append(s)

Now we can visualize the GPX file track data with three plots: Track, speed and elevation profile.

#PLOT TRACK
f,(track,speed,elevation)=plt.subplots(3,1)
f.set_figheight(8)
f.set_figwidth(5)
plt.subplots_adjust(hspace=0.5)
track.plot(lon_list,lat_list,'k')
track.set_ylabel("Latitude")
track.set_xlabel("Longitude")
track.set_title("Track Plot")

#PLOT SPEED
speed.bar(d_list,speed_list,30,color='w',edgecolor='w')
speed.set_title("Speed")
speed.set_xlabel("Distance(m)")
speed.set_ylabel("Speed(m/s)")

#PLOT ELEVATION PROFILE
base_reg=0
elevation.plot(d_list,h_list)
elevation.fill_between(d_list,h_list,base_reg,alpha=0.1)
elevation.set_title("Elevation Profile")
elevation.set_xlabel("Distance(m)")
elevation.set_ylabel("GPS Elevation(m)")
elevation.grid()

Lastly in a loop, dynamic plot is created in one second interval by plotting each tracking points over the track and elevation profile. The speed change is visualized in a bar chart.
 
#ANIMATION/DYNAMIC PLOT
for i in range(n_track):
    track.plot(lon_list[i],lat_list[i],'yo')
    speed.bar(d_list[i],speed_list[i],30,color='g',edgecolor='g')
    elevation.plot(d_list[i],h_list[i],'ro')
    display.display(plt.gcf())
    display.clear_output(wait=True)
    time.sleep(1)
    
plt.show()


That's all the tutorial how to create a GPX file viewer to visualize track, speed and elevation profile with Python in Jupyter Notebook. You can try to run the code by downloading the GPX file viewer Jupyter Notebook and a GPX file sample that used in this tutorial. Finally use your own GPX track file and try to visualize it with the code.