gdal python 图像金字塔,针对tif,tiff（生成图像金字塔，展示各层级图像）

已于 2022-09-30 10:56:02 修改
阅读量1.4k 收藏 9
点赞数 1
分类专栏： 图像处理GDAL 文章标签： python pytorch 深度学习 图像处理
于 2021-09-19 15:45:30 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/m0_46308537/article/details/120381150
版权

创建图像金字塔

原创不易
导包

import os
from osgeo import gdal
import cv2
import shutil
import skimage.io as io
from PIL import Image
Image.MAX_IMAGE_PIXELS = None

判断文件夹是否存在,存在就删除清空，不存在保留

if os.path.exists(self.dir_path):
       	shutil.rmtree(self.dir_path)
        os.makedirs(self.dir_path)
 else:
     	os.makedirs(self.dir_path)

检查输入和输出路径


class Image_pyramid:
    def __init__(self,inputfpath):
        self.inputtfpath=inputfpath//输入文件路劲
        # print(self.inputtfpath)
        dir_list = self.inputtfpath.strip(' ').split('.')[0].split("/")[0:-1]
        print(dir_list)
        dir_list_pic=self.inputtfpath.strip(' ').split('.')[0].split("/")[-1]
        dir_path2=""
        for i in range(0,len(dir_list)):
            if len(dir_list[i])==0:
                continue
            else:
                dir_list2=str(os.path.join('/',str(dir_list[i])))
                dir_path2=str(dir_path2)+dir_list2
        self.dir_path=os.path.join(dir_path2,dir_list_pic)
        print(self.dir_path)
        if os.path.exists(self.dir_path):
        	shutil.rmtree(self.dir_path)
             os.makedirs(self.dir_path)
    	else:
        	os.makedirs(self.dir_path)
    print(self.dir_path)
    self.outputtfpath = os.path.join(dir_path2,"10")
    print(self.outputtfpath)

生成新的tif图像的名称10.tif

def gdal_pymaid2(self):
      if os.path.exists(self.outputtfpath):
          shutil.rmtree(self.outputtfpath)
          os.makedirs(self.outputtfpath)
      else:
          os.makedirs(self.outputtfpath)
      outputfilename = os.path.join(self.outputtfpath, "10.tif")
      print(333333333333333)
      print(outputfilename)

获取tif图像的宽度和高度


 inputDataset = gdal.Open(self.inputtfpath)
 print(inputDataset==None)
  inputBand = inputDataset.GetRasterBand(1)
  img_arraytif = inputDataset.ReadAsArray(0, 0,
                                     inputDataset.RasterXSize, inputDataset.RasterYSize)

获取原始卫星图像获取仿射矩阵信息，获取投影信息，创建tiff图像驱动，获取波段数量

im_geotrans = inputDataset.GetGeoTransform()  # 获取仿射矩阵信息
im_proj = inputDataset.GetProjection()  # 获取投影信息
driver = gdal.GetDriverByName("GTiff")，#创建tiff图像驱动
im_bands=int(inputDataset.RasterCount)#获取波段数量
print(im_bands)

将原始卫星图像的矩阵信息宽度信息，高度信息，反射变化矩阵，投影信息，等写入到新的tif图像中

     if im_bands == 1:
            datasets = driver.Create(outputfilename,int((inputDataset.RasterXSize)),
            int((inputDataset.RasterYSize)), 1, inputBand.DataType)
            datasets.SetGeoTransform(im_geotrans)  # 写入仿射变换参数
            datasets.SetProjection(im_proj)  # 写入投影
      else:
            datasets = driver.Create(outputfilename, int((inputDataset.RasterXSize)), int((inputDataset.RasterYSize)), 3, inputBand.DataType)
            datasets.SetGeoTransform(im_geotrans)  # 写入仿射变换参数
            datasets.SetProjection(im_proj)  # 写入投影

卫星图像可能有一个波段，3个波段，甚至7个波段的信息，将每个波段的信息写入到对应的tif图像中，在这之前需要区分RGB通道，可以利用均值直方图信息，这里将图像一次缩小4倍，16倍，64倍…

 if im_bands >= 3:
            # for i in range(3):
            datasets.GetRasterBand(1).WriteArray(img_arraytif[2])
            datasets.GetRasterBand(2).WriteArray(img_arraytif[1])
            datasets.GetRasterBand(3).WriteArray(img_arraytif[0])
            for i in range(1, 4):
                datasets.GetRasterBand(i).ComputeStatistics(False)
                datasets.BuildOverviews('average', [1,2, 4, 8, 16, 32, 64])
 else:
            datasets.GetRasterBand(1).WriteArray(img_arraytif)
            datasets.GetRasterBand(1).ComputeStatistics(False)
            datasets.BuildOverviews('average', [1,2, 4, 8, 16, 32,64])
 del datasets

cv写入不同层数的图像金字塔所对应的图像

def show_pymaid(self):
    self.gdal_pymaid2()

    outputfilename = os.path.join(self.outputtfpath, "10.tif")

读取新创建的卫星图像基本信息，获取概率图深度，即金字塔层数和波段数量

 dataset = gdal.Open(outputfilename)
 # print(dataset)

 band = dataset.GetRasterBand(1)
 overviewNum = band.GetOverviewCount()

 print(dataset.RasterCount)
 print(overviewNum)
 im_bands = int(dataset.RasterCount)

不同波段的卫星图像，图像金字塔所对应的原始图像大小

if im_bands>1:
	rband = dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(6)
	r = roverviewBand.ReadAsArray()
	gband = dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(6)
	g = goverviewBand.ReadAsArray()
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(6)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([r, g, b])#正解
	# data = cv2.merge([b, g, r])
	print(data.shape)
	cv2.imwrite(f'{self.dir_path}/1.bmp', data)
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(5)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(5)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(5)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	print(data.shape)
	cv2.imwrite(f'{self.dir_path}/1.bmp',data)
	del data
	
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(4)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(4)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(4)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	print(data.shape)
	
	cv2.imwrite(f'{self.dir_path}/2.bmp',data)
	
	del data
	
	
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(3)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(3)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(3)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	print(data.shape)
	
	cv2.imwrite(f'{self.dir_path}/3.bmp',data)
	
	del data
	#
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(2)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(2)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(2)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	
	cv2.imwrite(f'{self.dir_path}/4.bmp',data)
	
	del data
	#
	#
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(1)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(1)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(1)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	print(data.shape)
	
	cv2.imwrite(f'{self.dir_path}/5.bmp',data)
	
	del data
	
	
	rband=dataset.GetRasterBand(3)
	roverviewBand = rband.GetOverview(0)
	r = roverviewBand.ReadAsArray()
	gband=dataset.GetRasterBand(2)
	goverviewBand = gband.GetOverview(0)
	g = goverviewBand.ReadAsArray()
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(0)
	b = broverviewBand.ReadAsArray()
	import numpy as np
	data=cv2.merge([r,g,b])
	# data = cv2.merge([b, g, r])
	print(data.shape)
	# io.imsave(f'{self.dir_path}/6.bmp',data)
	cv2.imwrite(f'{self.dir_path}/6.bmp',data)
	
	
	del data
else:
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(6)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	cv2.imwrite(f'{self.dir_path}/1.bmp', data)
	
	del data
	
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(5)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	
	cv2.imwrite(f'{self.dir_path}/1.bmp',data)
	
	del data
	
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(4)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	print(data.shape)
	
	cv2.imwrite(f'{self.dir_path}/2.bmp',data)
	
	del data
	
	
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(3)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	
	cv2.imwrite(f'{self.dir_path}/3.bmp',data)
	
	del data
	
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(2)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	
	cv2.imwrite(f'{self.dir_path}/4.bmp',data)
	
	del data
	
	bband = dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(1)
	b = broverviewBand.ReadAsArray()
	data = cv2.merge([b, b, b])
	
	cv2.imwrite(f'{self.dir_path}/5.bmp',data)
	
	del data
	
	bband=dataset.GetRasterBand(1)
	broverviewBand = bband.GetOverview(0)
	b = broverviewBand.ReadAsArray()
	data=cv2.merge([b,b,b])
	print(data.shape)
	
	cv2.imwrite(f'{self.dir_path}/6.bmp',data)
	
	# del data
m0_46308537
关注
  • 1
    点赞
  • 踩
  • 9
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
geotif、tif图像金字塔处理
09-14
利用C#编的关于geotif、tif图像利用GDAL处理成金字塔结构，代码可以直接运行。需安装GDAL
利用GDAL把TIF文件转换为bmp,png等格式图片并对图片进行压缩
12-15
资源声明：部分资源来自于互联网 代码利用GDAL转TIF文件为png、bmp等格式图片并根据屏幕分辨率对图片进行压缩。解决了因tif文件过大造成内存溢出等问题，欢迎相互学习
python+gdal+遥感图像拼接(mosaic)的实例
09-17
主要介绍了python+gdal+遥感图像拼接(mosaic)的实例，具有很好的参考价值，希望对大家有所帮助。一起跟随小编过来看看吧
gdal_opengl.rar_GDAL_gdal tiff_opengl gdal_tiff_读取tiff图像
09-20
同时使用gdal和opengl实现tiff图像的读取和现实
gdal裁剪tif_tif遥感图像_GDAL裁剪_遥感影像_GDAL_gdal批量裁剪_
10-01
基于gdal实现批量裁剪tif图像，对于遥感影像同样适用。
构建影像金字塔，成果无压缩TIF
03-02
读取影像，jpg，tif等格式，生成影像金字塔，并且以外部金字塔的方式写到硬盘。 包含opencv影像读取，降采样等技术； 包含GDAL数据集打开关闭，读写等技术。 可编译成exe可执行文件。
使用gdal对图像建立金字塔
liuzf1218的专栏
12-04 2410
对图像建立金字塔，参考了李民录的博客，给出了一个可以直接粘贴就能运行的代码 // CreatePyramids.cpp : 定义控制台应用程序的入口点。 #define _AFXDLL #pragma #include "stdafx.h" #include "gdal.h" #include "gdal_priv.h" #include "afx.h" #pragma comment(lib,
高分二号多光谱GEOTIFF影像批量转换为JPG图片
weixin_43490758的博客
03-13 1781
最近在做语义分割方面的研究，想对比一下高分二号4个波段的模型输入和JPG输入的分割结果有什么差别。使用了一些图片转换软件发现转换后的jpg图片一片漆黑，应该是geitif格式的特殊性吧！所有在网上找了一些代码，python调用gdal库，现在总结一下，分享给大家。
python/gdal处理遥感影像（读取、投影转换、裁剪、建立图像金字塔等）
热门推荐
Smilelulinlin的博客
06-25 1万+
python/gdal处理遥感影像（读取、投影转换、裁剪、建立图像金字塔等）gdal库简单介绍新的改变功能快捷键合理的创建标题，有助于目录的生成如何改变文本的样式插入链接与图片如何插入一段漂亮的代码片生成一个适合你的列表创建一个表格设定内容居中、居左、居右SmartyPants创建一个自定义列表如何创建一个注脚注释也是必不可少的KaTeX数学公式新的甘特图功能，丰富你的文章UML 图表FLowchart流程图导出与导入导出导入 gdal库简单介绍 gdal的全称含义是地理空间数据抽象库，是一个在X/MIT许
GDAL读取图像文件，建立金字塔
zou_ys88的博客
02-21 1294
GDAL是遥感数据的IO处理库，可以读写erdas img、GeoTiff等多种文件格式。本文通过建立金字塔读取大数据量的文件（测试数据1.45GB，文件格式 img）。 #include "gdal_priv.h" #include "math.h" unsigned char *uc_gray;//灰度数组 int iwidth;//图像宽度 int iheight;//图像高度 int CImageObject::GetSmpStep(fl...
Python空间数据处理之GDAL读写遥感图像
09-18
主要介绍了Python空间数据处理之GDAL读写遥感图像，文中通过示例代码介绍的非常详细，对大家的学习或者工作具有一定的参考学习价值，需要的朋友们下面随着小编来一起学习学习吧
利用python GDAL库读写geotiff格式的遥感影像方法
09-19
今天小编就为大家分享一篇利用python GDAL库读写geotiff格式的遥感影像方法，具有很好的参考价值，希望对大家有所帮助。一起跟随小编过来看看吧
tiff影像读取（主流库的对比，波段数等）python
大草坪
04-09 5623
一般来讲，读取时按照正常图片读取就可以，但是实际使用过程中，波段顺序和数量会因为不同的库表现不同。 涉及的库：numpy、cv2、matplot、tifffile
发布GeoTIFF时的Tile Caching属性和切片zoom：地图显示层级与切片
超频化石鱼的博客
10-07 2758
发布GeoTIFF图层时，查看Tile Caching： 最下方列出了所有可用的gridsets： 坐标系统gridsets 上面有两个gridsets，表示该地图将使用这两种坐标系统。下方有个Add grid subset的选项，可以添加更多的坐标系统。 一般来说，若非常肯定该地图将要使用的坐标系统，那么这里只需要留一个gridsets即可。因为后面切图的时候，会根据gridset...
转载：GeoServer发布多层级TIF地图大图数据
java_YHY的博客
08-25 2496
原文转载自：https://blog.csdn.net/gjysk/article/details/104442731 转载仅为方便自身查阅学习，请支持原创作者！ 一、概述 通常情况下，我们所拿到的数据是地图各层级构建好金字塔后的tif文件，每个层级对应一份tif文件，这篇文章主要是要解决多层级地图数据发布的问题。文章中以天地图的数据为例进行演示说明。 二、配置比例尺和坐标(Gridsets) 首先我们需要在Gridsets里添加下载的天地图各层级所对应的比例尺和坐标，按照如下的步骤进行添加。 1、如下图所
杂技01 Gdal 快速实现tif影像数据ovr金字塔文件建立方法
qq_24914861的博客
04-27 2269
Gdal 快速实现tif影像数据ovr金字塔文件建立方法。
如何将tiff转成金字塔格式tiff
cxsunn的博客
03-28 2193
因课题需要处理tiff格式的病理切片文件，在用openslide读取时发现tiff只有一个最大分辨率下的level dimension，不方便后续处理。 因此需要将tiff改成openslide支持的generic tiled tiff，具体步骤如下： 根据官方引导下载vips，博主用的MacOS，所以直接brew install vips就可以了嘿嘿 下载成功后在命令行输入 vips tiffsave input.tif output.tif —compression jpeg —tile ..
GDAL对TIF创建内建金字塔一个问题
weixin_34357928的博客
07-15 1187
gdalwarp输出tif图像的时候，默认如果没有使用BIGTIFF=YES选项，则会根据输出影像的大小进行判断，低于4G则不适用bigtiff格式。 对于非bigtiff图像，如果这时候使用gdaladdo创建金字塔（内建模式），如果出现了文件总大小超过4G的情况，则会导致图像文件出错。 gdaladdo会输出相关信息，但不会终止金字塔的创建，于是导致金字塔是黑色的。 C:\Program Fi...
【Python】Python中变量或者函数加下划线‘_‘的含义
最新发布
yangweipeng708的博客
05-14 670
理解这些约定可以帮助Python开发者编写更清晰、更Pythonic的代码，并理解其他Python代码中的命名习惯。）通常用来指示特定的编程约定或用法，这些约定主要关乎变量和函数的可访问性以及预期的使用范围。在Python中，变量或函数前面加一个下划线（
gdal 计算图像均值 python
08-09
在Python中使用GDAL计算图像均值的方法如下： 首先，需要导入必要的库文件。我们需要使用GDAL库来读取图像数据，以及NumPy库来处理数据。 ```python import gdal import numpy as np ``` 然后，我们需要打开图像文件，并读取其数据。 ```python dataset = gdal.Open("path/to/image.tif") # 替换为你的图像文件路径 band = dataset.GetRasterBand(1) # 获取第一个波段的数据 image_data = band.ReadAsArray() # 将波段数据读取为数组 ``` 接下来，我们可以使用NumPy库中的函数来计算图像的均值。 ```python mean_value = np.mean(image_data) ``` 最后，我们可以打印出计算得到的均值。 ```python print("图像均值为: ", mean_value) ``` 完整的代码如下： ```python import gdal import numpy as np dataset = gdal.Open("path/to/image.tif") # 替换为你的图像文件路径 band = dataset.GetRasterBand(1) # 获取第一个波段的数据 image_data = band.ReadAsArray() # 将波段数据读取为数组 mean_value = np.mean(image_data) print("图像均值为: ", mean_value) ``` 请注意，这只是计算了图像中的一个波段的均值。如果你想要计算每个波段的均值，可以使用循环来处理每个波段的数据。

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值