An example of QGis GeoAlgorithm: verify geometries with Oracle engine🔗

Posted by Médéric Ribreux 🗓 In blog/ Qgis/

Introduction

Since QGis 2.0, you can use GeoAlgorithms (aka Processing or geoprocessing) to manipulate your data. But you can also develop your own GeoAlgo to do things that are not (already) included in QGis. As GeoAlgo is relatively young, the development on this part of QGis is important.

Today I propose a very basic GeoAlgo that takes an Oracle layer, use internal Oracle Spatial procedure (SDO_GEOM.VALIDATE_GEOMETRY_WITH_CONTEXT) to find non valid geometries and returns a layer with those geometries and an explanation of the errors that have been found.

Install cx_Oracle

The code needs cx_Oracle Python library to work. This library is a low-level communication Python (DBAPIv2) library for inter-connections to Oracle databases. For the moment it is not included in QGis binary installation, so you have to install yourself first.

If you are running a GNU/Linux distribution and have installed QGis from the official repositories, you'll have to install Oracle Client and SDK and then install cx_Oracle. Read this…

If you are running QGis on MS-Windows, you have to download the good version of cx_Oracle binaries (32 or 64 bits) and extract and copy cx_Oracle.pyd file in C:\Program Files\QGis Chugiak/apps/Python2.7/DLL .

To check for a valid installation of cx_Oracle in the QGis Python environment:

import cx_Oracle

Source code for VerifyOracle.py

Here is the source of the GeoAlgo. Just save it in a file named VerifyOracle.py ans place it in the good directory. You'll have to find the directory that stores the official Processing plugin.

On MS-Windows systems, it is stored in C:\Program Files\QGIS Chugiak\apps\qgis\python\plugins\processing\algs\qgis.

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

#***************************************************************************
#    VerifyOracle.py
#    ---------------------
#    Date                 : September 2014
#    Copyright            : (C) 2014 by Médéric RIBREUX
#    Email                : mederic.ribreux@medspx.fr
#***************************************************************************
#*                                                                         *
#*   This program is free software; you can redistribute it and/or modify  *
#*   it under the terms of the GNU General Public License as published by  *
#*   the Free Software Foundation; either version 2 of the License, or     *
#*   (at your option) any later version.                                   *
#*                                                                         *
#***************************************************************************
#
#You need to install cx_Oracle under QGis Python directory.

__author__ = 'Médéric RIBREUX'
__date__ = 'September 2014'
__copyright__ = '(C) 2014, Médéric RIBREUX'

# This will get replaced with a git SHA1 when you do a git archive

__revision__ = '$Format:%H$'

from osgeo import gdal
import cx_Oracle
import re
from qgis.core import *
from PyQt4.QtCore import *
from processing.core.GeoAlgorithm import GeoAlgorithm
from processing.parameters.ParameterVector import ParameterVector
from processing.parameters.ParameterTableField import ParameterTableField
from processing.outputs.OutputVector import OutputVector
from processing.tools import dataobjects, vector
from processing.tools.general import *

class VerifyOracle(GeoAlgorithm):

INPUT_VECTOR = 'INPUT_VECTOR'
OUTPUT = 'OUTPUT'
FIELD = 'FIELD'

def defineCharacteristics(self):
	self.name = 'Verify layer geometries with Oracle Spatial engine'
	self.group = 'Vector analysis tools'

	self.addParameter(ParameterVector(self.INPUT_VECTOR, 'Vector layer',
					  [ParameterVector.VECTOR_TYPE_ANY]))
	self.addParameter(ParameterTableField(self.FIELD, 'UID Field for input vector layer',
										  self.INPUT_VECTOR))
	self.addOutput(OutputVector(self.OUTPUT, 'Result Vector layer'))

def processAlgorithm(self, progress):

	uri = self.getParameterValue(self.INPUT_VECTOR)
	layer = dataobjects.getObjectFromUri(uri)
	# Deals with fields
	fieldName = self.getParameterValue(self.FIELD)
	fieldIdx = layer.fieldNameIndex(fieldName)
	fields = layer.dataProvider().fields()

	# Add the Errors field
	fields.append(QgsField('Errors', QVariant.String))

	# Get connection parameters
	regexp = re.compile(".*dbname='([^']+)'.*user='([^']+).*password='([^']+)'.*table=([^ ]+).*\(([^\)]+)\)")
	dbname, user, password, table, geocol = regexp.match(uri).groups()

	query = (u"SELECT c.{0}, SDO_GEOM.VALIDATE_GEOMETRY_WITH_CONTEXT (c.{1},"
	         u"	0.001) FROM {2} c WHERE SDO_GEOM.VALIDATE_GEOMETRY_WITH_CONTEXT
	         u"	(c.{1},	0.001) <> 'TRUE'").format(self.getParameterValue(self.FIELD),geocol, table))

	# Make the connection and the query
	connection = cx_Oracle.connect( user, password, dbname)
	c = connection.cursor()
	c.execute(query)
rows = c.fetchall()
	c.close()

	# Open a writer to the output vector layer
	writer = self.getOutputFromName(
			self.OUTPUT).getVectorWriter(fields,
										 layer.dataProvider().geometryType(),
										 layer.dataProvider().crs())

	# We get all of the features from the input vector layer
	features = vector.features(layer)

	# And make some computations for the progress bar
	total = 100.0 / float(len(features))
	current = 0

	outFeat = QgsFeature()

	# Build a list of errors (at least as big as number of features of the layer)
	errors = []
	for row in rows:
		errors.append({'GID':row[0], 'ERROR':row[1]})

	# Main loop
	for feature in features:
		gid = feature.attributes()[fieldIdx]
		# if the feature has got an error
		if gid in [x['GID'] for x in errors]:
			error = (x['ERROR'] for x in errors if x['GID'] == gid).next()
			geom = feature.geometry()
			attrs = feature.attributes()

			# write the feature to the output layer
			outFeat.setGeometry(geom)
			attrs.append(error)
			outFeat.setAttributes(attrs)
			writer.addFeature(outFeat)

		current += 1
		progress.setPercentage(int(current * total))

	del writer

Register GeoAlgo

For QGis to display the GeoAlgo in the algorithms tree, you have to register the Python module first. Registering the module is done by editing the file named QGISAlgorithmProvider.py.

Just add the following lines:

# on the declaration of modules
from VerifyOracle import VerifyOracle
…
	self.alglist = [SumLines(), PointsInPolygon(), # just add a line to this listVerifyOracle(),

In action

Just open an Oracle layer in QGis, then launch the GeoAlgo. It is located in QGis GeoAlgorithms and in Vector analysis tools:

GeoAlgo menu
GeoAlgo menu

Then you'll find the following dialog box:

GeoAlgo dialog box
GeoAlgo dialog box

Once run and done, you can find a layer named "Result Vector Layer" that contains the non-valid geometries. Every attributes are taken from the original layer and one column has been added at the end, which is named Errorsand contains what errors have been found by Oracle. The ring and edges numbers are very useful to find where is the problem as QGis uses the same numbers in its nodes tool:

GeoAlgo results
GeoAlgo results

Conclusion

This code is provided as is ! There is plenty of room for improvement. You should consider it as something experimental… But it works !