www.pudn.com > vehiclehsv.rar > ImageProcessing.cpp, change:2002-04-02,size:41608b


#include <stdafx.h> 
#include <windows.h> 
#include <math.h> 
#include <malloc.h> 
#include "resource.h" 
#include "BmpTest.h" 
#include "BmpTestDoc.h" 
/********************************************************************   * 
*                                                                       * 
* ro_angle:rotated angle                                                * 
*	size :the size of the rotated bitmap                                  * 
* pav_cent:the rotate centre                                                                      * 
*   **********************************************************************/                                                                   
                                                                   
extern void Sh_Rotate(CBmpTestDoc* pDoc,double ro_angle,double size,POINT *pav_cent)  
{ 
	BYTE *lp=pDoc->m_Dib.m_Buffer,*ltemp; 
	int		 m,n;     
    unsigned int  a,b,d,j,k,stand=150,Height,Width,BPL,i; 
    DWORD    shift,y_value; 
    DWORD    large; 
    double   rotate,x_value,sinth,costh; 
    m=pav_cent->x; 
    n=pav_cent->y;  
	sinth=sin(-ro_angle); 
	Height=(int) pDoc->m_Dib.Height; 
	Width=(int) pDoc->m_Dib.Width; 
    BPL=(int) pDoc->m_Dib.Bpl; 
	costh=cos(-ro_angle); 
//loritemp 是本地的图像副本总指针,ltemp指向数据区,biSize是BITMAPINFORHEADER的大小为40 
	ltemp=new BYTE[9600]; 
    //initioalize the new field orign lpbi->biSizeImage 
    for(j=0;j<9600;j++) 
		*(ltemp+j)=0;  
   //if it is the bmppoint then enlage 
    for(k=0;k<Height;k++) 
    { 
	    for(i=0;i<=(int)(Width/8);i++) 
	    {  
		     if((i+1)*8<Width) 
		     	 d=8; 
		     else 
		     	 d=(unsigned)(Width-i*8); 
	    	 a=128; 
		     if(d!=0) 
		     { 
			     for(j=0;j<d;j++) 
			     { 
  	//if the bit is equal to 0 it is background 
					if((*(lp+(LONG)i+(LONG)k*BPL) & a) ==0) 
			     		a=a>>1; 
					else 
					{ 
		//i is the x coordinate so it is byte 
			     		x_value=(float)8*i+(float)j; 
			     		rotate=(x_value -m)*(costh)-(sinth) *((float)k-n); 
			    		rotate*=(double)110.0/size;   
			    		rotate+=(double)160.0; 
			     		large=(int)(rotate/8); 
			     		shift=(int)rotate-large*8; 
						b=128; 
	//shift>=0 it is in the next store,本身他已处于下一单元(i 从零开始)		     	  
			     		y_value=(LONG)(((x_value -m)*sinth+costh*((float)k-n)) 
			     	 	*110.0/size +(float)120 )*(LONG)(40); 
			     		large +=y_value; 
			     		b=b>>shift; 
			     		if(large < (9600-1)) 
				     		*(ltemp+large)|= b; 
			     		a=a>>1; 
					} 
				 } 
		     }     
	    } 
     } 
 // If the driver did not fill in the biSizeImage field, make one up  
     pDoc->m_Dib.Width=320; 
     pDoc->m_Dib.Height=300; 
     pDoc->m_Dib.Resize(9600); 
	 for(j=0;j<9600;j++) 
	 { 
		 *(pDoc->m_Dib.m_Buffer+j)=*(ltemp+j); 
		 *(pDoc->m_Dib.m_Buffer1+j)=*(ltemp+j); 
	 } 
 	 if(lp) 
	 { 
		delete [] lp; 
		lp=0; 
	 } 
 	 if(ltemp) 
	 { 
		delete [] ltemp; 
		ltemp=0; 
	 } 
} 
   
 
/************************************************************************** 
*  n                                      * 
*  gcc:the gravity contour centre	of the bitmap	                           * 
***************************************************************************/ 
extern BOOL  GCC(CBmpTestDoc* pDoc,double *n,POINT *gcc) 
{ 
//1:xnear 2:xfarer 3:ynear 4:yfarer 
	BYTE *lp=pDoc->m_Dib.m_Buffer; 
  	char name1[2]="1",name2[2]="2",name3[2]="3",name4[2]="4"; 
   	int	i,k,tempx,xImagecon[2000],BPL,black,white,Height,Width; 
    unsigned int		a=1,d,j,maxx,minx; 
    DWORD  x=0,y=0,centx=0,centy=0,counts=0;//counts represents the conunt of the contour 
    //cent represent the center_mass while av_cent represent the average center 
    DISTANCE   xdis[50]; 
    PDISTANCE     pdis; 
    //let x repreterative the near  and y represertive the far  
	POINT Imagecon[100],min,max,last,cent,av_cent,*pcent,*pav_cent; 
		//float distwo; 
	BOOL	maxfind=FALSE,minfind=FALSE,begin=TRUE,equal1,equal2,equal3; 
	double l,dis,midx,midy,d_centx,d_centy; 
	Height=(int) pDoc->m_Dib.Height; 
	Width=(int) pDoc->m_Dib.Width; 
	BPL=(int) pDoc->m_Dib.Bpl; 
    pcent=¢ 
    pav_cent=&av_cent; 
    pdis=xdis; 
    for(i=0;i<50;i++) 
	    (pdis+i)->number=0;     
    for(i=0;i<2000;i++) 
	    xImagecon[i]=0; 
    for(i=0;i<100;i++) 
    { 
	    Imagecon[i].x=0; 
		Imagecon[i].y=0; 
    } 
    for(k=0;k<Height;k++) 
	{ 
	      for(i=0;i<=(int)((Width)/8);i++) 
		  {   
	          a=255; 
			  if(((*(lp+i+k*BPL)) & a) ==0) 
				 continue; 	 
		      if((i+1)*8<(Width)) 
		     	 d=8; 
		      else 
		     	 d=(unsigned)((Width)-i*8); 
		      a=128; 
		      if(d!=0) 
			  { 
			     for(j=0;j<d;j++) 
			     { 
	 
	//if the bit is equal to 0 it is background 
			        if(((*(lp+i+k*BPL)) & a) ==0) 
			     	{ y++; 
			     	 a=a>>1;} 
			        else 
					{   x++; 
                        a=a>>1; 
					    centx=(WORD)i*8+j+centx; 
					    centy=(WORD)k+centy; 
 //分别求出从两个方向的投影点的轮廓 
 //从x方向求 k为其 坐标 
						if(xImagecon[k] ==0) 
							xImagecon[k] =i*8+j;     
						xImagecon[500+k] =max(xImagecon[500+k],i*8+(int)j);     
						xImagecon[k] =min(xImagecon[k],i*8+(int)j); 
 //从y方向求tempx为其坐标  
						tempx=i*8+j; 
						if(xImagecon[1000+tempx]==0)  
							xImagecon[1000+tempx] =k;     
			            xImagecon[1500+tempx] =max(xImagecon[1500+tempx],k); 
						xImagecon[1000+tempx] =min(xImagecon[1000+tempx],k);     
					} 
		         } 
		     }     
	    } 
    }                      
  
 //图的直方图 
    black=y; 
    white=x; 
 //求重心 
	d_centx=((double)centx/x); 
	d_centy=((double) centy/x); 
    for(k=0;k<4;k++) 
    { 
	    max.x=0; 
		max.y=0; 
		min.x=0; 
		min.y=0; 
		j=0; 
		maxx=0; 
		begin=TRUE; 
		maxfind=FALSE; 
		minfind=FALSE; 
	    for(i=1;i<500;i++) 
	    { 
		    //i  值为其 y方向的值 xImagecon 存的是X方向的最近和最远 
		    if(xImagecon[k*500+i] !=0)  
		    { 
			      if(max.x==0)  
			      { 
				       max.x=xImagecon[k*500+i]; 
//the first point is a contour point 
				     	if(j==0) 
						{   j++; 
							if(k<2) 
							{ 
								Imagecon[k*20+j].x=max.x; 
								Imagecon[k*20+j].y=i;  
								j++; 
							} 
					     else 
					     { 
							Imagecon[k*20+j].y=max.x; 
							Imagecon[k*20+j].x=i;  
							j++; 
					     } 
					} 
			} 
	        if(min.x==0)     min.x=xImagecon[k*500+i]; 
				if(!minfind || begin) 
				{ 
				     if(min.x>=xImagecon[k*500+i])   
				     { 
					     min.x=xImagecon[k*500+i]; 
					     min.y=i; 
						 minx=0; 
					 } 
	            else 
	            { 
	           		minx++;  
		          	if((minx>=2) && (min.x<=xImagecon[k*500+i]-2) && (min.x<=max.x+2)) 
		          	{ 
		          		minfind=TRUE; 
		          		begin=FALSE; 
		          	} 
	            } 
            } 
		    if(!maxfind || begin) 
			{ 
			     if(max.x<=xImagecon[k*500+i]) 
			     { 
				     max.x=xImagecon[k*500+i]; 
				     max.y=i;  
			         maxx=0; 
				  } 
				  else   
				  { 
				    	maxx++; 
						if((maxx>=2) && (max.x>=xImagecon[k*500+i]+2) && (max.x>=min.x+2)) 
						{ 
							maxfind=TRUE; 
							begin=FALSE; 
							minfind=FALSE; 
						} 
				   } 
			}  
	        last.x=xImagecon[k*500+i]; 
		    last.y=i; 
	 		if((maxfind && (maxx>=2)) || (minfind && (minx>=2))) 
	 		{		 
	 			if(minfind) 
				{ 
				     if(k==0) 
				     { 
					     Imagecon[k*20+j].x=min.x; 
					     Imagecon[k*20+j].y=min.y;  
					     j++; 
				     } 
			    	 if(k==2) 
				     { 
					     Imagecon[k*20+j].y=min.x; 
					     Imagecon[k*20+j].x=min.y;  
	                     j++; 
					 } 
				     max.x=min.x; 
				     max.y=min.y; 
	                 maxfind=FALSE;					     
               } 
			   if(maxfind) 
			   { 
				     if(k==1) 
				     { 
					     Imagecon[k*20+j].x=max.x; 
					     Imagecon[k*20+j].y=max.y;  
					     j++; 
				     } 
		 			 if(k==3) 
					 { 
					     Imagecon[k*20+j].y=max.x; 
					     Imagecon[k*20+j].x=max.y;  
					     j++; 
					 } 
					 min.x=max.x; 
					 min.y=max.y; 
			     } 
			     *(lp+Imagecon[k*20+j-1].x/8+Imagecon[k*20+j-1].y * BPL/4) ^ 255; 
			     *(lp+Imagecon[k*20+j-1].x/8+Imagecon[k*20+j-1].y * BPL/4) ^ 255; 
			     minx=0;			     
			     maxx=0; 
	    	} 
	    }  
	  } 
	  Imagecon[k*20].x=j; 
	  if(k<2) 
	  { 
		  Imagecon[k*20+j].x=last.x; 
		  Imagecon[k*20+j].y=last.y; 
	  } 
	  else 
	  { 
		  Imagecon[k*20+j].x=last.y; 
		  Imagecon[k*20+j].y=last.x; 
	  } 
	}   
 
//去掉相同点		 
	for(i=0;i<4;i++) 
	{	 
	    k=Imagecon[20*i].x; 
		Imagecon[20*i].x=0; 
		pav_cent->x =Imagecon[20*i+1].x; 
		pav_cent->y =Imagecon[20*i+1].y; 
		last.x=Imagecon[20*i+k].x; 
		last.y=Imagecon[20*i+k].y; 
		for(j=0;j<40;j++) 
		{ 
			if((i==0) && (Imagecon[40+20*(int)(j/20)+2].x<=pav_cent->x) && (Imagecon[40+j].y==pav_cent->y) ) 
				equal1=TRUE; 
			else 
         		equal1=FALSE; 
	        if((i==2) && (Imagecon[20*(int)(j/20)+2].y<=pav_cent->y) && (Imagecon[j].x==pav_cent->x) )		           
				equal2=TRUE; 
			else 
         		equal2=FALSE; 
       		tempx=40*(1-(int)(i/2))+j; 
		    if((Imagecon[tempx].x==pav_cent->x) && (Imagecon[tempx].y==pav_cent->y)) 
				 equal3=TRUE; 
			else 
			 	 equal3=FALSE; 
		    if(equal3 ||	equal1|| equal2) 
			{                               
				 Imagecon[20*i+1].x=0; 
				 Imagecon[20*i+1].y=0; 
			} 
	 	    if	 ((i==0) && (Imagecon[40+20*(int)(j/20)+2].x<=last.x) &&  
			 			 (Imagecon[40+j].y==last.y)) 
				equal1=TRUE; 
			else 
         		 equal1=FALSE; 
 
			if((i==2) &&(Imagecon[20*(int)(j/20)+2].y<=last.y) && (Imagecon[j].x==last.x)) 
				 equal2=TRUE; 
			else  
				 equal2=FALSE; 
		 	tempx=40*(1-i/2)+j; 
			if((Imagecon[tempx].x==last.x) && (Imagecon[tempx].y==last.y))		  
			 	 equal3=TRUE; 
			else 
				 equal3=FALSE; 
			if(equal3 ||	equal1 ||	 equal2)  
			{ 
				 Imagecon[20*i+k].x=0; 
				 Imagecon[20*i+k].y=0; 
			 }	 
       } 
    } 
//caculate the GCC   
	for(j=0;j<20;j++) 
	{ 
		 last.x=Imagecon[20+j].x; 
		 last.y=Imagecon[20+j].y; 
		 Imagecon[20+j].x=Imagecon[60+j].x; 
		 Imagecon[20+j].y=Imagecon[60+j].y; 
		 Imagecon[60+j].x=last.x; 
		 Imagecon[60+j].y=last.y; 
	}                                           
	l=0.; 
    midx=0; 
    midy=0; 
    max.x=0; 
    max.y=0; 
	for(i=0;i<80;i++) 
	{ 
		if(Imagecon[i].x !=0) 
		{ 
			min.x=Imagecon[i].x; 
			min.y=Imagecon[i].y;  
			if(max.x==0) 
			{ 
//max store the initional value				 
				max.x=min.x; 
				max.y=min.y; 
			} 
			j=i+1; 
			while((Imagecon[j].x==0) && (j<80)) 
				j++; 
			if((j>40) && (i<40)) 
			{ 
 
// pav_cent store the end of the y direction maxmun value 
				pav_cent->x=min.x; 
				pav_cent->y=min.y; 
				min.x=max.x; 
				min.y=max.y;  
			}            
			last.x=Imagecon[j].x; 
			last.y=Imagecon[j].y; 
			if(j==80)  
			{ 
			  last.x=pav_cent->x; 
			  last.y=pav_cent->y; 
			} 
			dis=sqrt((double)(min.x-last.x)*(min.x -last.x)+ 
			(double)(min.y-last.y)*(min.y-last.y)); 
			midx =(double)(min.x+last.x)/2 *dis +midx; 
			midy =(double)(min.y+last.y)/2 *dis +midy; 
			l+=dis; 
			i=j-1; 
		} 
 
	} 
    midx=midx/l; 
    midy=midy/l; 
    gcc->x=(int)midx; 
    gcc->y=(int)midy;          
    l=0; 
	for(i=0;i<80;i++) 
	{ 
		if(Imagecon[i].x !=0) 
		{ 
			dis=(Imagecon[i].x-midx)*(Imagecon[i].x-midx)+(Imagecon[i].y-midy)*(Imagecon[i].y-midy); 
			dis=sqrt(dis); 
			l=max(dis,l); 
        } 
	} 
    *n=l; 
	return TRUE; 
}                           
 
 
/*************************************************************************** 
*  pAnew_mom:the query bitmap's moments						                         * 
*  n:the size of the bitmap   					                                   * 
*  protate :the rotated angle            
* the moment name of the calculation below is the same as the maths define * 
***************************************************************************/ 
extern BOOL  Contour(CBmpTestDoc* pDoc,double *pAnew_mon,BOOL bmoment,double n,double *protate) 
{ 
 
//1:xnear 2:xfarer 3:ynear 4:yfarer 
	BYTE *lp=pDoc->m_Dib.m_Buffer; 
  	char	Outtext[300]="",ext1[15]="mom.txt",filename1[30]; 
   	int  i,k,tempx,black,white,BPL,Height,Width; 
    unsigned int		a=1,d,j; 
    DWORD  x=0,y=0,centx=0,centy=0,counts=0;//counts represents the conunt of the contour 
    //cent represent the center_mass while av_cent represent the average center 
    DISTANCE   xdis[50]; 
    PDISTANCE     pdis;  
    //let x repreterative the near  and y represertive the far  
	bool	maxfind=false,minfind=false,begin=true; 
	double  l,mid,d_centx,d_centy,A1,A2,A3,A4,A5,A6,A7,N,coefient, 
     		u02,u20,u03,u30,u11,u21,u12,a1,a2,a3,a4,a5,a6,a7, 
			m01,m02,m03,m20,m30,m11,m21,m12; 
	FILE	*ffp; 
	Height=(int) pDoc->m_Dib.Height; 
	Width=(int) pDoc->m_Dib.Width; 
    BPL=(int) pDoc->m_Dib.Bpl; 
	pdis=xdis;  
    N=n; 
    for(i=0;i<50;i++) 
	    (pdis+i)->number=0;     
    u02=0;u03=0;u20=0; 
//求直方图 
    for(k=0;k<(Height);k++) 
    { 
	    for(i=0;i<=Width/8;i++) 
	    {   
	       a=255; 
		   if(((int)*(lp+i+k*BPL) & a) ==0) 
			   continue; 	 
           if((i+1)*8<Width) 
	        	d=8; 
		   else 
		     	d=(unsigned)(Width-i*8); 
		   a=128; 
		   if(d!=0) 
		   { 
		      for(j=0;j<d;j++) 
			  { 
	//if the bit is equal to 0 it is background 
			     if(((int)*(lp+i+k*BPL) & a) ==0) 
				 {    
					 y++; 
			     	 a=a>>1; 
				 } 
			     else 
				 {   
					 x++; 
                    a=a>>1; 
		            centx=(WORD)i*8+j+centx; 
				    centy=(WORD)k+centy; 
				 } 
		      } 
		    }     
	    } 
    }                      
//图的面积=x+y 
    black=y; 
    white=x;  
//求重心 
    d_centx=((double)centx/x); 
	d_centy=((double) centy/x); 
//归一化系数     
    coefient=(double)x*x;   
    A7=sqrt((double)x); 
//A7=sqrt(a1);  
    mid=coefient*A7; 
    u11=0;u02=0;u03=0;u20=0;u30=0;u12=0;u21=0;  
    m11=0;m01=0;m02=0;m03=0;m20=0;m30=0;m12=0;m21=0; 
	A1=0.;A2=0.;A3=0.;A4=0.;A5=0.;A6=0.;A7=0.; 
    for(k=0;k<Height;k++) 
    { 
	    for(i=0;i<=(int)Width/8;i++) 
	    {   
	       a=255; 
		   if((((int)*(lp+i+k*BPL)) & a)==0) 
			   continue; 	 
    	   if((i+1)*8<Width) 
		       d=8; 
		   else 
		     	d=(unsigned)(Width-i*8); 
		   a=128; 
		   if(d!=0) 
		   { 
			     for(j=0;j<d;j++) 
			     { 
	 
	//if the bit is equal to 0 it is background 
				     if((((int)*(lp+i+k*BPL)) & a) ==0) 
				     	 a=a>>1; 
				     else 
				     {  
						a=a>>1;            
						tempx=i*8+j; 
						a1=((double)k-d_centy)*((double)k-d_centy); 
						u02+=a1; 
						a2=a1*(k- d_centy); 
						u03+=a2; 
						a3=(tempx- d_centx)*(tempx-d_centx); 
						u20+=a3; 
						u30+=a3*(tempx- d_centx);  
						a4=(tempx-d_centx)*(k-d_centy); 
						u11+=a4; 
						a5=a4*(tempx-d_centx); 
						u21+=a5; 
						a6=a4*(k-d_centy); 
						u12+=a6;   
					} 
				} 
			}     
		}    
	//to calculate the scale transform     
		a3=u02/coefient; 
		A1=A1+a3; 
		u02=0.;  
		a3=u20/coefient; 
		A2=A2+a3; 
		u20=0.; 
		a3=u11/coefient; 
		A3=a3+A3; 
		u11=0.; 
		a3=u03/mid; 
		A4=A4+a3; 
		u03=0.; 
		a3=u30/mid; 
		A5=A5+a3; 
		u30=0.; 
		a3=u21/mid; 
		A6=A6+a3; 
		u21=0.; 
		a3=u12/mid; 
		A7=A7+a3; 
		u12=0.; 
	}                      
    u02=A1; 
    u20=A2; 
    u11=A3; 
    u03=A4; 
    u30=A5; 
    u21=A6; 
    u12=A7;    
   
   
 //rotated invariant moment  
    a1=u02+u20;  
    a2=(u20-u02)*(u20-u02)+u11*u11*4; 
    a3=(u30- u12*3)*(u30-u12*3)+(u21*3-u03)*(u21*3-u03); 
    a4=(u30+u12)*(u30+u12)+(u03+u21)*(u03+u21); 
    a5=(u30-u12*3)*(u30+u12)*((u30+u12)*(u30+u12)-(u21+u03)*(u21+u03)*3) 
    	+(u21*3-u03)*(u21+u03)*((u30+u12)*(u30+u12)*3-(u21+u03)*(u21+u03)); 
    a6=(u20-u02)*((u30+u12)*(u30+u12)-(u21+u03)*(u21+u03)) 
    	+u11*(u30+u12)*(u21+u03)*4; 
      		 
    a7=(u21*3-u03)*(u30+u12)*((u30+u12)*(u30+u12)-(u21+u03)*(u21+u03)*3) 
    	-(u30-u12*3)*(u21+u03)*((u30+u12)*(u30+u12)*3 -(u21+u03)*(u21+u03)); 
//calculate the rotate angle  
    u11=2*u11;  
    u12=u20-u02; 
	l=atan2( u11, u12 ); 
	*protate=l/2; 
	*pAnew_mon=(double)n;                                                   
//to store the moment  
	A1=a1; 
	*(++pAnew_mon)=A1; 
	A2=a2;        
	*(++pAnew_mon)=A2; 
	A3=a3;  
	*(++pAnew_mon)=A3; 
	A4=a4;            
	*(++pAnew_mon)=A4; 
	A5=a5; 
	*(++pAnew_mon)=A5; 
	A6=a6;            
	*(++pAnew_mon)=A6; 
	A7=a7;            
	*(++pAnew_mon)=A7; 
    if(bmoment) 
    { 
        wsprintf(filename1,"%s",(LPSTR)ext1);  
		if((ffp=fopen(filename1,"ab"))==NULL) 
		{ 
			return false; 
		} 
		sprintf(Outtext,"\r\n %-20s%-7ld%-25le%-25le%-25le%-25le%-25le%-25le%-25le%-25le", 
			  pDoc->m_Dib.arcFileName,white,n,A1,A2,A3,A4,A5,A6,A7);                  
		fprintf (ffp,"%s",Outtext);	 
       	fclose(ffp);                      
     }  
     return true; 
}                           
 
void HorzMirror(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			lpOutput[x*3+y*nByteWidth]=lpInput[(nWidth-1-x)*3+y*nByteWidth]; 
			lpOutput[x*3+1+y*nByteWidth]=lpInput[(nWidth-1-x)*3+1+y*nByteWidth]; 
			lpOutput[x*3+2+y*nByteWidth]=lpInput[(nWidth-1-x)*3+2+y*nByteWidth]; 
		} 
		Progress(y); 
	} 
} 
 
void VertMirror(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nByteWidth;x++) 
		{ 
			lpOutput[x+y*nByteWidth]=lpInput[x+(nHeight-y-1)*nByteWidth]; 
		} 
		Progress(y); 
	} 
} 
 
void CornerMirror(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			lpOutput[x*3+y*nByteWidth]=lpInput[(nWidth-1-x)*3+(nHeight-y-1)*nByteWidth]; 
			lpOutput[x*3+1+y*nByteWidth]=lpInput[(nWidth-1-x)*3+1+(nHeight-y-1)*nByteWidth]; 
			lpOutput[x*3+2+y*nByteWidth]=lpInput[(nWidth-1-x)*3+2+(nHeight-y-1)*nByteWidth]; 
		} 
		Progress(y); 
	} 
} 
 
#define Point(x,y) lpPoints[(x)+(y)*nWidth] 
#define Point1(x,y) lpPoints1[(x)+(y)*nWidth] 
 
void GetPoints(int nWidth,int nHeight,BYTE *lpBits,BYTE *lpPoints) 
{ 
	int x,y,p; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			p=x*3+y*nByteWidth; 
			lpPoints[x+y*nWidth]=(BYTE)(0.299*(float)lpBits[p+2]+0.587*(float)lpBits[p+1]+0.114*(float)lpBits[p]+0.1); 
		} 
	} 
} 
 
void PutPoints(int nWidth,int nHeight,BYTE *lpBits,BYTE *lpPoints) 
{ 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	int x,y,p,p1; 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			p=x*3+y*nByteWidth; 
			p1=x+y*nWidth; 
			lpBits[p]=lpPoints[p1]; 
			lpBits[p+1]=lpPoints[p1]; 
			lpBits[p+2]=lpPoints[p1]; 
		} 
	} 
} 
 
void HistogramEq1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput) 
{ 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	 
	int r[256],s[256]; 
	ZeroMemory(r,1024); 
	ZeroMemory(s,1024); 
	for(y=0;y<nHeight;y++){ 
		for(x=0;x<nWidth;x++){ 
			r[Point(x,y)]++; 
		} 
	} 
	s[0]=r[0]; 
	for(y=1;y<256;y++) 
	{ 
		s[y]=s[y-1]; 
		s[y]+=r[y]; 
	} 
	for(y=0;y<nHeight;y++){ 
		for(x=0;x<nWidth;x++){ 
			Point(x,y)=s[Point(x,y)]*255/nWidth/nHeight; 
		} 
	} 
 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
	delete lpPoints; 
} 
 
//3x3中值滤波 
void Med1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	BYTE p[9],s; 
	int i,j; 
	for(y=1;y<nHeight-1;y++) 
	{ 
		for(x=3;x<nWidth*3-3;x++) 
		{ 
			p[0]=lpInput[x-3+(y-1)*nByteWidth]; 
			p[1]=lpInput[x+(y-1)*nByteWidth]; 
			p[2]=lpInput[x+3+(y-1)*nByteWidth]; 
			p[3]=lpInput[x-3+y*nByteWidth]; 
			p[4]=lpInput[x+y*nByteWidth]; 
			p[5]=lpInput[x+3+y*nByteWidth]; 
			p[6]=lpInput[x-3+(y+1)*nByteWidth]; 
			p[7]=lpInput[x+(y+1)*nByteWidth]; 
			p[8]=lpInput[x+3+(y+1)*nByteWidth]; 
			for(j=0;j<5;j++) 
			{ 
				for(i=j+1;i<9;i++) 
				{ 
					if (p[j]>p[i]) 
					{ 
						s=p[j]; 
						p[j]=p[i]; 
						p[i]=s; 
					} 
				} 
			} 
			lpOutput[x+y*nByteWidth]=p[4]; 
		} 
		Progress(y); 
	} 
} 
 
void Mean1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y,x1,y1,p; 
//	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
//	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	int sr,sg,sb; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	for(y=1;y<nHeight-1;y++) 
	{ 
		for(x=1;x<nWidth-1;x++) 
		{ 
			p=x*3+y*nByteWidth; 
			sb=0; 
			sg=0; 
			sr=0; 
			for(y1=-1;y1<=1;y1++) 
			for(x1=-1;x1<=1;x1++) 
			{ 
				sb+=lpInput[(y+y1)*nByteWidth+(x+x1)*3]; 
				sg+=lpInput[(y+y1)*nByteWidth+(x+x1)*3+1]; 
				sr+=lpInput[(y+y1)*nByteWidth+(x+x1)*3+2]; 
			} 
			lpOutput[p+2]=sr/9; 
			lpOutput[p+1]=sg/9; 
			lpOutput[p]=sb/9; 
		} 
		Progress(y); 
	} 
//	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
//	delete lpPoints; 
} 
 
 
#define pi (double)3.14159265359 
 
/*复数定义*/ 
typedef struct 
{ 
	double re; 
	double im; 
}COMPLEX; 
 
/*复数加运算*/ 
COMPLEX Add(COMPLEX c1, COMPLEX c2) 
{ 
	COMPLEX c; 
	c.re=c1.re+c2.re; 
	c.im=c1.im+c2.im; 
	return c; 
} 
 
/*复数减运算*/ 
COMPLEX Sub(COMPLEX c1, COMPLEX c2) 
{ 
	COMPLEX c; 
	c.re=c1.re-c2.re; 
	c.im=c1.im-c2.im; 
	return c; 
} 
 
/*复数乘运算*/ 
COMPLEX Mul(COMPLEX c1, COMPLEX c2) 
{ 
	COMPLEX c; 
	c.re=c1.re*c2.re-c1.im*c2.im; 
	c.im=c1.re*c2.im+c2.re*c1.im; 
	return c; 
} 
 
/*快速付里哀变换 
TD为时域值,FD为频域值,power为2的幂数*/ 
void FFT(COMPLEX * TD, COMPLEX * FD, int power) 
{ 
	int count; 
	int i,j,k,bfsize,p; 
	double angle; 
	COMPLEX *W,*X1,*X2,*X; 
 
	/*计算付里哀变换点数*/ 
	count=1<<power; 
	/*分配运算所需存储器*/ 
	W=(COMPLEX *)malloc(sizeof(COMPLEX)*count/2); 
	X1=(COMPLEX *)malloc(sizeof(COMPLEX)*count); 
	X2=(COMPLEX *)malloc(sizeof(COMPLEX)*count); 
	/*计算加权系数*/ 
	for(i=0;i<count/2;i++) 
	{ 
		angle=-i*pi*2/count; 
		W[i].re=cos(angle); 
		W[i].im=sin(angle); 
	} 
	/*将时域点写入存储器*/ 
	memcpy(X1,TD,sizeof(COMPLEX)*count); 
	/*蝶形运算*/ 
	for(k=0;k<power;k++) 
	{ 
		for(j=0;j<1<<k;j++) 
		{ 
			bfsize=1<<(power-k); 
			for(i=0;i<bfsize/2;i++) 
			{ 
				p=j*bfsize; 
				X2[i+p]=Add(X1[i+p],X1[i+p+bfsize/2]); 
				X2[i+p+bfsize/2]=Mul(Sub(X1[i+p],X1[i+p+bfsize/2]),W[i*(1<<k)]); 
			} 
		} 
		X=X1; 
		X1=X2; 
		X2=X; 
	} 
	/*重新排序*/ 
	for(j=0;j<count;j++) 
	{ 
		p=0; 
		for(i=0;i<power;i++) 
		{ 
			if (j&(1<<i)) p+=1<<(power-i-1); 
		} 
		FD[j]=X1[p]; 
	} 
	/*释放存储器*/ 
	free(W); 
	free(X1); 
	free(X2); 
} 
extern void str(int quotient,char addr[15]) 
{ 
	int remainder,i,j,k; 
	char  s1[15]; 
	i=0; 
	do  
	{ 
		remainder=quotient%10; 
		quotient=quotient/10; 
		remainder=remainder+48; 
        s1[i]=(char) remainder; 
		i++; 
	} while(quotient!=0); 
    j=0; 
	for(k=i-1;k>=0;k--) 
	{ 
		addr[j]=s1[k]; 
		j++; 
	} 
	addr[j]=(char) '\0'; 
    return; 
} 
 
/*快速离散余弦变换,利用快速付里哀变换 
f为时域值,F为变换域值,power为2的幂数*/ 
void DCT(double *f, double *F, int power) 
{ 
	int i,count; 
	COMPLEX *X; 
	double s; 
 
	/*计算离散余弦变换点数*/ 
	count=1<<power; 
	/*分配运算所需存储器*/ 
	X=(COMPLEX *)malloc(sizeof(COMPLEX)*count*2); 
	/*延拓*/ 
	memset(X,0,sizeof(COMPLEX)*count*2); 
	/*将时域点写入存储器*/ 
	for(i=0;i<count;i++) 
	{ 
		X[i].re=f[i]; 
	} 
	/*调用快速付里哀变换*/ 
	FFT(X,X,power+1); 
	/*调整系数*/ 
	s=1/sqrt(count); 
	F[0]=X[0].re*s; 
	s*=sqrt(2); 
	for(i=1;i<count;i++) 
	{ 
		F[i]=(X[i].re*cos(i*pi/(count*2))+X[i].im*sin(i*pi/(count*2)))*s; 
	} 
	/*释放存储器*/ 
	free(X); 
} 
 
/*快速沃尔什-哈达玛变换 
f为时域值,F为变换域值,power为2的幂数*/ 
void WALh(double *f, double *W, int power) 
{ 
	int count; 
	int i,j,k,bfsize,p; 
	double *X1,*X2,*X; 
 
	/*计算快速沃尔什变换点数*/ 
	count=1<<power; 
	/*分配运算所需存储器*/ 
	X1=(double *)malloc(sizeof(double)*count); 
	X2=(double *)malloc(sizeof(double)*count); 
	/*将时域点写入存储器*/ 
	memcpy(X1,f,sizeof(double)*count); 
	/*蝶形运算*/ 
	for(k=0;k<power;k++) 
	{ 
		for(j=0;j<1<<k;j++) 
		{ 
			bfsize=1<<(power-k); 
			for(i=0;i<bfsize/2;i++) 
			{ 
				p=j*bfsize; 
				X2[i+p]=X1[i+p]+X1[i+p+bfsize/2]; 
				X2[i+p+bfsize/2]=X1[i+p]-X1[i+p+bfsize/2]; 
			} 
		} 
		X=X1; 
		X1=X2; 
		X2=X; 
	} 
	/*调整系数*/ 
//	for(i=0;i<count;i++) 
//	{ 
//		W[i]=X1[i]/count; 
//	} 
	for(j=0;j<count;j++) 
	{ 
		p=0; 
		for(i=0;i<power;i++) 
		{ 
			if (j&(1<<i)) p+=1<<(power-i-1); 
		} 
		W[j]=X1[p]/count; 
	} 
	/*释放存储器*/ 
	free(X1); 
	free(X2); 
} 
 
 
void Fourier1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int w=1,h=1,wp=0,hp=0; 
	while(w*2<=nWidth) 
	{ 
		w*=2; 
		wp++; 
	} 
	while(h*2<=nHeight) 
	{ 
		h*=2; 
		hp++; 
	} 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	COMPLEX *TD=new COMPLEX[w*h]; 
	COMPLEX *FD=new COMPLEX[w*h]; 
 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			TD[x+w*y].re=Point(x,y); 
			TD[x+w*y].im=0; 
		} 
	} 
 
	for(y=0;y<h;y++) 
	{ 
		FFT(&TD[w*y],&FD[w*y],wp); 
//		Progress(y*nHeight/2/h); 
	} 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			TD[y+h*x]=FD[x+w*y]; 
//			TD[x+w*y]=FD[x*h+y]; 
		} 
	} 
	for(x=0;x<w;x++) 
	{ 
		FFT(&TD[x*h],&FD[x*h],hp); 
		Progress(x*nHeight/2/w+nHeight/2); 
	} 
 
	memset(lpPoints,0,nWidth*nHeight); 
	double m; 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			m=sqrt(FD[x*h+y].re*FD[x*h+y].re+FD[x*h+y].im*FD[x*h+y].im)/100; 
			if (m>255) m=255; 
			Point((x<w/2?x+w/2:x-w/2),nHeight-1-(y<h/2?y+h/2:y-h/2))=(BYTE)(m); 
		} 
	} 
	delete TD; 
	delete FD; 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
	delete lpPoints; 
} 
 
void Walsh1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int w=1,h=1,wp=0,hp=0; 
	while(w*2<=nWidth) 
	{ 
		w*=2; 
		wp++; 
	} 
	while(h*2<=nHeight) 
	{ 
		h*=2; 
		hp++; 
	} 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	double *f=new double[w*h]; 
	double *W=new double[w*h]; 
 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			f[x+y*w]=Point(x,y); 
		} 
	} 
 
	for(y=0;y<h;y++) 
	{ 
		WALh(f+w*y,W+w*y,wp); 
		Progress(y*nHeight/2/h); 
	} 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			f[x*h+y]=W[x+w*y]; 
		} 
	} 
	for(x=0;x<w;x++) 
	{ 
		WALh(f+x*h,W+x*h,hp); 
//		Progress(x*nHeight/2/w+nHeight/2); 
	} 
	double a; 
	memset(lpPoints,0,nWidth*nHeight); 
 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			a=fabs(W[x*h+y]*1000); 
			if (a>255) a=255; 
			Point(x,nHeight-y-1)=(BYTE)a; 
		} 
	} 
	delete f; 
	delete W; 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
	delete lpPoints; 
} 
 
void Dct1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int w=1,h=1,wp=0,hp=0; 
	while(w*2<=nWidth) 
	{ 
		w*=2; 
		wp++; 
	} 
	while(h*2<=nHeight) 
	{ 
		h*=2; 
		hp++; 
	} 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	double *f=new double[w*h]; 
	double *W=new double[w*h]; 
 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			f[x+y*w]=Point(x,y); 
		} 
	} 
 
	for(y=0;y<h;y++) 
	{ 
		DCT(&f[w*y],&W[w*y],wp); 
//		Progress(y*nHeight/2/h); 
	} 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			f[x*h+y]=W[x+w*y]; 
		} 
	} 
	for(x=0;x<w;x++) 
	{ 
		DCT(&f[x*h],&W[x*h],hp); 
		Progress(x*nHeight/2/w+nHeight/2); 
	} 
	double a; 
	memset(lpPoints,0,nWidth*nHeight); 
 
	for(y=0;y<h;y++) 
	{ 
		for(x=0;x<w;x++) 
		{ 
			a=fabs(W[x*h+y]); 
			if (a>255) a=255; 
			Point(x,nHeight-y-1)=(BYTE)(a); 
		} 
	} 
	delete f; 
	delete W; 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
	delete lpPoints; 
} 
 
int Conv2(int nWidth,int nHeight,BYTE *lpOutput,BYTE *lpInput,double mask[],int windowX,int windowY); 
int Laplacian(int nWidth,int nHeight,BYTE *lpOutput,double *lpInput); 
int Zerocross(int nWidth,int nHeight,BYTE *lpOutput,double *lpInput,double thresh); 
int LoG(int nWidth,int nHeight,BYTE *lpOutput,BYTE *lpInput,double thresh,double sigma); 
 
#define In(x,y) lpInput[(x)+(y)*nWidth] 
#define Out(x,y) lpOutput[(x)+(y)*nWidth] 
#define Mediate(x,y) lpMediate[(x)+(y)*nWidth] 
 
#include "NumberDlg.h" 
double GetNumber(const char *prompt,double number,double max) 
{ 
	CNumberDlg dlg; 
	dlg.m_Prompt=prompt; 
	dlg.m_Number=number; 
	dlg.DoModal(); 
	if (dlg.m_Number>max) dlg.m_Number=max; 
	return dlg.m_Number; 
}   
 
void Gauss1(int nWidth,int nHeight,BYTE *lpIn,BYTE *lpOut,void (* Progress)(int Pos)) 
{ 
	double sigma=GetNumber("输入高斯模糊半径:",1,10); 
	int x,y,xx,xk,yk; 
	BYTE *lpInput=new BYTE[nWidth*nHeight]; 
	BYTE *lpOutput=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpIn,lpInput); 
 
	int radius,window; //window of the converlutin kernel(window=2*radius+1) 
	double dev_inv=0.5/(sigma*sigma);  //1/(2*sigma^2) 
	radius = (int)ceil(3*sigma); 
	window = radius*2+1; 
//fast algorithm 
	memcpy(lpOutput,lpInput,nWidth*nHeight); 
	double* lpMediate=new double[nWidth*nHeight]; 
 
	double* mask=new double[window]; 
	double sum=0; 
	double temp; 
	double weight;//sum of weight of mask at edge of the image 
//generate mask 
	for(x=0;x<radius;x++) 
	{ 
		xx=(x-radius)*(x-radius); 
		mask[x]=exp(-xx*dev_inv); 
		mask[window-1-x]=mask[x]; 
		sum+=2*mask[x]; 
	} 
	mask[radius]=1; 
	sum+=1; 
	for(x=0;x<window;x++) 
	{ 
		mask[x]/=sum; 
	} 
//row process 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=radius;x<nWidth-radius;x++) 
		{ 
			temp=0; 
			for(xk=-radius;xk<radius+1;xk++) 
				temp+=In(x+xk,y)*mask[xk+radius]; 
			Mediate(x,y)=temp; 
		} 
		for(x=0;x<radius;x++) 
		{ 
			temp=0; 
			weight=0; 
			for(xk=-x;xk<radius+1;xk++) 
			{ 
				temp+=In(x+xk,y)*mask[xk+radius]; 
				weight+=mask[xk+radius]; 
			} 
			Mediate(x,y)=temp/weight; 
		} 
		for(x=nWidth-radius;x<nWidth;x++) 
		{ 
			temp=0; 
			weight=0; 
			for(xk=-radius;xk<nWidth-x;xk++) 
			{ 
				temp+=In(x+xk,y)*mask[xk+radius]; 
				weight+=mask[xk+radius]; 
			} 
			Mediate(x,y)=temp/weight; 
		} 
		Progress(y/2); 
	} 
//column process 
	double* Column = new double[nHeight]; 
	for(x=0;x<nWidth;x++) 
	{ 
		for(y=radius;y<nHeight-radius;y++) 
		{ 
			temp=0; 
			for(yk=-radius;yk<radius+1;yk++) 
				temp+=Mediate(x,y+yk)*mask[yk+radius]; 
			Column[y]=temp; 
		} 
		for(y=0;y<radius;y++) 
		{ 
			temp=0; 
			weight=0; 
			for(yk=-y;yk<radius+1;yk++) 
			{ 
				temp+=Mediate(x,y+yk)*mask[yk+radius]; 
				weight+=mask[yk+radius]; 
			} 
			Column[y]=temp/weight; 
		} 
        for(y=nHeight-radius;y<nHeight;y++) 
		{ 
			temp=0; 
			weight=0; 
			for(yk=-radius;yk<nHeight-y;yk++) 
			{ 
				temp+=Mediate(x,y+yk)*mask[yk+radius]; 
				weight+=mask[yk+radius]; 
			} 
			Column[y]=temp/weight; 
		} 
 
		for(y=0;y<nHeight;y++) 
		{ 
			Out(x,y)=(int)Column[y]; 
		} 
		Progress(nHeight/2+x*nHeight/2/nWidth); 
	} 
	delete []Column; 
	delete []mask; 
//    Laplacian(nWidth,nHeight,lpOutput,lpMediate); 
	delete []lpMediate; 
//general method(hide) 
	PutPoints(nWidth,nHeight,lpOut,lpOutput); 
	delete lpInput; 
	delete lpOutput; 
} 
 
int Conv2(int nWidth,int nHeight,double *lpOutput,BYTE *lpInput,double mask[],int windowX,int windowY) 
{ 
	int x,y,xk,yk; 
	int radiusX=windowX/2; 
	int radiusY=windowY/2; 
	double temp; 
    for(y=radiusY;y<nHeight-radiusY;y++) 
	{ 
		for(x=radiusX;x<nWidth-radiusX;x++) 
		{ 
			temp=0; 
			for(xk=-radiusX;xk<radiusX+1;xk++) 
				for(yk=-radiusY;yk<radiusY+1;yk++) 
					temp+=In(x+xk,y+yk)*mask[xk+radiusX+(yk+radiusY)*windowX]; 
			Out(x,y)=temp; 
		} 
	}     
	return 0; 
} 
 
int Zerocross(int nWidth,int nHeight,BYTE *lpOutput,double *lpInput,double thresh) 
{ 
	int x,y; 
//zerocross cal 
	for(y=1;y<nHeight-1;y++) 
	{ 
		for(x=1;x<nWidth-1;x++) 
		{ 
			if(In(x,y)<0 && In(x+1,y)>0 && 
				fabs(In(x,y)-In(x+1,y))>thresh) 
			{ 
				Out(x,y)=255; 
			} 
			else if(In(x-1,y)>0 && In(x,y)<0 && 
				fabs(In(x-1,y)-In(x,y))>thresh) 
			{ 
				Out(x,y)=255; 
			} 
			else if(In(x,y)<0 && In(x,y+1)>0 && 
				fabs(In(x,y)-In(x,y+1))>thresh) 
			{ 
				Out(x,y)=255; 
			} 
			else if(In(x,y-1)>0 && In(x,y)<0 && 
				fabs(In(x,y)-In(x,y-1))>thresh) 
			{ 
				Out(x,y)=255; 
			} 
		/*	if((In(x,y)<0 && In(x+1,y)>0)|| 
			   (In(x,y)<0 && In(x,y+1)>0)|| 
			   (In(x,y)>0 && In(x+1,y)<0)|| 
			   (In(x,y)>0 && In(x,y+1)<0)) 
			{ 
				Out(x,y)=255; 
			} 
		*/	else if(In(x,y)==0) 
			{ 
				if(In(x-1,y)<0 && In(x+1,y)>0 && 
				   fabs(In(x-1,y)-In(x+1,y))>2*thresh) 
				{ 
					Out(x,y)=255; 
				} 
			    else if(In(x-1,y)>0 && In(x+1,y)<0 && 
				    fabs(In(x-1,y)-In(x+1,y))>2*thresh) 
				{ 
				    Out(x,y)=255; 
				} 
			    else if(In(x,y-1)>0 && In(x,y+1)<0 && 
			    	fabs(In(x,y-1)-In(x,y+1))>2*thresh) 
				{ 
				    Out(x,y)=255; 
				} 
			    else if(In(x,y+1)>0 && In(x,y-1)<0 && 
				    fabs(In(x,y+1)-In(x,y-1))>2*thresh) 
				{ 
				    Out(x,y)=255; 
				} 
				else  
				{ 
				    Out(x,y)=0; 
				} 
			} 
			else  
			{ 
				Out(x,y)=0; 
			} 
		} 
	} 
	return 0; 
} 
 
void LoG(int nWidth,int nHeight,BYTE *lpIn,BYTE *lpOut,void (* Progress)(int Pos)) 
{ 
 
	BYTE *lpInput=new BYTE[nWidth*nHeight]; 
	BYTE *lpOutput=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpIn,lpInput); 
	double thresh=0,sigma=1; 
 
    int x,y,xx,yy; 
	int radius,window; //window of the converlutin kernel(window=2*radius+1) 
	double dev_inv=0.5/(sigma*sigma);  //1/(2*sigma^2) 
	radius = (int)ceil(3*sigma); 
	window = radius*2+1; 
    //generate mask of convolution (matlab code) 
	//std2 = std*std; 
    //h1 = exp(-(x.*x + y.*y)/(2*std2)); 
     
	double* mask = new double [window*window]; 
	double sum=0; 
	double std2=2*sigma*sigma; 
    for(y=0;y<radius;y++) 
	{ 
		for(x=0;x<radius;x++) 
		{ 
			yy=(y-radius)*(y-radius); 
			xx=(x-radius)*(x-radius); 
			mask[x+y*window]=exp(-(xx+yy)*dev_inv); 
//			mask[window-1-x+y*window]=mask[x+y*window]; 
//			mask[window-1-x+(window-1-y)*window]=mask[x+y*window]; 
//			mask[x+(window-1-y)*window]=mask[x+y*window]; 
			sum=sum+4*mask[x+y*window]; 
		} 
	} 
	x=radius; 
	for(y=0;y<radius;y++) 
	{ 
		yy=(y-radius)*(y-radius); 
		mask[x+y*window]=exp(-yy*dev_inv); 
//		mask[x+(window-1-y)*window]=mask[x+y*window]; 
		sum=sum+2*mask[x+y*window]; 
	} 
	y=radius; 
	for(x=0;x<radius;x++) 
	{ 
		xx=(x-radius)*(x-radius); 
		mask[x+y*window]=exp(-xx*dev_inv); 
//		mask[window-1-x+y*window]=mask[x+y*window]; 
		sum=sum+2*mask[x+y*window]; 
	} 
	mask[radius+radius*window]=1; 
	sum=sum+mask[radius+radius*window]; 
//h = h1.*(x.*x + y.*y - 2*std2)/(2*pi*(std^6));	 
	double denominator=2*3.14159265*pow(sigma,6)*sum; 
	sum=0; 
	for(x=0;x<radius;x++) 
	{ 
		for(y=0;y<radius;y++) 
		{ 
			yy=(y-radius)*(y-radius); 
			xx=(x-radius)*(x-radius); 
			mask[x+y*window]=mask[x+y*window]*(xx+yy-2*std2)/denominator; 
//			mask[window-1-x+y*window]=mask[x+y*window]; 
//			mask[window-1-x+(window-1-y)*window]=mask[x+y*window]; 
//			mask[x+(window-1-y)*window]=mask[x+y*window]; 
			sum=sum+4*mask[x+y*window]; 
		} 
	} 
	x=radius; 
	for(y=0;y<radius;y++) 
	{ 
		yy=(y-radius)*(y-radius); 
		mask[x+y*window]=mask[x+y*window]*(yy-2*std2)/denominator; 
//		mask[x+(window-1-y)*window]=mask[x+y*window]; 
		sum=sum+2*mask[x+y*window]; 
	} 
	y=radius; 
	for(x=0;x<radius;x++) 
	{ 
		xx=(x-radius)*(x-radius); 
		mask[x+y*window]=mask[x+y*window]*(xx-2*std2)/denominator; 
//		mask[window-1-x+y*window]=mask[x+y*window]; 
		sum=sum+2*mask[x+y*window]; 
	} 
	mask[radius+radius*window]=mask[radius+radius*window]*(-2*std2)/denominator; 
	sum=sum+mask[radius+radius*window]; 
//h = h - sum(h(:))/prod(size(h)); % make the filter sum to zero 
	double average=sum/(window*window); 
	for(x=0;x<radius;x++) 
	{ 
		for(y=0;y<radius;y++) 
		{ 
			mask[x+y*window]=mask[x+y*window]-average; 
			mask[window-1-x+y*window]=mask[x+y*window]; 
			mask[window-1-x+(window-1-y)*window]=mask[x+y*window]; 
			mask[x+(window-1-y)*window]=mask[x+y*window]; 
		} 
	} 
	x=radius; 
	for(y=0;y<radius;y++) 
	{ 
		mask[x+y*window]=mask[x+y*window]-average; 
		mask[x+(window-1-y)*window]=mask[x+y*window]; 
	} 
	y=radius; 
	for(x=0;x<radius;x++) 
	{ 
		mask[x+y*window]=mask[x+y*window]-average; 
		mask[window-1-x+y*window]=mask[x+y*window]; 
	} 
	mask[radius+radius*window]=mask[radius+radius*window]-average; 
 
    double* lpMediate = new double[nWidth*nHeight]; 
	//2D convolution 
	if(Conv2(nWidth,nHeight,lpMediate,lpInput,mask,window,window)) 
	{ 
		TRACE("error in Conv2"); 
		return; 
	}; 
//calcaulate the thresh for default 
//  thresh = .75*mean2(abs(b(rr,cc))); 
	sum = 0; 
    for(y=radius;y<nHeight-radius;y++) 
	{ 
		for(x=radius;x<nWidth-radius;x++) 
		{ 
			sum = sum + fabs(Mediate(x,y)); 
		} 
	} 
	thresh=0.5*sum/((nWidth-2*radius)*(nHeight-2*radius)); 
//zerocross cal 
    if(Zerocross(nWidth,nHeight,lpOutput,lpMediate,thresh)) 
	{ 
		TRACE("error in Zerocross"); 
		return; 
	} 
	delete mask; 
	PutPoints(nWidth,nHeight,lpOut,lpOutput); 
	delete lpInput; 
	delete lpOutput; 
} 
 
void BiValue1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput) 
{ 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
	double t=GetNumber("请输入阈值",128,255); 
	BYTE threshold=(BYTE)t; 
	int x,y; 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			if (Point(x,y)>threshold) Point(x,y)=255; 
			else Point(x,y)=0; 
		} 
	} 
 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints); 
	delete lpPoints; 
} 
 
void Sobel(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y,x1,y1,i; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	BYTE *lpPoints1=new BYTE[nWidth*nHeight]; 
	memset(lpPoints1,0,nWidth*nHeight); 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	int d,max; 
	static s[8][9]={ 
		{-1,-2,-1,0,0,0,1,2,1}, 
		{0,-1,-2,1,0,-1,2,1,0}, 
		{1,0,-1,2,0,-2,1,0,-1}, 
		{2,1,0,1,0,-1,0,-1,-2}, 
		{1,2,1,0,0,0,-1,-2,-1}, 
		{0,1,2,-1,0,1,-2,-1,0}, 
		{-1,0,1,-2,0,2,-1,0,1}, 
		{-2,-1,0,-1,0,1,0,1,2} 
	}; 
	for(y=1;y<nHeight-1;y++) 
	{ 
		for(x=1;x<nWidth-1;x++) 
		{ 
			max=0; 
			for(i=0;i<8;i++) 
			{ 
				d=0; 
				for(y1=0;y1<3;y1++) 
				for(x1=0;x1<3;x1++) 
				{ 
					d+=s[i][x1+y1*3]*Point(x+x1-1,y+y1-1); 
				} 
				if (d>max) max=d; 
			} 
			if (max>255) max=255; 
			Point1(x,y)=(BYTE)max; 
		} 
		Progress(y); 
	} 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints1); 
	delete lpPoints; 
	delete lpPoints1; 
} 
 
void Prewitte(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y,x1,y1,i; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	BYTE *lpPoints1=new BYTE[nWidth*nHeight]; 
	memset(lpPoints1,0,nWidth*nHeight); 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	int d,max; 
	static s[8][9]={ 
		{-1,-1,-1,0,0,0,1,1,1}, 
		{0,-1,-1,1,0,-1,1,1,0}, 
		{1,0,-1,1,0,-1,1,0,-1}, 
		{1,1,0,1,0,-1,0,-1,-1}, 
		{1,1,1,0,0,0,-1,-1,-1}, 
		{0,1,1,-1,0,1,-1,-1,0}, 
		{-1,0,1,-1,0,1,-1,0,1}, 
		{-1,-1,0,-1,0,1,0,1,1} 
	}; 
	for(y=1;y<nHeight-1;y++) 
	{ 
		for(x=1;x<nWidth-1;x++) 
		{ 
			max=0; 
			for(i=0;i<8;i++) 
			{ 
				d=0; 
				for(y1=0;y1<3;y1++) 
				for(x1=0;x1<3;x1++) 
				{ 
					d+=s[i][x1+y1*3]*Point(x+x1-1,y+y1-1); 
				} 
				if (d>max) max=d; 
			} 
			if (max>255) max=255; 
			Point1(x,y)=(BYTE)max; 
		} 
		Progress(y); 
	} 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints1); 
	delete lpPoints; 
	delete lpPoints1; 
} 
 
void Roberts(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	BYTE *lpPoints1=new BYTE[nWidth*nHeight]; 
	memset(lpPoints1,0,nWidth*nHeight); 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
 
	int d; 
	for(y=0;y<nHeight-1;y++) 
	{ 
		for(x=0;x<nWidth-1;x++) 
		{ 
			d=2*(abs(Point(x,y)-Point(x+1,y+1))+abs(Point(x+1,y)-Point(x,y+1))); 
			if (d>255) d=255; 
			Point1(x,y)=(BYTE)d; 
		} 
		Progress(y); 
	} 
	PutPoints(nWidth,nHeight,lpOutput,lpPoints1); 
	delete lpPoints; 
	delete lpPoints1; 
} 
 
void PseudoColor1(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
	BYTE b,R,G,B; 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			b=Point(x,y); 
			if (b<64) 
			{ 
				B=255; 
				G=4*b; 
				R=0; 
			} 
			else if (b<128) 
			{ 
				B=(127-b)*4; 
				G=255; 
				R=0; 
			} 
			else if (b<192) 
			{ 
				B=0; 
				G=255; 
				R=(b-128)*4; 
			} 
			else 
			{ 
				B=0; 
				G=(255-b)*4; 
				R=255; 
			} 
 
			lpOutput[x*3+nByteWidth*y]=B; 
			lpOutput[x*3+1+nByteWidth*y]=G; 
			lpOutput[x*3+2+nByteWidth*y]=R; 
		} 
		Progress(y); 
	} 
	delete lpPoints; 
} 
 
BYTE psf2(BYTE x) 
{ 
	double x1=x; 
	return (BYTE)((5-x1/42.5)*x1); 
} 
 
void PseudoColor2(int nWidth,int nHeight,BYTE *lpInput,BYTE *lpOutput,void (* Progress)(int Pos)) 
{ 
	int x,y; 
	BYTE *lpPoints=new BYTE[nWidth*nHeight]; 
	int nByteWidth=nWidth*3; 
	if (nByteWidth%4) nByteWidth+=4-(nByteWidth%4); 
	GetPoints(nWidth,nHeight,lpInput,lpPoints); 
	BYTE b,R,G,B; 
	for(y=0;y<nHeight;y++) 
	{ 
		for(x=0;x<nWidth;x++) 
		{ 
			b=Point(x,y); 
			if (b<85) 
			{ 
				B=psf2(84-b); 
				G=psf2(b); 
				R=0; 
			} 
			else if (b<170) 
			{ 
				B=0; 
				G=psf2(169-b); 
				R=psf2(b-85); 
			} 
			else 
			{ 
				B=psf2(b-170); 
				G=0; 
				R=psf2(255-b); 
			} 
 
			lpOutput[x*3+nByteWidth*y]=B; 
			lpOutput[x*3+1+nByteWidth*y]=G; 
			lpOutput[x*3+2+nByteWidth*y]=R; 
		} 
		Progress(y); 
	} 
	delete lpPoints; 
}