7

I'm trying to make a line pattern fill in QGIS using the geometry generator and a Python function.

The code is inspired from this answer : Save line geometry as points geometry in QGIS (pyqgis)

from qgis.core import *
from qgis.gui import *

@qgsfunction(args='auto', group='Custom')
def fill_line_pattern(x, feature, parent):

#create a base line that is 45 degrees
base_xmin = 2285865.0
base_xmax = 2669093.0
base_ymin = 1017594.0
base_ymax = 1399851.0

#calculating interval between parrallels
countX = int((2800000.0-2285865.0) / xInterval)
#empty list
pt=[]
pts = []
intersected_lines = []
#raw lines creation
for xOff in range(countX+1):
    ptXmin = base_xmin + xOff*(xInterval)
    ptXmax = base_xmax + xOff*(xInterval)
    ptmin = QgsPoint(ptXmin,base_ymin)
    ptmax = QgsPoint(ptXmax,base_ymax)
    pt = ([ptmin, ptmax])
    pts.append(pt)
#intersect lines with polygons
#for f in feat :
for p in pts :
    rawlines = QgsGeometry.fromPolyline(p)
    if feature.geometry().intersects(rawlines):
        geom = feature.geometry().intersection(rawlines)
        intersected_lines.append(geom.asPolyline())
#create lines
for line_ok in intersected_lines:
    return QgsGeometry.fromPolylineXY(line_ok)

The issue is that the first line only is crossing the polygon...

enter image description here

How can I have all the lines crossing the polygons?

In comparison to the inspirational code above, there's no loop for the polygon features.

In the function, the feature called is a QgsFeature which can't use getQgsFeatures().

If a for loop is used :

for f in feature :
    for p in pts :
    rawlines = QgsGeometry.fromPolyline(p)
        if f.geometry().intersects(rawlines):
            geom = f.geometry().intersection(rawlines)
            intersected_lines.append(geom.asPolyline())

This error message appears:

Evaluation error: 'int' object has no attribute 'geometry'

EDIT

The code below does almost what I want, it seems a little bit buggy because not all the intersections appear... The other disadvantage is the long time QGIS takes to render it on the canvas and in the composer.

On the other hand, the pdf output size is so light and is competitive with other software.

Maps are made to be shared: QGIS biggest issue is PDF output file size

The inspiration of the source : https://www.lutraconsulting.co.uk/blog/2015/06/05/qgis-function-editor/

Lutra Consulting uses a cache but I didn't find the way to do it...

from qgis.core import *
from qgis.gui import *

def _layer_intersection(layer, xInterval):
#make 45 degrees lines from layer extent
    ext = layer.extent()
    xmax = ext.xMaximum()
    ymax = ext.yMaximum()
    xmin = ext.xMinimum()
    ymin = ext.yMinimum()

    pts_hash_box = []
    line_in_polygon = []

    x_base_min = xmin - (ymax-ymin)
    countX = int(((xmax-xmin)+(ymax-ymin)) / xInterval)

    for xOff in range(countX+1):
        ptXmin = x_base_min + xOff*(xInterval)
        ptXmax = xmin + xOff*(xInterval)

        ptmin_hash_box = QgsPointXY(ptXmin,ymin)
        ptmax_hash_box = QgsPointXY(ptXmax,ymax)

        pt_hash_box = ([ptmin_hash_box, ptmax_hash_box])
        pts_hash_box.append(pt_hash_box)

    multipolylines = QgsGeometry.fromMultiPolylineXY(pts_hash_box)

    #intersection between 45 degrees virtual lines and the polygon features
    for layer_feature in layer.getFeatures() :
        if layer_feature.geometry().intersects(multipolylines):
            pts = layer_feature.geometry().intersection(multipolylines).asMultiPolyline()
            for pt in pts :
                line_in_polygon.append(y)

return QgsGeometry.fromMultiPolylineXY(line_in_polygon)

@qgsfunction(args='auto', group='Custom')
def line_pattern_fill(layer_name, xInterval, feature, parent):
for layer in QgsProject.instance().mapLayers().values():
    if layer.name() == layer_name :
        return _layer_intersection(layer, xInterval)
6
  • Must the lines be true geometry or could you use the QgsLinePatternFillSymbolLayer() instead? – Matt Needle Apr 26 '19 at 5:52
  • The lines doesn't have to be true geometry. QgsLinePatternFillSymbolLayer() is rasterized when exported in PDF. The consequence is that the PDF size increase considerably... So I was thinking about a new methodology to keep vectorized line pattern to preserve the size of my PDFs with QGIS. – geavcon Apr 29 '19 at 13:23
  • @geavcon did you find a solution to this? – TeddyTedTed Mar 16 '20 at 15:36
  • @geavcon: are you really sure that printing with line pattern is rasterized ? For me, it's exported to vectors (pdf). – J. Monticolo Mar 19 '20 at 14:46
  • @J.Monticolo I get the same rasterized result as geavcon reported. Which version of QGIS are you using? Are you exporting a PDF with a polygon that has a Line Pattern Fill? – TeddyTedTed Mar 19 '20 at 15:53
9
+100

Here is my working solution. I took a great part of your code but rearranged and improved it.

Create a custom function in the function editor and paste the code below.

Use it as : line_pattern_fill(@layer_id, 45, 1000) where 45 is the angle in degrees and 1000 is the X interval in layer units.

I create, as Lutra Consulting blog post, a cache. I reinitialize the variable if the angle or x_interval change so that QGIS compute new pattern fills.

from qgis.core import *
from qgis.gui import *
from math import ceil, cos, radians, sin

cache = {}


@qgsfunction(args='auto', group='geom_generator')
def line_pattern_fill(layer_id, angle, x_interval, feature, parent):
    """
    Returns a MultiLineString geometry with line pattern fill, contained
    in the current feature geometry.
    <h4>Syntax</h4>
    <div class="syntax">
    <code><span class="functionname">line_pattern_fill</span>
    (<span class="argument">@layer_id</span>,
    <span class="argument">angle</span>,
    <span class="argument">x interval</span>)</code>
    <h4>Arguments</h4>
    <div class="arguments">
    <table><tr>
    <td class="argument">@layer_id</td><td>ID of current layer</td></tr><tr>
    <td class="argument">angle</td><td>angle in degrees, from 0 to 179</td></tr><tr>
    <td class="argument">x interval</td><td>spacing between lines, in layer units</td></tr>
    </table>
    </div>
    <h4>Examples</h4>
    <div class="examples">
    <ul>
    <li><code>line_pattern_fill(@layer_id, 45, 1000)</code> &rarr; 
    <code>a line pattern fill at 45 degrees, spaced by 1000 layer units</code></li>
    <li><code>line_pattern_fill(@layer_id, 0, 60)</code> &rarr; 
    <code>a vertical line pattern fill, spaced by 60 layer units</code></li>
    <li><code>line_pattern_fill(@layer_id, 90, 150)</code> &rarr; 
    <code>a horizontal line pattern fill, spaced by 150 layer units</code></li>
    </ul>
    </div>
    """
    global cache
    key1 = f"{angle}_{x_interval}"
    key2 = feature.id()
    if layer_id in cache:
        if not key1 in cache[layer_id]:
            cache[layer_id] = {key1: {}}
    else:
        cache[layer_id] = {key1: {}}

    if key2 not in cache[layer_id][key1]:
        pts_hash_box = []
        bbox = feature.geometry().boundingBox()
        xcoords = [bbox.xMinimum(), bbox.xMaximum()]
        ycoords = [bbox.yMinimum(), bbox.yMaximum()]

        xcoef = round(cos(radians(angle)), 5)
        ycoef = round(sin(radians(angle)), 5)
        fn = lambda x, y: ycoef * y - xcoef * x

        ccorners = [fn(x, y) for x in xcoords for y in ycoords]
        cmin = min(ccorners)
        cmax = max(ccorners)
        cvar = cmax - cmin

        for i in range(0, ceil(cvar), x_interval):
            c = cmin + i
            if xcoef == 0:
                xs = []
            else:
                xs = [((ycoef * y - c) / xcoef, y) for y in ycoords]

            if ycoef == 0:
                ys = []
            else:
                ys = [(x, (xcoef * x + c) / ycoef) for x in xcoords]

            pts = [
                (x, y)
                for x, y in xs + ys
                if xcoords[0] <= x <= xcoords[1] and ycoords[0] <= y <= ycoords[1]
            ]
            pt_hash_box = [QgsPointXY(x, y) for x, y in sorted(pts)]
            pts_hash_box.append(pt_hash_box)

        multipolylines = QgsGeometry.fromMultiPolylineXY(pts_hash_box)
        geom = multipolylines.intersection(feature.geometry())
        if geom.type() == QgsWkbTypes.LineGeometry:
            cache[layer_id][key1][key2] = geom
        elif geom.wkbType() == QgsWkbTypes.GeometryCollection:
            if geom.convertGeometryCollectionToSubclass(QgsWkbTypes.LineGeometry):
                cache[layer_id][key1][key2] = geom
            else:
                cache[layer_id][key1][key2] = QgsGeometry()
        else:
            cache[layer_id][key1][key2] = QgsGeometry()

    return cache[layer_id][key1][key2]

EDIT

  • Code improved : QgsGeometry.intersection() produce sometimes a GeometryCollection with lines and points. If the result of intersection is a collection, the code try to convert it to Line geometry type.
  • Angle choice added : the angle argument is now in the code. The angle is in degrees. Rewritten code for this.
  • Change wkbType to type : LineGeometry is more generic and works with MultiLineStringZ and all other Line types.
  • Formatting help : adapt function help to QGIS expression function help format
8
  • I don't know why, but it doesn't work for some features in layer. I apply to different layers, but the result is the same. – Kadir Şahbaz Mar 20 '20 at 11:56
  • Is the geometry of the feature is valid ? Can you give a geometry as an example and test ? – J. Monticolo Mar 20 '20 at 12:15
  • Sample data – Kadir Şahbaz Mar 20 '20 at 12:58
  • When I try this the "Output previous" says "<geometry; Point>". And when I apply the styling (as a point geometry type) it just shows the first point on the first feature. EDIT: Nevermind! That was with a scratch layer, working fine with Spatialite so far – TeddyTedTed Mar 20 '20 at 13:11
  • 1
    @TeddyTedTed: angle choice done. You can play with all angles as you want ! – J. Monticolo Mar 20 '20 at 20:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.