python读取tif格式的遥感影像(灰度图)并用迭代法阈值选择算法进行二值化处理

在python中,一般的图像处理库或者绘图库无法处理或显示遥感影像,例如Matplotlib,opencv,scipy等等,用这些库读取遥感影像数值会发生一定的错误,从而导致后续工作无法展开。因此我们需要用到GDAL这个库。

GDAL的安装

GDAL轮子安装网址
在这里插入图片描述

找到与自己python相对应的版本即可
轮子安装完成后,在终端cd到轮子的安装位置,输入pip install “文件名”
这样就完成GDAL的安装了

批量读取遥感影像

import os
from osgeo import gdal

filedir = "D:/img"		# 图片所在路径
for image in os.listdir(filedir):
	driver = gdal.GetDriverByName("GTiff")
	driver.Register()
	img = gdal.Open(filedir + '/' + image)
	
	im_width = img.RasterXSize	# 栅格矩阵的列数
	print("im_width:", im_width)
	
	im_height = img.RasterYSize	# 栅格矩阵的行数
	print("im_height:", im_height)
	
	band = img.GetRasterBand(1)
	im_data = img.ReadAsArray(0, 0, im_width, im_height)	# 获取数据
	print(im_data)

迭代法阈值选择算法

迭代法阈值选择算法是对双峰法的改进,首先选择一个近似的阈值T,将图像分割成两份,R1和R2,分别计算出R1和R2的均值u1和u2,再选择新的阈值T=(u1+u2)/2,重复以上步骤,直到u1和u2不再变化。

迭代法是基于逼近的思想,其步骤如下:
(1)求出图像的最大灰度值和最小灰度值,分别记为ZMAX和ZMIN,令初始阈值为T=(ZMAX+ZMIN)/2
(2)根据阈值T将图像分割为前景和背景,分别求出两者的平均灰度值ZO和ZB
(3)求出新阈值T=(ZO+ZB)/2
(4)若两个平均灰度值ZO和ZB不再发生变化(或T不再变化),则T为阈值,否则转(2)迭代计算

完整代码

import cv2
import os
import tiff
import numpy as np
import pandas as pd
from osgeo import gdal

filedir = "D:/image"		# 遥感图片所在路径
savedir = "D:/image_save"	# 保存路径

for image in os.listdir(filedir):
	driver = gdal.GetDriverByName("GTiff")
	driver.Register()
	img = gdal.Open(filedir + '/' + image)
	
	im_width = img.RasterXSize	# 栅格矩阵的列数
	print("im_width:", im_width)
	
	im_height = img.RasterYSize	# 栅格矩阵的行数
	print("im_height:", im_height)
	
	band = img.GetRasterBand(1)
	im_data = img.ReadAsArray(0, 0, im_width, im_height)	# 获取数据
	print(im_data)
	
	img = np.array(im_data)
	ZMAX = 0; ZMIN = 255
	ZO = 0; ZB = 255
	
	for i in img:
		if max(i) > ZMAX:
			ZMAX = max(i)
		if min(i) < ZMIN:
			ZMIN = min(i)
			
	T = (ZMAX + ZMIN) / 2
	img = pd.DataFrame(img)
	
	while True:
		ZO_ = img[img>T].mean().mean()
		ZB_ = img[img<=T].mean().mean()
		if (ZO_ == ZO) and (ZB_ == ZB):
			break
		else:
			ZO = ZO_; ZB = ZB_
			T = (ZO + ZB) / 2
		print(T)
		
	img[img>T] = 255
	img[img<=T] = 0
	tiff.imsave(savedir + '/' + image, np.array(img))	# 保存二值化后的图像
图像的二值化(需要是灰度图像) import cv2 as cv import numpy as np def threshold_demo(image): gray = cv.cvtColor(image,cv.COLOR_BGR2GRAY) #cv.THRESH_TRIANGLE与cv.THRESH_OTSU是自动寻找阈值,这个时候threshold值必须是0 #如果threshold不为0,即是自己指定的(如127),那么type就直接THRESH_BINARY就够了 ret,binary = cv.threshold(gray,0,255,cv.THRESH_BINARY|cv.THRESH_OTSU) print(ret) cv.imshow("binary",binary) def local_threshold(image): gray = cv.cvtColor(image, cv.COLOR_BGR2GRAY) binary =cv.adaptiveThreshold(gray,255,cv.ADAPTIVE_THRESH_GAUSSIAN_C,cv.THRESH_BINARY,25,10) cv.imshow("binary1", binary) src =cv.imread("E:/opencv/picture/dog.jpg") cv.imshow("initial_window",src) threshold_demo(src) local_threshold(src) cv.waitKey(0) cv.destroyAllWindows() 分析: 1. 全局阈值 def threshold_demo(image): gray = cv.cvtColor(image,cv.COLOR_BGR2GRAY) #cv.THRESH_TRIANGLE与cv.THRESH_OTSU是自动寻找阈值,这个时候threshold值必须是0 #如果threshold不为0,即是自己指定的(如127),那么type就直接THRESH_BINARY就够了 ret,binary = cv.threshold(gray,0,255,cv.THRESH_BINARY|cv.THRESH_OTSU) print(ret) cv.imshow("binary",binary) threshold其函数原型为:threshold(src, thresh, maxval, type[, dst]) -> retval, dst src参数表示输入图像(多通道,8位或32位浮点)。 thresh参数表示阈值。(如果type是自动二值化就设置为0) maxval参数表示与THRESH_BINARY和THRESH_BINARY_INV阈值类型一起使用设置的最大值。 type参数表示阈值类型。 retval参数表示返回的阈值。若是全局固定阈值算法,则返回thresh参数值。若是全局自适应阈值算法,则返回自适应计算得出的合适阈值。 dst参数表示输出与src相同大小和类型以及相同通道数的图像。 type阈值类型说明 cv.THRESH_BINARY | cv.THRESH_OTSU)#大律法,全局自适应阈值 参数0可改为任意数字但不起作用 cv.THRESH_BINARY | cv.THRESH_TRIANGLE)#TRIANGLE法,,全局自适应阈值, 参数0可改为任意数字但不起作用,适用于单个波峰 cv.THRESH_BINARY)# 自定义阈值为127,大于127的是白色 小于的是黑色 cv.THRESH_BINARY_INV)# 自定义阈值为127,大于127的是黑色 小于的是白色 cv.THRESH_TRUNC)# 截断 大于127的是改为127 小于127的保留 cv.THRESH_TOZERO)# 截断 小于127的是改为127 大于127的保留 2. 局部阈值法 3. def local_threshold(image): gray = cv.cvtColor(image, cv.COLOR_BGR2GRAY) binary =cv.adaptiveThreshold(gray,255,cv.ADAPTIVE_THRESH_GAUSSIAN_C,cv.THRESH_BINARY,25,10) cv.imshow("binary1", binary) adaptiveThreshold函数进行局部阈值 函数原型为:adaptiveThreshold(src, maxValue, adaptiveMethod, thresholdType, blockSize, C[, dst]) -> dst src参数表示输入图像(8位单通道图像)。 maxValue参数表示使用 THRESH_BINARY 和 THRESH_BINARY_INV 的最大值. adaptiveMethod参数表示自适应阈值算法, 平均— cv2.ADAPTIVE_THRESH_MEAN_C :领域内均值 高斯—cv2.ADAPTIVE_THRESH_GAUSSIAN_C :领域内像素点加权和,权 重为一个高斯窗口 thresholdType参数表示阈值类型,必须为THRESH_BINARY或THRESH_BINARY_INV的阈值类型。 blockSize参数表示块大小,规定领域大小(奇数且大于1,比如3,5,7........ )。 C参数是常数,表示从平均值或加权平均值中减去的数。 通常情况下,这是正值,但也可能为零或负值。 二.对超大图像进行二值化 如果这个时候只是单纯的用二值化api,图像上会出现很多噪声,所以我们特地介绍了一种对大图像进行二值化的方法import cv2 as cv import numpy as np from matplotlib import pyplot as plt #对超大图像进行二值化 def big_image_threshold(image): cw = 256 ch = 256 h,w = image.shape[:2] gray = cv.cvtColor(image,cv.COLOR_BGR2GRAY) for row in range(0, h, ch): for col in range(0, w, cw): #gray[0:3] 从索引0开始取,直到索引3为止,但不包括索引3。即索引0,1,2,正好是3个元素。 roi = gray[row:row+ch,col:col+ch]#表示从[row,col]到[row+ch,col+ch]的所有元素所组成的矩阵 #推荐使用局部阈值二值化 dst = cv.adaptiveThreshold(roi,255,cv.ADAPTIVE_THRESH_GAUSSIAN_C,cv.THRESH_BINARY,127,20) gray[row:row+ch,col:col+ch] =dst cv.imwrite("E:/opencv/picture/bigpicture1.jpg",gray) src = cv.imread("E:/opencv/picture/bigpicture.jpg") t1 = cv.getTickCount() big_image_threshold(src) t2 = cv.getTickCount() time= (t2-t1)/cv.getTickFrequency() #print("time =%sms\n"%(time)*1000) cv.waitKey(0) cv.destroyAllWindows()
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

余磬

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

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

余额充值