很多遙感衛星資料使用的量化層級都要比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 charDT_8U;/*! 16U */typedef unsigned shortDT_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;elseeDT = 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;elsepDstData[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的範圍。