Python core.QgsDistanceArea类代码示例

 def nearest_neighbour_analysis(self, vlayer):
     vprovider = vlayer.dataProvider()
     sumDist = 0.00
     distance = QgsDistanceArea()
     A = vlayer.extent()
     A = float(A.width() * A.height())
     index = ftools_utils.createIndex(vprovider)
     nFeat = vprovider.featureCount()
     nElement = 0
     if nFeat > 0:
         self.emit(SIGNAL("runStatus(PyQt_PyObject)"), 0)
         self.emit(SIGNAL("runRange(PyQt_PyObject)"), (0, nFeat))
     feat = QgsFeature()
     neighbour = QgsFeature()
     fit = vprovider.getFeatures()
     while fit.nextFeature(feat):
         neighbourID = index.nearestNeighbor(feat.geometry().asPoint(), 2)[1]
         vprovider.getFeatures(QgsFeatureRequest().setFilterFid(neighbourID).setSubsetOfAttributes([])).nextFeature(neighbour)
         nearDist = distance.measureLine(neighbour.geometry().asPoint(), feat.geometry().asPoint())
         sumDist += nearDist
         nElement += 1
         self.emit(SIGNAL("runStatus(PyQt_PyObject)"), nElement)
     nVal = vprovider.featureCount()
     do = float(sumDist) / nVal
     de = float(0.5 / math.sqrt(nVal / A))
     d = float(do / de)
     SE = float(0.26136 / math.sqrt((nVal * nVal) / A))
     zscore = float((do - de) / SE)
     lstStats = []
     lstStats.append(self.tr("Observed mean distance:") + unicode(do))
     lstStats.append(self.tr("Expected mean distance:") + unicode(de))
     lstStats.append(self.tr("Nearest neighbour index:") + unicode(d))
     lstStats.append(self.tr("N:") + unicode(nVal))
     lstStats.append(self.tr("Z-Score:") + unicode(zscore))
     return (lstStats, [])

    def testMeasurePolygon(self):
        # +-+-+
        # |   |
        # + +-+
        # | |
        # +-+
        polygon = QgsGeometry.fromPolygon(
            [
                [
                    QgsPoint(0, 0),
                    QgsPoint(1, 0),
                    QgsPoint(1, 1),
                    QgsPoint(2, 1),
                    QgsPoint(2, 2),
                    QgsPoint(0, 2),
                    QgsPoint(0, 0),
                ]
            ]
        )

        da = QgsDistanceArea()
        area = da.measure(polygon)
        assert area == 3, "Expected:\n%f\nGot:\n%f\n" % (3, area)

        perimeter = da.measurePerimeter(polygon)
        assert perimeter == 8, "Expected:\n%f\nGot:\n%f\n" % (8, perimeter)

    def regularMatrix(self, inLayer, inField, targetLayer, targetField,
                      nPoints, progress):
        features = vector.features(inLayer)
        total = 100.0 / len(features) if len(features) > 0 else 1

        targetIdx = targetLayer.fieldNameIndex(targetField)

        first = True
        distArea = QgsDistanceArea()
        index = vector.spatialindex(targetLayer)

        for current, inFeat in enumerate(features):
            inGeom = inFeat.geometry()
            if first:
                featList = index.nearestNeighbor(inGeom.asPoint(), nPoints)
                first = False
                data = ['ID']
                request = QgsFeatureRequest().setFilterFids(featList).setSubsetOfAttributes([targetIdx])
                for f in targetLayer.getFeatures(request):
                    data.append(unicode(f[targetField]))
                self.writer.addRecord(data)

            data = [unicode(inFeat[inField])]
            for f in targetLayer.getFeatures(request):
                outGeom = f.geometry()
                dist = distArea.measureLine(inGeom.asPoint(), outGeom.asPoint())
                data.append(unicode(float(dist)))
            self.writer.addRecord(data)

            progress.setPercentage(int(current * total))

    def prepareAlgorithm(self, parameters, context, feedback):
        source = self.parameterAsSource(parameters, 'INPUT', context)
        mapping = self.parameterAsFieldsMapping(parameters, self.FIELDS_MAPPING, context)

        self.fields = QgsFields()
        self.expressions = []

        da = QgsDistanceArea()
        da.setSourceCrs(source.sourceCrs())
        da.setEllipsoid(context.project().ellipsoid())

        # create an expression context using thread safe processing context
        self.expr_context = self.createExpressionContext(parameters, context, source)

        for field_def in mapping:
            self.fields.append(QgsField(name=field_def['name'],
                                        type=field_def['type'],
                                        typeName="",
                                        len=field_def.get('length', 0),
                                        prec=field_def.get('precision', 0)))
            expression = QgsExpression(field_def['expression'])
            expression.setGeomCalculator(da)
            expression.setDistanceUnits(context.project().distanceUnits())
            expression.setAreaUnits(context.project().areaUnits())
            if expression.hasParserError():
                raise QgsProcessingException(
                    self.tr(u'Parser error in expression "{}": {}')
                    .format(str(expression.expression()),
                            str(expression.parserErrorString())))
            self.expressions.append(expression)
        return True

    def regularMatrix(self, inLayer, inField, targetLayer, targetField,
                      nPoints, progress):
        index = vector.spatialindex(targetLayer)

        inIdx = inLayer.fields().lookupField(inField)

        distArea = QgsDistanceArea()

        first = True
        features = vector.features(inLayer)
        total = 100.0 / len(features)
        for current, inFeat in enumerate(features):
            inGeom = inFeat.geometry()
            inID = str(inFeat.attributes()[inIdx])
            featList = index.nearestNeighbor(inGeom.asPoint(), nPoints)
            if first:
                first = False
                data = ['ID']
                for i in range(len(featList)):
                    data.append('DIST_{0}'.format(i + 1))
                self.writer.addRecord(data)

            data = [inID]
            for i in featList:
                request = QgsFeatureRequest().setFilterFid(i)
                outFeat = next(targetLayer.getFeatures(request))
                outGeom = outFeat.geometry()
                dist = distArea.measureLine(inGeom.asPoint(),
                                            outGeom.asPoint())
                data.append(str(float(dist)))
            self.writer.addRecord(data)

            progress.setPercentage(int(current * total))

def measurePerimeter(geom):
    measure = QgsDistanceArea()
    value = 0.00
    polygon = geom.asPolygon()
    for line in polygon:
        value += measure.measureLine(line)
    return value

 def loopThruPolygons(self, inLayer, numRand, design):
     sProvider = inLayer.dataProvider()
     sFeat = QgsFeature()
     sGeom = QgsGeometry()
     sPoints = []
     if design == self.tr("field"):
       index = sProvider.fieldNameIndex(numRand)
     count = 10.00
     add = 60.00 / sProvider.featureCount()
     sFit = sProvider.getFeatures()
     featureErrors = []
     while sFit.nextFeature(sFeat):
         sGeom = sFeat.geometry()
         if design == self.tr("density"):
             sDistArea = QgsDistanceArea()
             value = int(round(numRand * sDistArea.measure(sGeom)))
         elif design == self.tr("field"):
             sAtMap = sFeat.attributes()
             try:
               value = int(sAtMap[index])
             except (ValueError,TypeError):
               featureErrors.append(sFeat)
               continue
         else:
             value = numRand
         sExt = sGeom.boundingBox()
         sPoints.extend(self.simpleRandom(value, sGeom, sExt.xMinimum(), sExt.xMaximum(), sExt.yMinimum(), sExt.yMaximum()))
         count = count + add
         self.progressBar.setValue(count)
     return sPoints, featureErrors

    def processAlgorithm(self, parameters, context, feedback):
        source = self.parameterAsSource(parameters, self.INPUT, context)
        layer = self.parameterAsVectorLayer(parameters, self.INPUT, context)
        field_name = self.parameterAsString(parameters, self.FIELD_NAME, context)
        field_type = self.TYPES[self.parameterAsEnum(parameters, self.FIELD_TYPE, context)]
        width = self.parameterAsInt(parameters, self.FIELD_LENGTH, context)
        precision = self.parameterAsInt(parameters, self.FIELD_PRECISION, context)
        new_field = self.parameterAsBool(parameters, self.NEW_FIELD, context)
        formula = self.parameterAsString(parameters, self.FORMULA, context)

        expression = QgsExpression(formula)
        da = QgsDistanceArea()
        da.setSourceCrs(source.sourceCrs())
        da.setEllipsoid(context.project().ellipsoid())
        expression.setGeomCalculator(da)

        expression.setDistanceUnits(context.project().distanceUnits())
        expression.setAreaUnits(context.project().areaUnits())

        fields = source.fields()
        field_index = fields.lookupField(field_name)
        if new_field or field_index < 0:
            fields.append(QgsField(field_name, field_type, '', width, precision))

        (sink, dest_id) = self.parameterAsSink(parameters, self.OUTPUT, context,
                                               fields, source.wkbType(), source.sourceCrs())

        exp_context = self.createExpressionContext(parameters, context)
        if layer is not None:
            exp_context.appendScope(QgsExpressionContextUtils.layerScope(layer))

        if not expression.prepare(exp_context):
            raise QgsProcessingException(
                self.tr('Evaluation error: {0}').format(expression.parserErrorString()))

        features = source.getFeatures()
        total = 100.0 / source.featureCount() if source.featureCount() else 0

        for current, f in enumerate(features):
            if feedback.isCanceled():
                break

            rownum = current + 1
            exp_context.setFeature(f)
            exp_context.lastScope().setVariable("row_number", rownum)
            value = expression.evaluate(exp_context)
            if expression.hasEvalError():
                feedback.reportError(expression.evalErrorString())
            else:
                attrs = f.attributes()
                if new_field or field_index < 0:
                    attrs.append(value)
                else:
                    attrs[field_index] = value
                f.setAttributes(attrs)
                sink.addFeature(f, QgsFeatureSink.FastInsert)
            feedback.setProgress(int(current * total))

        return {self.OUTPUT: dest_id}

    def testCrs(self):
        # test setting/getting the source CRS
        da = QgsDistanceArea()

        # try setting using a CRS object
        crs = QgsCoordinateReferenceSystem(3111, QgsCoordinateReferenceSystem.EpsgCrsId)
        da.setSourceCrs(crs, QgsProject.instance().transformContext())
        self.assertEqual(da.sourceCrs().srsid(), crs.srsid())

 def getSemiMajorAndSemiMinorAxis(self):
     """Obtains the semi major axis and semi minor axis from the used ellipsoid
     """
     currentLayer = self.iface.mapCanvas().currentLayer()
     distanceArea = QgsDistanceArea()
     distanceArea.setEllipsoid(currentLayer.crs().ellipsoidAcronym())
     a = distanceArea.ellipsoidSemiMajor()
     b = distanceArea.ellipsoidSemiMinor()
     
     return (a,b)

def simpleMeasure(geom, method=0, ellips=None, crs=None):
    # Method defines calculation type:
    # 0 - layer CRS
    # 1 - project CRS
    # 2 - ellipsoidal

    if geom.wkbType() in [QGis.WKBPoint, QGis.WKBPoint25D]:
        pt = geom.asPoint()
        attr1 = pt.x()
        attr2 = pt.y()
    elif geom.wkbType() in [QGis.WKBMultiPoint, QGis.WKBMultiPoint25D]:
        pt = geom.asMultiPoint()
        attr1 = pt[0].x()
        attr2 = pt[0].y()
    else:
        measure = QgsDistanceArea()

        if method == 2:
            measure.setSourceCrs(crs)
            measure.setEllipsoid(ellips)
            measure.setEllipsoidalMode(True)

        attr1 = measure.measure(geom)
        if geom.type() == QGis.Polygon:
            attr2 = measure.measurePerimeter(geom)
        else:
            attr2 = None

    return (attr1, attr2)

 def testMeasureLine(self):
     #   +-+
     #   | |
     # +-+ +
     linestring = QgsGeometry.fromPolyline(
         [QgsPoint(0, 0), QgsPoint(1, 0), QgsPoint(1, 1), QgsPoint(2, 1), QgsPoint(2, 0)]
     )
     da = QgsDistanceArea()
     length = da.measure(linestring)
     myMessage = "Expected:\n%f\nGot:\n%f\n" % (4, length)
     assert length == 4, myMessage

 def compute(self, inPoly, inLns, inField, outPath, progressBar):
     polyLayer = ftools_utils.getVectorLayerByName(inPoly)
     lineLayer = ftools_utils.getVectorLayerByName(inLns)
     polyProvider = polyLayer.dataProvider()
     lineProvider = lineLayer.dataProvider()
     if polyProvider.crs() != lineProvider.crs():
         QMessageBox.warning(self, self.tr("CRS warning!"), self.tr("Warning: Input layers have non-matching CRS.\nThis may cause unexpected results."))
     fieldList = ftools_utils.getFieldList(polyLayer)
     index = polyProvider.fieldNameIndex(unicode(inField))
     if index == -1:
         index = polyProvider.fields().count()
         fieldList.append(QgsField(unicode(inField), QVariant.Double, "real", 24, 15, self.tr("length field")))
     sRs = polyProvider.crs()
     inFeat = QgsFeature()
     inFeatB = QgsFeature()
     outFeat = QgsFeature()
     inGeom = QgsGeometry()
     outGeom = QgsGeometry()
     distArea = QgsDistanceArea()
     start = 0.00
     add = 100.00 / polyProvider.featureCount()
     check = QFile(self.shapefileName)
     if check.exists():
         if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
             return
     writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fieldList, polyProvider.geometryType(), sRs)
     spatialIndex = ftools_utils.createIndex(lineProvider)
     polyFit = polyProvider.getFeatures()
     while polyFit.nextFeature(inFeat):
         inGeom = QgsGeometry(inFeat.geometry())
         atMap = inFeat.attributes()
         lineList = []
         length = 0
         lineList = spatialIndex.intersects(inGeom.boundingBox())
         if len(lineList) > 0:
             check = 0
         else:
             check = 1
         if check == 0:
             for i in lineList:
                 lineProvider.getFeatures(QgsFeatureRequest().setFilterFid(int(i))).nextFeature(inFeatB)
                 tmpGeom = QgsGeometry(inFeatB.geometry())
                 if inGeom.intersects(tmpGeom):
                     outGeom = inGeom.intersection(tmpGeom)
                     length = length + distArea.measure(outGeom)
         outFeat.setGeometry(inGeom)
         atMap.append(length)
         outFeat.setAttributes(atMap)
         writer.addFeature(outFeat)
         start = start + 1
         progressBar.setValue(start * (add))
     del writer

 def testMeasureMultiLine(self):
     #   +-+ +-+-+
     #   | | |   |
     # +-+ + +   +-+
     linestring = QgsGeometry.fromMultiPolyline(
         [
             [QgsPoint(0, 0), QgsPoint(1, 0), QgsPoint(1, 1), QgsPoint(2, 1), QgsPoint(2, 0)],
             [QgsPoint(3, 0), QgsPoint(3, 1), QgsPoint(5, 1), QgsPoint(5, 0), QgsPoint(6, 0)],
         ]
     )
     da = QgsDistanceArea()
     length = da.measure(linestring)
     myMessage = "Expected:\n%f\nGot:\n%f\n" % (9, length)
     assert length == 9, myMessage

    def _calc_north(self):
        extent = self.canvas.extent()
        if self.canvas.layerCount() == 0 or extent.isEmpty():
            print "No layers or extent"
            return 0

        outcrs = self.canvas.mapSettings().destinationCrs()

        if outcrs.isValid() and not outcrs.geographicFlag():
            crs = QgsCoordinateReferenceSystem()
            crs.createFromOgcWmsCrs("EPSG:4326")

            transform = QgsCoordinateTransform(outcrs, crs)

            p1 = QgsPoint(extent.center())
            p2 = QgsPoint(p1.x(), p1.y() + extent.height() * 0.25)

            try:
                pp1 = transform.transform(p1)
                pp2 = transform.transform(p2)
            except QgsCsException:
                roam.utils.warning("North arrow. Error transforming.")
                return None

            area = QgsDistanceArea()
            area.setEllipsoid(crs.ellipsoidAcronym())
            area.setEllipsoidalMode(True)
            area.setSourceCrs(crs)
            distance, angle, _ = area.computeDistanceBearing(pp1, pp2)
            angle = math.degrees(angle)
            return angle
        else:
            return 0