# coding UTF-8
# author:huangZengli
from osgeo import gdal
from pylab import * # 支持中文
import os
import pandas as pd
mpl.rcParams['font.sans-serif'] = ['SimHei']
#遍历指定类型文件
def getFilePathList(path, filetype):
pathList = []
for root, dirs, files in os.walk(path):
for file in files:
if file.endswith(filetype):
pathList.append(os.path.join(root, file))
return sorted(pathList) # 输出以filetype为后缀的列表,已按名称排序,注意这个地方可能因为命名方式出错,可以输出list看下
if __name__ == "__main__":
for year in range(2001,2021):
filePath = '/Users/sunergao/Desktop/VM/NDVI2/geo/NDVIresample/' + str(year) # tif文件路径
filePathList = getFilePathList(filePath,'tif')
for i in range(0,12):
dataset = gdal.Open(filePathList[i]) # 打开tif
# 获取行数列数和地理信息
# geo_information(0):左上像素左上角的x坐标。
# geo_information(1):w - e像素分辨率 / 像素宽度。
# geo_information(2):行旋转(通常为零)。
# geo_information(3):左上像素左上角的y坐标。
# geo_information(4):列旋转(通常为零)。
# geo_information(5):n - s像素分辨率 / 像素高度(北半球上图像为负值)
geo_information = dataset.GetGeoTransform()
col = dataset.RasterXSize # 438
row = dataset.RasterYSize # 671
band = dataset.RasterCount
dem = dataset.GetRasterBand(1).ReadAsArray()
# 获取行列数,对应其经纬度,j对于x坐标
lon = []
lat = []
ndvi = []
for y in range(row): # 行
for x in range(col): # 列
# 有效高程
# if dem[y][x] > 0:
# 输出经纬度
lon.append(geo_information[0] + x * geo_information[1] + y * geo_information[2])
lat.append(geo_information[3] + x * geo_information[4] + y * geo_information[5])
ndvi.append(dem[y][x])
outfilePath = '/Users/sunergao/Desktop/VM/NDVI2/geo/NDVItable3/' + str(year) + '/' + str(i+1) + '_NDVI.csv'
name = ['lon','lat','NDVI']
data = pd.DataFrame({name[0]:lon,name[1]:lat,name[2]:ndvi})
save_name = outfilePath
data.to_csv(save_name,index =False ,sep = ',')
print(year,'-',i+1,'NDVI表已经生成')
原文:使用python读取tiff文件中的经纬度,并将数据以excel表的形式输出(详细步骤)_迷沉的博客-CSDN博客_tif的经纬度坐标