Python遥感开发之快速判断TIF数据为空
前言:介绍一下如何使用python下的gdal读取tif数据的时候,快速判断该tif数据是否为空,如果为空的话就把当前的tif删掉。
如图所示,通过arcgis查看箭头指向的为空值。
仅通过文件的大小无法判断
实现思路:我们可以借助python下的gdal读取tif,以第一个像元为基准像元,如果后面的像元有不一样的值就认定该tif不为空。
代码实现:
import numpy as np
from osgeo import gdal
import os
def read_tif01(filepath):
dataset = gdal.Open(filepath)
col = dataset.RasterXSize # 图像长度
row = dataset.RasterYSize # 图像宽度
geotrans = dataset.GetGeoTransform() # 读取仿射变换
proj = dataset.GetProjection() # 读取投影
data = dataset.ReadAsArray() # 转为numpy格式
data = data.astype(np.float32) # 转为float类型
# data[(data < 0) | (data > 15000)] = np.nan
# a = data[0][0]
# data[data == a] = np.nan #原因:读取某一个行政区的影像图的时候,往往第一行的第一列值为空值
# band = dataset.GetRasterBand(1)
# print(band)
return [col, row, geotrans, proj, data]
def get_data_list(file_path, out=""):
list1 = [] # 文件的完整路径
if os.path.isdir(file_path):
fileList = os.listdir(file_path)
if out != "":
for f in fileList:
out_data = out + "\\" + f
out_data = out_data.replace(".HDF", "_ndvi.tif")
list1.append(out_data)
else:
for f in fileList:
pre_data = file_path + '\\' + f # 文件的完整路径
list1.append(pre_data)
return list1
if __name__ == '__main__':
inputfile_path = r"F:\风云数据\MERSI-II陆表反射比1KM段产品\b19\01原始"
infile_list = get_data_list(inputfile_path)
for file_path in infile_list:
col, row, geotrans, proj, data = read_tif01(file_path)
print(data[0, 0])
are_all_elements_equal = np.all(data == data[0, 0])
if are_all_elements_equal:
print("所有元素都相同")
os.remove(file_path)
else:
print("数组中有不同的元素")
print("已去除空数据")