本文整理匯總了Python中qgis.analysis.QgsRasterCalculator類的典型用法代碼示例。如果您正苦於以下問題：Python QgsRasterCalculator類的具體用法？Python QgsRasterCalculator怎麽用？Python QgsRasterCalculator使用的例子？那麽, 這裏精選的類代碼示例或許可以為您提供幫助。
在下文中一共展示了QgsRasterCalculator類的14個代碼示例，這些例子默認根據受歡迎程度排序。您可以為喜歡或者感覺有用的代碼點讚，您的評價將有助於係統推薦出更棒的Python代碼示例。
示例1: split_bands
def split_bands(pathIn,pathOut):
# Recibe un path de entrada (raster multibanda) y devuelve las bandas por separado
fileInfo=QFileInfo(pathIn)
baseName=fileInfo.baseName()
layer=QgsRasterLayer(pathIn, baseName)
if not layer.isValid():
print "Error importing Micasense Mosaic to split"
print "Splitting bands from " + baseName
numBands=layer.bandCount()
i=1
entries=[]
output=[]
while(i<=numBands):
band = QgsRasterCalculatorEntry()
band.ref = "[email protected]"+str(i)
band.raster=layer
band.bandNumber=i
entries.append(band)
# Saves the current band as a separate file
calc=QgsRasterCalculator(band.ref, pathOut+ "/" +baseName+"_band_"+str(i)+".tif","GTiff",layer.extent(),layer.width(),layer.height(), entries)
calc.processCalculation()
output.append(pathOut+"/"+baseName+"_band_"+str(i)+".tif")
i=i+1
return output示例2: processAlgorithm
def processAlgorithm(self, parameters, context, feedback):
expression = self.getParameterValue(self.EXPRESSION)
layersValue = self.getParameterValue(self.LAYERS)
layersDict = {}
if layersValue:
layers = [QgsProcessingUtils.mapLayerFromString(f, context) for f in layersValue.split(";")]
layersDict = {os.path.basename(lyr.source().split(".")[0]): lyr for lyr in layers}
for lyr in QgsProcessingUtils.compatibleRasterLayers(context.project()):
name = lyr.name()
if (name + "@") in expression:
layersDict[name] = lyr
entries = []
for name, lyr in layersDict.items():
for n in range(lyr.bandCount()):
entry = QgsRasterCalculatorEntry()
entry.ref = '{:s}@{:d}'.format(name, n + 1)
entry.raster = lyr
entry.bandNumber = n + 1
entries.append(entry)
output = self.getOutputValue(self.OUTPUT)
extentValue = self.getParameterValue(self.EXTENT)
if not extentValue:
extentValue = QgsProcessingUtils.combineLayerExtents(layersValue)
if extentValue:
extent = extentValue.split(',')
bbox = QgsRectangle(float(extent[0]), float(extent[2]),
float(extent[1]), float(extent[3]))
else:
if layersDict:
bbox = list(layersDict.values())[0].extent()
for lyr in layersDict.values():
bbox.combineExtentWith(lyr.extent())
else:
raise GeoAlgorithmExecutionException(self.tr("No layers selected"))
def _cellsize(layer):
return (layer.extent().xMaximum() - layer.extent().xMinimum()) / layer.width()
cellsize = self.getParameterValue(self.CELLSIZE) or min([_cellsize(lyr) for lyr in layersDict.values()])
width = math.floor((bbox.xMaximum() - bbox.xMinimum()) / cellsize)
height = math.floor((bbox.yMaximum() - bbox.yMinimum()) / cellsize)
driverName = GdalUtils.getFormatShortNameFromFilename(output)
calc = QgsRasterCalculator(expression,
output,
driverName,
bbox,
width,
height,
entries)
res = calc.processCalculation()
if res == QgsRasterCalculator.ParserError:
raise GeoAlgorithmExecutionException(self.tr("Error parsing formula"))示例3: calculo
def calculo(expresion,capa):
calc = QgsRasterCalculator(expresion,
os.path.join(carpeta,troncoresumido+'_'+capa+'.tif'),
'GTiff',
layerglobal.extent(),
layerglobal.width(),
layerglobal.height(),
entries )
calc.processCalculation()
del(calc)示例4: processAlgorithm
def processAlgorithm(self, parameters, context, feedback):
expression = self.parameterAsString(parameters, self.EXPRESSION, context)
layers = self.parameterAsLayerList(parameters, self.LAYERS, context)
layersDict = {}
if layers:
layersDict = {os.path.basename(lyr.source().split(".")[0]): lyr for lyr in layers}
for lyr in QgsProcessingUtils.compatibleRasterLayers(context.project()):
name = lyr.name()
if (name + "@") in expression:
layersDict[name] = lyr
entries = []
for name, lyr in layersDict.items():
for n in range(lyr.bandCount()):
entry = QgsRasterCalculatorEntry()
entry.ref = '{:s}@{:d}'.format(name, n + 1)
entry.raster = lyr
entry.bandNumber = n + 1
entries.append(entry)
output = self.parameterAsOutputLayer(parameters, self.OUTPUT, context)
bbox = self.parameterAsExtent(parameters, self.EXTENT, context)
if bbox.isNull():
bbox = QgsProcessingUtils.combineLayerExtents(layers)
if bbox.isNull():
if layersDict:
bbox = list(layersDict.values())[0].extent()
for lyr in layersDict.values():
bbox.combineExtentWith(lyr.extent())
else:
raise QgsProcessingException(self.tr("No layers selected"))
def _cellsize(layer):
return (layer.extent().xMaximum() - layer.extent().xMinimum()) / layer.width()
cellsize = self.parameterAsDouble(parameters, self.CELLSIZE, context) or min([_cellsize(lyr) for lyr in layersDict.values()])
width = math.floor((bbox.xMaximum() - bbox.xMinimum()) / cellsize)
height = math.floor((bbox.yMaximum() - bbox.yMinimum()) / cellsize)
driverName = GdalUtils.getFormatShortNameFromFilename(output)
calc = QgsRasterCalculator(expression,
output,
driverName,
bbox,
width,
height,
entries)
res = calc.processCalculation()
if res == QgsRasterCalculator.ParserError:
raise QgsProcessingException(self.tr("Error parsing formula"))
return {self.OUTPUT: output}示例5: calculate_PCD
def calculate_PCD(red,nir,pcdPath):
# Obtain file information and create the layers
redInfo=QFileInfo(red)
nirInfo=QFileInfo(nir)
redBaseName=redInfo.baseName()
nirBaseName=nirInfo.baseName()
folderPath = redInfo.absolutePath()
redReflectancePath = folderPath + "/red_reflectance.tif"
nirReflectancePath = folderPath + "/nir_reflectance.tif"
redLayer = QgsRasterLayer(red,redBaseName)
if not redLayer.isValid():
print "Error importing red band to calculate reflectances"
nirLayer = QgsRasterLayer(nir,nirBaseName)
if not nirLayer.isValid():
print "Error importing NIR band to calculate reflectances"
# The images are transformed into reflectances by dividing by 32768
entries=[]
redReflectance = QgsRasterCalculatorEntry()
redReflectance.ref = "[email protected]"
redReflectance.raster=redLayer
redReflectance.bandNumber = 1
entries.append(redReflectance)
# Converts the DN raster into a reflectance raster
calc=QgsRasterCalculator('float(' + redReflectance.ref + ')/32768', redReflectancePath,"GTiff",redLayer.extent(),redLayer.width(),redLayer.height(), entries)
calc.processCalculation()
nirReflectance = QgsRasterCalculatorEntry()
nirReflectance.ref = "[email protected]"
nirReflectance.raster=nirLayer
nirReflectance.bandNumber = 1
entries.append(nirReflectance)
# Converts the DN raster into a reflectance raster
calc=QgsRasterCalculator('float(' + nirReflectance.ref + ')/32768', nirReflectancePath,"GTiff",nirLayer.extent(),nirLayer.width(),nirLayer.height(), entries)
calc.processCalculation()
# Calculate the PCD index
calc=QgsRasterCalculator("float(" + nirReflectance.ref + ")/float(" + redReflectance.ref + ")", pcdPath,"GTiff",nirLayer.extent(),nirLayer.width(),nirLayer.height(), entries)
calc.processCalculation()
print "PCD calculated"示例6: reclassifyRaster
def reclassifyRaster(prjpath, inRasterName, bandnum, minValue, tupleUpValues,
outRasterName):
""" Reclassify a raster to groups defined by tupleUpValues."""
# Get raster
inRaster=getRasterLayerByName(inRasterName)
if not inRaster:
message=inRasterName+ " not loaded or not a raster!"
QMessageBox.critical(None,'Error',message, QMessageBox.Ok)
return False
# Check prjpath exists
if not os.path.isdir(prjpath):
message= prjpath + " does not exist!"
QMessageBox.critical(None,'Error',message, QMessageBox.Ok)
return False
# Define band
boh = QgsRasterCalculatorEntry()
bandName=inRasterName+'@'+str(bandnum)
boh.ref = bandName
boh.raster = inRaster
boh.bandNumber =bandnum
# Prepare entries
entries = []
entries.append( boh )
# Prepare calculation command
bandNameAddStr= '<='+ bandName + ' AND ' + bandName + '<'
i = 1
lowerVal=0
calcCommand=""
for upValue in tupleUpValues:
calcCommand = calcCommand + '( ' + str(minValue) + bandNameAddStr
calcCommand = calcCommand + str(upValue) + ')' + '*' + str(i)
if i!=len(tupleUpValues):
calcCommand = calcCommand + ' + '
minValue = upValue
i = i + 1
# Process calculation with input extent and resolution
pathFilename=os.path.join( prjpath, outRasterName) + '.tif'
calc = QgsRasterCalculator(calcCommand, pathFilename, 'GTiff',
inRaster.extent(), inRaster.width(),
inRaster.height(), entries )
if not calc: return False
ok= (calc.processCalculation() == 0)
return ok示例7: raster_subtract
def raster_subtract(self):
# if the layer does not exist it has to be created
rlayer1 = QgsMapLayerRegistry.instance().mapLayersByName( self.cb_input1.currentText() )[0]
rlayer2 = QgsMapLayerRegistry.instance().mapLayersByName( self.cb_input2.currentText() )[0]
fileName = self.lineEdit.text()
entries = []
# Define band1
boh1 = QgsRasterCalculatorEntry()
boh1.ref = '[email protected]'
boh1.raster = rlayer1
boh1.bandNumber = 1
entries.append( boh1 )
# Define band2
boh2 = QgsRasterCalculatorEntry()
boh2.ref = '[email protected]'
boh2.raster = rlayer2
boh2.bandNumber = 1
entries.append( boh2 )
# Process calculation with input extent and resolution
calc = QgsRasterCalculator( '[email protected] - [email protected]', fileName, \
'GTiff', rlayer1.extent(), \
rlayer1.width(), rlayer1.height(), entries )
calc.processCalculation()
# Load the file into the map
fileInfo = QFileInfo(fileName)
baseName = fileInfo.baseName()
root = QgsProject.instance().layerTreeRoot()
node_group1 = root.insertGroup(0, "Group 1")
node_subgroup1 = node_group1.addGroup("Sub-group 1")
# Check out signals from nodes section
# http://www.lutraconsulting.co.uk/blog/2014/07/25/qgis-layer-tree-api-part-2/
# if the layer does not exist it has to be created
if not QgsMapLayerRegistry.instance().mapLayersByName(baseName):
rOutput = QgsRasterLayer(fileName, baseName)
QgsMapLayerRegistry.instance().addMapLayer(rOutput, False)
setRamp(rOutput, self.iface)
node_layer1 = node_subgroup1.addLayer(rOutput)
# if the layer already exists trigger a refresh
else:
rOutput = QgsMapLayerRegistry.instance().mapLayersByName(baseName)[0]
rOutput.triggerRepaint()示例8: adjustRasterToBaseRaster
def adjustRasterToBaseRaster (baseRaster, adjustingRaster, outputPath):
entries = []
rCalcObj = QgsRasterCalculatorEntry()
rCalcObj.ref = '[email protected]'
rCalcObj.raster = adjustingRaster
rCalcObj.bandNumber = 1
entries.append(rCalcObj)
#print '---'
#print baseRaster.extent()
#print baseRaster.width()
#print baseRaster.height()
#print '---'
calc = QgsRasterCalculator('[email protected]', outputPath, 'GTiff',
baseRaster.extent(),
baseRaster.width(), baseRaster.height(), entries)
#print 'calc'
calc.processCalculation()示例9: generateOneValueRasterQGisCalc
def generateOneValueRasterQGisCalc (baseRaster, value, outputPath):
"""
Generates raster with extent and resolution of baseRaster, each pixel is value.
QgsRasterCalculator used
:param baseRaster:
:param value: output value
:param outputPath: path for output file
"""
entries = []
rCalcObj = QgsRasterCalculatorEntry()
rCalcObj.ref = '[email protected]'
rCalcObj.raster = baseRaster
rCalcObj.bandNumber = 1
entries.append(rCalcObj)
calc = QgsRasterCalculator(str(value), outputPath, 'GTiff',
baseRaster.extent(),
baseRaster.width(), baseRaster.height(), entries)
calc.processCalculation()示例10: processAlgorithm
def processAlgorithm(self, parameters, context, feedback):
expression = self.parameterAsString(parameters, self.EXPRESSION, context)
layers = self.parameterAsLayerList(parameters, self.LAYERS, context)
layersDict = {}
if layers:
layersDict = {lyr.source(): lyr for lyr in layers}
crs = self.parameterAsCrs(parameters, self.CRS, context)
if crs is None or not crs.isValid():
if not layers:
raise QgsProcessingException(self.tr("No reference layer selected nor CRS provided"))
else:
crs = list(layersDict.values())[0].crs()
bbox = self.parameterAsExtent(parameters, self.EXTENT, context)
if bbox.isNull() and not layers:
raise QgsProcessingException(self.tr("No reference layer selected nor extent box provided"))
if not bbox.isNull():
bboxCrs = self.parameterAsExtentCrs(parameters, self.EXTENT, context)
if bboxCrs != crs:
transform = QgsCoordinateTransform(bboxCrs, crs, context.project())
bbox = transform.transformBoundingBox(bbox)
if bbox.isNull() and layers:
bbox = QgsProcessingUtils.combineLayerExtents(layers, crs, context)
cellsize = self.parameterAsDouble(parameters, self.CELLSIZE, context)
if cellsize == 0 and not layers:
raise QgsProcessingException(self.tr("No reference layer selected nor cellsize value provided"))
def _cellsize(layer):
ext = layer.extent()
if layer.crs() != crs:
transform = QgsCoordinateTransform(layer.crs(), crs, context.project())
ext = transform.transformBoundingBox(ext)
return (ext.xMaximum() - ext.xMinimum()) / layer.width()
if cellsize == 0:
cellsize = min([_cellsize(lyr) for lyr in layersDict.values()])
# check for layers available in the model
layersDictCopy = layersDict.copy() # need a shallow copy because next calls invalidate iterator
for lyr in layersDictCopy.values():
expression = self.mappedNameToLayer(lyr, expression, layersDict, context)
# check for layers available in the project
for lyr in QgsProcessingUtils.compatibleRasterLayers(context.project()):
expression = self.mappedNameToLayer(lyr, expression, layersDict, context)
# create the list of layers to be passed as inputs to RasterCalculaltor
# at this phase expression has been modified to match available layers
# in the current scope
entries = []
for name, lyr in layersDict.items():
for n in range(lyr.bandCount()):
ref = '{:s}@{:d}'.format(name, n + 1)
if ref in expression:
entry = QgsRasterCalculatorEntry()
entry.ref = ref
entry.raster = lyr
entry.bandNumber = n + 1
entries.append(entry)
# Append any missing entry from the current project
for entry in QgsRasterCalculatorEntry.rasterEntries():
if not [e for e in entries if e.ref == entry.ref]:
entries.append(entry)
output = self.parameterAsOutputLayer(parameters, self.OUTPUT, context)
width = round((bbox.xMaximum() - bbox.xMinimum()) / cellsize)
height = round((bbox.yMaximum() - bbox.yMinimum()) / cellsize)
driverName = GdalUtils.getFormatShortNameFromFilename(output)
calc = QgsRasterCalculator(expression,
output,
driverName,
bbox,
crs,
width,
height,
entries,
context.transformContext())
res = calc.processCalculation(feedback)
if res == QgsRasterCalculator.ParserError:
raise QgsProcessingException(self.tr("Error parsing formula"))
return {self.OUTPUT: output}示例11: agregado
def agregado(rasterdeentrada):
#filtro gausian para dar valor en funcion de los vecinos
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'1.tif')
sigma=0.2
mode=1
radius=2
result=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2.tif')
processing.runalg('saga:gaussianfilter', input, sigma, mode, radius, result)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2.tif'),rasterdeentrada+str("g2"))
#filtro y me quedo con lo mayor de 0,40
calc = QgsRasterCalculator("'"+rasterdeentrada+'[email protected] > 0.40',
os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2s.tif'),
'GTiff',
layerglobal.extent(),
layerglobal.width(),
layerglobal.height(),
entries )
calc.processCalculation()
del(calc)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2s.tif'),rasterdeentrada+str("g2s"))
#convierto en nodata lo que no me interesa
calc = QgsRasterCalculator(("'"+rasterdeentrada+'[email protected]'>0)*"'"+rasterdeentrada+'[email protected]',
os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2sn.tif'),
'GTiff',
layerglobal.extent(),
layerglobal.width(),
layerglobal.height(),
entries )
calc.processCalculation()
del(calc)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2sn.tif'),rasterdeentrada+str("g2sn"))
#filtro filter clums eliminar los huecos menores de 1300 m2
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'g2sn.tif')
threshold=13
result=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'3.tif')
processing.runalg('saga:filterclumps', input, threshold, result)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'3.tif'),rasterdeentrada+str("3"))
#filtro mayorityffilter para dar valor en funcion de los vecinos
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'3.tif')
mode=0
radius=1
threshold=4
result=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'4.tif')
try:
processing.runalg('saga:majorityfilter', input, mode, radius, threshold, result)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'4.tif'),rasterdeentrada+str("4"))
#filtro para rellenar huecos pequenos
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'4.tif')
distance=3
iterations=0
band=1
mask=None
no_default_mask='True'
output=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'5.tif')
processing.runalg('gdalogr:fillnodata', input, distance, iterations, band,mask,no_default_mask, output)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'5.tif'),rasterdeentrada+str("5"))
#filtro filter clums eliminar los huecos
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'5.tif')
threshold=5
result=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'6.tif')
processing.runalg('saga:filterclumps', input, threshold, result)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'6.tif'),rasterdeentrada+str("6"))
#filtro para rellenar huecos pequenos
input=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'6.tif')
distance=3
iterations=0
band=1
mask=None
no_default_mask='True'
output=os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'7.tif')
processing.runalg('gdalogr:fillnodata', input, distance, iterations, band,mask,no_default_mask, output)
StringToRaster(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'7.tif'),rasterdeentrada+str("7"))
#lo vectorizo
processing.runalg("gdalogr:polygonize",os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'7.tif'),"DN",os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'.shp'))
#seleciono lo que me interesa
lyr=QgsVectorLayer(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'.shp'),rasterdeentrada,"ogr")
QgsMapLayerRegistry.instance().addMapLayers([lyr])
selection = lyr.getFeatures(QgsFeatureRequest().setFilterExpression(u'"DN" = 1'))
selecionado = lyr.setSelectedFeatures([s.id() for s in selection])
nbrSelected=lyr.selectedFeatureCount()
if nbrSelected > 0:
#guardo lo selecionado
processing.runalg("qgis:saveselectedfeatures",os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'.shp'),os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'2.shp'))
#calcula la superficie de esta capa pero no en todos los registros
layer=QgsVectorLayer(os.path.join(carpeta,troncoresumido+'_'+rasterdeentrada+'2.shp'),rasterdeentrada+str("2"),"ogr")
provider = layer.dataProvider()
areas = [ feat.geometry().area() for feat in layer.getFeatures() ]
indice = [ feat.id() for feat in layer.getFeatures() ]
#.........這裏部分代碼省略.........
示例12: processAlgorithm
def processAlgorithm(self, parameters, context, feedback):
expression = self.parameterAsString(parameters, self.EXPRESSION, context)
layers = self.parameterAsLayerList(parameters, self.LAYERS, context)
layersDict = {}
if layers:
layersDict = {os.path.basename(lyr.source().split(".")[0]): lyr for lyr in layers}
crs = self.parameterAsCrs(parameters, self.CRS, context)
if not layers and not crs.isValid():
raise QgsProcessingException(self.tr("No reference layer selected nor CRS provided"))
if not crs.isValid() and layers:
crs = list(layersDict.values())[0].crs()
bbox = self.parameterAsExtent(parameters, self.EXTENT, context)
if not layers and bbox.isNull():
raise QgsProcessingException(self.tr("No reference layer selected nor extent box provided"))
if not bbox.isNull():
bboxCrs = self.parameterAsExtentCrs(parameters, self.EXTENT, context)
if bboxCrs != crs:
transform = QgsCoordinateTransform(bboxCrs, crs, context.project())
bbox = transform.transformBoundingBox(bbox)
if bbox.isNull() and layers:
bbox = QgsProcessingUtils.combineLayerExtents(layers, crs)
cellsize = self.parameterAsDouble(parameters, self.CELLSIZE, context)
if not layers and cellsize == 0:
raise QgsProcessingException(self.tr("No reference layer selected nor cellsize value provided"))
def _cellsize(layer):
ext = layer.extent()
if layer.crs() != crs:
transform = QgsCoordinateTransform(layer.crs(), crs, context.project())
ext = transform.transformBoundingBox(ext)
return (ext.xMaximum() - ext.xMinimum()) / layer.width()
if cellsize == 0:
cellsize = min([_cellsize(lyr) for lyr in layersDict.values()])
for lyr in QgsProcessingUtils.compatibleRasterLayers(context.project()):
name = lyr.name()
if (name + "@") in expression:
layersDict[name] = lyr
entries = []
for name, lyr in layersDict.items():
for n in range(lyr.bandCount()):
ref = '{:s}@{:d}'.format(name, n + 1)
if ref in expression:
entry = QgsRasterCalculatorEntry()
entry.ref = ref
entry.raster = lyr
entry.bandNumber = n + 1
entries.append(entry)
output = self.parameterAsOutputLayer(parameters, self.OUTPUT, context)
width = math.floor((bbox.xMaximum() - bbox.xMinimum()) / cellsize)
height = math.floor((bbox.yMaximum() - bbox.yMinimum()) / cellsize)
driverName = GdalUtils.getFormatShortNameFromFilename(output)
calc = QgsRasterCalculator(expression,
output,
driverName,
bbox,
crs,
width,
height,
entries)
res = calc.processCalculation(feedback)
if res == QgsRasterCalculator.ParserError:
raise QgsProcessingException(self.tr("Error parsing formula"))
return {self.OUTPUT: output}示例13: StringToRaster
def StringToRaster(raster):
# Check if string is provided
fileInfo = QFileInfo(raster)
path = fileInfo.filePath()
baseName = fileInfo.baseName()
layer = QgsRasterLayer(path, baseName)
#QgsMapLayerRegistry.instance().addMapLayer(layer)
entries = []
# Define band1
boh1 = QgsRasterCalculatorEntry()
boh1.ref = '[email protected]'
boh1.raster = layer
boh1.bandNumber = 5
entries.append( boh1 )
# Process calculation with input extent and resolution
calc = QgsRasterCalculator( '([email protected]) * 166.67 + 111.67', 'C:/Hackathon Farmhack data/Output/outputfile.tif', 'GTiff', layer.extent(), layer.width(), layer.height(), entries )
calc.processCalculation()
fileInfo = QFileInfo('C:/Hackathon Farmhack data/Output/outputfile.tif')
path = fileInfo.filePath()
baseName = fileInfo.baseName()
layer = QgsRasterLayer(path, baseName)
#QgsMapLayerRegistry.instance().addMapLayer(layer)
if layer.isValid() is True:
print "Layer was loaded successfully!"
else:
print "Unable to read basename and file path - Your string is probably invalid"
shape = QgsVectorLayer('C:/Hackathon Farmhack data/perceel-hier-rdnew.geojson', 'perceel', 'ogr')
#QgsMapLayerRegistry.instance().addMapLayer(shape)
xmin = (shape.extent().xMinimum()) #extract the minimum x coord from our layer
xmax = (shape.extent().xMaximum()) #extract our maximum x coord from our layer
ymin = (shape.extent().yMinimum()) #extract our minimum y coord from our layer
ymax = (shape.extent().yMaximum()) #extract our maximum y coord from our layer
#prepare the extent in a format the VectorGrid tool can interpret (xmin,xmax,ymin,ymax)
extent = str(xmin)+ ',' + str(xmax)+ ',' +str(ymin)+ ',' +str(ymax)
# raster the given shape
processing.runalg('qgis:vectorgrid', extent, 20, 20, 0, 'C:/Hackathon Farmhack data/Output/rasterShapes.geojson')
shapeRaster = QgsVectorLayer('C:/Hackathon Farmhack data/Output/rasterShapes.geojson', 'perceelRaster', 'ogr')
shapeRaster.setCrs(QgsCoordinateReferenceSystem(28992, QgsCoordinateReferenceSystem.EpsgCrsId))
#QgsMapLayerRegistry.instance().addMapLayer(shapeRaster)
#clip the raster returned
processing.runalg('qgis:clip', shapeRaster, shape, 'C:/Hackathon Farmhack data/Output/clippedRaster.shp')
#define oldPath and newPath
ogr2ogr.main(["","-f", "ESRI Shapefile", "-s_srs", "epsg:28992", "-t_srs", "epsg:32632", "C:/Hackathon Farmhack data/Output/clippedRasterNew.shp", "C:/Hackathon Farmhack data/Output/clippedRaster.shp"])
clippedRaster = QgsVectorLayer('C:/Hackathon Farmhack data/Output/clippedRasterNew.shp', 'clippedRaster', 'ogr')
clippedRaster.setCrs(QgsCoordinateReferenceSystem(32632, QgsCoordinateReferenceSystem.EpsgCrsId))
#QgsMapLayerRegistry.instance().addMapLayer(clippedRaster)
#zonalstatistics
processing.runalg('qgis:zonalstatistics', layer, 1, clippedRaster, '', False, 'C:/Hackathon Farmhack data/Output/filledRaster.geojson')
filledRaster = QgsVectorLayer('C:/Hackathon Farmhack data/Output/filledRaster.geojson', 'filledRaster', 'ogr')
filledRaster.setCrs(QgsCoordinateReferenceSystem(32632, QgsCoordinateReferenceSystem.EpsgCrsId))
#QgsMapLayerRegistry.instance().addMapLayer(filledRaster)
ogr2ogr.main(["","-f", "GeoJSON", "-s_srs", "epsg:32632", "-t_srs", "epsg:4326", "C:/Hackathon Farmhack data/Output/taakkaart.geojson", "C:/Hackathon Farmhack data/Output/filledRaster.geojson"])
taakkaart = QgsVectorLayer('C:/Hackathon Farmhack data/Output/taakkaart.geojson', 'taakkaart', 'ogr')
QgsMapLayerRegistry.instance().addMapLayer(taakkaart) 示例14:
##reflectance_mult=number 0
##reflectance_add=number 0
##sun_elevation=number 0
##inImage=raster
##outImage=output raster
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
lyr = processing.getObject(inImage)
entries=[]
rasterCalcEntry=QgsRasterCalculatorEntry()
rasterCalcEntry.ref='[email protected]'
rasterCalcEntry.raster=lyr
rasterCalcEntry.bandNumber=1
entries.append(rasterCalcEntry)
if not ".tif" in outImage:
outImage=outImage+".tif"
noData=-3.4028234663852886e+38
radElev = sun_elevation*3.14159/180
calc=QgsRasterCalculator('([email protected] != 0)*('+str(reflectance_mult)+' * [email protected] + '+str(reflectance_add)+')/sin('+str(radElev)+') + ([email protected] = 0) * '+str(noData), outImage, "GTiff", lyr.extent(), lyr.crs(), lyr.width(), lyr.height(), entries)
calc.processCalculation()