Python 循环读取影像 计算NDVI

from osgeo import gdal
import numpy as np
from osgeo import gdal
import os
#循环读取
FilePath='D:\\影像练习'
tif_name=os.listdir(FilePath) #读取FilePath下的目录名,包含扩展名
SavePath='D:\\影像练习\\NDVI输出'
print(tif_name)

for i in range(len(tif_name)):  #i从0开始到len(tif_name)-1
    qianzui=os.path.splitext(tif_name[i])[0]  #将文件名和扩展名分隔开
    houzui=os.path.splitext(tif_name[i])[1]
    if houzui=='.tif':
        ReadName=FilePath+'\\'+qianzui+houzui
        SaveName=SavePath+'\\'+qianzui+'_NDVI.tif'
        dataset = gdal.Open(ReadName)
        cols = dataset.RasterXSize  # 列数
        rows = dataset.RasterYSize  # 行数  空间坐标 可以理解为行列号
        # 计算 NDVI=（NIR-RED）/（NIR+RED）
        # 提取红波段
        band_red = dataset.GetRasterBand(3)
        data_red = band_red.ReadAsArray(0, 0, cols, rows)  # 根据位置获取像素值
        # 提取近红外波段
        band_nir = dataset.GetRasterBand(4)
        data_nir = band_nir.ReadAsArray(0, 0, cols, rows)
        ndvi = (data_nir - data_red) / (data_nir + data_red)

        #创建栅格数据集
        driver = gdal.GetDriverByName('GTiff')  # 将GTiff图纸实体化为driver，括号内为文件类型
        dataset_out = driver.Create(SaveName, cols, rows, 1, gdal.GDT_Float32)  # 这里的1指的是创建的栅格数据文件只能容纳一个波段
        # 定义仿射矩阵和投影信息
        dataset_out.SetGeoTransform(dataset.GetGeoTransform())
        dataset_out.SetProjection(dataset.GetProjection())
        # 写入像素值
        band_out = dataset_out.GetRasterBand(1)
        band_out.WriteArray(ndvi)  # 把数组写入到栅格数据相应的波段中
        # 清空缓存，关闭数据集
        dataset.FlushCache()
        dataset_out.FlushCache()  # 把内存中的数据写入到文件中
        del dataset
        del dataset_out