我试图创建一个脚本,它读取一个光栅数组并计算一个参数,但是我得到了一个错误:操作数不能与形状一起广播(27184310)(70,86)。在
# -*- coding: cp1252 -*-
from osgeo import gdal
import numpy as np
import sys ,os
from osgeo.gdalconst import *
output_dir='C:\\wamp\\www\\Donnees_CRTS\\irrisat_dessimination-smiej'
etp = os.path.join(output_dir,"L8_L8_ETpot_24_30m_2016_223.tif")
rain = os.path.join(output_dir,"GLDAS_NOAH025SUBP_3H.A2016223.sumday.rain.tif")
eff = os.path.join(output_dir,"eff_v11.tif")
driver=gdal.GetDriverByName('GTIff')
driver.Register()
paths = []
paths.append(etp)
paths.append(rain)
paths.append(eff)
raster_px = []
rasters_px = []
bands = []
def open_raster(raster):
for i in range(len(paths)) :
raster = gdal.Open(paths[i])
columns = raster.RasterXSize
lines = raster.RasterYSize
band = raster.GetRasterBand(1).ReadAsArray(0,0,columns,lines)
bands.append(band)
gt=raster.GetGeoTransform()
raster_px.append(gt)
band=None
if raster is None :
print ("Erreur : Impossible d'ouvrir le raster: ")
try :
smallest = min(raster_px)
raster.SetGeoTransform(smallest)
except:
print ("Erreur : Impossible d'extraire la bande")
return smallest, columns, lines
def parameters(parameter) :
diff = bands[0]-bands[1]
if parameter == "IWR" :
iwr=diff/bands[2]
return iwr
def createImage(new_image,columns,lines,smallest):
new_image=driver.Create("iwr9.tif", columns,lines, 1, gdal.GDT_Byte)
new_image.SetProjection("EPSG:4326")
new_image.SetGeoTransform(smallest)
new_image.GetRasterBand(1).WriteRaster( 0, 0, columns, lines, output_dir
)
new_image=None
return columns, lines
def main() :
smallest, columns, lines = open_raster(paths)
p = "IWR"
parameters(p)
createImage("iwr9.tif",columns, lines,smallest)
main()
谢谢你的帮助
据我所知,你要做以下事情:
您有三个不同的.tif文件,每个文件包含一个波段。将这些图像的内容读入(numpy)数组,并将它们存储在一个数组列表中,称为
bands
。所以band[0]
是etp的数组,band[1]
是雨的数组,band[2]
是eff的数组。这三个数组具有x&y维(即图像的维数/“picture”)。在bands[1] - bands[0]
是一个数组除法,它意味着从另一个数组的项中减去每个数组项的值。在例如
上面的例子处理的是1D数组,但是2D数组也是如此。现在,如果python告诉您它不能将一个数组“广播”到另一个数组中,那么这意味着逐项计算失败,因为数组的形状不同。在上面的示例中,如果
x
的长度为5,而y
的长度为4,则numpy/python无法知道如何处理丢失的项并抛出一个错误。在你可以通过查看你的带子项目的形状来检查(确保你
^{pr2}$import numpy as np
)可能他们会有所不同,这会导致你的错误。因此,在你执行像素代数之前,首先你必须确保你所有的图像都有相同的大小。在
相关问题 更多 >
编程相关推荐