Python QgsGeometry.buffer方法代码示例

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
def buffering(progress, writer, distance, field, useField, layer, dissolve,
              segments):

    if useField:
        field = layer.fieldNameIndex(field)

    outFeat = QgsFeature()
    inFeat = QgsFeature()
    inGeom = QgsGeometry()
    outGeom = QgsGeometry()

    current = 0
    features = vector.features(layer)
    total = 100.0 / float(len(features))

    # With dissolve
    if dissolve:
        first = True
        for inFeat in features:
            attrs = inFeat.attributes()
            if useField:
                value = attrs[field]
            else:
                value = distance

            inGeom = QgsGeometry(inFeat.geometry())
            outGeom = inGeom.buffer(float(value), segments)
            if first:
                tempGeom = QgsGeometry(outGeom)
                first = False
            else:
                tempGeom = tempGeom.combine(outGeom)

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

        outFeat.setGeometry(tempGeom)
        outFeat.setAttributes(attrs)
        writer.addFeature(outFeat)
    else:
        # Without dissolve
        for inFeat in features:
            attrs = inFeat.attributes()
            if useField:
                value = attrs[field]
            else:
                value = distance
            inGeom = QgsGeometry(inFeat.geometry())
            outGeom = inGeom.buffer(float(value), segments)
            outFeat.setGeometry(outGeom)
            outFeat.setAttributes(attrs)
            writer.addFeature(outFeat)
            current += 1
            progress.setPercentage(int(current * total))

    del writer

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def mask_geometry(self):
        if not self.parameters.geometry:
            geom = QgsGeometry()
            return geom, QgsRectangle()

        geom = QgsGeometry(self.parameters.geometry)  # COPY !!

        if self.parameters.do_simplify:
            if hasattr(self.canvas, 'mapSettings'):
                tol = self.parameters.simplify_tolerance * \
                    self.canvas.mapSettings().mapUnitsPerPixel()
            else:
                tol = self.parameters.simplify_tolerance * \
                    self.canvas.mapRenderer().mapUnitsPerPixel()

            if tol in list(self.simplified_geometries.keys()):
                geom, bbox = self.simplified_geometries[tol]
            else:
                if self.has_simplifier:
                    simplifier = QgsMapToPixelSimplifier(QgsMapToPixelSimplifier.SimplifyGeometry,
                                                         tol)
                    geom = simplifier.simplify(geom)
                    if not geom.isGeosValid():
                        # make valid
                        geom = geom.buffer(0.0, 1)
                bbox = geom.boundingBox()
                self.simplified_geometries[tol] = (QgsGeometry(geom), QgsRectangle(bbox))
        else:
            bbox = geom.boundingBox()

        return geom, bbox

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def compute_mask_geometries(self, parameters, poly):
        geom = None
        for g in poly:
            if geom is None:
                geom = QgsGeometry(g)
            else:
                # do an union here
                geom = geom.combine(g)

        if parameters.do_buffer:
            geom = geom.buffer(parameters.buffer_units, parameters.buffer_segments)

        # reset the simplified geometries dict
        self.simplified_geometries = {}

        return geom

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
 def vectorRandom(self, n, layer, xmin, xmax, ymin, ymax):
     provider = layer.dataProvider()
     index = ftools_utils.createIndex(provider)
     seed()
     points = []
     feat = QgsFeature()
     i = 1
     count = 40.00
     add = ( 70.00 - 40.00 ) / n
     while i <= n:
         point = QgsPoint(xmin + (xmax-xmin) * random(), ymin + (ymax-ymin) * random())
         pGeom = QgsGeometry().fromPoint(point)
         ids = index.intersects(pGeom.buffer(5,5).boundingBox())
         for id in ids:
             provider.getFeatures( QgsFeatureRequest().setFilterFid( int(id) ) ).nextFeature( feat )
             tGeom = QgsGeometry(feat.geometry())
             if pGeom.intersects(tGeom):
                 points.append(pGeom)
                 i = i + 1
                 count = count + add
                 self.progressBar.setValue(count)
                 break
     return points

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def processAlgorithm(self, progress):
        target = dataobjects.getObjectFromUri(
            self.getParameterValue(self.TARGET))
        join = dataobjects.getObjectFromUri(
            self.getParameterValue(self.JOIN))
        predicates = self.getParameterValue(self.PREDICATE)
        precision = self.getParameterValue(self.PRECISION)

        summary = self.getParameterValue(self.SUMMARY) == 1
        keep = self.getParameterValue(self.KEEP) == 1

        sumList = self.getParameterValue(self.STATS).lower().split(',')

        targetFields = target.fields()
        joinFields = join.fields()

        fieldList = QgsFields()

        if not summary:
            joinFields = vector.testForUniqueness(targetFields, joinFields)
            seq = list(range(len(targetFields) + len(joinFields)))
            targetFields.extend(joinFields)
            targetFields = dict(list(zip(seq, targetFields)))
        else:
            numFields = {}
            for j in range(len(joinFields)):
                if joinFields[j].type() in [QVariant.Int, QVariant.Double, QVariant.LongLong, QVariant.UInt, QVariant.ULongLong]:
                    numFields[j] = []
                    for i in sumList:
                        field = QgsField(i + str(joinFields[j].name()), QVariant.Double, '', 24, 16)
                        fieldList.append(field)
            field = QgsField('count', QVariant.Double, '', 24, 16)
            fieldList.append(field)
            joinFields = vector.testForUniqueness(targetFields, fieldList)
            targetFields.extend(fieldList)
            seq = list(range(len(targetFields)))
            targetFields = dict(list(zip(seq, targetFields)))

        fields = QgsFields()
        for f in list(targetFields.values()):
            fields.append(f)

        writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(
            fields, target.wkbType(), target.crs())

        outFeat = QgsFeature()
        inFeatB = QgsFeature()
        inGeom = QgsGeometry()

        index = vector.spatialindex(join)

        mapP2 = dict()
        features = vector.features(join)
        for f in features:
            mapP2[f.id()] = QgsFeature(f)

        features = vector.features(target)
        total = 100.0 / len(features)
        for c, f in enumerate(features):
            atMap1 = f.attributes()
            outFeat.setGeometry(f.geometry())
            inGeom = vector.snapToPrecision(f.geometry(), precision)
            none = True
            joinList = []
            if inGeom.type() == QgsWkbTypes.PointGeometry:
                bbox = inGeom.buffer(10, 2).boundingBox()
            else:
                bbox = inGeom.boundingBox()
            bufferedBox = vector.bufferedBoundingBox(bbox, 0.51 * precision)
            joinList = index.intersects(bufferedBox)
            if len(joinList) > 0:
                count = 0
                for i in joinList:
                    inFeatB = mapP2[i]
                    inGeomB = vector.snapToPrecision(inFeatB.geometry(), precision)

                    res = False
                    for predicate in predicates:
                        res = getattr(inGeom, predicate)(inGeomB)
                        if res:
                            break

                    if res:
                        count = count + 1
                        none = False
                        atMap2 = inFeatB.attributes()
                        if not summary:
                            atMap = atMap1
                            atMap2 = atMap2
                            atMap.extend(atMap2)
                            atMap = dict(list(zip(seq, atMap)))
                            break
                        else:
                            for j in list(numFields.keys()):
                                numFields[j].append(atMap2[j])

                if summary and not none:
                    atMap = atMap1
                    for j in list(numFields.keys()):
                        for k in sumList:
#.........这里部分代码省略.........

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def processAlgorithm(self, progress):
        target = dataobjects.getObjectFromUri(
            self.getParameterValue(self.TARGET))
        join = dataobjects.getObjectFromUri(
            self.getParameterValue(self.JOIN))
        predicates = self.getParameterValue(self.PREDICATE)

        summary = self.getParameterValue(self.SUMMARY) == 1
        keep = self.getParameterValue(self.KEEP) == 1

        sumList = self.getParameterValue(self.STATS).lower().split(',')

        targetProvider = target.dataProvider()
        joinProvider = join.dataProvider()

        targetFields = targetProvider.fields()
        joinFields = joinProvider.fields()

        fieldList = QgsFields()

        if not summary:
            joinFields = vector.testForUniqueness(targetFields, joinFields)
            seq = range(0, len(targetFields) + len(joinFields))
            targetFields.extend(joinFields)
            targetFields = dict(zip(seq, targetFields))
        else:
            numFields = {}
            for j in xrange(len(joinFields)):
                if joinFields[j].type() in [QVariant.Int, QVariant.Double]:
                    numFields[j] = []
                    for i in sumList:
                        field = QgsField(i + unicode(joinFields[j].name()), QVariant.Double, '', 24, 16)
                        fieldList.append(field)
            field = QgsField('count', QVariant.Double, '', 24, 16)
            fieldList.append(field)
            joinFields = vector.testForUniqueness(targetFields, fieldList)
            targetFields.extend(fieldList)
            seq = range(0, len(targetFields))
            targetFields = dict(zip(seq, targetFields))

        fields = QgsFields()
        for f in targetFields.values():
            fields.append(f)

        writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(
            fields, targetProvider.geometryType(), targetProvider.crs())

        outFeat = QgsFeature()
        inFeatB = QgsFeature()
        inGeom = QgsGeometry()

        index = vector.spatialindex(join)

        mapP2 = dict()
        features = vector.features(join)
        for f in features:
            mapP2[f.id()] = QgsFeature(f)

        features = vector.features(target)
        total = 100.0 / len(features)
        for c, f in enumerate(features):
            inGeom = f.geometry()
            atMap1 = f.attributes()
            outFeat.setGeometry(inGeom)
            none = True
            joinList = []
            if inGeom.type() == QGis.Point:
                joinList = index.intersects(inGeom.buffer(10, 2).boundingBox())
                if len(joinList) > 0:
                    check = 0
                else:
                    check = 1
            else:
                joinList = index.intersects(inGeom.boundingBox())
                if len(joinList) > 0:
                    check = 0
                else:
                    check = 1

            if check == 0:
                count = 0
                for i in joinList:
                    inFeatB = mapP2[i]
                    inGeomB = inFeatB.geometry()

                    res = False
                    for predicate in predicates:
                        if predicate == 'intersects':
                            res = inGeom.intersects(inGeomB)
                        elif predicate == 'contains':
                            res = inGeom.contains(inGeomB)
                        elif predicate == 'equals':
                            res = inGeom.equals(inGeomB)
                        elif predicate == 'touches':
                            res = inGeom.touches(inGeomB)
                        elif predicate == 'overlaps':
                            res = inGeom.overlaps(inGeomB)
                        elif predicate == 'within':
                            res = inGeom.within(inGeomB)
                        elif predicate == 'crosses':
#.........这里部分代码省略.........

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
def multi_buffering(layer, radii, callback=None):
    """Buffer a vector layer using many buffers (for volcanoes or rivers).

    This processing algorithm will keep the original attribute table and
    will add a new one for the hazard class name according to
    safe.definitions.fields.hazard_value_field.

    radii = OrderedDict()
    radii[500] = 'high'
    radii[1000] = 'medium'
    radii[2000] = 'low'

    Issue https://github.com/inasafe/inasafe/issues/3185

    :param layer: The layer to polygonize.
    :type layer: QgsVectorLayer

    :param radii: A dictionary of radius.
    :type radii: OrderedDict

    :param callback: A function to all to indicate progress. The function
        should accept params 'current' (int), 'maximum' (int) and 'step' (str).
        Defaults to None.
    :type callback: function

    :return: The buffered vector layer.
    :rtype: QgsVectorLayer
    """
    # Layer output
    output_layer_name = buffer_steps['output_layer_name']
    processing_step = buffer_steps['step_name']

    input_crs = layer.crs()
    feature_count = layer.featureCount()

    fields = layer.fields()
    # Set the new hazard class field.
    new_field = create_field_from_definition(hazard_class_field)
    fields.append(new_field)
    # Set the new buffer distances field.
    new_field = create_field_from_definition(buffer_distance_field)
    fields.append(new_field)

    buffered = create_memory_layer(
        output_layer_name, QGis.Polygon, input_crs, fields)
    data_provider = buffered.dataProvider()

    # Reproject features if needed into UTM if the layer is in 4326.
    if layer.crs().authid() == 'EPSG:4326':
        center = layer.extent().center()
        utm = QgsCoordinateReferenceSystem(
            get_utm_epsg(center.x(), center.y(), input_crs))
        transform = QgsCoordinateTransform(layer.crs(), utm)
        reverse_transform = QgsCoordinateTransform(utm, layer.crs())
    else:
        transform = None
        reverse_transform = None

    for i, feature in enumerate(layer.getFeatures()):
        geom = QgsGeometry(feature.geometry())

        if transform:
            geom.transform(transform)

        inner_ring = None

        for radius in radii:
            attributes = feature.attributes()

            # We add the hazard value name to the attribute table.
            attributes.append(radii[radius])
            # We add the value of buffer distance to the attribute table.
            attributes.append(radius)

            circle = geom.buffer(radius, 30)

            if inner_ring:
                circle.addRing(inner_ring)

            inner_ring = circle.asPolygon()[0]

            new_feature = QgsFeature()
            if reverse_transform:
                circle.transform(reverse_transform)

            new_feature.setGeometry(circle)
            new_feature.setAttributes(attributes)

            data_provider.addFeatures([new_feature])

        if callback:
            callback(current=i, maximum=feature_count, step=processing_step)

    # We transfer keywords to the output.
    buffered.keywords = layer.keywords
    buffered.keywords['layer_geometry'] = 'polygon'
    buffered.keywords['layer_purpose'] = layer_purpose_hazard['key']
    buffered.keywords['inasafe_fields'][hazard_class_field['key']] = (
        hazard_class_field['field_name'])

#.........这里部分代码省略.........

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def compute(self, inName, joinName, outName, summary, sumList, keep, progressBar):
        layer1 = ftools_utils.getVectorLayerByName(inName)
        provider1 = layer1.dataProvider()
        fieldList1 = ftools_utils.getFieldList(layer1)

        layer2 = ftools_utils.getVectorLayerByName(joinName)
        provider2 = layer2.dataProvider()

        fieldList2 = ftools_utils.getFieldList(layer2)
        fieldList = QgsFields()
        if provider1.crs() != provider2.crs():
            QMessageBox.warning(self, self.tr("CRS warning!"), self.tr("Warning: Input layers have non-matching CRS.\nThis may cause unexpected results."))
        if not summary:
            fieldList2 = ftools_utils.testForUniqueness(fieldList1, fieldList2)
            seq = range(0, len(fieldList1) + len(fieldList2))
            fieldList1.extend(fieldList2)
            fieldList1 = dict(zip(seq, fieldList1))
        else:
            numFields = {}
            for j in xrange(len(fieldList2)):
                if fieldList2[j].type() == QVariant.Int or fieldList2[j].type() == QVariant.Double:
                    numFields[j] = []
                    for i in sumList:
                        field = QgsField(i + unicode(fieldList2[j].name()), QVariant.Double, "real", 24, 16, self.tr("Summary field"))
                        fieldList.append(field)
            field = QgsField("COUNT", QVariant.Double, "real", 24, 16, self.tr("Summary field"))
            fieldList.append(field)
            fieldList2 = ftools_utils.testForUniqueness(fieldList1, fieldList)
            fieldList1.extend(fieldList)
            seq = range(0, len(fieldList1))
            fieldList1 = dict(zip(seq, fieldList1))

        sRs = provider1.crs()
        progressBar.setValue(13)
        check = QFile(self.shapefileName)
        if check.exists():
            if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
                QMessageBox.warning(
                    self, self.tr('Error deleting shapefile'),
                    self.tr("Can't delete existing shapefile\n%s") % (self.shapefileName))
                return False
        fields = QgsFields()
        for f in fieldList1.values():
            fields.append(f)
        writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fields, provider1.geometryType(), sRs)
        #writer = QgsVectorFileWriter(outName, "UTF-8", fieldList1, provider1.geometryType(), sRs)
        inFeat = QgsFeature()
        outFeat = QgsFeature()
        inFeatB = QgsFeature()
        inGeom = QgsGeometry()
        progressBar.setValue(15)
        start = 15.00
        add = 85.00 / provider1.featureCount()

        index = ftools_utils.createIndex(provider2)

        # cache all features from provider2 to avoid huge number of feature requests in the inner loop
        mapP2 = {}
        for f in provider2.getFeatures():
            mapP2[f.id()] = QgsFeature(f)

        fit1 = provider1.getFeatures()
        while fit1.nextFeature(inFeat):
            inGeom = inFeat.geometry()
            atMap1 = inFeat.attributes()
            outFeat.setGeometry(inGeom)
            none = True
            joinList = []
            if inGeom.type() == QGis.Point:
                #(check, joinList) = layer2.featuresInRectangle(inGeom.buffer(10,2).boundingBox(), True, True)
                #layer2.select(inGeom.buffer(10,2).boundingBox(), False)
                #joinList = layer2.selectedFeatures()
                joinList = index.intersects(inGeom.buffer(10, 2).boundingBox())
                if len(joinList) > 0:
                    check = 0
                else:
                    check = 1
            else:
                #(check, joinList) = layer2.featuresInRectangle(inGeom.boundingBox(), True, True)
                #layer2.select(inGeom.boundingBox(), False)
                #joinList = layer2.selectedFeatures()
                joinList = index.intersects(inGeom.boundingBox())
                if len(joinList) > 0:
                    check = 0
                else:
                    check = 1
            if check == 0:
                count = 0
                for i in joinList:
                    inFeatB = mapP2[i]  # cached feature from provider2
                    if inGeom.intersects(inFeatB.geometry()):
                        count = count + 1
                        none = False
                        atMap2 = inFeatB.attributes()
                        if not summary:
                            atMap = atMap1
                            atMap2 = atMap2
                            atMap.extend(atMap2)
                            atMap = dict(zip(seq, atMap))
                            break
#.........这里部分代码省略.........

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
def buffering(progress, writer, distance, field, useField, layer, dissolve,
              segments):

    if useField:
        field = layer.fieldNameIndex(field)

    outFeat = QgsFeature()
    inFeat = QgsFeature()
    inGeom = QgsGeometry()
    outGeom = QgsGeometry()

    current = 0
    features = vector.features(layer)
    total = 100.0 / float(len(features))

    # With dissolve
    if dissolve:
        first = True
        for inFeat in features:
            attrs = inFeat.attributes()
            if useField:
                value = attrs[field]
            else:
                value = distance

            inGeom = QgsGeometry(inFeat.geometry())
            if inGeom.isGeosEmpty() or not inGeom.isGeosValid():
                ProcessingLog.addToLog(ProcessingLog.LOG_WARNING, 'Feature {} has empty or invalid geometry. Skipping...'.format(inFeat.id()))
                continue
            outGeom = inGeom.buffer(float(value), segments)
            if first:
                tempGeom = QgsGeometry(outGeom)
                first = False
            else:
                tempGeom = tempGeom.combine(outGeom)

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

        outFeat.setGeometry(tempGeom)
        outFeat.setAttributes(attrs)
        writer.addFeature(outFeat)
    else:
        # Without dissolve
        for inFeat in features:
            attrs = inFeat.attributes()
            if useField:
                value = attrs[field]
            else:
                value = distance
            inGeom = QgsGeometry(inFeat.geometry())
            if inGeom.isGeosEmpty() or not inGeom.isGeosValid():
                ProcessingLog.addToLog(ProcessingLog.LOG_WARNING, 'Feature {} has empty or invalid geometry. Skipping...'.format(inFeat.id()))
                continue

            outGeom = inGeom.buffer(float(value), segments)
            outFeat.setGeometry(outGeom)
            outFeat.setAttributes(attrs)
            writer.addFeature(outFeat)
            current += 1
            progress.setPercentage(int(current * total))

    del writer

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def in_mask(self, feature, srid=None):
        if feature is None:  # expression overview
            return False

        if self.layer is None:
            return False

        try:
            # layer is not None but destroyed ?
            self.layer.id()
        except:
            self.reset_mask_layer()
            return False

        # mask layer empty due to unloaded memlayersaver plugin > no filtering
        if self.layer.featureCount() == 0:
            return True

        mask_geom, bbox = self.mask_geometry()
        geom = QgsGeometry(feature.geometry())
        if not geom.isGeosValid():
            geom = geom.buffer(0.0, 1)

        if geom is None:
            return False

        if srid is not None and self.layer.crs().postgisSrid() != srid:
            src_crs = QgsCoordinateReferenceSystem(srid)
            dest_crs = self.layer.crs()
            xform = QgsCoordinateTransform(src_crs, dest_crs, QgsProject.instance())
            try:
                geom.transform(xform)
            except:
                # transformation error. Check layer projection.
                pass

        if geom.type() == QgsWkbTypes.PolygonGeometry:
            if self.parameters.polygon_mask_method == 2 and not self.has_point_on_surface:
                self.parameters.polygon_mask_method = 1

            if self.parameters.polygon_mask_method == 0:
                # this method can only work when no geometry simplification is involved
                return (mask_geom.overlaps(geom) or mask_geom.contains(geom))
            elif self.parameters.polygon_mask_method == 1:
                # the fastest method, but with possible inaccuracies
                pt = geom.vertexAt(0)
                return bbox.contains(QgsPointXY(pt)) and mask_geom.contains(geom.centroid())
            elif self.parameters.polygon_mask_method == 2:
                # will always work
                pt = geom.vertexAt(0)
                return bbox.contains(QgsPointXY(pt)) and mask_geom.contains(geom.pointOnSurface())
            else:
                return False
        elif geom.type() == QgsWkbTypes.LineGeometry:
            if self.parameters.line_mask_method == 0:
                return mask_geom.intersects(geom)
            elif self.parameters.line_mask_method == 1:
                return mask_geom.contains(geom)
            else:
                return False
        elif geom.type() == QgsWkbTypes.PointGeometry:
            return mask_geom.intersects(geom)
        else:
            return False

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import buffer [as 别名]
    def traf(self,geometry,trafic,spatial,lines,echelle,angle_max,dist_min,double_sens):
        lignes=lines
        conn = db.connect(':memory:')
        conn.enable_load_extension(True)
        conn.execute("select load_extension('mod_spatialite')")
        c = conn.cursor()
        proj=str(self.lines.crs().postgisSrid())
        t=geometry[trafic]
        
        texte="select astext(st_buffer(st_geomfromtext('"+geometry.geometry().asWkt()+"',"+proj+"),"+str(echelle*t)+"))"
        rs = c.execute(texte)
        resultat=c.fetchall()
        conn.commit()
        texte_buf= resultat[0][0]
        texte_buf=texte_buf.replace("Polygon","MultiLineString")
        buf=QgsGeometry.fromWkt(texte_buf)
        #buf=buf.convertToType(QgsWkbTypes.LineGeometry,False)
        texte="select astext(st_union(st_exteriorring(st_geomfromtext('"+buf.asWkt()+"',"+proj+"))))"
        rs = c.execute(texte)
        resultat=c.fetchall()
        conn.commit()
        buf= QgsGeometry.fromWkt(resultat[0][0])
        #QgsMessageLog.logMessage(texte)
        #QgsMessageLog.logMessage(buf.asWkt())
        buf_poly=buf
        geom=geometry.geometry()
        buf_sav=buf
        buf_l=buf.length()
        l2=QgsLineString()
        #print(geom.asWkt())
        geom.convertToSingleType()
        #print(geom.asWkt())
        l2.fromWkt(geom.asWkt())
        #print(l2.asWkt())
        pt1=QgsGeometry(l2.startPoint())
        pt2=QgsGeometry(l2.endPoint())
        ###point debut
        debut=spatial.intersects(QgsGeometry.buffer(pt1,10,3).boundingBox())
        fe= [f for f in debut]

        feats=[]
        for f in fe:
            ff=QgsLineString()
            #QgsMessageLog.logMessage(lignes[f].geometry().asWkt())
            geom=lignes[f].geometry()
            geom.convertToSingleType()
            ff.fromWkt(geom.asWkt())
            #QgsMessageLog.logMessage(ff.asWkt())
            #print(pt1.distance(QgsGeometry(ff.startPoint())))
            if pt1.distance(QgsGeometry(ff.startPoint()))<10:
                if lignes[f] not in feats:
                    feats.append(f)
            elif pt1.distance(QgsGeometry(ff.endPoint()))<10:
                if lignes[f] not in feats:
                    feats.append(f)
        #QgsMessageLog.logMessage(str(fe))
        #QgsMessageLog.logMessage(str(feats))

                
            
        distances={}
        angles={}
        for i in feats:
            longueur=lignes[i].geometry().length()
            if not(geometry.id()==i):
                distances[i]=lignes[i].geometry().lineLocatePoint(pt1)
                if distances[i]<dist_min:
                    angles[i]=((lignes[i].geometry().interpolateAngle(min(dist_min,longueur))*180/math.pi)+180)%360
                else:
                    angles[i]=lignes[i].geometry().interpolateAngle(longueur-min(dist_min,longueur))*180/math.pi
            else:
                angle1=lignes[i].geometry().interpolateAngle(min(dist_min,longueur))*180/math.pi
        angle_maxi=1e38
        voisin=None
        angle_voisin=None
        angle2=None
        if len(distances)==0:
            angle=(angle1)%360
            angle2=angle1
            angle_voisin=angle2
        for i in distances:
            if distances[i]<dist_min:
                angle=(angles[i])%360
                min_angle=min(abs((angle+180)%360-(angle1+180)%360),abs(angle-angle1))
                if min_angle<angle_maxi:
                    angle_maxi=min_angle
                    angle_voisin=angle
                    voisin=i
            else:
                angle=angles[i]
                min_angle=min(abs((angle+180)%360-(angle1+180)%360),abs(angle-angle1))
                if min_angle<angle_maxi:
                    angle_maxi=min_angle
                    angle_voisin=angle
                    voisin=i

        if  min(abs((angle_voisin+180)%360-(angle1+180)%360),abs(angle_voisin-angle1))<angle_max:
            if abs((angle_voisin+180)%360-(angle1+180)%360)<abs(angle_voisin-angle1):
                angle2=(0.5*(((angle_voisin+180)%360+(angle1+180)%360))+180)%360
            else:
#.........这里部分代码省略.........