Selecting and clipping layers within user defined area using PyQGIS

Question

I am new to QGIS and am trying to write a script that will allow the user to define a region on the map (either polygon or freehand) and clip all displayed map layers within this region to the user-defined selection area so that line lengths can be calculated and totaled within the region for each layer.

Using the example below, suppose there are two vector layers of polylines. The long vertical line is on one layer, and the pink lines are on the 2nd layer.

Instead of selecting the entire intersecting lines from endpoint to endpoint, I need to select/clip only the part of the lines that are within the polygon boundary.

The output I am looking for would be something like:

Layer Name 1: Total Lines Length: 55
Layer Name 2: Total Lines Length: 145
Layer Name 3: etc.

I would store this output in a table or csv file.

I am not sure how to get the user defined polygon region as input, how to get the displayed layers that fall within that boundary, and how to clip their extents to the region so that I can get their clipped lengths.

d = QgsDistanceArea()
d.setEllipsoid('WGS84')

#this is just selecting the active layer, 
#but I need to get layers from boundary

layer = qgis.utils.iface.activeLayer()

"""
what I need is something like this <pseudocode>:

selection_layers = interface.CurrentSelectionPolygon()
for layer in selection_layers:
    new_temp_layer = clip_to_selection(layer)
    for feature in new_temp_layer.selectedFeatures():
       line_length = d.measureLength(feature.geometry())
       #add line length to total for layer etc.

"""



#getting the user selected features from active layer just to test
features = layer.selectedFeatures()


#
lengths = []
if len(features) != 0:
    for f in features:
        geom = f.geometry()
        #use measureLength to get accurate projected values
        length = d.measureLength(geom)
        #convert this to feet
        c_length = d.convertLengthMeasurement(length, 2)             
        #0 is meters, 2 is feet           
        lengths.append(c_length)

# this is outputting the sum of the single layer's line lengths        
total_length = sum(lengths)
print(f"total length: {total_length}")

I cannot figure out how to get the user's map drawn selection polygon to use as a clipping mask.

Answer 1

For this purpose, I think you will want to use the QgsRubberBand class inside a custom map tool which will allow the user to digitize a polygon on the canvas. You can use or adapt my script below which should do more or less what you want. To test this script, paste it into a new editor in the Python console and click Run. Then use the left mouse button to digitize a polygon and right-click to finish. This will launch a message box with the sum of line feature intersections with the rubber band polygon for all line layers in the current project which you can choose to save as a simple text file.

class MapToolSumLineLengths(QgsMapToolEmitPoint):
    def __init__(self, canvas):
        self.canvas = canvas
        QgsMapToolEmitPoint.__init__(self, self.canvas)
        self.rubber_band_line = None
        self.rubber_band_fixed = None
        self.rubber_band_feature = None
        self.points = []
        self.result_text = ''
        self.msg = QMessageBox()
        self.msg.btn_save_txt = self.msg.addButton('Save to text file', QMessageBox.AcceptRole)
        self.msg.btn_save_txt.clicked.connect(lambda: self.save_txt(self.result_text))
        self.msg.finished.connect(self.remove_feature)
        
        self.deactivated.connect(self.clean_up)
        
        self.project = QgsProject.instance()
        self.da = QgsDistanceArea()
        self.da.setEllipsoid('WGS84')
        
    def canvasPressEvent(self, e):
        if self.result_text:
            self.result_text = ''
        set_point = self.toMapCoordinates(e.pos())
        if self.rubber_band_fixed is not None:
            self.rubber_band_fixed.reset()
        if e.button() == Qt.LeftButton:
            self.points.append(set_point)
            if len(self.points) == 2:
                self.rubber_band_fixed = QgsRubberBand(self.canvas, QgsWkbTypes.LineGeometry)
                self.rubber_band_fixed.setColor(QColor('Red'))
                geom = QgsGeometry().fromPolylineXY(self.points)
                self.rubber_band_fixed.setToGeometry(geom)
                self.rubber_band_fixed.show()
            elif len(self.points) > 2:
                self.rubber_band_fixed = QgsRubberBand(self.canvas, QgsWkbTypes.PolygonGeometry)
                self.rubber_band_fixed.setColor(QColor('Red'))
                self.rubber_band_fixed.setFillColor(QColor(210, 153, 90, 125))
                pnt_ring = [p for p in self.points]
                pnt_ring.append(pnt_ring[0])
                geom = QgsGeometry().fromPolygonXY([pnt_ring])
                self.rubber_band_fixed.setToGeometry(geom)
                self.rubber_band_fixed.show()
        elif e.button() == Qt.RightButton:
            self.points.append(set_point)
            if len(self.points) == 2:
                self.clean_up()
            elif len(self.points) > 2:
                self.rubber_band_feature = QgsRubberBand(self.canvas, QgsWkbTypes.PolygonGeometry)
                self.rubber_band_feature.setColor(QColor('Black'))
                self.rubber_band_feature.setFillColor(QColor(65, 220, 245, 125))
                pnt_ring = [p for p in self.points]
                pnt_ring.append(pnt_ring[0])
                temp_geom = QgsGeometry().fromPolygonXY([pnt_ring])
                self.rubber_band_feature.setToGeometry(temp_geom)
                self.rubber_band_feature.show()
                self.sum_lines(temp_geom)
                self.clean_up()
                
    def remove_feature(self, result):
        self.rubber_band_feature.reset()
        self.rubber_band_feature = None
        
    def canvasMoveEvent(self, e):
        point = self.toMapCoordinates(e.pos())
        if self.rubber_band_line is not None:
            self.rubber_band_line.reset()
        if self.points:
            self.rubber_band_line = QgsRubberBand(self.canvas, QgsWkbTypes.LineGeometry)
            self.rubber_band_line.setBrushStyle(Qt.Dense7Pattern)
            self.rubber_band_line.setColor(QColor('Green'))
            if len(self.points) == 1:
                geom = QgsGeometry().fromPolylineXY([self.points[0], point])
                self.rubber_band_line.setToGeometry(geom)
                self.rubber_band_line.show()
            elif len(self.points) > 1:
                geom = QgsGeometry().fromPolylineXY([self.points[0], point, self.points[-1]])
                self.rubber_band_line.setToGeometry(geom)
                self.rubber_band_line.show()
    
    def clean_up(self):
        if self.rubber_band_line is not None:
            self.rubber_band_line.reset()
            self.rubber_band_line = None
        if self.rubber_band_fixed is not None:
            self.rubber_band_fixed.reset()
            self.rubber_band_fixed = None
        self.points.clear()
        
    def sum_lines(self, temp_polygon):
        all_line_layers = [l for l in self.project.mapLayers().values() if l.type() == QgsMapLayer.VectorLayer and l.geometryType() == 1]
        for layer in all_line_layers:
            total = 0
            features = [f for f in layer.getFeatures() if f.geometry().intersects(temp_polygon)]
            for f in features:
                intersection = self.da.measureLength(f.geometry().intersection(temp_polygon))
                metres = self.da.convertLengthMeasurement(intersection, QgsUnitTypes.DistanceFeet)
                total += metres
            if total:
                self.result_text += f'{layer.name()}: Total line lengths: {round(total, 3)} feet\n'
        self.msg.setText(self.result_text)
        self.msg.show()
        
    def save_txt(self, results):
        file_name = QFileDialog.getSaveFileName(None, 'Save Results', filter='*.txt')
        if file_name:
            f = open(file_name[0], 'w')
            f.write(results)
            f.close()
            
    
canvas = iface.mapCanvas()
T = MapToolSumLineLengths(canvas)
canvas.setMapTool(T)

See gif below for example of use: