Solved: QGIS 3.20 Python Console shows an error when running Clip Vector by Extent tool

Fellow researchers and open-source GIS enthusiasts,

Welcome to my blog!

I’d like to start with a disclaimer – I may be a researcher of this very area but that doesn’t mean everything I do or write here will work for you, in your own desktop configurations and package versions. I have no responsibility if you lose data or mess up your installation. I also do not authorize any copies of my content.

Today, I am writing about an error I obtained while using QGIS 3.20 Python Console. I was trying to clip a vector layer populated with points I had just created, by an extent of the type $xmin,xmax,ymin,ymax$, using the tool Clip Vector by Extent, from GDAL.

The error

First, I ran the following line, to set the tool parameters:

paramgdal= {'INPUT': points_layer, 'EXTENT':extent_final, 'OUTPUT':output_path}

The variables “extent_final” and “output_path” were previously set. “points_layer” is the vector layer of the type shapefile in which the points I want to clip by the extent are stored. It was loaded previously, using the following snippet:

points_layer = QgsVectorLayer(path_to_file, 'points_loaded', 'ogr')

Well, so far, so good.

But, after setting the parameters, I ran:

processing.run("gdal:clipvectorbyextent", paramgdal)

and then, this error appeared:

ERROR 1: Attempt to write non-multipoint (POINT) geometry to multipoint shapefile. 
ERROR 1: Unable to write feature 0 from layer points_layer. 
ERROR 1: Terminating translation prematurely after failed translation of layer points_layer (use -skipfailures to skip errors)

And the resulting layer was not correctly generated.

Interpretation of the issue

Re-reading the error message, I noticed that GDAL, for some reason, tries to save the generated layer as a Multipoint shapefile, instead of, simply, a Point shapefile. For this reason, it is unable to write any of the features from the original layer on the new one, because their types do not match.

Solution

Checking on the QGIS documentation we can notice that there are four possible parameters to run gdal:clipvectorbyextent. One of them, “OPTIONS”, is Optional, and refers to GDAL creation options. Because the algorithm Clip Vector by Extent is based on the GDAL utility ogr2ogr, the options available for the original algorithm can be used to build the new layer through the Python Console as well.

Especially, one of the options is useful in this scenario. The “nlt” option of ogr2ogr defines the geometry of the generated file, which is just what we were looking for. In the case here shown, I simply added

'OPTIONS': '-nlt point'

to the $gdalparams$ Python dictionary. This tells GDAL that the layer to be created is a point layer, and it solved my issue completely. The code snippet shown at the start of this post was adapted by me as:

paramgdal= {'INPUT': points_layer, 'EXTENT':extent_final, 'OPTIONS': '-nlt point', 'OUTPUT':output_path}
processing.run("gdal:clipvectorbyextent", paramgdal)

And when I ran the snippet above, the clipped point layer was correctly generated, without any errors or warnings.

A simple fix!

Extras

• Why don’t you simply use GDAL from the command line, or click on the tool through the Processing Toolbox and use it directly from the GUI? Why do you want to use the Python Console on QGIS?

Because sometimes the Python Console is the most practical way to perform a calculation, and this code snippet shown here is a small part of a large code. For this motive, it is not worth it to change the whole interface or language in which we are developing code.

• I loaded my vector file using the QGIS Python Console, but it did not appear on the layers list on QGIS! What happened?

This is not an issue. The code snippet used to load the vector layer does not list the layer on QGIS Layers panel, nor shows it in the map visualization. To do that, you should run

QgsProject.instance().addMapLayers([points_layer])

• I’m not sure what I should import before using the tool “gdal:clipvectorbyextent” on QGIS Python Console.

You should import the QGIS algorithms and the processing tool. See the snippet below:

from qgis.core import *
import processing