GDAL ReadAsArray() 不断返回 None
GDAL ReadAsArray() keeps returning None
我无法阻止 ReadAsArray() 方法返回 None。我不明白为什么因为我用完全相同的方式编写了另一个运行良好的脚本。
我想提早指出我确实启用了 gdal.UseExceptions()。
这似乎是搜索中的一个常见问题,但我还没有找到使其工作的解决方案,这很奇怪,我在一个现有的脚本中做了完全相同的事情,只是找到了。我在其他脚本中使用过这个文件,所以我知道它没有损坏。另外,将另一个图像加载到其中会导致同样的问题。
有人知道这是怎么回事吗?
这是当前的脚本:
import gdal
from gdal import Open, GetDriverByName
def qa_mask(in_qa_band, in_cols, in_rows, in_geotransform, mask_tiff):
qa_np = in_qa_band.ReadAsArray(0, 0, in_cols, in_rows)
print(qa_np)
qa_np_str = qa_np.astype('str')
mask = qa_np_str[(qa_np_str[0:2] == '10')] = -99
geotiff = GetDriverByName('GTiff')
output = geotiff.Create(mask_tiff, in_cols, in_rows, 1, gdal.GDT_Byte)
output_band = output.GetRasterBand(1)
output_band.WriteArray(mask)
output.SetGeoTransform(in_geotransform)
return None
gdal.UseExceptions()
tiff = Open(r'D:\landsat_data160710_landsat8qa_test\LC80410222016190LGN00_BQA.TIF')
band = tiff.GetRasterBand(1)
rows, cols, geotransform = tiff.RasterYSize, tiff.RasterXSize, tiff.GetGeoTransform()
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
我总是得到错误:
Traceback (most recent call last):
File "D:/gis_image/landsat8qa_gdal.py", line 27, in <module>
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
File "D:/gis_image/landsat8qa_gdal.py", line 8, in qa_mask
qa_np_str = qa_np.astype('str')
AttributeError: 'NoneType' object has no attribute 'astype'
然而,这个脚本工作得很好:
import numpy as np
from numpy import subtract, add, divide, multiply
import gdal
from gdal import GetDriverByName
def ndvi(in_nir_band, in_colour_band, in_rows, in_cols, in_geotransform, out_tiff, data_type=gdal.GDT_Float32):
"""
Performs an NDVI calculation given two input bands, as well as other information that can be retrieved from the
original image.
@param in_nir_band A GDAL band object representing the near-infrared image data.
@type in_nir_band GDALRasterBand
@param in_colour_band A GDAL band object representing the colour image data.
@type: in_colour_band GDALRasterBand
@param in_rows The number of rows in both input bands.
@type: in_rows int
@param in_cols The number of columns in both input bands.
@type: in_cols int
@param in_geotransform The geographic transformation to be applied to the output image.
@type in_geotransform Tuple (as returned by GetGeoTransform())
@param out_tiff Path to the desired output .tif file.
@type: out_tiff String (should end in ".tif")
@param data_type Data type of output image. Valid values are gdal.UInt16 and gdal.Float32. Default is
gdal.Float32
@type data_type GDALDataType
@return None
"""
# Read the input bands as numpy arrays.
np_nir = in_nir_band.ReadAsArray(0, 0, in_cols, in_rows)
np_colour = in_colour_band.ReadAsArray(0, 0, in_cols, in_rows)
# Convert the np arrays to 32-bit floating point to make sure division will occur properly.
np_nir_as32 = np_nir.astype(np.float32)
np_colour_as32 = np_colour.astype(np.float32)
# Calculate the NDVI formula.
numerator = subtract(np_nir_as32, np_colour_as32)
denominator = add(np_nir_as32, np_colour_as32)
result = divide(numerator, denominator)
# Remove any out-of-bounds areas
result[result == -0] = -99
# Initialize a geotiff driver.
geotiff = GetDriverByName('GTiff')
# If the desired output is an int16, map the domain [-1,1] to [0,255], create an int16 geotiff with one band and
# write the contents of the int16 NDVI calculation to it. Otherwise, create a float32 geotiff with one band and
# write the contents of the float32 NDVI calculation to it.
if data_type == gdal.GDT_UInt16:
ndvi_int8 = multiply((result + 1), (2**7 - 1))
output = geotiff.Create(out_tiff, in_cols, in_rows, 1, gdal.GDT_Byte)
output_band = output.GetRasterBand(1)
output_band.SetNoDataValue(-99)
output_band.WriteArray(ndvi_int8)
elif data_type == gdal.GDT_Float32:
output = geotiff.Create(out_tiff, in_cols, in_rows, 1, gdal.GDT_Float32)
output_band = output.GetRasterBand(1)
output_band.SetNoDataValue(-99)
output_band.WriteArray(result)
else:
raise ValueError('Invalid output data type. Valid types are gdal.UInt16 or gdal.Float32.')
# Set the geographic transformation as the input.
output.SetGeoTransform(in_geotransform)
return None
# Open NIR image and get its only band.
nir_tiff = Open(r'D:\landsat_data160621_fort_mac\nir.tif')
nir_band = nir_tiff.GetRasterBand(1)
# Open red image and get its only band.
red_tiff = Open(r'D:\landsat_data160621_fort_mac\red.tif')
red_band = red_tiff.GetRasterBand(1)
# Get the rows and cols from one of the images (both should always be the same)
rows, cols, geotransform = nir_tiff.RasterYSize, nir_tiff.RasterXSize, nir_tiff.GetGeoTransform()
#print(geotransform)
# Set an output for a 16-bit unsigned integer (0-255)
#out_tiff_int16 = r'NDVI_INT16.tif'
# Set the output for a 32-bit floating point (-1 to 1)
out_tiff_float32 = r'D:\landsat_data160621_fort_mac\out3.tif'
# Run the function for unsigned 16-bit integer
#ndvi(nir_band, red_band, rows, cols, geotransform, out_tiff_int16, gdal.GDT_UInt16)
# Run the function for 32-bit floating point
ndvi(nir_band, red_band, rows, cols, geotransform, out_tiff_float32, gdal.GDT_Float32)
print('done')
这两个脚本都接受栅格波段作为输入并尝试将它们读取为 numpy 数组
嗯,我不觉得自己很蠢吗...
函数定义为:
def qa_mask(in_qa_band, in_cols, in_rows, in_geotransform, mask_tiff):
而我称它为:
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
所以,这里要吸取的教训是,显然 GDAL 不会告诉您交换 rows/cols 是否有问题。它不会定义 numpy 数组...
我无法阻止 ReadAsArray() 方法返回 None。我不明白为什么因为我用完全相同的方式编写了另一个运行良好的脚本。
我想提早指出我确实启用了 gdal.UseExceptions()。
这似乎是搜索中的一个常见问题,但我还没有找到使其工作的解决方案,这很奇怪,我在一个现有的脚本中做了完全相同的事情,只是找到了。我在其他脚本中使用过这个文件,所以我知道它没有损坏。另外,将另一个图像加载到其中会导致同样的问题。
有人知道这是怎么回事吗?
这是当前的脚本:
import gdal
from gdal import Open, GetDriverByName
def qa_mask(in_qa_band, in_cols, in_rows, in_geotransform, mask_tiff):
qa_np = in_qa_band.ReadAsArray(0, 0, in_cols, in_rows)
print(qa_np)
qa_np_str = qa_np.astype('str')
mask = qa_np_str[(qa_np_str[0:2] == '10')] = -99
geotiff = GetDriverByName('GTiff')
output = geotiff.Create(mask_tiff, in_cols, in_rows, 1, gdal.GDT_Byte)
output_band = output.GetRasterBand(1)
output_band.WriteArray(mask)
output.SetGeoTransform(in_geotransform)
return None
gdal.UseExceptions()
tiff = Open(r'D:\landsat_data160710_landsat8qa_test\LC80410222016190LGN00_BQA.TIF')
band = tiff.GetRasterBand(1)
rows, cols, geotransform = tiff.RasterYSize, tiff.RasterXSize, tiff.GetGeoTransform()
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
我总是得到错误:
Traceback (most recent call last):
File "D:/gis_image/landsat8qa_gdal.py", line 27, in <module>
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
File "D:/gis_image/landsat8qa_gdal.py", line 8, in qa_mask
qa_np_str = qa_np.astype('str')
AttributeError: 'NoneType' object has no attribute 'astype'
然而,这个脚本工作得很好:
import numpy as np
from numpy import subtract, add, divide, multiply
import gdal
from gdal import GetDriverByName
def ndvi(in_nir_band, in_colour_band, in_rows, in_cols, in_geotransform, out_tiff, data_type=gdal.GDT_Float32):
"""
Performs an NDVI calculation given two input bands, as well as other information that can be retrieved from the
original image.
@param in_nir_band A GDAL band object representing the near-infrared image data.
@type in_nir_band GDALRasterBand
@param in_colour_band A GDAL band object representing the colour image data.
@type: in_colour_band GDALRasterBand
@param in_rows The number of rows in both input bands.
@type: in_rows int
@param in_cols The number of columns in both input bands.
@type: in_cols int
@param in_geotransform The geographic transformation to be applied to the output image.
@type in_geotransform Tuple (as returned by GetGeoTransform())
@param out_tiff Path to the desired output .tif file.
@type: out_tiff String (should end in ".tif")
@param data_type Data type of output image. Valid values are gdal.UInt16 and gdal.Float32. Default is
gdal.Float32
@type data_type GDALDataType
@return None
"""
# Read the input bands as numpy arrays.
np_nir = in_nir_band.ReadAsArray(0, 0, in_cols, in_rows)
np_colour = in_colour_band.ReadAsArray(0, 0, in_cols, in_rows)
# Convert the np arrays to 32-bit floating point to make sure division will occur properly.
np_nir_as32 = np_nir.astype(np.float32)
np_colour_as32 = np_colour.astype(np.float32)
# Calculate the NDVI formula.
numerator = subtract(np_nir_as32, np_colour_as32)
denominator = add(np_nir_as32, np_colour_as32)
result = divide(numerator, denominator)
# Remove any out-of-bounds areas
result[result == -0] = -99
# Initialize a geotiff driver.
geotiff = GetDriverByName('GTiff')
# If the desired output is an int16, map the domain [-1,1] to [0,255], create an int16 geotiff with one band and
# write the contents of the int16 NDVI calculation to it. Otherwise, create a float32 geotiff with one band and
# write the contents of the float32 NDVI calculation to it.
if data_type == gdal.GDT_UInt16:
ndvi_int8 = multiply((result + 1), (2**7 - 1))
output = geotiff.Create(out_tiff, in_cols, in_rows, 1, gdal.GDT_Byte)
output_band = output.GetRasterBand(1)
output_band.SetNoDataValue(-99)
output_band.WriteArray(ndvi_int8)
elif data_type == gdal.GDT_Float32:
output = geotiff.Create(out_tiff, in_cols, in_rows, 1, gdal.GDT_Float32)
output_band = output.GetRasterBand(1)
output_band.SetNoDataValue(-99)
output_band.WriteArray(result)
else:
raise ValueError('Invalid output data type. Valid types are gdal.UInt16 or gdal.Float32.')
# Set the geographic transformation as the input.
output.SetGeoTransform(in_geotransform)
return None
# Open NIR image and get its only band.
nir_tiff = Open(r'D:\landsat_data160621_fort_mac\nir.tif')
nir_band = nir_tiff.GetRasterBand(1)
# Open red image and get its only band.
red_tiff = Open(r'D:\landsat_data160621_fort_mac\red.tif')
red_band = red_tiff.GetRasterBand(1)
# Get the rows and cols from one of the images (both should always be the same)
rows, cols, geotransform = nir_tiff.RasterYSize, nir_tiff.RasterXSize, nir_tiff.GetGeoTransform()
#print(geotransform)
# Set an output for a 16-bit unsigned integer (0-255)
#out_tiff_int16 = r'NDVI_INT16.tif'
# Set the output for a 32-bit floating point (-1 to 1)
out_tiff_float32 = r'D:\landsat_data160621_fort_mac\out3.tif'
# Run the function for unsigned 16-bit integer
#ndvi(nir_band, red_band, rows, cols, geotransform, out_tiff_int16, gdal.GDT_UInt16)
# Run the function for 32-bit floating point
ndvi(nir_band, red_band, rows, cols, geotransform, out_tiff_float32, gdal.GDT_Float32)
print('done')
这两个脚本都接受栅格波段作为输入并尝试将它们读取为 numpy 数组
嗯,我不觉得自己很蠢吗...
函数定义为:
def qa_mask(in_qa_band, in_cols, in_rows, in_geotransform, mask_tiff):
而我称它为:
qa_mask(band, rows, cols, geotransform, 'D:\landsat_data160710_landsat8qa_test\cloudmask.TIF')
所以,这里要吸取的教训是,显然 GDAL 不会告诉您交换 rows/cols 是否有问题。它不会定义 numpy 数组...