hyperspectralimagcode.cpp

改进的JPEG-LS算法

#include "stdafx.h"
#include <math.h>
#include ".\AdvanceJPEG-LS\AdvanceJPEG-LSCode.h"
#include "HyperspectralImagCode.h"
//高光谱图像编码器
CHyperspectralImageCode::CHyperspectralImageCode()
{
	HyperspectralDataFormat=BFILE;
	BandNumber=Height=Width=OnePixelBytes=OnePixelBits=LinesPerCompressBlock=AllowMaxError=0;
	lpHead=NULL;
	lpImageData=NULL;
	HImageWidth=HImageHeight=HImageBand=HBitPerPixel=HDataStartPosition=HImageFormat=OneLineAuxDataWidth=0;
	lpIMF=NULL;
	lpImageIndex=NULL;
	lpCompressUnit=NULL;
	lpPredictTable=NULL;
	lpCodeStreamBuf=NULL;
}
CHyperspectralImageCode::~CHyperspectralImageCode()
{
	EndHyperspectralCode();
}
void CHyperspectralImageCode::EndHyperspectralCode()
{
	int i;
	if(HMF.IsOpen()==TRUE)HMF.Close();
	if(lpImageIndex!=NULL){delete lpImageIndex;lpImageIndex=NULL;}
	if(lpIMF!=NULL)
	{
		for(i=0;i<BandNumber;i++)
		{
			if(lpIMF[i]->IsOpen()==TRUE)lpIMF[i]->Close();
			delete lpIMF[i];
		}
		delete lpIMF;
		lpIMF=NULL;
	}
	if(lpPredictTable!=NULL)delete lpPredictTable;lpPredictTable=NULL;
	if(CMF.IsOpen()==TRUE)
	{
		CMF.Close();
		lpCodeStreamBuf=NULL;
	}
	else
	{ 
		if(lpCodeStreamBuf!=NULL)delete lpCodeStreamBuf;
		lpCodeStreamBuf=NULL;
	}
	if(lpCompressUnit!=NULL)
	{
		for(i=0;i<BandNumber;i++)
		{
			lpCompressUnit[i]->EndCode();
			delete lpCompressUnit[i];
		}
		delete lpCompressUnit;
		lpCompressUnit=NULL;
	}
	BandNumber=Height=Width=OnePixelBytes=OnePixelBits=LinesPerCompressBlock=AllowMaxError=0;
}
//OMIS相关图像格式时
BOOL CHyperspectralImageCode::LoadHyperspectralImageFile(LPCTSTR pHFile,
														 int OneLineImageAuxDataL/*每扫描行图像辅助数据数*/)
{
	LPBYTE lpHpic;
	int i=::GetFileLength(pHFile);//取得完整高光谱数据长度
	if(HMF.CreateFileMap(pHFile,//创建文件映射内存
						PAGE_READONLY,
						i,
						NULL,
						OPEN_ALWAYS)!=NULL)
	{
		lpHpic=(LPBYTE)HMF.MapOfFile(0,i);
		lpHead=(stAuxMsg *)lpHpic;
		if(lpHead->iImageBand<=MAXBANDNUMBER&&lpHead->iBitPerPixel<=MAXBITPERPIXEL)
		{
			lpImageData=lpHpic+lpHead->iDataStartPosition;//高光谱数据指针
			BandNumber=HImageBand=lpHead->iImageBand;//波段数
			HImageWidth=lpHead->iImageWidth;//高光谱数据宽度
			Height=HImageHeight=lpHead->iImageHeight;//高光谱数据高度
			OnePixelBits=HBitPerPixel=lpHead->iBitPerPixel;//每象素比特数
			HDataStartPosition=lpHead->iDataStartPosition;//高光谱数据在文件中的实际位置
			HyperspectralDataFormat=HImageFormat=lpHead->iImageFormat;//高光谱数据格式
			OneLineAuxDataWidth=OneLineImageAuxDataL;
			Width=HImageWidth-OneLineImageAuxDataL;//每行图像的宽度要减出辅助数据项
			if(OnePixelBits<=8)OnePixelBytes=1;//由象素比特数产生单象素字节数
			else OnePixelBytes=2;

			return TRUE;
		}
		else
		{
			HMF.Close();
			return FALSE;
		}
	}
	return FALSE;
}
//纯图像格式时,任何高光谱数据都可以先转化为单波段纯图像格式再行压缩.
BOOL CHyperspectralImageCode::LoadHyperspectralImageFile(LPCTSTR *lpIFileList,
														int BandN,
														int H,
														int W,
														int Bytes,
														int Bits)
{//这是从单波段纯图像列表中加载高光谱数据
	int i,l;
	if(BandN<=MAXBANDNUMBER&&Bits<=MAXBITPERPIXEL)
	{
		HyperspectralDataFormat=BFILE;//数据格式标志
		OneLineAuxDataWidth=0;//没有辅助数据
		BandNumber=HImageBand=BandN;//波段数
		Height=HImageHeight=H;//图像高度
		Width=HImageWidth=W;//图像宽度
		OnePixelBits=HBitPerPixel=Bits;//单象素比特数
		OnePixelBytes=Bytes;//单象素字节数
		l=H*W*Bytes;//最小文件长度
		for(i=0;i<BandN;i++)
		{//检查每个文件是否有足够的长度
			if((int)::GetFileLength(lpIFileList[i])<l)
			{
				break;
			}
		}
		if(i<BandN)
		{//中间有文件不足最小长度时,加载高光谱数据失败
			return FALSE;
		}
		lpIMF=new CMapFile *[BandN];//创建映射文件表
		for(i=0;i<BandN;i++)
		{
			lpIMF[i]=new CMapFile;
		}
		lpImageIndex=new LPBYTE[BandN];//创建各波段数据指针表
		for(i=0;i<BandN;i++)
		{
			if(lpIMF[i]->CreateFileMap(lpIFileList[i],
									PAGE_READONLY,
									l,
									NULL,
									OPEN_ALWAYS)!=NULL)
			{
				lpImageIndex[i]=(LPBYTE)lpIMF[i]->MapOfFile(0,l);//从映射文件中获得数据指针表
			}
			else break;
		}
		if(i<BandN)
		{//若其中有文件映射失败,清理已创建的对象.
			for(i=0;i<BandN;i++)
			{
				delete lpIMF[i];
			}
			delete lpIMF;lpIMF=NULL;
			delete lpImageIndex;lpImageIndex=NULL;
			return FALSE;
		}
		return TRUE;
	}
	return FALSE;
}

LPBYTE CHyperspectralImageCode::GetOneLineImageData(int Band,int Line)//直接从原始数据中取得指针
{//这种方式支持纯图像格式及BSQ_H、BIL_H格式
	if(HyperspectralDataFormat==BFILE)//纯图像格式
	{
		return lpImageIndex[Band]+Line*Width*OnePixelBytes;
	}
	else if(HyperspectralDataFormat==BSQ_H)//BSQ_H图像格式
	{
		return lpImageData+(Band*HImageWidth*HImageHeight*OnePixelBytes)+(Line*HImageWidth*OnePixelBytes);
	}
	else if(HyperspectralDataFormat==BSQ_H)//BIL_H图像格式
	{
		return lpImageData+(Line*HImageWidth*HImageBand*OnePixelBytes)+(Band*HImageWidth*OnePixelBytes);
	}
	return NULL;
}

void CHyperspectralImageCode::GetOneLineImageData(LPBYTE lpOneLineBuf,int Band,int Line)//从原始数据中取得数据装入到指定缓冲
{//这种方式支持所有格式，尤其是BIP_H格式。
	LPBYTE lpLineData;
	int i,j;
	TCHAR Andcod1;
	short int Andcod2;
	if(HyperspectralDataFormat==BIP_H)//BIP_H图像格式
	{
		lpLineData=lpImageData+(Line*HImageWidth*HImageBand*OnePixelBytes);
		//须按波段间隔采样才能得到单一波段连续一行的数据。
		if(OnePixelBytes==1)
		{
			Andcod1=~((-1)<<OnePixelBits);
			for(i=0,j=Band;i<Width;i++,j+=HImageBand)
			{
				lpOneLineBuf[i]=lpLineData[j]&Andcod1;
			}
		}
		else
		{
			Andcod2=~((-1)<<OnePixelBits);
			for(i=0,j=Band;i<Width;i++,j+=HImageBand)
			{
				*(((short int *)lpOneLineBuf)+i)=(*(((short int *)lpLineData)+j))&Andcod2;
			}
		}
	}
	else if(HyperspectralDataFormat==BIL_H||HyperspectralDataFormat==BSQ_H||HyperspectralDataFormat==BFILE)
	{
		MoveMemory(lpOneLineBuf,GetOneLineImageData(Band,Line),Width*OnePixelBytes);
		if(OnePixelBytes==1)
		{
			Andcod1=~((-1)<<OnePixelBits);
			for(i=0;i<Width;i++)
			{
				lpOneLineBuf[i]&=Andcod1;
			}
		}
		else
		{
			Andcod2=~((-1)<<OnePixelBits);
			for(i=0;i<Width;i++)
			{
				*(((short int *)lpOneLineBuf)+i)&=Andcod2;
			}
		}
	}
}
//从原始数据中取得每行辅助数据装入到指定缓冲
int CHyperspectralImageCode::GetOneLineAuxData(LPBYTE lpOneLineBuf,int Band,int Line)
{//这种方式支持BSQ_H、BIL_H格式。
	if(HyperspectralDataFormat==BIL_H||HyperspectralDataFormat==BSQ_H)
	{
		MoveMemory(lpOneLineBuf,
					GetOneLineImageData(Band,Line)+Width*OnePixelBytes,
					OneLineAuxDataWidth*OnePixelBytes);
		return OneLineAuxDataWidth*OnePixelBytes;
	}
	return 0;
}
//为输出码流准备缓冲
void CHyperspectralImageCode::PrepareCodeStreamBuf(int OneBlockLines,LPCTSTR CodeFileName,int AllowError)
{
	CompressFile=CodeFileName;
	LinesPerCompressBlock=OneBlockLines;
	StreamBufLength=OneBlockLines*Width*BandNumber*OnePixelBits/8;
	if(AllowError>=4)StreamBufLength/=2;
	if(AllowError>=24)StreamBufLength/=2;
	if(AllowError>=48)StreamBufLength/=2;
	if(lpCodeStreamBuf!=NULL)delete lpCodeStreamBuf;
	lpCodeStreamBuf=new BYTE[StreamBufLength];
	ClearCodeStreamBuf();
}
void CHyperspectralImageCode::ClearCodeStreamBuf()
{
	ZeroMemory(lpCodeStreamBuf,StreamBufLength);
}
//从输入的字符串中装入预测波段的配置
void CHyperspectralImageCode::LoadPredictTable(LPCTSTR lpAsc)
{
	TCHAR *lpasc;
	int i,j,k,l;
	if(lpPredictTable!=NULL)delete lpPredictTable;
	lpPredictTable=new short int[BandNumber];
	l=strlen(lpAsc)+2;
	lpasc=new TCHAR[l];
	strcpy(lpasc,lpAsc);
	for(i=0,j=0;i<BandNumber&&j<l;i++)
	{
		for(;lpasc[j]!='-'&&(lpasc[j]>'9'||lpasc[j]<'0')&&j<l;j++);
		k=j;
		for(;lpasc[j]=='-'||(lpasc[j]<='9'&&lpasc[j]>='0')&&j<l;j++);
		lpasc[j]='\0';
		lpPredictTable[i]=atoi(lpasc+k);
	}
	if(i<BandNumber)for(;i<BandNumber;i++)lpPredictTable[i]=-1;
	delete lpasc;
}
//一般在调用预测波段配置后才可调用
void CHyperspectralImageCode::CreatCompressObject()//根据波段数创建相应数量的压缩对象
{
	int i;
	lpCompressUnit=new CAdvanceJepgLSOneLineCode *[BandNumber];
	for(i=0;i<BandNumber;i++)
	{
		lpCompressUnit[i]=new CAdvanceJepgLSOneLineCode;
		if(lpPredictTable[i]>=0)lpCompressUnit[i]->PrepareCode(Width,lpCompressUnit[lpPredictTable[i]]);
		else lpCompressUnit[i]->PrepareCode(Width);
	}
}
//创建高光谱数据压缩信息头，返回头字节数
int CHyperspectralImageCode::CreatHyperspectralCompressHead(LPBYTE lpBuf,LPBYTE &lpAuxBuf,int &AuxBufL)
{
	lpBuf[0]=HyperspectralDataFormat;//高光谱数据格式
	*((short int *)(lpBuf+1))=BandNumber;//波段数
	*((short int *)(lpBuf+3))=Width;//图像宽
	*((short int *)(lpBuf+5))=Height;//图像高
	if(HyperspectralDataFormat==BFILE)
	{//纯图像格式时
		*((char *)(lpBuf+7))=(char)OnePixelBits;//单象素比特数
		*((char *)(lpBuf+8))=(char)OnePixelBytes;//单象素字节数
		MoveMemory(lpBuf+9,lpPredictTable,BandNumber*sizeof(short int));//以双字节存放预测波段号
		lpAuxBuf=NULL;AuxBufL=0;
		return 9+BandNumber*sizeof(short int);
	}
	else 
	{//高光谱数据格式时
		//以双字节存放预测波段号
		MoveMemory(lpBuf+7,lpPredictTable,BandNumber*sizeof(short int));
		lpAuxBuf=(LPBYTE)lpHead;
		AuxBufL=sizeof(stAuxMsg)+lpHead->stLineAuxDataInfo.ladHeight*lpHead->stLineAuxDataInfo.ladWidth;
		//高光谱数据格式头
//		MoveMemory(lpBuf+7+BandNumber*sizeof(short int),lpHead,sizeof(stAuxMsg));
		//高光谱数据格式头中的辅助数据
//		MoveMemory(lpBuf+7+BandNumber*sizeof(short int)+sizeof(stAuxMsg),
//					lpHead+sizeof(stAuxMsg),
//					lpHead->stLineAuxDataInfo.ladHeight*lpHead->stLineAuxDataInfo.ladWidth);
		return 7+BandNumber*sizeof(short int)+sizeof(stAuxMsg)+lpHead->stLineAuxDataInfo.ladHeight*lpHead->stLineAuxDataInfo.ladWidth;
	}
	return 0;
}
//创建高光谱数据块压缩信息头，返回头字节数
int CHyperspectralImageCode::CreatBlockCompressHead(LPBYTE lpBuf,int Maxerr,int LineBG,int LineN)
{
	//前四字节用于保存块的总字节数
	lpBuf[4]=Maxerr;//当前块误差界
	*((short int *)(lpBuf+5))=LineBG;//当前块起始行
	*((short int *)(lpBuf+7))=LineN;//当前块所包含的图像行数
	return 9+OneLineAuxDataWidth*OnePixelBytes*BandNumber*LineN;//返回块信息头总字节数
}
int CHyperspectralImageCode::LoadBlockCompressHead(LPBYTE lpCodeStream,int &AllowMaxErr,int &LineBG,int &LineN)
{
	AllowMaxErr=lpCodeStream[4];
	LineBG=*((short int *)(lpCodeStream+5));
	LineN=*((short int *)(lpCodeStream+7));
	return 9+OneLineAuxDataWidth*OnePixelBytes*BandNumber*LineN;

}
//执行压缩操作
void CHyperspectralImageCode::DoCompress(int AllowMaxErr,double *lpErrMsePerPixel,double *lpBitPerPixel)
{
	/*
			a、建立压缩块循环。
				b、创建压缩块信息头。
				c、建立图像行循环
					d、建立波段循环
						e、取某行某波段数据
						f、用相应的波段编码器编码该图像行数据，码流汇入总码流。
					g、下一个波段。
				h、下一个图像行
				i、导出压缩块码流。
			j、下一个压缩图像块
	*/
	int Line,i,j,LineBG,LineN,AuxData_cp,CodeBit_cp;
	LPBYTE lpAuxData=new BYTE[9+OneLineAuxDataWidth*OnePixelBytes*BandNumber*LinesPerCompressBlock];
	LPBYTE lpImageData=new BYTE[Width*OnePixelBytes+4];
	LPBYTE lpAuxD;
	int AuxL;
	double ErrMse,BitPerPixel;
	AllowMaxError=AllowMaxErr;
	//创建并保存高光谱数据压缩头信息
	i=CreatHyperspectralCompressHead(lpCodeStreamBuf,lpAuxD,AuxL);
	WriteDataToFile(lpCodeStreamBuf,i-AuxL);
	if(AuxL>0)AppendDataToFile(lpAuxD,AuxL);
	//执行压缩循环
	for(Line=0;Line<Height;)
	{//这是块循环层
		LineBG=Line;
		LineN=((Height-LineBG)>=LinesPerCompressBlock)?LinesPerCompressBlock:(Height-LineBG);
		AuxData_cp=9;
		CreatBlockCompressHead(lpAuxData,AllowMaxError,LineBG,LineN);
		for(i=0;i<BandNumber;i++)
		{//各波段编码器复位
			lpCompressUnit[i]->InitialisationsCode(OnePixelBits,AllowMaxError);
		}
		CodeBit_cp=0;
		ErrMse=0;
		BitPerPixel=0;
		for(j=0;j<LineN;j++)
		{
			for(i=0;i<BandNumber;i++)
			{
				//取当前行图像数据
				GetOneLineImageData(lpImageData,i,j+LineBG);
				//取当前行辅助数据
				AuxData_cp+=GetOneLineAuxData(lpAuxData+AuxData_cp,i,j+LineBG);
				//压缩当前行图像,码流以比特为单位存放
				CodeBit_cp+=lpCompressUnit[i]->DoEncodeLine(lpImageData,
															OnePixelBytes,
															lpCodeStreamBuf,
															CodeBit_cp,
															&ErrMse);
			}
		}
		Line+=LineN;
		//累计总的码流比特数
		BitPerPixel+=CodeBit_cp;
		//计算当前压缩块总字节数
		*((int *)lpAuxData)=AuxData_cp+(CodeBit_cp>>3)+((CodeBit_cp&7)?1:0);
		//保存当前压缩块信息头
		AppendDataToFile(lpAuxData,AuxData_cp);
		//保存当前压缩块码流