7

I have a CSV which contains the geometry information as a GeoJSON string, e.g.

{ "type": "MultiLineString",
    "coordinates": [
        [[10, 10], [20, 20], [10, 40]],
        [[40, 40], [30, 30], [40, 20], [30, 10]]
    ]
}

which is just a text field.

I would like to display this geometry in QGIS 3. I could go via PostGIS and use ST_GeomFromGeoJSON, but would prefer a direct way, without adding my data to a PostGIS database first.

Unfortunately, virtual layers do not recognize this function and as expressions there only seem to be geom_from_wkt() and geom_from_gml().

Is there a possibilty to read this GeoJSON-String as Geometry directly from a CSV column?

MrXsquared
  • 34,292
  • 21
  • 67
  • 117

3 Answers3

9

If you can retrieve the GeoJSON text in a Python object, you can execute the code below :

# geoj is the string object that contains the GeoJSON
geoj = """
{ "type": "MultiLineString",
    "coordinates": [
        [[10, 10], [20, 20], [10, 40]],
        [[40, 40], [30, 30], [40, 20], [30, 10]]
    ]
}
"""
# PyQGIS has a parser class for JSON and GeoJSON
feats = QgsJsonUtils.stringToFeatureList(geoj, QgsFields(), None)
# if there are features in the list
if len(feats) > 0:
    # define the geometry type of the layer
    geom_type = feats[0].geometry().type()
    if geom_type == QgsWkbTypes.PointGeometry:
        geotype = "MultiPoint"
    elif geom_type == QgsWkbTypes.LineGeometry:
        geotype = "MultiLineString"
    elif geom_type == QgsWkbTypes.PolygonGeometry:
        geotype = "MultiLinePolygon"
    else:
        print("geometry type not defined in the script")
    # create a new memory layer with the features
    vl = QgsVectorLayer(geotype, "geojson_layer", "memory")
    with edit(vl):
        vl.addFeatures(feats)
        vl.updateExtents()
    # add this brand new layer
    QgsProject.instance().addMapLayer(vl)
else:
    print("no features found in the geoJSON")
J. Monticolo
  • 15,695
  • 1
  • 29
  • 64
  • 1
    I compared 3.4 and latest API, and they added a QgsFields() default value, so I made a little edit in my code. – J. Monticolo Dec 21 '20 at 15:11
  • I'm getting an AssertionError (assert self.layer.startEditing() AssertionError). This happens even if I'm applying layer.startEditing() before the 'with edit' block, and layer.commitChanges() after it. – Casian Pavel Jan 11 '21 at 11:27
4

Based on J. Monticolo's answer I have modified his script. The first one is for python console and the second one a processing algorithm. Both will ignore features, that have an invalid GeoJSON-String as content.

Script for Python-Console:

# Define settings:
source_layer = QgsProject.instance().mapLayersByName("mylayer")[0] # Define Source-layer
new_layer = QgsVectorLayer('MultiLineString', "new_geojson_layer", "memory") # Define new layer and its type
source_geojsonfield = 'geojsonstring' # Fieldname containing the GeoJSON

No Changes needed below

Copy fields from old layer to new layer

source_fields = source_layer.fields() new_layer_pr = new_layer.dataProvider() new_layer.startEditing() new_layer_pr.addAttributes(source_fields) new_layer.commitChanges()

for feature in source_layer.getFeatures(): # geoj is the string object that contains the GeoJSON geoj = feature.attributes()[source_fields.indexFromName(source_geojsonfield)] # PyQGIS has a parser class for JSON and GeoJSON geojfeats = QgsJsonUtils.stringToFeatureList(geoj, QgsFields(), None) # if there are features in the list if len(geojfeats) > 0: new_geom = geojfeats[0].geometry() with edit(new_layer): new_feat = QgsFeature(feature) new_feat.setGeometry(new_geom) new_layer.addFeature(new_feat) new_layer.updateExtents() else: print("No features found in the GeoJSON or no valid GeoJSON-String")

add this brand new layer

QgsProject.instance().addMapLayer(new_layer)

Script for Processing-Tool:

The selection of the target-wkb-type in the processing tool is not an optimal implementation... if you know a better solution, feel free to edit it

from PyQt5.QtCore import QCoreApplication, QVariant
from qgis.core import (QgsJsonUtils, QgsField, QgsFeature, QgsProcessing, QgsExpression, QgsGeometry, QgsPoint, QgsFields, QgsWkbTypes,
                       QgsFeatureSink, QgsFeatureRequest, QgsProcessingAlgorithm,
                       QgsProcessingParameterFeatureSink, QgsProcessingParameterCrs, QgsProcessingParameterField, QgsProcessingParameterFeatureSource, QgsProcessingParameterEnum, QgsProcessingParameterString, QgsProcessingParameterNumber)

class GeometryLayerFromGeojsonStringField(QgsProcessingAlgorithm): SOURCE_LYR = 'SOURCE_LYR' GEOJSON_FIELD = 'GEOJSON_FIELD' #GEOMETRYTYPE_STRING = 'GEOMETRYTYPE_STRING' GEOMETRYTYPE_ENUM = 'GEOMETRYTYPE_ENUM' CRS = 'CRS' OUTPUT = 'OUTPUT'

def initAlgorithm(self, config=None):  
    self.addParameter(
        QgsProcessingParameterFeatureSource(
            self.SOURCE_LYR, self.tr('Source'), [QgsProcessing.TypeMapLayer])) # Take any source layer, unfortunately no-geometry layers will not be available...
    self.addParameter(
        QgsProcessingParameterField(
            self.GEOJSON_FIELD, self.tr('Field containing the GeoJSON as string'),'GeoJSON','SOURCE_LYR', 1)) # Choose the field containing the GeoJSON as string
    #self.addParameter(
    #    QgsProcessingParameterNumber(
    #        self.GEOMETRYTYPE_STRING, self.tr('Geometry type of the target layer / of the GeoJSON content as number (lookup at https://qgis.org/pyqgis/3.0/core/Wkb/QgsWkbTypes.html)'),0,5)) # Unfortunately there is no WKB-Type-Input available...
    self.addParameter(
        QgsProcessingParameterEnum(
            self.GEOMETRYTYPE_ENUM, self.tr('Geometry type of the target layer / of the GeoJSON content'),
            ['Unknown','Point','LineString','Polygon','MultiPoint','MultiLineString','MultiPolygon','GeometryCollection','CircularString','CompoundCurve','CurvePolygon'],defaultValue=5)) # Only Works because these are ascending numerated in QGIS... NOT A GOOD SOLUTION!! But better than typing in a number by hand... see https://qgis.org/api/classQgsWkbTypes.html
    self.addParameter(
        QgsProcessingParameterCrs(
            self.CRS, self.tr('CRS of the target layer / of the GeoJSON content'),'EPSG:4326')) # CRS of the targetlayer
    self.addParameter(
        QgsProcessingParameterFeatureSink(
            self.OUTPUT, self.tr('new_geojson_layer'))) # Output

def processAlgorithm(self, parameters, context, feedback):
    # Get Parameters and assign to variable to work with
    source_layer = self.parameterAsLayer(parameters, self.SOURCE_LYR, context)
    source_geojsonfield = self.parameterAsString(parameters, self.GEOJSON_FIELD, context)
    #wkbgeometrytype = self.parameterAsInt(parameters, self.GEOMETRYTYPE_STRING, context)
    wkbgeometrytype_fromenum = self.parameterAsInt(parameters, self.GEOMETRYTYPE_ENUM, context)
    wkbgeometrytype = wkbgeometrytype_fromenum # testing assignment        
    crsgeometry = self.parameterAsCrs(parameters, self.CRS, context)

    total = 100.0 / source_layer.featureCount() if source_layer.featureCount() else 0 # Initialize progress for progressbar

    source_fields = source_layer.fields() # get all fields of the sourcelayer

    (sink, dest_id) = self.parameterAsSink(parameters, self.OUTPUT, context, source_fields, wkbgeometrytype, crsgeometry)

    for current, feature in enumerate(source_layer.getFeatures()): # iterate over source 
        # geoj is the string object that contains the GeoJSON
        geoj = feature.attributes()[source_fields.indexFromName(source_geojsonfield)]
        # PyQGIS has a parser class for JSON and GeoJSON
        geojfeats = QgsJsonUtils.stringToFeatureList(geoj, QgsFields(), None)
        # if there are features in the list
        if len(geojfeats) > 0:
            new_geom = geojfeats[0].geometry()
            new_feat = QgsFeature(feature)
            new_feat.setGeometry(new_geom)
            sink.addFeature(new_feat, QgsFeatureSink.FastInsert) # add feature to the output

        if feedback.isCanceled(): # Cancel algorithm if button is pressed
            break

        feedback.setProgress(int(current * total)) # Set Progress in Progressbar

    return {self.OUTPUT: dest_id} # Return result of algorithm



def tr(self, string):
    return QCoreApplication.translate('Processing', string)

def createInstance(self):
    return GeometryLayerFromGeojsonStringField()

def name(self):
    return 'GeometryLayerFromGeojsonStringField'

def displayName(self):
    return self.tr('New Layer from GeoJSON String')

def group(self):
    return self.tr('FROM GISSE')

def groupId(self):
    return 'from_gisse'

def shortHelpString(self):
    return self.tr('This Algorithm takes a source layer containing a GeoJSON as a String in a field and creates a copy of this layer with the geometry of this GeoJSON field')

MrXsquared
  • 34,292
  • 21
  • 67
  • 117
  • 1
    For the Wkb type, why don't you use a Enum parameter ? https://qgis.org/pyqgis/3.2/core/Processing/QgsProcessingParameterEnum.html?highlight=qgsprocessingparameterenum – J. Monticolo Dec 21 '20 at 20:46
  • @J.Monticolo just done that, but its not a good solution as I need to rewrite the entire list. Was trying to solve it as described here https://stackoverflow.com/questions/29503339/how-to-get-all-values-from-python-enum-class but getting TypeError: 'sip.enumtype' object is not iterable when trying to access QgsWkbTypes.Type like this. – MrXsquared Dec 21 '20 at 20:52
  • Why don't you delete this parameter and guess the geometry type after geojfeats as I wrote it in my code ? – J. Monticolo Dec 21 '20 at 21:40
  • @J.Monticolo because I need to define the wkb-type when creating the output layer. If I would guess it, I'd either have to iterate over the source layer twice (one time before creating the output, and one time after creating it), or pick a random example in hope it is a valid one. – MrXsquared Dec 21 '20 at 21:55
  • I think it's not the case, look at my code in my post, I just take the first feature after geojfeats and extract the geometry type. – J. Monticolo Dec 21 '20 at 22:00
  • 1
    And they are dirty methods to get types list with values : [(eval(f"QgsWkbTypes.{i}"), i) for i in dir(QgsWkbTypes) if type(eval(f"QgsWkbTypes.{i}")) == QgsWkbTypes.Type] – J. Monticolo Dec 21 '20 at 22:29
1

As I am not allowed to comment I have to make an annotation answer and just because it took me a while testing the answer and finding the problem:

There is no geotype = "MultiLinePolygon" rather should it be geotype = "MultiPolygon"

mach0
  • 31
  • 2