使用GDAL将12bit量化图像转为16bit或者8bit

    很多遥感卫星数据使用的量化级别都要比8bit高，比如常用的WorldView用的是12bit的量化，对于一些图像处理软件，是不能直接处理12bit量化的图像，所以需要将12bit的数据转为16bit数据或者8bit数据来进行处理。

    下面写了一个简单的函数来进行处理，具体原理很简单，就是使用GDAL将12bit的数据读进来，然后再使用线性拉伸为8bit存出去，或者直接保存为16bit数据。注意12bit的数据在GDAL中读取的时候会显示为16bit数据，就好比2bit的数据在GDAL中是8bit一样，因为在C或者C++中很难找到一个类型来表示2bit或者12bit的东西，最小的char是8bit，short是16bit。代码如下，首先是头文件：

/***************************************************************************
*
* Time: 2012-02-23
* Project: 遥感平台 
* Purpose: 将12bit数据转换为8bit或者16bit
* Author:  李民录
* Copyright (c) 2012, liminlu0314@gmail.com
* Describe:将12bit数据转换为8bit或者16bit
*
****************************************************************************/
#ifndef DATARESCALE_H
#define DATARESCALE_H

/*! 8U */
typedef unsigned char	DT_8U;
/*! 16U */
typedef unsigned short	DT_16U;
/**
* @brief 释放数组
*/
#define RELEASE(x)	if(x!=NULL) {delete []x; x = NULL;}

/**
* @brief 图像转换，将图像存为16bit，前提确保输入的数据是12bit的
* @param pszSrcFile		输入文件路径
* @param pszDstFile			输出文件路径
* @param bTo8			是否转为8位，false为专为16bit数据，true表示转为8bit数据
* @param pszFormat			输出文件格式，详细参考GDAL支持数据类型
* @return 返回值，表示计算过程中出现的各种错误信息
*/
int ImageDataRescale216(const char* pszSrcFile, const char* pszDstFile, bool bTo8 = true, const char* pszFormat = "GTiff");

#endif /* DATARESCALE_H */

下面是函数实现代码：

/***************************************************************************
*
* Time: 2012-02-23
* Project: 遥感平台 
* Purpose: 将12bit数据转换为8bit或者16bit
* Author:  李民录
* Copyright (c) 2012, liminlu0314@gmail.com
* Describe:将12bit数据转换为8bit或者16bit
*
****************************************************************************/
#include "DataRescale.h"
#include "gdal_priv.h"

/**
* @brief 图像转换，将图像存为16bit，前提确保输入的数据是12bit的
* @param pszSrcFile		输入文件路径
* @param pszDstFile			输出文件路径
* @param bTo8			是否转为8位，false为专为16bit数据，true表示转为8bit数据
* @param pszFormat			输出文件格式，详细参考GDAL支持数据类型
* @return 返回值，表示计算过程中出现的各种错误信息
*/
int ImageDataRescale(const char* pszSrcFile, const char* pszDstFile, bool bTo8 = true, const char* pszFormat = "GTiff")
{
	//判断输入路径是否为空
	if(pszSrcFile == NULL || pszDstFile == NULL)
		return -1;

	GDALAllRegister();

	GDALDataset *poSrcDS = (GDALDataset *) GDALOpen( pszSrcFile, GA_ReadOnly );
	if( poSrcDS == NULL )
	{
		//图像打开失败
		return -2;
	}

	GDALDriver *poDriver = GetGDALDriverManager()->GetDriverByName(pszFormat);
	if( poDriver == NULL )
	{
		//不能创建制定类型的文件，请检查该文件类型GDAL是否支持创建
		GDALClose( (GDALDatasetH) poSrcDS );
		return -3;
	}

	//获取图像宽高和波段数
	int iXSize = poSrcDS->GetRasterXSize();
	int iYSize = poSrcDS->GetRasterYSize();
	int iBandCount = poSrcDS->GetRasterCount();

	//确定输出图像的位数
	GDALDataType eDT = GDT_UInt16;
	if(bTo8)
		eDT = GDT_Byte;
	else
		eDT = GDT_UInt16;
		
	//创建16bit的数据
	GDALDataset *poDstDS = poDriver->Create(pszDstFile, iXSize, iYSize, iBandCount, eDT, NULL);
	double dGeoTrans[6] = {0};
	
	//设置仿射变换参数
	poSrcDS->GetGeoTransform(dGeoTrans);
	poDstDS->SetGeoTransform(dGeoTrans);
	//设置图像投影信息
	poDstDS->SetProjection(poSrcDS->GetProjectionRef());

	//用于保存读取的12bit数据
	DT_16U *pSrcData = new DT_16U[iXSize];

	if(bTo8)	//转换为8bit
	{
		//定义结果数据存储空间
		DT_8U *pDstData = new DT_8U[iXSize];
		
		//循环波段
		for(int iBand=1; iBand<=iBandCount; iBand++)
		{
			GDALRasterBand *pSrcBand = poSrcDS->GetRasterBand(iBand);
			GDALRasterBand *pDstBand = poDstDS->GetRasterBand(iBand);
			
			for(int i=0; i<iYSize; i++)	//循环图像高
			{
				//将数据读出来
				pSrcBand->RasterIO(GF_Read, 0, i, iXSize, 1, pSrcData, iXSize, 1, GDT_UInt16, 0, 0);
				
				//循环，将12bit数据专为8bit数据，使用线性拉伸方式
				for(int j=0; j<iXSize; j++)
				{
					double dTemp = pSrcData[j] * 255.0 / 4095.0;
					if(dTemp > 255.0)
						pDstData[j] = 255;
					else
						pDstData[j] = (DT_8U)dTemp;
				}
				
				pDstBand->RasterIO(GF_Write, 0, i, iXSize, 1, pDstData, iXSize, 1, GDT_Byte, 0, 0);
			}
		}
		
		RELEASE(pDstData);	//释放内存
	}
	else	//转换为16bit数据
	{
		//循环波段
		for(int iBand=1; iBand<=iBandCount; iBand++)
		{
			GDALRasterBand *pSrcBand = poSrcDS->GetRasterBand(iBand);
			GDALRasterBand *pDstBand = poDstDS->GetRasterBand(iBand);
			
			for(int i=0; i<iYSize; i++)	//循环图像高
			{
				//将数据读出来，然后写入结果数据
				pSrcBand->RasterIO(GF_Read, 0, i, iXSize, 1, pSrcData, iXSize, 1, GDT_UInt16, 0, 0);
				pDstBand->RasterIO(GF_Write, 0, i, iXSize, 1, pSrcData, iXSize, 1, GDT_UInt16, 0, 0);
			}
		}
	}

	RELEASE(pSrcData);

	//关闭原始图像和结果图像
	GDALClose( (GDALDatasetH) poDstDS );
	GDALClose( (GDALDatasetH) poSrcDS );

	return 0;
}

    对于上面的实现做一个简单的说明，12bit的数据读进来，对于16bit的直接写到结果图像里面，没有拉伸到16bit的范围，这样就是完全保留了原始数据的所有图像信息，对于保存为8bit的数据，肯定会造成部分信息的丢失，使用最简单的线性方程进行拉伸，将0～4095的范围拉伸到0～255的范围。