影像拉伸优化

已于 2024-01-08 00:11:22 修改
阅读量439 收藏 8
点赞数 10
分类专栏： 影像处理 文章标签： opencv 计算机视觉 人工智能
于 2024-01-08 00:05:11 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/yangfahe1/article/details/135446971
版权

Qt、GDAL遥感影像显示

影像拉伸优化

使用GDAL读取影像数据之后，需要进行拉伸显示。比如读取landsat数据时，由于数据类型时ushort，如果直接显示，那么显示的结果就是一片空白。因此我们需要对数据进行拉伸，将ushort拉伸到0-255。因此出现以下三种拉伸方法。

方案1（最慢）

	for(int index = 0; index < pImpObject->bandMap.length(); ++index) {
		int bandMap = pImpObject->bandMap[index] + 1;
		cv::Mat tmp = cv::Mat(buf.y, buf.x, GisEntsHelper::GDALType2GCType(dType, 1));
		GDALRasterBand *pBand = pDataset->GetRasterBand(bandMap);
		CPLErr err = pBand->RasterIO(GF_Read, read.x, read.y, size.x, size.y, tmp.data, tmp.cols, tmp.rows, dType, 0, 0);
		GisEntsHelper::GisLineStretchImage(tmp, pImpObject->histArray[index], sType, GisEntsHelper::GDALTypeMaxValue(dType));
		mats.push_back(tmp.clone());
	}

void GisEntsHelper::GisLineStretchImage(cv::Mat &mat, const HistogramData &histData, OGis::StretchType sType, uint maskValue)
{
	cv::Mat mask1;
	cv::compare(mat, histData._noDataValue, mask1, cv::CMP_EQ);
	mat.setTo(maskValue, mask1);
	if(sType != OGis::kStretch_No) {
		switch(sType) {
			case OGis::kStretch_Line:
			case OGis::kStretch_Line1:
			case OGis::kStretch_Line2:
			case OGis::kStretch_Line5:
				cv::threshold(mat, mat, histData.minmax[sType - 1].x, 0, cv::THRESH_TOZERO);
				cv::threshold(mat, mat, histData.minmax[sType - 1].y, 255, cv::THRESH_TRUNC);
				mat = (mat - histData.minmax[sType - 1].x) * histData.minmax[sType - 1].z;
				mat.convertTo(mat, CV_8U);
				break;
			case OGis::kStretch_Equalize:
				cv::normalize(mat, mat, 255, 0, CV_MINMAX, CV_8U);
				cv::equalizeHist(mat, mat); break;
		}
	}
}

方案2（中等）

for(int index = 0; index < pImpObject->bandMap.length(); ++index) {
	int bandMap = pImpObject->bandMap[index] + 1;
	cv::Mat tmp = cv::Mat(buf.y, buf.x, GisEntsHelper::GDALType2GCType(dType, 1));
	GDALRasterBand *pBand = pDataset->GetRasterBand(bandMap);
	CPLErr err = pBand->RasterIO(GF_Read, read.x, read.y, size.x, size.y, tmp.data, tmp.cols, tmp.rows, dType, 0, 0);
	mats.push_back(tmp.clone());
}
GisEntsHelper::GisLineStretchImage(mats, pImpObject->histArray, sType, 0);

void GisEntsHelper::GisLineStretchImage(std::vector<cv::Mat> &mats, const HistogramDataArray &hists, OGis::StretchType sType, uint trans)
{
	if(sType == OGis::kStretch_No) {
		return;
	}
	std::vector<cv::Mat> temps = mats;
	cv::Mat mask; mats.clear();
	if(temps.size() == 3) {
		mask = ((temps[0] == temps[0]) & (temps[1] == temps[0]) & (temps[2] == temps[0]));
		temps[0].setTo(hists[0]._max, mask); temps[1].setTo(hists[1]._max, mask); temps[2].setTo(hists[2]._max, mask);
	}
	else {
		mask = ~cv::Mat(temps[0] == temps[0]);
		temps[0].setTo(hists[0]._max, mask);
	}
	for(int idx = 0; idx < temps.size(); ++idx) {
		cv::Mat mat = temps[idx]; HistogramData hist = hists[idx];
		cv::Mat result = cv::Mat(temps[0].rows, temps[0].cols, CV_8U);
		if(sType <= OGis::kStretch_Line5) {
			Double3 minmax = hist.minmax[sType - 1];
			for(int row = 0; row < mat.rows; ++row) {
				for(int col = 0; col < mat.cols; ++col) {
					double val = GetMatValue(mat, row, col);
					if(isnan(val)) {
						result.at<byte>(row, col) = 255;
						continue;
					}
					if(val < minmax.x) {
						val = 0;
					}
					else if(val > minmax.y) {
						val = 255;
					}
					else {
						val = (val - minmax.x) * minmax.z;
					}
					result.at<byte>(row, col) = val;
				}
			}
		}
		else if(OGis::kStretch_SquareRoot == sType) {
			if(hist._min < 0) {
				mat -= hist._min;
			}
			mat.convertTo(mat, CV_32F);
			mat = (mat - hist._min) * hist.minmax[sType - 1].x + hist.minmax[sType - 1].y;
			mat.convertTo(result, CV_8U);
		}
		else if(OGis::kStretch_Logaithmic == sType) {
			if(hist._min < 1) {
				mat -= (hist._min - 1);
			}
			mat.convertTo(mat, CV_32F);
			mat = (mat - hist._min) * hist.minmax[sType - 1].x + hist.minmax[sType - 1].y;
			mat.convertTo(result, CV_8U);
		}
		else if(OGis::kStretch_Gaussian == sType) {
		}
		else if(OGis::kStretch_Equalize == sType) {
			cv::normalize(mat, result, 255, 0, CV_MINMAX, CV_8U);
			cv::equalizeHist(result, result);
		}

		mats.push_back(result.clone());
	}
}

方案2的代码量虽然比较多，但是速度比方案1确实快了一点点。

方案3（最快）

for(int index = 0; index < pImpObject->bandMap.length(); ++index) {
	int bandMap = pImpObject->bandMap[index] + 1;
	cv::Mat tmp = cv::Mat(buf.y, buf.x, GisEntsHelper::GDALType2GCType(dType, 1));
	GDALRasterBand *pBand = pDataset->GetRasterBand(bandMap);
	CPLErr err = pBand->RasterIO(GF_Read, read.x, read.y, size.x, size.y, tmp.data, tmp.cols, tmp.rows, dType, 0, 0);
	mats.push_back(tmp.clone());
}
GisEntsHelper::GisLineStretchImageEx(mats, pImpObject->histArray, sType, 0);

void GisEntsHelper::GisLineStretchImageEx(std::vector<cv::Mat> &mats, const HistogramDataArray &hists, OGis::StretchType sType, uint trans)
{
	if(sType == OGis::kStretch_No) {
		return;
	}
	std::vector<cv::Mat> temps;
	cv::Mat mask;

	if(mats.size() == 3) {
		mask = ((mats[0] == hists[0]._noDataValue) & (mats[1] == hists[1]._noDataValue) & (mats[2] == hists[2]._noDataValue));
	}
	else {
		mask = ~cv::Mat(mats[0] == mats[0]);
	}
	for(int idxMat = 0; idxMat < mats.size(); ++idxMat) {
		mats[idxMat].setTo(hists[idxMat]._max, mask);
		HistogramData hist = hists[idxMat]; cv::Mat result;
		mats[idxMat].convertTo(result, CV_8U, 255.0 / (hist._max - hist._min), -255.0 * hist._min / (hist._max - hist._min));
		if(OGis::kStretch_Equalize == sType) {
			cv::equalizeHist(result, mats[idxMat]);
		}
		temps.push_back(result);
	}
	if(OGis::kStretch_Equalize == sType) {
		return;
	}
	mats.clear();
	for(int idxMat = 0; idxMat < temps.size(); ++idxMat) {
		HistogramData hist = hists[idxMat];
		cv::Mat lookTableMat = cv::Mat(1, 256, CV_8U);
		uchar *pData = lookTableMat.data;
		double valMin, valMax, factor;
		Double3 minmax = hist.minmax[sType - 1];
		double dMin = hist._min, dMax = hist._max;
		if(sType <= OGis::kStretch_Line5) {
			valMin = 255.0 * (minmax.x - hist._min) / (hist._max - hist._min);
			valMax = 255.0 * (minmax.y - hist._min) / (hist._max - hist._min);
		}
		else if(sType == OGis::kStretch_SquareRoot) {
			valMin = sqrt(0);
			valMax = sqrt(255);
		}
		else if(sType == OGis::kStretch_Logaithmic) {
			valMin = log(1);
			valMax = log(256);
		}
		factor = 255.0 / (valMax - valMin);
		for(int idx = 0; idx < 256; ++idx) {
			double val = idx;
			if(OGis::kStretch_SquareRoot == sType)  val = sqrt(idx);
			else if(OGis::kStretch_Logaithmic == sType)  val = log(idx + 1);
			if(val < valMin) pData[idx] = 0;
			else if(val > valMax) pData[idx] = 255;
			else pData[idx] = (val - valMin) * factor;
		}
		cv::Mat result = cv::Mat(temps[idxMat].rows, temps[idxMat].cols, CV_8U);
		cv::LUT(temps[idxMat], lookTableMat, result);
		mats.push_back(result.clone());
	}
}

方案3使用OpenCV的查找表进行数据的拉伸。

分析原因

方案原因
方案1对mat进行了6次遍历
方案2对mat进行了2次遍历
方案3对mat进行了2次遍历，第二次使用的是cv::LUT，这个效率不是遍历像素可比的
杨发荷
关注
  • 10
    点赞
  • 踩
  • 8
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
GEE：栅格计算（GoogleEarthEngine平台）
养乐多的博客
04-01 1万+
ArcGIS中常用到的栅格计算器公式 https://blog.csdn.net/ccnucb/article/details/81010308 在GEE中怎么对影像进行栅格计算呢？ 我们可以用.where()和.clamp() 具体语法可以参考https://zhuanlan.zhihu.com/p/76343298 这里抽取他的例子。 var DEM = ee.Image("USGS/SRTMGL1_003"); var Land_Cover = ee.Image("ESA/GLOBCOVER_L4_
基于GDAL读写影像并做线性拉伸，调整影像亮度
04-14
线性拉伸影像提升亮度，使用GDAL读取影像，做基本的线性拉伸，可以将暗色的影像变的更加清晰
【Python】GDAL基本操作/遥感大图显示
兴趣使然的创作者
02-08 3479
遥感图像往往尺寸较大，无法用默认的图像浏览器加载。GDAL是空间数据处理的开源包，支持多种数据格式的读写。遥感图像是一种带大地坐标的栅格数据，因此，可以借用GDAL对遥感图像进行读写，本文就来记录一些相关操作。
GDAL 2.1.0工具开发之gdal_translate（图像裁剪、缩放、拉伸、赋坐标投影、格式转换等）
热门推荐
杂货铺
06-12 1万+
模拟gdal_translate功能
Python遥感图像处理应用篇(二十九)：遥感图像拉伸处理
soderayer的博客
02-13 1572
将原始遥感图像各个波段拉伸到0-255范围显示。运行环境：windows10 pycharm python3.7.7包等。
C# GDAL 数字图像处理Part2 Band合成BMP
Im_Nobody33的博客
10-23 614
关于如何在Gdal中取出波段合成对应的BMP图像
GDAL实现标准差拉伸渲染影像
TJLCY的博客
11-15 2678
import os from sys import path from osgeo import gdal from osgeo import osr import numpy as np import matplotlib.pyplot as plt import matplotlib.colors as mpc import argparse def stdStretch(imgFile, colorList, flag, outPath): # imgFile 影像 # color
遥感影像百分比线性拉伸
AIEarth的博客
12-28 1254
AI Earth地球科学云平台开发者模式提供了丰富的遥感数据和函数计算能力，下面介绍结合AIE Notebook，实现遥感数据的百分比线性灰度拉伸。
对遥感图像进行2%线性拉伸
sunnygirl's house
05-28 9855
参考链接：http://www.voidcn.com/article/p-ethxsjor-du.html
C++实现遥感图像2%线性拉伸
闫小亘的博客
12-25 1429
1.接口函数 cv::Mat ColorImageStretch(const cv::Mat img, const float ratio) { //影像分RGB计算灰度直方图 cv::Mat mimg = img.clone(); cv::Mat bgr[3]; split(mimg, bgr); std::vector<cv::Point2f> hb = IS_Histogram(bgr[0]); std::vector<cv::Point2f> hg = IS_Hi
c++&python实现遥感影像1、2、5...%线性拉伸
qq_41824159的博客
06-29 686
遥感影像处理
GDAL百分比截断拉伸调亮影像
09-13
GDAL调整图像亮度之百分比截断拉伸。简介百分比截断拉伸原理，以及验证过的核心代码，相关注释，效果图
遥感图像拉伸代码
04-13
使用C++语言+OpenCV，实现了对数拉伸、高斯平滑以及直方图均衡化这三种图像拉伸的方法。
遥感图像的显示和拉伸
05-22
遥感图像的显示和拉伸！！！！！！！！！！！！
影像亮度拉伸程序
08-16
IDL 影像亮度拉伸程序
基于ArcEngine的影像拉伸算法
12-04
基于ArcEngine的影像拉伸算法,单波段和多波段组合的拉伸显示。
python2%线性拉伸
01-20
核心要点 本方法是针对遥感数据的，因此使用了Gdal，如果你针对的是普通的图片，就直接使用cv2.imread()来读取即可，就不需要gdal了。 np.percentile()函数 np.clip()函数 cv2.split()函数 cv2.merge()函数 ...
Qt、GDAL遥感影像显示
小刀砍流氓的博客
04-15 5617
影像显示原理 显示遥感影像在业界做得比较牛逼的有ENVI，Erdas，以及ArcGIS软件。其实我并不想列举ArcGIS软件，但是ArcGIS毕竟是行业的大佬，所以。。 遥感影像一般小到几百兆大到几个G，像这个么大的数据ENVI以及Erdas等软件是如何将它们显示出来呢！其实我一开始也比较好奇，同时也比较纳闷。这么大的数据，加载进来，软件居然一点都不卡，响应效率非常...
C++搭建框架，利用OpenGL、GDAL、Qt进行分块显示遥感影像
小刀砍流氓的博客
11-11 3479
  主要是利用C++搭建的框架，利用OpenGL、GDAL及Qt进行影像分块显示遥感影像，目前测试显示600M的资源3号卫星影像，仅仅需要15秒左右的时间。 此文章不对OpenGL以及GDAL做解释，如果对OpenGL和GDAL不熟悉，请自行查阅相应的文档。 利用OpenGL显示影像，那么第一步首先是要对OpenGL进行初始化。 void GeoGraphicsView::initializ...
python遥感影像线性拉伸
最新发布
12-22
Python遥感影像线性拉伸是一种用于增强遥感影像对比度的方法。在遥感影像处理中，由于不同地物的反射率和亮度不同，使得影像中的细节不易被直观地识别。线性拉伸就是通过调整影像的亮度范围，将原有的像素值分布拉伸...

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

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

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值