import osgeo
from osgeo import ogr
from osgeo import gdal
from osgeo import osr
import os
import numpy as np
gdal.AllRegister()
osgeo.gdal.SetConfigOption('GDAL_FILENAME_IS_UTF8', 'NO')
osgeo.gdal.SetConfigOption('SHAPE_ENCODING', 'gb2312')
filePath = r'H:\task\mask_14.tif'#365个波段,25列,26行的tif
dataset = gdal.Open(filePath)
adfGeoTransform = dataset.GetGeoTransform()
im_proj = dataset.GetProjection()
nXSize = dataset.RasterXSize #列数
nYSize = dataset.RasterYSize #行数
im_bands = dataset.RasterCount #波段数
#print(im_bands)
data = dataset.ReadAsArray(0, 0, nXSize, nYSize)
#print(data)
##########以下操作属博主瞎搞,可略过
result=np.maximum(data,0)#判断小于0的像素值都赋值为0
print(result.shape)
test=result[10:80,0:24,0:25]
dataarray=np.mean(test,axis=0)对0轴方向取均值
print(dataarray.shape)
##########
driver = gdal.GetDriverByName("GTiff")
##一开始没有仔细看数据类型(gdal.GDT_Float32),导致出错,望下次注意
datasetnew = driver.Create( r'H:\task\new_80.tif', nXSize, nYSize, 1, gdal.GDT_Float32)
datasetnew.SetGeoTransform(adfGeoTransform)
datasetnew.SetProjection(im_proj)
band = datasetnew.GetRasterBand(1)
band.WriteArray(dataarray)
练习代码,勿喷,个人学习使用