问题描述

现需要对11景某一个区域的NDVI数据进行出图，且需要使用同一个拉伸的色带，但是拉伸色带的间断值是根据影像的直方图确定的，意味着11景影像会有11个不同的拉伸色带，不符合需求。

解决方法

目前想到的解决方法就是将11景影像通过改变地理位置的方法，排列为一个4*3的格网，再将11景影像镶嵌为一张影像。原理可以参考下面这张图片：
在这里插入图片描述

代码如下


import gdal
import os,glob

file = r'G:\temp\0817\hjx\1'
outfile = r'G:\temp\0817\hjx\2'
TifFile = glob.glob(file + '//*.tif')
index = 1
for tif in TifFile:
    print(tif)
    raster = gdal.Open(tif)
    x = raster.RasterYSize
    y = raster.RasterXSize
    data = raster.GetRasterBand(1).ReadAsArray()
    geotransform = raster.GetGeoTransform()
    proj = raster.GetProjection()
    # 通过index计算新的geotransform
    if index < 4:
        print('index:',index)
        new_x = geotransform[0]+geotransform[1]*y*index
        newgeotransform = (new_x, geotransform[1], 0.0, geotransform[3], 0.0, geotransform[5])
        print('-------0--------')
        print(new_x, geotransform[3])
    elif 4 <= index < 8:
        new_x = geotransform[0] + geotransform[1] * y * (index - 4)
        new_y = geotransform[3] + geotransform[5] * x
        newgeotransform = (new_x, geotransform[1], 0.0, new_y, 0.0, geotransform[5])
        print('-------4-------')
        print(new_x,new_y)
    else:
        new_x = geotransform[0] + geotransform[1] * y * (index - 8)
        new_y = geotransform[3] + geotransform[5] * x * 2
        newgeotransform = (new_x, geotransform[1], 0.0, new_y, 0.0, geotransform[5])
        print('-------8-------')
        print(new_x,new_y)
	# 影像输出
    driver = gdal.GetDriverByName('GTiff')
    out_name = outfile + os.sep + os.path.basename(tif)
    out_ds = driver.Create(out_name,y,x,1,gdal.GDT_Float32)
    out_ds.GetRasterBand(1).WriteArray(data)
    # out_ds.GetRasterBand(2).WriteArray(landdata)
    out_ds.SetGeoTransform(newgeotransform)
    out_ds.SetProjection(proj)
    out_ds.FlushCache()
    index+=1