Python QgsGeometry.distance方法代码示例

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import distance [as 别名]
 def hdist(self, feata, lDist):
     geoma = QgsGeometry(feata.geometry()) 
     feat = QgsFeature()
     provider = self.ml.dataProvider()
     feats = provider.getFeatures()
     self.emit(SIGNAL("runStatus(PyQt_PyObject)"), 0)
     self.emit(SIGNAL("runRange(PyQt_PyObject)"), (0, provider.featureCount())) 
     ne = 0              
     neighbours = ""
     while feats.nextFeature(feat):
         ne += 1
         self.emit(SIGNAL("runStatus(PyQt_PyObject)"), ne)                       
         geomb = QgsGeometry(feat.geometry())            
         if feata.id()!=feat.id():
             if geoma.distance(geomb)<=lDist:
                 neighbours = neighbours + '%s,' % (feat.id()+self.p)
     return neighbours[:-1]            

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import distance [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:
#.........这里部分代码省略.........

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import distance [as 别名]
def setSelectFeatures(canvas, selectGeometry, doContains, doDifference, singleSelect=None):
    """
    QgsMapCanvas* canvas,
    QgsGeometry* selectGeometry,
    bool doContains,
    bool doDifference,
    bool singleSelect 
    """
    if selectGeometry.type() != QGis.Polygon:
        return

    vlayer = getCurrentVectorLayer(canvas)

    if vlayer == None:
        return

    #toLayerCoordinates will throw an exception for any 'invalid' points in
    #the rubber band.
    #For example, if you project a world map onto a globe using EPSG 2163
    #and then click somewhere off the globe, an exception will be thrown.
    selectGeomTrans = QgsGeometry(selectGeometry)

    if canvas.mapSettings().hasCrsTransformEnabled():
        try:
            ct = QgsCoordinateTransform(canvas.mapSettings().destinationCrs(), vlayer.crs())
            selectGeomTrans.transform( ct )
        except QgsCsException as cse:
            Q_UNUSED(cse)
            #catch exception for 'invalid' point and leave existing selection unchanged
            """
            QgsLogger::warning( "Caught CRS exception " + QString( __FILE__ ) + ": " + QString::number( __LINE__ ) );
            QgisApp::instance()->messageBar()->pushMessage(
            QObject::tr( "CRS Exception" ),
            QObject::tr( "Selection extends beyond layer's coordinate system" ),
            QgsMessageBar::WARNING,
            QgisApp::instance()->messageTimeout() );
            """
            return

    QApplication.setOverrideCursor(Qt.WaitCursor)
    """
    QgsDebugMsg( "Selection layer: " + vlayer->name() );
    QgsDebugMsg( "Selection polygon: " + selectGeomTrans.exportToWkt() );
    QgsDebugMsg( "doContains: " + QString( doContains ? "T" : "F" ) );
    QgsDebugMsg( "doDifference: " + QString( doDifference ? "T" : "F" ) );
    """

    context = QgsRenderContext().fromMapSettings(canvas.mapSettings())
    r = vlayer.rendererV2()

    if r:
        r.startRender(context, vlayer.pendingFields())

    request = QgsFeatureRequest()
    request.setFilterRect(selectGeomTrans.boundingBox())
    request.setFlags(QgsFeatureRequest.ExactIntersect)

    if r:
        request.setSubsetOfAttributes(r.usedAttributes(), vlayer.pendingFields())
    else:
        request.setSubsetOfAttributes(QgsAttributeList)

    fit = vlayer.getFeatures(request)

    newSelectedFeatures = [] #QgsFeatureIds
    f = QgsFeature() 
    closestFeatureId = 0 #QgsFeatureId 
    foundSingleFeature = False
    #double closestFeatureDist = std::numeric_limits<double>::max();
    closestFeatureDist = sys.float_info.max

    while fit.nextFeature(f):
        # make sure to only use features that are visible
        if r and not r.willRenderFeature( f ):
            continue;
        g = QgsGeometry(f.geometry())
        if doContains:
            if not selectGeomTrans.contains(g):
                continue
        else:
            if not selectGeomTrans.intersects(g):
                continue
        if singleSelect:
            foundSingleFeature = True
            distance = float(g.distance(selectGeomTrans))
            if ( distance <= closestFeatureDist ):
                closestFeatureDist = distance
                closestFeatureId = f.id()
        else:
            newSelectedFeatures.insert(0, f.id())

    if singleSelect and foundSingleFeature:
        newSelectedFeatures.insert(0, closestFeatureId)
    if r:
        r.stopRender(context)
    #QgsDebugMsg( "Number of new selected features: " + QString::number( newSelectedFeatures.size() ) 

    if doDifference:
        layerSelectedFeatures = vlayer.selectedFeaturesIds()

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

# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import distance [as 别名]
    def do_indexjoin(self, feat):
        '''Find the nearest neigbour using an index, if possible

        Parameter: feat -- The feature for which a neighbour is
                           sought
        '''
        infeature = feat
        infeatureid = infeature.id()
        inputgeom = QgsGeometry(infeature.geometry())
        # Shall approximate input geometries be used?
        if self.approximateinputgeom:
            # Use the centroid as the input geometry
            inputgeom = QgsGeometry(infeature.geometry()).centroid()
        # Check if the coordinate systems are equal, if not,
        # transform the input feature!
        if (self.inpvl.crs() != self.joinvl.crs()):
            try:
                inputgeom.transform(QgsCoordinateTransform(
                    self.inpvl.crs(), self.joinvl.crs()))
            except:
                import traceback
                self.error.emit(self.tr('CRS Transformation error!') +
                                ' - ' + traceback.format_exc())
                self.abort = True
                return
        nnfeature = None
        mindist = float("inf")
        ## Find the closest feature!
        if (self.approximateinputgeom or
                self.inpvl.wkbType() == QGis.WKBPoint or
                self.inpvl.wkbType() == QGis.WKBPoint25D):
            # The input layer's geometry type is point, or shall be
            # approximated to point (centroid).
            # Then a join index will always be used.
            if (self.usejoinlayerapprox or
                    self.joinvl.wkbType() == QGis.WKBPoint or
                    self.joinvl.wkbType() == QGis.WKBPoint25D):
                # The join layer's geometry type is point, or the
                # user wants approximate join geometries to be used.
                # Then the join index nearest neighbour function can
                # be used without refinement.
                if self.selfjoin:
                    # Self join!
                    # Have to get the two nearest neighbours
                    nearestids = self.joinlind.nearestNeighbor(
                                             inputgeom.asPoint(), 2)
                    if nearestids[0] == infeatureid and len(nearestids) > 1:
                        # The first feature is the same as the input
                        # feature, so choose the second one
                        nnfeature = self.joinvl.getFeatures(
                            QgsFeatureRequest(nearestids[1])).next()
                    else:
                        # The first feature is not the same as the
                        # input feature, so choose it
                        nnfeature = self.joinvl.getFeatures(
                            QgsFeatureRequest(nearestids[0])).next()
                    ## Pick the second closest neighbour!
                    ## (the first is supposed to be the point itself)
                    ## Should we check for coinciding points?
                    #nearestid = self.joinlind.nearestNeighbor(
                    #    inputgeom.asPoint(), 2)[1]
                    #nnfeature = self.joinvl.getFeatures(
                    #    QgsFeatureRequest(nearestid)).next()
                else:
                    # Not a self join, so we can search for only the
                    # nearest neighbour (1)
                    nearestid = self.joinlind.nearestNeighbor(
                                           inputgeom.asPoint(), 1)[0]
                    nnfeature = self.joinvl.getFeatures(
                                 QgsFeatureRequest(nearestid)).next()
                mindist = inputgeom.distance(nnfeature.geometry())
            elif (self.joinvl.wkbType() == QGis.WKBPolygon or
                  self.joinvl.wkbType() == QGis.WKBPolygon25D or
                  self.joinvl.wkbType() == QGis.WKBLineString or
                  self.joinvl.wkbType() == QGis.WKBLineString25D):
                # Use the join layer index to speed up the join when
                # the join layer geometry type is polygon or line
                # and the input layer geometry type is point or an
                # approximation (point)
                nearestindexid = self.joinlind.nearestNeighbor(
                    inputgeom.asPoint(), 1)[0]
                # Check for self join
                if self.selfjoin and nearestindexid == infeatureid:
                    # Self join and same feature, so get the two
                    # first two neighbours
                    nearestindexes = self.joinlind.nearestNeighbor(
                                             inputgeom.asPoint(), 2)
                    nearestindexid = nearestindexes[0]
                    if (nearestindexid == infeatureid and
                                  len(nearestindexes) > 1):
                        nearestindexid = nearestindexes[1]
                nnfeature = self.joinvl.getFeatures(
                    QgsFeatureRequest(nearestindexid)).next()
                mindist = inputgeom.distance(nnfeature.geometry())
                px = inputgeom.asPoint().x()
                py = inputgeom.asPoint().y()
                closefids = self.joinlind.intersects(QgsRectangle(
                    px - mindist,
                    py - mindist,
                    px + mindist,
#.........这里部分代码省略.........