本文整理汇总了Python中qgis.core.QgsGeometry.transform方法的典型用法代码示例。如果您正苦于以下问题:Python QgsGeometry.transform方法的具体用法?Python QgsGeometry.transform怎么用?Python QgsGeometry.transform使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类qgis.core.QgsGeometry
的用法示例。
在下文中一共展示了QgsGeometry.transform方法的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: transformGeom
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def transformGeom(self, item):
src_crs = QgsCoordinateReferenceSystem()
src_crs.createFromSrid(item.srid)
dest_crs = self.mapCanvas.mapRenderer().destinationCrs()
geom = QgsGeometry(item.geometry)
geom.transform(QgsCoordinateTransform(src_crs, dest_crs))
return geom
示例2: processAlgorithm
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def processAlgorithm(self, parameters, context, feedback):
source = self.parameterAsSource(parameters, self.INPUT, context)
extent = self.parameterAsExtent(parameters, self.TARGET_AREA, context)
target_crs = self.parameterAsCrs(parameters, self.TARGET_AREA_CRS, context)
target_geom = QgsGeometry.fromRect(extent)
fields = QgsFields()
fields.append(QgsField('auth_id', QVariant.String, '', 20))
(sink, dest_id) = self.parameterAsSink(parameters, self.OUTPUT, context,
fields, QgsWkbTypes.NoGeometry, QgsCoordinateReferenceSystem())
# make intersection tests nice and fast
engine = QgsGeometry.createGeometryEngine(target_geom.constGet())
engine.prepareGeometry()
layer_bounds = QgsGeometry.fromRect(source.sourceExtent())
crses_to_check = QgsCoordinateReferenceSystem.validSrsIds()
total = 100.0 / len(crses_to_check)
found_results = 0
transform_context = QgsCoordinateTransformContext()
for current, srs_id in enumerate(crses_to_check):
if feedback.isCanceled():
break
candidate_crs = QgsCoordinateReferenceSystem.fromSrsId(srs_id)
if not candidate_crs.isValid():
continue
transform_candidate = QgsCoordinateTransform(candidate_crs, target_crs, transform_context)
transformed_bounds = QgsGeometry(layer_bounds)
try:
if not transformed_bounds.transform(transform_candidate) == 0:
continue
except:
continue
try:
if engine.intersects(transformed_bounds.constGet()):
feedback.pushInfo(self.tr('Found candidate CRS: {}').format(candidate_crs.authid()))
f = QgsFeature(fields)
f.setAttributes([candidate_crs.authid()])
sink.addFeature(f, QgsFeatureSink.FastInsert)
found_results += 1
except:
continue
feedback.setProgress(int(current * total))
if found_results == 0:
feedback.reportError(self.tr('No matching projections found'))
return {self.OUTPUT: dest_id}
示例3: set_next_analysis_extent
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def set_next_analysis_extent(self, extent, crs):
"""Setter for the next analysis extent.
This function will redraw the rubberband if needed.
:param extent: The next analysis extent.
:type extent: QgsGeometry
:param crs: The CRS of the extent.
:type crs: QgsCoordinateReferenceSystem
"""
extent = QgsGeometry(extent)
transform = QgsCoordinateTransform(crs, self.crs)
extent.transform(transform)
self._next_analysis_extent = extent
if self._show_rubber_bands:
self.display_next_analysis()
示例4: set_user_extent
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def set_user_extent(self, extent, crs):
"""Setter for the user requested extent.
This function will redraw the rubberband if needed.
:param extent: The user extent.
:type extent: QgsGeometry
:param crs: The CRS of the extent.
:type crs: QgsCoordinateReferenceSystem
"""
extent = QgsGeometry(extent)
transform = QgsCoordinateTransform(crs, self.crs)
extent.transform(transform)
self._user_extent = extent
set_setting('user_extent', extent.exportToWkt())
set_setting('user_extent_crs', crs.authid())
if self._show_rubber_bands:
self.display_user_extent()
示例5: _addHighlightGeom
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def _addHighlightGeom(self):
def highlight():
rb = QgsRubberBand( self.canvas, QGis.Polygon)
rb.setBorderColor( QColor( 255, 0, 0 ) )
rb.setWidth( 2 )
rb.setToGeometry( geomRB, None )
return rb
crsCanvas = self.canvas.mapSettings().destinationCrs()
crsLayer = self.layerSrc.crs()
if not crsCanvas == crsLayer:
geomRB = QgsGeometry( self.geom )
ct = QgsCoordinateTransform( crsLayer, crsCanvas )
geomRB.transform( ct )
else:
geomRB = self.geom
return highlight()
示例6: dissolvePolygonsOnCanvas
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def dissolvePolygonsOnCanvas(writer, layer):
"""dissolve polygons of the layer and clip the dissolution with base extent"""
settings = writer.settings
baseExtent = settings.baseExtent
baseExtentGeom = baseExtent.geometry()
rotation = baseExtent.rotation()
transform = QgsCoordinateTransform(layer.crs(), settings.crs)
combi = None
request = QgsFeatureRequest()
request.setFilterRect(transform.transformBoundingBox(baseExtent.boundingBox(), QgsCoordinateTransform.ReverseTransform))
for f in layer.getFeatures(request):
geometry = f.geometry()
if geometry is None:
logMessage("null geometry skipped")
continue
# coordinate transformation - layer crs to project crs
geom = QgsGeometry(geometry)
if geom.transform(transform) != 0:
logMessage("Failed to transform geometry")
continue
# check if geometry intersects with the base extent (rotated rect)
if rotation and not baseExtentGeom.intersects(geom):
continue
if combi:
combi = combi.combine(geom)
else:
combi = geom
# clip geom with slightly smaller extent than base extent
# to make sure that the clipped polygon stays within the base extent
geom = combi.intersection(baseExtent.clone().scale(0.999999).geometry())
if geom is None:
return None
# check if geometry is empty
if geom.isGeosEmpty():
logMessage("empty geometry")
return None
return geom
示例7: processAlgorithm
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT_LAYER), context)
extent = self.getParameterValue(self.TARGET_AREA).split(',')
if not extent:
extent = QgsProcessingUtils.combineLayerExtents([layer])
target_crs = QgsCoordinateReferenceSystem(self.getParameterValue(self.TARGET_AREA_CRS))
target_geom = QgsGeometry.fromRect(QgsRectangle(float(extent[0]), float(extent[2]),
float(extent[1]), float(extent[3])))
output_file = self.getOutputValue(self.OUTPUT_HTML_FILE)
# make intersection tests nice and fast
engine = QgsGeometry.createGeometryEngine(target_geom.geometry())
engine.prepareGeometry()
layer_bounds = QgsGeometry.fromRect(layer.extent())
results = []
for srs_id in QgsCoordinateReferenceSystem.validSrsIds():
candidate_crs = QgsCoordinateReferenceSystem.fromSrsId(srs_id)
if not candidate_crs.isValid():
continue
transform_candidate = QgsCoordinateTransform(candidate_crs, target_crs)
transformed_bounds = QgsGeometry(layer_bounds)
try:
if not transformed_bounds.transform(transform_candidate) == 0:
continue
except:
continue
if engine.intersects(transformed_bounds.geometry()):
results.append(candidate_crs.authid())
self.createHTML(output_file, results)
示例8: geomTransform
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def geomTransform(self, geom, crs_orig, crs_dest):
g = QgsGeometry(geom)
crsTransform = QgsCoordinateTransform(
crs_orig, crs_dest, QgsCoordinateTransformContext()) # which context ?
g.transform(crsTransform)
return g
示例9: multi_buffering
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [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'])
#.........这里部分代码省略.........
示例10: setSelectFeatures
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [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()
#.........这里部分代码省略.........
示例11: export_geometry_info
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def export_geometry_info( self ):
ellips = None
crs = None
coordTransform = None
# calculate with:
# 0 - layer CRS
# 1 - project CRS
# 2 - ellipsoidal
if self.myCalcType == 2:
settings = QSettings()
ellips = settings.value( "/qgis/measure/ellipsoid", "WGS84" )
crs = self.vlayer.crs().srsid()
elif self.myCalcType == 1:
mapCRS = self.parent.iface.mapCanvas().mapRenderer().destinationCrs()
layCRS = self.vlayer.crs()
coordTransform = QgsCoordinateTransform( layCRS, mapCRS )
inFeat = QgsFeature()
outFeat = QgsFeature()
inGeom = QgsGeometry()
nElement = 0
vprovider = self.vlayer.dataProvider()
self.emit( SIGNAL( "runStatus( PyQt_PyObject )" ), 0)
self.emit( SIGNAL( "runRange( PyQt_PyObject )" ), ( 0, vprovider.featureCount() ) )
( fields, index1, index2 ) = self.checkMeasurementFields( self.vlayer, not self.writeShape )
if self.writeShape:
writer = QgsVectorFileWriter( self.myName, self.myEncoding, fields,
vprovider.geometryType(), vprovider.crs() )
fit = vprovider.getFeatures()
while fit.nextFeature(inFeat):
self.emit( SIGNAL( "runStatus( PyQt_PyObject )" ), nElement )
nElement += 1
inGeom = inFeat.geometry()
if self.myCalcType == 1:
inGeom.transform( coordTransform )
( attr1, attr2 ) = self.simpleMeasure( inGeom, self.myCalcType, ellips, crs )
if self.writeShape:
outFeat.setGeometry( inGeom )
atMap = inFeat.attributes()
maxIndex = index1 if index1 > index2 else index2
if maxIndex >= len(atMap):
atMap += [ "" ] * ( index2+1 - len(atMap) )
atMap[ index1 ] = attr1
if index1 != index2:
atMap[ index2 ] = attr2
outFeat.setAttributes( atMap )
writer.addFeature( outFeat )
else:
changeMap = {}
changeMap[ inFeat.id() ] = {}
changeMap[ inFeat.id() ][ index1 ] = attr1
if index1!=index2:
changeMap[ inFeat.id() ][ index2 ] = attr2
vprovider.changeAttributeValues( changeMap )
self.vlayer.updateFields()
if self.writeShape:
del writer
return True
示例12: processAlgorithm
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def processAlgorithm(self, progress):
layer = dataobjects.getObjectFromUri(
self.getParameterValue(self.INPUT))
method = self.getParameterValue(self.METHOD)
geometryType = layer.geometryType()
idx1 = -1
idx2 = -1
fields = layer.pendingFields()
if geometryType == QGis.Polygon:
(idx1, fields) = vector.findOrCreateField(layer, fields, 'area',
21, 6)
(idx2, fields) = vector.findOrCreateField(layer, fields,
'perimeter', 21, 6)
elif geometryType == QGis.Line:
(idx1, fields) = vector.findOrCreateField(layer, fields, 'length',
21, 6)
idx2 = idx1
else:
(idx1, fields) = vector.findOrCreateField(layer, fields, 'xcoord',
21, 6)
(idx2, fields) = vector.findOrCreateField(layer, fields, 'ycoord',
21, 6)
writer = self.getOutputFromName(
self.OUTPUT).getVectorWriter(fields.toList(),
layer.dataProvider().geometryType(), layer.crs())
ellips = None
crs = None
coordTransform = None
# Calculate with:
# 0 - layer CRS
# 1 - project CRS
# 2 - ellipsoidal
if method == 2:
ellips = QgsProject.instance().readEntry('Measure', '/Ellipsoid',
'NONE')[0]
crs = layer.crs().srsid()
elif method == 1:
mapCRS = iface.mapCanvas().mapRenderer().destinationCrs()
layCRS = layer.crs()
coordTransform = QgsCoordinateTransform(layCRS, mapCRS)
outFeat = QgsFeature()
inGeom = QgsGeometry()
outFeat.initAttributes(len(fields))
outFeat.setFields(fields)
current = 0
features = vector.features(layer)
total = 100.0 / float(len(features))
for f in features:
inGeom = f.geometry()
if method == 1:
inGeom.transform(coordTransform)
(attr1, attr2) = vector.simpleMeasure(inGeom, method, ellips, crs)
outFeat.setGeometry(inGeom)
attrs = f.attributes()
attrs.insert(idx1, attr1)
if attr2 is not None:
attrs.insert(idx2, attr2)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat)
current += 1
progress.setPercentage(int(current * total))
del writer
示例13: features
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [as 别名]
def features(self, request=None, clipGeom=None):
settings = self.writer.settings
mapTo3d = settings.mapTo3d()
baseExtent = settings.baseExtent
baseExtentGeom = baseExtent.geometry()
rotation = baseExtent.rotation()
prop = self.prop
useZ = prop.useZ()
if useZ:
srs_from = osr.SpatialReference()
srs_from.ImportFromProj4(str(self.layer.crs().toProj4()))
srs_to = osr.SpatialReference()
srs_to.ImportFromProj4(str(self.writer.settings.crs.toProj4()))
ogr_transform = osr.CreateCoordinateTransformation(srs_from, srs_to)
clipGeomWkb = clipGeom.asWkb() if clipGeom else None
ogr_clipGeom = ogr.CreateGeometryFromWkb(clipGeomWkb) if clipGeomWkb else None
else:
# z_func: function to get elevation at given point (x, y) on surface
if prop.isHeightRelativeToDEM():
if self.geomType == QGis.Polygon and prop.type_index == 1: # Overlay
z_func = lambda x, y: 0
else:
# get elevation from DEM
z_func = lambda x, y: self.writer.demProvider.readValue(x, y)
else:
z_func = lambda x, y: 0
request = request or QgsFeatureRequest()
for f in self.layer.getFeatures(request):
geometry = f.geometry()
if geometry is None:
logMessage("null geometry skipped")
continue
# coordinate transformation - layer crs to project crs
geom = QgsGeometry(geometry)
if geom.transform(self.transform) != 0:
logMessage("Failed to transform geometry")
continue
# check if geometry intersects with the base extent (rotated rect)
if rotation and not baseExtentGeom.intersects(geom):
continue
# create feature
feat = Feature(self.writer, self, f)
# transform_func: function to transform the map coordinates to 3d coordinates
relativeHeight = prop.relativeHeight(f)
def transform_func(x, y, z):
return mapTo3d.transform(x, y, z + relativeHeight)
if useZ:
ogr_geom = ogr.CreateGeometryFromWkb(geometry.asWkb())
# transform geometry from layer CRS to project CRS
if ogr_geom.Transform(ogr_transform) != 0:
logMessage("Failed to transform geometry")
continue
# clip geometry
if ogr_clipGeom and self.geomType == QGis.Line:
ogr_geom = ogr_geom.Intersection(ogr_clipGeom)
if ogr_geom is None:
continue
# check if geometry is empty
if ogr_geom.IsEmpty():
logMessage("empty geometry skipped")
continue
feat.geom = self.geomClass.fromOgrGeometry25D(ogr_geom, transform_func)
else:
# clip geometry
if clipGeom and self.geomType in [QGis.Line, QGis.Polygon]:
geom = geom.intersection(clipGeom)
if geom is None:
continue
# check if geometry is empty
if geom.isGeosEmpty():
logMessage("empty geometry skipped")
continue
if self.geomType == QGis.Polygon:
feat.geom = self.geomClass.fromQgsGeometry(geom, z_func, transform_func, self.hasLabel())
if prop.type_index == 1 and prop.isHeightRelativeToDEM(): # Overlay and relative to DEM
feat.geom.splitPolygon(self.writer.triangleMesh())
else:
feat.geom = self.geomClass.fromQgsGeometry(geom, z_func, transform_func)
if feat.geom is None:
continue
#.........这里部分代码省略.........
示例14: in_mask
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [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
示例15: do_indexjoin
# 需要导入模块: from qgis.core import QgsGeometry [as 别名]
# 或者: from qgis.core.QgsGeometry import transform [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,
#.........这里部分代码省略.........