使用python批处理读取tiff文件中的经纬度和值,并将数据以excel表的形式输出(详细步骤)
目标:把文件夹内的tif文件中的经纬度和值提取并写入excel中
为了实现这一目标,分为以下几步:
1、提取一个tif内栅格的经纬度和值,将其导入表中
2、提取一个文件内所有tif的值,将其导入表中,
3、提取多个文件内所有tif的值,将其导入一个excel的多个sheet(一个文件夹对应一个sheet),或者将其导入一个表
以下图片是文件的存放路径:
co2_tif文件夹内有三个文件夹来存放tif,各个文件内有一个月每天四个时间段,分别为0点,6点、12点、18点的co2浓度的tif文件。
1、实现第一步(提取一个tif内栅格的经纬度和值,将其导入表中)使用的代码
# coding UTF-8
from osgeo import gdal
from pylab import * # 支持中文
mpl.rcParams['font.sans-serif'] = ['SimHei']
from openpyxl import Workbook
# 创建一个Workbook对象
work = Workbook()
def out(data, name):
ws = work.active
ws['A1'] = '经度'
ws['B1'] = '纬度'
ws['C1'] = '二氧化碳'
ws['D1'] = '所在栅格行'
ws['E1'] = '所在栅格列'
for i in range(len(data)):
rows = []
row_length = len(data[i])
if row_length != 0:
for j in range(row_length):
rows.append(data[i][j])
ws.append(rows[j])
work.save(name)
if __name__ == "__main__":
filePath = 'F:/data/co2_tif/201301_975/c_cty_Band_1.tif' # tif文件路径
dataset = gdal.Open(filePath) # 打开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 # 25
row = dataset.RasterYSize # 20
band = dataset.RasterCount #波段数,该tif文件中波段为 1
conc = dataset.GetRasterBand(1).ReadAsArray()
# 获取行列数,对应其经纬度,j对于x坐标
co2 = []
for y in range(row): # 行
rows = []
for x in range(col): # 列
# 有效co2值
if conc[y][x] > 0:
# 输出经纬度
lon = geo_information[0] + x * geo_information[1] + y * geo_information[2]
lat = geo_information[3] + x * geo_information[4] + y * geo_information[5]
child = [lon, lat, conc[y][x], y, x]#每个栅格的经纬度、值、所在行、列 存放在一个列表内
rows.append(child)
co2.append(rows)
out(co2, 'F:/data/co2_tif/test11.xlsx')
print('表已经生成')
这段代码来源:
使用python读取tiff文件中的经纬度,并将数据以excel表的形式输出(详细步骤)_迷沉的博客-CSDN博客_tif的经纬度坐标
运行结束,所得到的表格
2、实现第二步(提取一个文件内所有tif文件的栅格值,将其导入表中)代码
由于每个栅格的经纬度及其所在行、列都一致,之后的代码就不再提取经纬和行列
import xlwt
from pylab import * # 支持中文
from osgeo import gdal
mpl.rcParams['font.sans-serif'] = ['SimHei']
from openpyxl import Workbook
work_book = xlwt.Workbook(encoding='utf-8')
sheet = work_book.add_sheet(sheetname='co2值')
def out(data, name,row):
for i in range(len(data)):
sheet.write(row,i+4, data[i])#row=写入的值所在的行,i+4=写入值所在列,data是存放所在值的列表
work_book.save(name)
if __name__ == "__main__":
for z in range(1,121):#z代表文件中tif文件名称的波段数,例如:c_cty_Band_1、c_cty_Band_2等等,范围设置(1,121),是因为该文件中有120个命名为该形式tif文件
filePath = r'F:\data\co2_tif\201401_975\c_cty_Band_'+str(z)+'.tif' # tif文件路径
dataset = gdal.Open(filePath) # 打开tif
col = dataset.RasterXSize # 25
row = dataset.RasterYSize # 20
conc = dataset.GetRasterBand(1).ReadAsArray()
co2=[]
for y in range(row): # 行
cols = []
for x in range(col): # 列
if conc[y][x] > 0:
co2.append(float(conc[y][x]))#co2这个列表用来存放所有栅格的值
out(co2,'F:/data/test12.xlsx',int(z)-1)#int(z)-1 被用来控制表格的行,每读取一个栅格图,则行数增加一行
print('表已经生成')
运行完成的结果
3、实现第三步(提取多个文件内所有tif的值,将其导入一个excel的多个sheet,或者将其导入一个表)
先看一下读取文件名称的语句
import os
path='F:/data/co2_tif'
datanames=os.listdir(path)
print(datanames)
txtnumber=len(os.listdir('F:/data/co2_tif/201301_975'))
print(txtnumber)
这是运行的结果 (co2_tif)文件下的三个文件名都被读取进datanames这个列表
os.listdir 读取该文件夹内的文件名
len(os.listdir())读取一个文件夹内存储文件的个数
提取多个文件内所有tif的值,将其导入一个excel的多个sheet(一个文件对应一个sheet)的代码
import os
import xlwt
from pylab import * # 支持中文
from osgeo import gdal
mpl.rcParams['font.sans-serif'] = ['SimHei']
work_book = xlwt.Workbook(encoding='utf-8')
#sheet = work_book.add_sheet(sheetname='co2值')
def out(data, name,col):
for i in range(len(data)):
sheet.write(col+4,i, data[i])
work_book.save(name)
if __name__ == "__main__":
path='F:/data/co2_tif'
datanames=os.listdir(path)#datanames用来存放co2_tif文件内的文件名称
for u in datanames:
file='F:/data/co2_tif/'+u #u=datanames列表内存放的文件名
txtnumber=len(os.listdir('F:/data/co2_tif/'+u)) #txtnumber=文件夹内所存储的文件个数
sheet = work_book.add_sheet(sheetname=u) # 一个文件内的所有tif栅格值,导入一个sheet中(一个文件对应一个sheet),因此,就把第九行创建新sheet的语句放在在了循环里面
for z in range(1,int(txtnumber/5+1)):#txtnumber/5 看图三就会理解。 文件内总的文件数除5=文件内tif文件的个数,z还被用来读取路径,因此加txtnumber/5+1
filePath = 'F:/data/co2_tif/'+u+'/c_cty_Band_' + str(z) + '.tif' # tif文件路径
dataset = gdal.Open(filePath) # 打开tif
col = dataset.RasterXSize # 25
row = dataset.RasterYSize # 20
conc = dataset.GetRasterBand(1).ReadAsArray()
co2=[]
for y in range(row): # 行
cols = []
for x in range(col): # 列
if conc[y][x] > 0:
co2.append(float(conc[y][x]))
out(co2, 'F:/data/test13.xlsx',int(z)-1)
print('表已经生成')
运行完成结果(注意看下面sheet)
提取多个文件内所有tif的值,将其导入一个表的代码
import xlwt
from pylab import * # 支持中文
from osgeo import gdal
import os
mpl.rcParams['font.sans-serif'] = ['SimHei']
work_book = xlwt.Workbook(encoding='utf-8')
sheet = work_book.add_sheet(sheetname='co值')
list=['0点','6点','12点','18点']
def out(data, name,row):
sheet.write(row,0,u[:4])#第一列表头:哪一年
sheet.write(row,1,u[4:6])#第二列表头:月份
sheet.write(row,2,int((z+3)/4))#第三列表头 日
sheet.write(row,3,list[int(row%4)])# 第四列表头 时间点
for i in range(len(data)):
sheet.write(row,i+4, data[i])
work_book.save(name)
work_book = xlwt.Workbook(encoding='utf-8')
sheet = work_book.add_sheet(sheetname='co值')
if __name__ == "__main__":
n=4# n 用来控制行数,第一个tif文件的栅格值写入表格的四行 每循环一次n就+1,另外,n=4主要是为了便于,第四列表头 时间点的添加
path='F:/data/co2_tif'
datanames=os.listdir(path)#datanames用来存放co2_tif文件内的文件名称
for u in datanames:
file='F:/data/co2_tif/'+u #u=datanames列表内存放的文件名
txtnumber=len(os.listdir('F:/data/co2_tif/'+u))#txtnumber=文件夹内所存储的文件个数
for z in range(1,int(txtnumber/5+1)):# 一个文件内的所有tif栅格值,导入一个sheet中(一个文件对应一个sheet),因此,就把第九行创建新sheet的语句放在在了循环里面
filePath = 'F:/data/co2_tif/'+u+'/c_cty_Band_' + str(z) + '.tif' #txtnumber/5 看图三就会理解。 文件内总的文件数除5=文件内tif文件的个数,z还被用来读取文件路径的改变,因此加txtnumber/5+1
dataset = gdal.Open(filePath) # 打开tif
col = dataset.RasterXSize # 栅格的列数
row = dataset.RasterYSize # 栅格的行数
conc= dataset.GetRasterBand(1).ReadAsArray()
co2=[]
for y in range(row): # 行
for x in range(col): # 列
if conc[y][x] > 0:
co2.append(float(conc[y][x]))
out(co2, 'F:/data/test14.xlsx',n)
n=n+1
运行结果