【影像配准】目标影像在参考影像中的自动定位与裁剪（附有完整代码）

       在影像配准时通常要求输入的参考影像和待配准影像是同一地区且尺寸相等，但是在实际应用

时如果要求输入的参考影像和待配准影像是同一地区且尺寸相等的话，则需人工手动裁剪，无形中

增加了工作量！！！

       所以在配准前需要增加一个自动裁剪功能，即输入的参考影像很大，然后根据输入的待配准影

像从参考影像中裁剪出目标区域，然后进行配准！！！

预备知识：了解像素坐标与投影坐标间的转换关系！！！

转换关系见：

像素坐标转投影坐标、投影坐标转像素坐标（附有完整代码及测试结果）_一米九零小胖子的博客-CSDN博客_tfw文件

一、利用待配准影像的投影信息从参考影像中裁剪出目标区域原理 

1、获取待配准影像的左上角投影坐标和右下角投影坐标；

2、该投影坐标即为参考影像中相应区域的左上角投影坐标和右下角投影坐标；

3、把投影坐标转化为像素坐标；

4、利用gdal从原参考影像中裁剪出目标区域(.tif)； 

二、利用待配准影像的投影信息从参考影像中裁剪出目标区域代码

#include<iostream>
#include<gdal_priv.h>
#include<gdal.h>
#include<windows.h>
#include<opencv2\imgproc\types_c.h>
#include<opencv2\imgproc\imgproc.hpp>
#include<opencv2\highgui\highgui.hpp>

using namespace std;


//获取tif影像的投影信息
void getProj(double* geo, const char* srcTifPath)
{
	GDALAllRegister();

	GDALDataset* podataset = (GDALDataset*)GDALOpen(srcTifPath, GA_ReadOnly);

	//保存 .tif 文件的投影信息' x=geo[0]、y=geo[3]、像素大小=geo[1]
	//double geo[6];

	podataset->GetGeoTransform(geo);
}


//裁剪参考影像中的目标区域
void imageCut(string refpath, string senpath, string dstpath)
{
	//获取两幅影像的投影信息
	double refgeo[6], sengeo[6];

	getProj(refgeo, refpath.c_str());
	getProj(sengeo, senpath.c_str());

	/******使用投影信息在参考影像中截取目标影像******/
	cv::Mat refImage = cv::imread(refpath, -1);
	cv::Mat senImage = cv::imread(senpath, -1);

	int senRow = senImage.rows;
	int senCol = senImage.cols;

	//找到目标影像左上角和右下角坐标：（右下x-左上x）/size=像素宽度，（右下y-左上y）/size=像素宽度
	double senX = senCol * sengeo[1] + sengeo[0];
	double senY = senRow * sengeo[5] + sengeo[3];

	//计算目标影像坐标在参考影像中对应的像素坐标
	int refLeftX = (sengeo[0] - refgeo[0]) / refgeo[1];
	int refLeftY = (refgeo[3] - sengeo[3]) / refgeo[1];
	int refRightX = (senX - refgeo[0]) / refgeo[1];
	int refRightY = (refgeo[3] - senY) / refgeo[1];

	/********创建裁剪后的tif********/

	//裁剪区域的投影信息
	refgeo[0] = sengeo[0];
	refgeo[3] = sengeo[3];

	GDALDataset* poDataset1; //待创建的GDAL数据集 
	GDALDataset* poDataset2; //参考影像影像
	GDALDriver* poDriver;    //驱动，用于创建新的文件

	//gdal打开参考图像
	GDALAllRegister();
	poDataset2 = (GDALDataset*)GDALOpen(refpath.c_str(), GA_ReadOnly);

	//判断是否打开文件
	if (poDataset2 == NULL) { printf("文件打开失败"); }
	else printf("文件打开成功");

	//获取待加入坐标的图像的通道长宽
	int nBandCount = poDataset2->GetRasterCount();
	int nImgSizeX = refRightX - refLeftX;
	int nImgSizeY = refRightY - refLeftY;

	//创建新文件
	char fomat[] = "GTiff";
	poDriver = GetGDALDriverManager()->GetDriverByName(fomat);

	//获取格式类型
	char** papszMetadata = poDriver->GetMetadata();

	//根据文件路径文件名，图像宽，高，波段数，数据类型，文件类型，创建新的数据集
	poDataset1 = poDriver->Create(dstpath.c_str(), nImgSizeX, nImgSizeY, nBandCount, GDT_Byte, papszMetadata);

	//坐标赋值
	poDataset1->SetGeoTransform(refgeo);

	//将原图像数据读出，进行相应处理后，写入新文件
	BYTE* ppafScan = new BYTE[nImgSizeX * nImgSizeY * nBandCount];

	//读取参考影像中的目标区域信息
	poDataset2->RasterIO(GF_Read, refLeftX, refLeftY, nImgSizeX, nImgSizeY, ppafScan, nImgSizeX, nImgSizeY, GDT_Byte, nBandCount, 0, 0, 0, 0);

	//将图像数据写入新图像中
	poDataset1->RasterIO(GF_Write, 0, 0, nImgSizeX, nImgSizeY, ppafScan, nImgSizeX, nImgSizeY, GDT_Byte, nBandCount, 0, 0, 0, 0);

	delete poDataset1;
	delete poDataset2;
}


int main()
{
	string refpath = "..\\testImage\\ref.tif";
	string senpath = "..\\testImage\\sen.tif";
	string dstpath = ".\\imageSave\\dst.tif";

	imageCut(refpath, senpath, dstpath);

	return 0;
}