imgprocessing.cpp

直方图的相关操作

	if(file.Open(str,CFile::modeRead)==NULL)
	{
		AfxMessageBox("读向量错误");
		return false;
	}
	CArchive ar(&file,CArchive::load);
	cf.Serialize(ar);
	return true;
}
void Img_Mirror(CByteArray& img,int height,int width)
{
	int i,j;
	int buf;
	for(i=0;i<height;i++)
		for(j=0;j<width/2;j++)
		{
			buf=img[i*width+j];
			img[i*width+j]=img[i*width+width-j-1];
			img[i*width+width-j-1]=buf;
		}
}
void Img_ToDbFormat(DoubleArray& klc,CString& str)
{
	int size=klc.GetSize();
	str="";
	if(size!=0)
    {
	     CString strbuf="";
	     strbuf.Format("%.3d",size);
	     str+=strbuf;
		 for(int i=0;i<size;i++)
		 {
			 strbuf.Format("%8.3f",klc[i]);
			 str+=strbuf;
		 }
	}
}
void Img_FromDbFormat(DoubleArray& klc,unsigned short* str)
{
	 int index=0;
	 char* strbuf;
	 strbuf=new char[9];
     strbuf[0]=(char)str[0];
     strbuf[1]=(char)str[1];
     strbuf[2]=(char)str[2];
     strbuf[3]='\0';
	 int size=atoi(strbuf);
	 klc.SetSize(size);
	 for(int i=0;i<size;i++)
	 {
         strbuf[0]=(char)str[i*8+3];
         strbuf[1]=(char)str[i*8+4];
         strbuf[2]=(char)str[i*8+5];
         strbuf[3]=(char)str[i*8+6];
         strbuf[4]=(char)str[i*8+7];
         strbuf[5]=(char)str[i*8+8];
         strbuf[6]=(char)str[i*8+9];
         strbuf[7]=(char)str[i*8+10];
		 strbuf[8]='\0';
		 klc[i]=atof(strbuf);
	 }
}
void Img_FromDbFormat(CByteArray& img,unsigned short* str)
{
	 int index=0;
	 char strbuf[4];
     strbuf[0]=(char)str[index++];
     strbuf[1]=(char)str[index++];
     strbuf[2]=(char)str[index++];
     strbuf[3]='\n';
	 int width=atoi(strbuf);
     strbuf[0]=(char)str[index++];
     strbuf[1]=(char)str[index++];
     strbuf[2]=(char)str[index++];
     strbuf[3]='\n';
	 int height=atoi(strbuf);
	 img.SetSize(width*height);
	 for(int i=0;i<width*height;i++)
	 {
		 strbuf[0]=(char)str[index++];
		 strbuf[1]=(char)str[index++];
		 strbuf[2]=(char)str[index++];
		 strbuf[3]='\n';
		 img[i]=atoi(strbuf);
	 }
}
/*void SetDBList(COleSafeArray& vFieldList,COleSafeArray& vValueList,
						long No,CString& TypeStr,CString& MemoStr,int flag)
{
		
	    vFieldList.CreateOneDim(VT_VARIANT,3);
		vValueList.CreateOneDim(VT_VARIANT,3);
        long lArrayIndex[1];
		//ImgNo域
		lArrayIndex[0] = 0;
		if(flag==0)
		{
			vFieldList.PutElement(lArrayIndex,&(_variant_t("PersonNo")));
			vValueList.PutElement(lArrayIndex,&(_variant_t(No)));
		}
		else
		{
            vFieldList.PutElement(lArrayIndex,
                  &(_variant_t("ImgNo")));
            vValueList.PutElement(lArrayIndex,
                  &(_variant_t(No)));
		}
        //RecType域
		lArrayIndex[0] = 1;
		if(flag==0)
		{
			vFieldList.PutElement(lArrayIndex,&(_variant_t("ViewType")));
			vValueList.PutElement(lArrayIndex,&(_variant_t(TypeStr)));
		}
		else
		{
		    vFieldList.PutElement(lArrayIndex,
                       &(_variant_t("RecInfoType")));
		    vValueList.PutElement(lArrayIndex,
			     &(_variant_t(TypeStr)));
		}
		//RecInfo域
		lArrayIndex[0] = 2;
		if(flag==0)
		{
            vFieldList.PutElement(lArrayIndex,
                        &(_variant_t("ImgBody")));
			vValueList.PutElement(lArrayIndex,
			        &(_variant_t(MemoStr)));
		}
		else
		{
		    vFieldList.PutElement(lArrayIndex,
                       &(_variant_t("RecInfo")));
		    vValueList.PutElement(lArrayIndex,
			     &(_variant_t(MemoStr)));
		}
}
*/
/*
///////////////////////////////////////////////////////////////////////////// 
// CFFT 
// Calculate the frequency function by Fast Forier Transformation 
// 
//              Input:  CFloatArray& pData 
//              Output: CFloatArray& 
//                              the first element is the direct current portion 
//                              the second element is the one time wave portion 
//                              the third element is the two times wave portion 
//                              the forth element is the three times wave portion 
//                              ... 
//                              the unit is percent 
///////////////////////////////////////////////////////////////////////////// 
void fft( DoubleArray& pData ) 
{ 
        int Count = pData.GetSize(); 
        double *TR = new double[Count+1]; 
        double *TI = new double[Count+1]; 
        int K,G,D,Q,C; 
        for( int I = 0; I <= Count; I ++ ) 
                TR[I] = TI[I] = 0.0f; 
 
        for( I = 0; I < Count; I ++ ) 
                TR[I] = pData[I]; 
 
        int M=int(log((double)Count)/log((double)2) + 0.5); 
        int T=1; 
        Count = 2; 
 
        for( I = 1; I < M; I ++ ) 
                Count *= 2; 
         
         
        for(int L=1; L<= M; L ++ ) 
        { 
                int N1 = 2; 
                for( I = 1; I < M-L; I ++ ) 
                        N1 *= 2; 
                if( M-L == 0 ) 
                        N1 = 1; 
                 
                int P=0; 
entry:  G=P; 
 
                D=0; 
                for( Q=1; Q <= M; Q ++ ) 
                { 
                        C=G/2; 
                        D=2*D+G-2*C; 
                        G=C; 
                } 
 
                double A=6.283f*D/Count; 
                double A1=cos(A); 
                double A2=sin(A); 
 
                for( K=P*N1; K <= (P+1)*N1-1; K ++ ) 
                { 
                        double TP=TR[K+N1]*A1+TI[K+N1]*A2; 
                        double TQ=TI[K+N1]*A1-TR[K+N1]*A2; 
                        TR[K+N1]=TR[K]-TP; 
                        TI[K+N1]=TI[K]-TQ; 
                        TR[K]+=TP; 
                        TI[K]+=TQ; 
                } 
 
                P=P+2; 
                if( (P+1)*N1<Count ) 
                        goto entry; 
        } 
        for(K=0; K < Count; K ++ ) 
        { 
                G=K; 
 
                D=0; 
                for( Q=1; Q <= M; Q ++ ) 
                { 
                        C=G/2; 
                        D=2*D+G-2*C; 
                        G=C; 
                } 
 
                if( K>=D ) 
                        continue; 
                double XP=TR[K]; 
                double XQ=TI[K]; 
                TR[K]=TR[D]; 
                TI[K]=TI[D]; 
                TR[D]=XP; 
                TI[D]=XQ; 
        } 
 
        int OutCount = Count/2; 
        pData.SetSize( OutCount + 1 ); 
        float TOTAL = 0.0f; 
        for( I=0; I <= OutCount; I ++ ) 
        { 
                pData[I] = sqrt((TR[I]*T/Count)*(TR[I]*T/Count)+(TI[I]*T/Count)*(TI[I]*T/Count)); 
                TOTAL+=pData[I]; 
        } 
        if( TOTAL == 0 ) 
                TOTAL = 1; 
 
        for( K=0; K <= OutCount; K ++ ) 
                pData[K]=pData[K]/TOTAL*100; 
 
        delete [] TR; 
        delete [] TI; 
} */
#define Pi          3.1415926
void fft(DoubleArray& Real ,DoubleArray& Imag,DoubleArray& NewReal,DoubleArray& NewImag)
{
	register int	loop, loop1, loop2;
	int	    i,j,i1;			/* going to right shift this */
	int		        i2, i3, i4, y;
	double		    a1, a2, b1, b2, z1, z2;
    //样点数目
	int temp,Count;
	temp=Count=Real.GetSize();
	i1 = Count >> 1;
	i2 = 1;
    //2的幂数
	int power=0;
	while(temp>=2)
	{
	    temp=temp/2;
		power++;
	}
    //计算sine值和Br_table
	DoubleArray Sine,Cos;
	Sine.SetSize(Count);
	CUIntArray Br_table;
    Br_table.SetSize(Count);
	char* buf=new char[power+1];
    CString str;
	for(i=0;i<Count;i++)
	{
//        Real[i]/=(double)Count;
		Sine[i]=sin(i*2*Pi/Count); 
		_itoa(i,buf,2);
		str.Format("%s",buf);
		while(str.GetLength()<power)
		{
			str="0"+str;
		}
		str.MakeReverse();
		j=0;
		int x=str.GetAt(j)-'0';
		j++;
		while(j<str.GetLength())
		{
			x*=2;
			if(str.GetAt(j)=='1') x++;
			j++;
		}
		Br_table[i]=x;
	}
	/* perform the butterfly's */
	for (loop = 0; loop < power; loop++)  {
		i3 = 0;
		i4 = i1;

		for (loop1 = 0; loop1 < i2; loop1++)  {
			y = Br_table[(i3 / (int)i1)];
			z1 = Sine[(y+Count/4)%Count];
			z2 = -Sine[y];

			for (loop2 = i3; loop2 < i4; loop2++)  {

				a1 = Real[loop2];
				a2 = Imag[loop2];

				b1  = z1*Real[loop2+i1] - z2*Imag[loop2+i1];
				b2  = z2*Real[loop2+i1] + z1*Imag[loop2+i1];

				Real[loop2] = a1 + b1;
				Imag[loop2] = a2 + b2;

				Real[loop2+i1] = a1 - b1;
				Imag[loop2+i1] = a2 - b2;
			}

			i3 += (i1 << 1);
			i4 += (i1 << 1);
		}

		i1 >>= 1;
		i2 <<= 1;
	}
	NewReal.SetSize(Count);
	NewImag.SetSize(Count);
	for(i=0;i<Count;i++)
	{
		NewReal[i]=Real[Br_table[i]];
		NewImag[i]=Imag[Br_table[i]];
	}
	
	delete[] buf;
}
//fft2:  2D Fourie transform 
void fft2(CByteArray& input,DoubleArray& output,int Height,int Width)
{
    DoubleArray Real,Imag;
    DoubleArray RealBuf,ImagBuf,RealBuf1,ImagBuf1;
	int i,j;
	int count=input.GetSize();
	double x;
	//初始化
	Real.SetSize(count);
	Imag.SetSize(count);
	for(i=0;i<count;i++)
	{
		Real[i]=input[i];
		Imag[i]=0;
	}
	//对行做一维fourie变换
	RealBuf.SetSize(Width);
	ImagBuf.SetSize(Width);
	RealBuf1.SetSize(Width);
	ImagBuf1.SetSize(Width);
	for(i=0;i<Height;i++)
	{
	    for(j=0;j<Width;j++)
		{
			ImagBuf[j]=0;
			RealBuf[j]=Real[i*Width+j];
			x=Real[i*Width+j];
		}
		fft(RealBuf,ImagBuf,RealBuf1,ImagBuf1);
		for(j=0;j<Width;j++)
		{
			Real[i*Width+j]=RealBuf1[j];
			Imag[i*Width+j]=ImagBuf1[j];
			x=Real[i*Width+j];
		}
	}
	//对列做一维Fourie变换
	RealBuf.SetSize(Height);
	ImagBuf.SetSize(Height);
	RealBuf1.SetSize(Height);
	ImagBuf1.SetSize(Height);
	for(i=0;i<Width;i++)
	{
		for(j=0;j<Height;j++)
		{
			RealBuf[j]=Real[j*Width+i];
			ImagBuf[j]=Imag[j*Width+i];
			x=Real[i*Width+j];
		}
		fft(RealBuf,ImagBuf,RealBuf1,ImagBuf1);
		for(j=0;j<Height;j++)
		{
			Real[j*Width+i]=RealBuf1[j];
			Imag[j*Width+i]=ImagBuf1[j];
			x=Real[i*Width+j];
		}
	}
	output.SetSize(count);
	for(i=0;i<count;i++)
		x=output[i]=sqrt(Real[i]*Real[i]+Imag[i]*Imag[i]);
}
void dinv(DoubleArray& mat,int width,int& flag)
{
	flag=1;
	int i,j,k;
	double temp,d=0.0;
	CWordArray IS,JS;
	IS.SetSize(width);
	JS.SetSize(width);
	for(k=0;k<width;k++)
	{
		d=0.0;
		for(i=k;i<width;i++)
			for(j=k;j<width;j++)
			{
				temp=mat[i*width+j];
				if(temp<0.0) temp=-temp;
				if(temp>d)
				{
					d=temp;
					IS[k]=i;
					JS[k]=j;
				}
			}
		if(d==0.0)//奇异
		{
			flag=0;
			return;
		}
		if(IS[k]!=k)//行交换
			for(j=0;j<width;j++)
			{
				temp=mat[k*width+j];
				mat[k*width+j]=mat[IS[k]*width+j];
				mat[IS[k]*width+j]=temp;
			}
		if(JS[k]!=k)//列交换
			for(i=0;i<width;i++)
			{
				temp=mat[i*width+k];
				mat[i*width+k]=mat[i*width+JS[k]];
				mat[i*width+JS[k]]=temp;
			}
		mat[k*width+k]=1.0/mat[k*width+k];
		for(j=0;j<width;j++)
           if(j!=k)
			   mat[k*width+j]=mat[k*width+j]*mat[k*width+k];
		for(i=0;i<width;i++)
			if(i!=k)
				for(j=0;j<width;j++)
					if(j!=k)
					  mat[i*width+j]=mat[i*width+j]-mat[i*width+k]*mat[k*width+j];
		for(i=0;i<width;i++)
			if(i!=k) mat[i*width+k]=-1*mat[i*width+k]*mat[k*width+k];
	}
	for(k=width-1;k>=0;k--)
	{
		for(j=0;j<width;j++)//行交换
			if(JS[k]!=k)
			{
				temp=mat[k*width+j];
				mat[k*width+j]=mat[JS[k]*width+j];
				mat[JS[k]*width+j]=temp;
			}
		for(i=0;i<width;i++)//列交换
			if(IS[k]!=k)
			{
				temp=mat[i*width+k];
				mat[i*width+k]=mat[i*width+IS[k]];
				mat[i*width+IS[k]]=temp;
			}
	}
}