//-----------------------------------------------------------------
// ipcommon.cpp:
//     画像処理一般
//                 Last Update: <2006/10/24 12:05:46 A.Murakami>
//-----------------------------------------------------------------
#include <stdlib.h>
#include <time.h>
#include <math.h>
#include <windows.h>
#include "wingui.h"
#include "ipcommon.h"
//-----------------------------------------------------------------
// ユークリッドの互除法による最小公倍数
//-----------------------------------------------------------------
// ユークリッドの互除法による最小公倍数
int lcm(int a,int b)
{
	int x=a,y=b,c;
	if(a<=0||b<=0) return 1;
	if(a<b) swapv(a,b);
	if(a<b){ c=a; a=b; b=c; }
	while(b!=0){ c=a%b; a=b; b=c; }
	return (x*y/a);
}
//-----------------------------------------------------------------
// 乱数発生 [0..n-1 までの乱数]
//-----------------------------------------------------------------
void init_random() // 乱数の初期化
{
    int seed = time(NULL);
    if(seed<1) seed=1;
    srand(seed);
}
int random_int(int n)
{
	return rint((double)rand()/((double)RAND_MAX+1.0)*n);
}
//-----------------------------------------------------------------
// データの入れ替え
//-----------------------------------------------------------------
void swapv(int c1,int c2)
{
	int c=c1; c1=c2; c2=c;
}
// 点の交換
void swap_point(POINT p1,POINT p2)
{
	POINT p={0,0};
	p.x  = p1.x; p.x  = p1.y;
	p1.x = p2.x; p1.y = p2.y;
	p2.x = p.x;  p2.y = p.y;
}
//-----------------------------------------------------------------
// 色の比較
//-----------------------------------------------------------------
int colcmp(LPBYTE inb,CPIXEL c)
{
	if(*(inb+2)==c.r && *(inb+1)==c.g && *(inb+0)==c.b)
		return 1;
	return 0;
}
//-----------------------------------------------------------------
// 色の設定
//-----------------------------------------------------------------
void setcol(LPBYTE inb,CPIXEL c)
{
	*(inb+2) = c.r;
	*(inb+1) = c.g;
	*(inb+0) = c.b;
}
//-----------------------------------------------------------------
// 距離の計算
//-----------------------------------------------------------------
double calc_distance(POINT p1,POINT p2)
{
	return sqrt((p2.x-p1.x)*(p2.x-p1.x)+(p2.y-p1.y)*(p2.y-p1.y));
}
double calc_distance(int x1,int y1,int x2,int y2)
{
	return sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
}
//-----------------------------------------------------------------
//  平均値 (mean)
//-----------------------------------------------------------------
double calc_mean(int n,double* a)
{
    int i; double mean=0;
    for(i=0;i<n;i++) mean += a[i];
    return mean/n;
}
//-----------------------------------------------------------------
// 標準偏差 (standard deviation)
//-----------------------------------------------------------------
double calc_sd(int n,double* a)
{
    int i;
    double sd=0,mean=calc_mean(n,a);
    for(i=0;i<n;i++) sd+=(a[i]*a[i]);
    sd -= n*mean*mean;
    sd /= (double)n;
    return sqrt(sd);
}
//-----------------------------------------------------------------
// 標準化[平均,標準偏差=0]
//-----------------------------------------------------------------
void standalize_data(int n,double* a)
{
    int i;
    double mv=calc_mean(n,a);
    double sd=calc_sd(n,a);
    for(i=0;i<n;i++)
        a[i] = (a[i]-mv)/sd;
}
//-----------------------------------------------------------------
// 正規化[0..1]
//-----------------------------------------------------------------
void normalize_data(int n,double* a)
{
    int i; double maxv,minv;
	maxv = -10000.; minv = 10000.;
	for(i=0;i<n;i++){
		if(a[i]-maxv > 0) maxv = a[i];
		if(a[i]-minv < 0) minv = a[i];
    }
    // 正規化
    for(i=0;i<n;i++)
		a[i] = maxv*(a[i]-minv) / fabs(maxv-minv);
}
//-----------------------------------------------------------------
// 表色系の変換
//-----------------------------------------------------------------
// RGB から 濃淡値の計算
BYTE rgb2gray(BYTE r,BYTE g,BYTE b)
{
	int gray=rint(r*0.299+g*0.587+b*0.144);
	if(gray>iWHITE) gray=iWHITE;
	return gray;
}
//-----------------------------------------------------------------
// RGB から YUV 色調への計算
//-----------------------------------------------------------------
void rgb2yuv(BYTE r,BYTE g,BYTE b,BYTE &y,BYTE &u,BYTE &v)
{
	double yv,uv,vv;
	yv = (0.2989 *(r+175.45)+0.5871*(g-132.40)+0.1140*(b+221.74));
	uv = (0.5115 *(r+175.45)-0.4283*(g-132.40)-0.0832*(b+221.74));
	vv = (-0.1726*(r+175.45)-0.3389*(g-132.40)+0.5114*(b+221.74));
	if(yv>iWHITE) yv=iWHITE; if(yv<0) yv=0;
	if(uv>iWHITE) uv=iWHITE; if(uv<0) uv=0;
	if(vv>iWHITE) vv=iWHITE; if(vv<0) vv=0;
	y = (BYTE)yv; u = (BYTE)uv; v = (BYTE)vv;
}
BYTE rgb2y(BYTE r,BYTE g,BYTE b)
{
	double yv;
	yv = (0.2989 *(r+175.45)+0.5871*(g-132.40)+0.1140*(b+221.74));
	if(yv>iWHITE) yv=iWHITE; if(yv<0) yv=0;
	return (BYTE)yv;
}
BYTE rgb2u(BYTE r,BYTE g,BYTE b)
{
	double uv;
	uv = (0.5115 *(r+175.45)-0.4283*(g-132.40)-0.0832*(b+221.74));
	if(uv>iWHITE) uv=iWHITE; if(uv<0) uv=0;
	return (BYTE)uv;
}
BYTE rgb2v(BYTE r,BYTE g,BYTE b)
{
	double vv;
	vv = (-0.1726*(r+175.45)-0.3389*(g-132.40)+0.5114*(b+221.74));
	if(vv>iWHITE) vv=iWHITE; if(vv<0) vv=0;
	return (BYTE)vv;
}
//-----------------------------------------------------------------
// RGB から ρ/θ の計算
//-----------------------------------------------------------------
void rgb2rth(BYTE r,BYTE g,BYTE b,BYTE &rho,BYTE &th)
{
	// YUV 色調への変換
	BYTE y,u,v;
	rgb2yuv(r,g,b,y,u,v);
	// 色彩距離(ρ),色相(θ) の計算
	double tmp1,tmp2,rhov,thv;
	tmp1=u-128.; tmp2=v-128.;
	rhov = sqrt(tmp1*tmp1+tmp2*tmp2);
	thv  = rad2deg(atan2(tmp1,tmp2));
	if(rhov<0) rhov=0;   if(rhov>180) rhov=180;
	if(thv<0)  thv+=360; if(thv>360)  thv-=360;
	rhov = (double)iWHITE*(rhov)/(180.);
	thv  = (double)iWHITE*(thv) /(360.);
	rho = (BYTE)rhov;
	th  = (BYTE)thv;
}
BYTE rgb2rho(BYTE r,BYTE g,BYTE b)
{
	// YUV 色調への変換
	BYTE y,u,v;
	rgb2yuv(r,g,b,y,u,v);
	// 色彩距離(ρ)の計算
	double tmp1,tmp2,rhov;
	tmp1=u-128.; tmp2=v-128.;
	rhov  = sqrt(tmp1*tmp1+tmp2*tmp2);
	if(rhov<0) rhov=0;   if(rhov>180) rhov=180;
	rhov = (double)iWHITE*(rhov)/(180.);
	return (BYTE)rhov;
}
BYTE rgb2theta(BYTE r,BYTE g,BYTE b)
{
	// YUV 色調への変換
	BYTE y,u,v;
	rgb2yuv(r,g,b,y,u,v);
	// 色相(θ) の計算
	double tmp1,tmp2,thv;
	tmp1=u-128.; tmp2=v-128.;
	thv = rad2deg(atan2(tmp1,tmp2));
	if(thv<0)  thv+=360; if(thv>360)  thv-=360;
	thv  = (double)iWHITE*(thv) /(360.);
	return (BYTE)thv;
}
//-----------------------------------------------------------------
// 描画関数
//-----------------------------------------------------------------
// 中心:p(x,y),半径:r の円を描く
void draw_circle(LPBYTE oBuf,POINT cp,int radius)
{
	int xs,ys,th,th_max;
	double cirstep;
	// １ピクセル毎の計算となるように
	th_max  = rint(2*M_PI*radius);
	cirstep = 360./th_max;
	for(th=0;th<th_max/2.;th++){
		xs = rint(cp.x + radius * cos(deg2rad(cirstep*th)));
		ys = rint(cp.y + radius * sin(deg2rad(cirstep*th)));
		if(ip_safe(xs,ys,iWidth,iHeight))
			oBuf[XY2PI(xs,ys,iWidth)] = iWHITE;
		xs = rint(cp.x + radius * cos(deg2rad(360-cirstep*th)));
		ys = rint(cp.y + radius * sin(deg2rad(360-cirstep*th)));
		if(ip_safe(xs,ys,iWidth,iHeight))
			oBuf[XY2PI(xs,ys,iWidth)] = iWHITE;
	}
}
//-----------------------------------------------------------------
// p1:(x1,y1), p2:(x2,y2) 間の直線を描く
//-----------------------------------------------------------------
void draw_line(LPBYTE oBuf,POINT p1,POINT p2,int col)
{
	draw_line(oBuf,p1.x,p1.y,p2.x,p2.y,col);
}
void draw_line(LPBYTE oBuf,int x1,int y1,int x2,int y2,int col)
{
	int xp,yp;
    double distance,step,t;
	int dx1,dx2,dy1,dy2;
	dx1 = x1; dx2 = x2; dy1 = y1; dy2 = y2;
	// 描画位置のクリッピング処理
	if(scrossln(&dx1,&y1,&x2,&y2)==0) return;
	// 描画開始
	t = 0.0;
    distance = calc_distance(dx1,dy1,dx2,dy2);
    step = 1.0/(1.5*distance);
    while(t<0.5){
        xp = rint(t*dx2 + (1.0-t)*dx1);
        yp = rint(t*dy2 + (1.0-t)*dy1);
		if(ip_safe(xp,yp,iWidth,iHeight))
			oBuf[XY2PI(xp,yp,iWidth)] = col;
		xp = rint((1.0-t)*dx2 + t*dx1);
        yp = rint((1.0-t)*dy2 + t*dy1);
		if(ip_safe(xp,yp,iWidth,iHeight))
			oBuf[XY2PI(xp,yp,iWidth)] = col;
        t = t + step;
    }
}
//-----------------------------------------------------------------
// 直線と(x=xc)の直線との交点のy座標を求める
//-----------------------------------------------------------------
int croslnxc(int x1,int y1,int x2,int y2,int xc,int *yc)
{
    double wk;
    if(y1 == y2){ *yc = y1; return  1; }
    if(x1 == x2){ *yc = 0;  return -1; }
    if(x1 != x2){
		wk  = (y2-y1)/(x2-x1);
		*yc = rint(wk*(xc-x1)+y1);
		return 1;
    }
    return 0;
}
//-----------------------------------------------------------------
// 直線と(y=yc)の直線との交点のx座標を求める
//-----------------------------------------------------------------
int croslnyc(int x1,int y1,int x2,int y2,int yc,int *xc)
{
    double wk;
    if(x1 == x2){ *xc = x1; return  1; }
    if(y1 == y2){ *xc = 0;  return -1; }
    if(x1 != x2){
		wk  = (y2-y1)/(x2-x1);
		*xc = rint(x1+(yc-y1)/wk);
		return 1;
    }
    return 0;
}
//-----------------------------------------------------------------
// 直線描画領域の２次元クリッピング処理:
//        直線の両端点の座標から実際に描く線分の端点を求める
//-----------------------------------------------------------------
int scrossln(int* x1,int* y1,int* x2,int* y2)
{
	int dx1,dx2,dy1,dy2;
    int i,is,ip,mkp1,mkp2;
    int wk,xc,yc,xw1,xw2,yw1,yw2;
    int xx[2],yy[2],ds1[2],ds2[2];
	int gr_wnx[2],gr_wny[2];
	//--------------------------------------------------
	// 初期化
	//--------------------------------------------------
	dx1 = *x1; dx2 = *x2; dy1 = *y1; dy2 = *y2;
    for(i=0;i<2;i++){
		xx[i] =0; yy[i] =0; ds1[i]=0; ds2[i]=0;
    }
    xw1=gr_wnx[0]=0; xw2=gr_wnx[1]=iWidth;
    yw1=gr_wny[0]=0; yw2=gr_wny[1]=iHeight;
    mkp1=mkp2=0;
	//--------------------------------------------------
    // 直線の端点がウィンドウの枠の内部かどうか
	//--------------------------------------------------
    if(dx1>=xw1 && dx1<=xw2 && dy1>=yw1 && dy1<=yw2){
		*x1  = dx1; *y1 = dy1;
		mkp1 = 1;
    }
    if(dx2>=xw1 && dx2<=xw2 && dy2>=yw1 && dy2<=yw2){
		*x2  = dx2; *y2 = dy2;
		mkp2 = 1;
    }
	//--------------------------------------------------
    // ウィンドウ内であれば終了
	//--------------------------------------------------
    if(mkp1 == 1 && mkp2 == 1) return 1;
	ip = 0;
	//--------------------------------------------------
    // 直線とウィンドウの上下の枠と交差する点を求める
	//--------------------------------------------------
    for(i=0;i<2;i++){
		wk=gr_wnx[i];
		is=croslnxc(dx1,dy1,dx2,dy2,wk,&yc); 
		if(is == 1 && yc>=yw1 && yc<=yw2){
			xx[ip] = wk; yy[ip] = yc;
			ip++;
		}
    }
    //--------------------------------------------------
    // 直線とウィンドウの左右の枠と交差する点を求める
	//--------------------------------------------------
    if(ip<2){
		for(i=0;i<2;i++){
			wk=gr_wny[i];
			is=croslnyc(dx1,dy1,dx2,dy2,wk,&xc); 
			if(is == 1 && xc>=xw1 && xc<=xw2){
				xx[ip] = xc; yy[ip] = wk;
				ip++;
			}
		}	
    }
	//--------------------------------------------------
    // 描く直線が枠と交差しない場合は直線の全てが
	// 枠の外側にあるので、描く必要がないことをしめす
	//--------------------------------------------------
    if(ip==0) return 0;
	//--------------------------------------------------
    // 枠と交差する点を描く直線の対応する端点として代入
	//--------------------------------------------------
    for(i=0;i<2;i++){
		ds1[i] = ((xx[i]-dx1)*(xx[i]-dx1)+(yy[i]-dy1)*(yy[i]-dy1));
		ds2[i] = ((xx[i]-dx2)*(xx[i]-dx2)+(yy[i]-dy2)*(yy[i]-dy2));
    }
    if(mkp1 == 0){
		if(ds1[0]<ds1[1]){
			*x1 = xx[0]; *y1 = yy[0];
		}
		if(ds1[0]>ds1[1]){
			*x1 = xx[1]; *y1 = yy[1];
		}
    }
    if(mkp2 == 0){
		if(ds2[0]<ds2[1]){
			*x2 = xx[0]; *y2 = yy[0];
		}
		if(ds2[0]>ds2[1]){
			*x2 = xx[1]; *y2 = yy[1];
		}
    }
    return 1;
}
//-----------------------------------------------------------------
// テンプレート位置の描画[網掛け]
//-----------------------------------------------------------------
void mark_template_area(LPBYTE inBuf,POINT tp,int txsize,int tysize)
{
	int x,y,xstep=5,ystep=5;
	for(y=tp.y;y<tp.y+tysize;y++) for(x=tp.x;x<tp.x+txsize;x++){
		if(ip_safe(x,y,iWidth,iHeight))
			inBuf[XY2PI(x,y,iWidth)] = iWHITE-inBuf[XY2PI(x,y,iWidth)];
	}
	// 網掛け
	for(y=tp.y;y<tp.y+tysize;y+=ystep){
		draw_line(inBuf,tp.x,y,tp.x+txsize,y,175);
	}
	for(x=tp.x;x<tp.x+txsize;x+=xstep){
		draw_line(inBuf,x,tp.y,x,tp.y+tysize,175);
	}
#if 0
	// 囲む
	for(y=tp.y;y<=tp.y+tysize;y+=tysize){
		draw_line(inBuf,tp.x,y,tp.x+txsize,y,0);
	}
	for(x=tp.x;x<=tp.x+txsize;x+=txsize){
		draw_line(inBuf,x,tp.y,x,tp.y+tysize,0);
	}
#endif
}
//-----------------------------------------------------------------
// 画像処理関数
//-----------------------------------------------------------------
//-----------------------------------------------------------------
// 画像のヒストグラムの作成-->iGHistへ
//-----------------------------------------------------------------
void make_histogram(LPBYTE inBuf)
{
	int x,y;
	FillMemory(iGHist,sizeof(UINT)*iMAX_GRAY,0);
    for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++)
		iGHist[inBuf[XY2PI(x,y,iWidth)]]++;
}
//-----------------------------------------------------------------
// 表示用
//-----------------------------------------------------------------
void GrayToColor(LPBYTE iGray,LPBYTE iColor)
{
	int i,x,y;
	int iadd=(4-(iWidth*3)%4);
	// カラー画像へ
	for(y=0;y<iHeight;y++){
		for(x=0;x<iWidth;x++) {
			FillMemory(iColor+x*3+y*iLength,
					   sizeof(BYTE)*3,iGray[XY2PI(x,y,iWidth)]);
		}
		for(i=0;i<iadd;i++)
			iColor[x*3+y*iLength+i] = 0;
	}
}
//-----------------------------------------------------------------
// 濃淡化
//-----------------------------------------------------------------
void toGray(LPBYTE inBuf)
{
	int i,j; BYTE r,g,b,gray;
	FillMemory(iGHist,sizeof(UINT)*512,0);
	// 濃淡化
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++) {
		b=inBuf[XY2PI(j*3,i,iLength)+0];
		g=inBuf[XY2PI(j*3,i,iLength)+1];
		r=inBuf[XY2PI(j*3,i,iLength)+2];
		gray=rgb2gray(r,g,b);
		iGHist[gray]++;
		inBuf[XY2PI(j,i,iWidth)] = gray;
	}
	FillMemory(inBuf+iSize,sizeof(BYTE)*(biSize-iSize),0);
}
void toCGray(LPBYTE inBuf,int csel)
{
	int i,j; BYTE gray,r,g,b;
	FillMemory(iGHist,sizeof(UINT)*512,0);
	// 濃淡化
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++) {
		gray=inBuf[XY2PI(j*3,i,iLength)+csel];
		// RGB成分
		if(csel <= COLOR_B){
			switch(csel){
			case COLOR_R: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 2]; break;
			case COLOR_G: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 1]; break;
			case COLOR_B: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 0]; break;
			}
		} else {
			r = lpOrgBMP[XY2PI(j*3,i,iLength) + 2];
			g = lpOrgBMP[XY2PI(j*3,i,iLength) + 1];
			b = lpOrgBMP[XY2PI(j*3,i,iLength) + 0];
			switch(csel){
			case COLOR_Y:     gray = rgb2y(r,g,b);     break;
			case COLOR_U:     gray = rgb2u(r,g,b);     break;
			case COLOR_V:     gray = rgb2v(r,g,b);     break;
			case COLOR_RHO:   gray = rgb2rho(r,g,b);   break;
			case COLOR_THETA: gray = rgb2theta(r,g,b); break;
			}
		}
		iGHist[gray]++;
		inBuf[j+i*iWidth] = gray;
	}
	FillMemory(inBuf+iSize,sizeof(BYTE)*(biSize-iSize),0);
}
//-----------------------------------------------------------------
// 濃淡化
//-----------------------------------------------------------------
void toGray()
{
	int i,j; BYTE r,g,b,gray;
	FillMemory(iGHist,sizeof(UINT)*512,0);
	// 濃淡化
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++) {
		b=lpOrgBMP[XY2PI(j*3,i,iLength)+0];
		g=lpOrgBMP[XY2PI(j*3,i,iLength)+1];
		r=lpOrgBMP[XY2PI(j*3,i,iLength)+2];
		gray=rgb2gray(r,g,b);
		iGHist[gray]++;
		FillMemory(lpBMP+XY2PI(j*3,i,iLength),sizeof(BYTE)*3,gray);
	}
}
//-----------------------------------------------------------------
// 濃淡[R,G,B]画像の取得
//-----------------------------------------------------------------
LPBYTE GetGray()
{
	int i,j; BYTE r,g,b,gray;
	LPBYTE iGray=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	FillMemory(iGHist,sizeof(UINT)*512,0);
	// 濃淡化
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++) {
		b=lpOrgBMP[XY2PI(j*3,i,iLength)+0];
		g=lpOrgBMP[XY2PI(j*3,i,iLength)+1];
		r=lpOrgBMP[XY2PI(j*3,i,iLength)+2];
		gray=rgb2gray(r,g,b);
		if(gray>iWHITE) gray=iWHITE;
		iGHist[gray]++;
		iGray[j+i*iWidth] = gray;
	}
	return iGray;
}
LPBYTE GetColor(int csel)
{
	int i,j; BYTE r,g,b,gray;
	LPBYTE iColor=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	FillMemory(iGHist,sizeof(UINT)*512,0);
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++) {
		// RGB成分
		if(csel <= COLOR_B){
			switch(csel){
			case COLOR_R: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 2]; break;
			case COLOR_G: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 1]; break;
			case COLOR_B: gray = lpOrgBMP[XY2PI(j*3,i,iLength) + 0]; break;
			}
		} else {
			r = lpOrgBMP[XY2PI(j*3,i,iLength) + 2];
			g = lpOrgBMP[XY2PI(j*3,i,iLength) + 1];
			b = lpOrgBMP[XY2PI(j*3,i,iLength) + 0];
			switch(csel){
			case COLOR_Y:     gray = rgb2y(r,g,b);     break;
			case COLOR_U:     gray = rgb2u(r,g,b);     break;
			case COLOR_V:     gray = rgb2v(r,g,b);     break;
			case COLOR_RHO:   gray = rgb2rho(r,g,b);   break;
			case COLOR_THETA: gray = rgb2theta(r,g,b); break;
			}
		}
		iColor[j+i*iWidth] = gray;
		iGHist[gray]++;
	}
	return iColor;
}
//-----------------------------------------------------------------
// 線形補間法による画像サイズの縮小・拡大
//-----------------------------------------------------------------
void affine_transform(LPBYTE inBuf,int intp_opt,
					  double zx,double zy,
					  double deg,double px,double py)
// zx,zy: 倍率x,y
// deg  : 回転角(degree)
// px,py: 移動量x,y
{
	int i,j,xc=iWidth/2,yc=iHeight/2;
	double x,y,u,v,r,c,s;
	LPBYTE ib;
	// 画像確保
	ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	FillMemory(inBuf,sizeof(BYTE)*iSize,0);
	// パラメータ計算
	r = deg2rad(deg); c = cos(r); s = sin(r);
	// アフィン変換
	for(i=-yc;i<yc;i++) for(j=-xc;j<xc;j++){
		v = i-py;         u = j-px;
		y = (u*s+v*c)/zy; x = (u*c-v*s)/zx;
		// 補間形式の選択
		inBuf[XY2PI(j+xc,i+yc,iWidth)] =
			image_interpolation(ib,intp_opt,x,y,xc,yc);
	}
	GlobalFree(ib);
}
void c_affine_transform(LPBYTE inCBuf,int intp_opt,
						double zx,double zy,
						double deg,double px,double py)
// zx,zy: 倍率x,y
// deg  : 回転角(degree)
// px,py: 移動量x,y
{
	int i,j,xc=iWidth/2,yc=iHeight/2;
	double x,y,u,v,r,c,s;
	CPIXEL rgb;
	LPBYTE icb;
	// 画像確保
	icb=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*biSize);
	CopyMemory(icb,inCBuf,sizeof(BYTE)*biSize);
	FillMemory(inCBuf,sizeof(BYTE)*biSize,0);
	// パラメータ計算
	r = deg2rad(deg); c = cos(r); s = sin(r);
	// アフィン変換
	for(i=-yc;i<yc;i++) for(j=-xc;j<xc;j++){
		v = i-py;         u = j-px;
		y = (u*s+v*c)/zy; x = (u*c-v*s)/zx;
		// 補間形式の選択
		rgb = c_image_interpolation(icb,intp_opt,x,y,xc,yc);
		inCBuf[XY2PI((j+xc)*3,i+yc,iLength) + 0] = rgb.b;
		inCBuf[XY2PI((j+xc)*3,i+yc,iLength) + 1] = rgb.g;
		inCBuf[XY2PI((j+xc)*3,i+yc,iLength) + 2] = rgb.r;
	}
	GlobalFree(icb);
}
//-----------------------------------------------------------------
// 画像の補間
//-----------------------------------------------------------------
BYTE image_interpolation(LPBYTE ib,int intp_opt,
						 double x,double y,int xc,int yc)
{
	int i,m,n,gray;
	// for linear
	BYTE p1,p2,p3,p4; double p,q;
	// for bicubic
	double c[8],sum;
	// point of interpolation
	if(intp_opt>INTP_NN){
		if(y>0) m=(int)y; else m=(int)(y-1);
		if(x>0) n=(int)x; else n=(int)(x-1);
	} else {
		m = rint(y); n = rint(x);
	}
	// out of area
	if(!ip_safe(n+xc,m+yc,iWidth,iHeight))
		return 0;
	// interpolation
	switch(intp_opt){
	case INTP_NN: // nearest neighbour
		gray = ib[XY2PI(n+xc,m+yc,iWidth)];
		break;
	case INTP_LN: // linear
		if(!ip_safe(n+xc+1,m+yc+1,iWidth,iHeight))
			return 0;
		// 線形補間
		q = y-m; p = x-n;
		p1 = ib[XY2PI(n+xc  ,m+yc  ,iWidth)];
		p2 = ib[XY2PI(n+xc+1,m+yc  ,iWidth)];
		p3 = ib[XY2PI(n+xc  ,m+yc+1,iWidth)];
		p4 = ib[XY2PI(n+xc+1,m+yc+1,iWidth)];
		gray = rint((1.0-q)*((1.0-p)*p1+p*p2) + q*((1.0-p)*p3+p*p4));
		break;
	case INTP_BCC: // bicubic
		if(!ip_safe(n+xc+2,m+yc+2,iWidth,iHeight) ||
		   !ip_safe(n+xc-1,m+yc-1,iWidth,iHeight))
			return 0;
		// 16個の格子点を用いた補間(3次補間)
		c[0]=(m-1)-y; c[1]=m-y;
		c[2]=(m+1)-y; c[3]=(m+2)-y;
		c[4]=(n-1)-x; c[5]=n-x;
		c[6]=(n+1)-x; c[7]=(n+2)-x;
		// 補間関数
		for(i=0;i<8;i++){
			if(c[i]<0) c[i]=-c[i];
			if(0<=c[i] && c[i]<1)
				c[i]=1-2*c[i]*c[i]+c[i]*c[i]*c[i];
			else if(1<=c[i] && c[i]<2)
				c[i]=4-8*c[i]+5*c[i]*c[i]-c[i]*c[i]*c[i];
			else if(2<=c[i]) c[i]=0;
		}
		// 補間濃度
		sum = 0;
		sum += c[0]*c[4]*ib[XY2PI(n+xc-1,m+yc-1,iWidth)];
		sum += c[0]*c[5]*ib[XY2PI(n+xc  ,m+yc-1,iWidth)];
		sum += c[0]*c[6]*ib[XY2PI(n+xc+1,m+yc-1,iWidth)];
		sum += c[0]*c[7]*ib[XY2PI(n+xc+2,m+yc-1,iWidth)];
		sum += c[1]*c[4]*ib[XY2PI(n+xc-1,m+yc  ,iWidth)];
		sum += c[1]*c[5]*ib[XY2PI(n+xc  ,m+yc  ,iWidth)];
		sum += c[1]*c[6]*ib[XY2PI(n+xc+1,m+yc  ,iWidth)];
		sum += c[1]*c[7]*ib[XY2PI(n+xc+2,m+yc  ,iWidth)];
		sum += c[2]*c[4]*ib[XY2PI(n+xc-1,m+yc+1,iWidth)];
		sum += c[2]*c[5]*ib[XY2PI(n+xc  ,m+yc+1,iWidth)];
		sum += c[2]*c[6]*ib[XY2PI(n+xc+1,m+yc+1,iWidth)];
		sum += c[2]*c[7]*ib[XY2PI(n+xc+2,m+yc+1,iWidth)];
		sum += c[3]*c[4]*ib[XY2PI(n+xc-1,m+yc+2,iWidth)];
		sum += c[3]*c[5]*ib[XY2PI(n+xc  ,m+yc+2,iWidth)];
		sum += c[3]*c[6]*ib[XY2PI(n+xc+1,m+yc+2,iWidth)];
		sum += c[3]*c[7]*ib[XY2PI(n+xc+2,m+yc+2,iWidth)];
		if(sum<0) sum = 0; if(sum>iWHITE) sum = iWHITE;
		gray = rint(sum);
		break;
	}
	if(gray<0) gray=0; if(gray>iWHITE) gray=iWHITE;
	return (BYTE)gray;
}
CPIXEL c_image_interpolation(LPBYTE icb,int intp_opt,
						   double x,double y,int xc,int yc)
{
	CPIXEL rgb={0,0,0};
	int i,m,n,gray,r,g,b;
	// for linear
	BYTE p1,p2,p3,p4; double p,q;
	// for bicubic
	double c[8],sum;
	// point of interpolation
	if(intp_opt>INTP_NN){
		if(y>0) m=(int)y; else m=(int)(y-1);
		if(x>0) n=(int)x; else n=(int)(x-1);
	} else {
		m = rint(y); n = rint(x);
	}
	// out of area
	if(!ip_safe(n+xc,m+yc,iWidth,iHeight))
		return rgb;
	// interpolation
	switch(intp_opt){
	case INTP_NN: // nearest neighbour
		r = icb[XY2PI((n+xc)*3,m+yc,iLength)+2];
		g = icb[XY2PI((n+xc)*3,m+yc,iLength)+1];
		b = icb[XY2PI((n+xc)*3,m+yc,iLength)+0];
		break;
	case INTP_LN: // linear
		if(!ip_safe(n+xc+1,m+yc+1,iWidth,iHeight))
			return rgb;
		// 線形補間
		for(i=0;i<3;i++){
			q = y-m; p = x-n;
			p1 = icb[XY2PI((n+xc  )*3,m+yc  ,iLength)+i];
			p2 = icb[XY2PI((n+xc+1)*3,m+yc  ,iLength)+i];
			p3 = icb[XY2PI((n+xc  )*3,m+yc+1,iLength)+i];
			p4 = icb[XY2PI((n+xc+1)*3,m+yc+1,iLength)+i];
			gray = rint((1.0-q)*((1.0-p)*p1+p*p2) + q*((1.0-p)*p3+p*p4));
			switch(i){
			case 0: b = gray; break;
			case 1: g = gray; break;
			case 2: r = gray; break;
			}
		}
		break;
	case INTP_BCC: // bicubic
		if(!ip_safe(n+xc+2,m+yc+2,iWidth,iHeight) ||
		   !ip_safe(n+xc-1,m+yc-1,iWidth,iHeight))
			return rgb;
		// 16個の格子点を用いた補間(3次補間)
		c[0]=(m-1)-y; c[1]=m-y;
		c[2]=(m+1)-y; c[3]=(m+2)-y;
		c[4]=(n-1)-x; c[5]=n-x;
		c[6]=(n+1)-x; c[7]=(n+2)-x;
		// 補間関数
		for(i=0;i<8;i++){
			if(c[i]<0) c[i]=-c[i];
			if(0<=c[i] && c[i]<1)
				c[i]=1-2*c[i]*c[i]+c[i]*c[i]*c[i];
			else if(1<=c[i] && c[i]<2)
				c[i]=4-8*c[i]+5*c[i]*c[i]-c[i]*c[i]*c[i];
			else if(2<=c[i]) c[i]=0;
		}
		// 補間濃度
		for(i=0;i<3;i++){
			sum = 0;
			sum += c[0]*c[4]*icb[XY2PI((n+xc-1)*3,m+yc-1,iLength)+i];
			sum += c[0]*c[5]*icb[XY2PI((n+xc  )*3,m+yc-1,iLength)+i];
			sum += c[0]*c[6]*icb[XY2PI((n+xc+1)*3,m+yc-1,iLength)+i];
			sum += c[0]*c[7]*icb[XY2PI((n+xc+2)*3,m+yc-1,iLength)+i];
			sum += c[1]*c[4]*icb[XY2PI((n+xc-1)*3,m+yc  ,iLength)+i];
			sum += c[1]*c[5]*icb[XY2PI((n+xc  )*3,m+yc  ,iLength)+i];
			sum += c[1]*c[6]*icb[XY2PI((n+xc+1)*3,m+yc  ,iLength)+i];
			sum += c[1]*c[7]*icb[XY2PI((n+xc+2)*3,m+yc  ,iLength)+i];
			sum += c[2]*c[4]*icb[XY2PI((n+xc-1)*3,m+yc+1,iLength)+i];
			sum += c[2]*c[5]*icb[XY2PI((n+xc  )*3,m+yc+1,iLength)+i];
			sum += c[2]*c[6]*icb[XY2PI((n+xc+1)*3,m+yc+1,iLength)+i];
			sum += c[2]*c[7]*icb[XY2PI((n+xc+2)*3,m+yc+1,iLength)+i];
			sum += c[3]*c[4]*icb[XY2PI((n+xc-1)*3,m+yc+2,iLength)+i];
			sum += c[3]*c[5]*icb[XY2PI((n+xc  )*3,m+yc+2,iLength)+i];
			sum += c[3]*c[6]*icb[XY2PI((n+xc+1)*3,m+yc+2,iLength)+i];
			sum += c[3]*c[7]*icb[XY2PI((n+xc+2)*3,m+yc+2,iLength)+i];
			if(sum<0) sum = 0; if(sum>iWHITE) sum = iWHITE;
			gray = rint(sum);
			switch(i){
			case 0: b = gray; break;
			case 1: g = gray; break;
			case 2: r = gray; break;
			}
		}
		break;
	}
	if(r<0) r=0; if(r>iWHITE) r=iWHITE;
	if(g<0) g=0; if(g>iWHITE) g=iWHITE;
	if(b<0) b=0; if(b>iWHITE) b=iWHITE;
	rgb.r = r; rgb.g = g; rgb.b = b;
	return rgb;
}
//-----------------------------------------------------------------
// 最近傍法による画像サイズの縮小・拡大[nearest neighbour]
//-----------------------------------------------------------------
void nn_scale(LPBYTE inBuf,double zx,double zy)
{
	int i,j,m,n;
	int xsize=(zx>1.0)?iWidth : rint(iWidth*zx);
	LPBYTE ib;
	// 画像確保
	ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	FillMemory(inBuf,sizeof(BYTE)*iSize,0);
	// 変換
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++){
		m = rint(i/zy); n = rint(j/zx);
		if(ip_safe(n,m,iWidth,iHeight))
			inBuf[XY2PI(j,i,xsize)] = ib[XY2PI(n,m,iWidth)];
		else
			inBuf[XY2PI(j,i,xsize)] = 0;
	}
	GlobalFree(ib);
}
//-----------------------------------------------------------------
// 線形補間法による画像サイズの縮小・拡大
//-----------------------------------------------------------------
void linear_scale(LPBYTE inBuf,double zx,double zy)
{
	int i,j,m,n;
	int p1,p2,p3,p4,d;
	int xsize=(zx>1.0)?iWidth:rint(iWidth*zx);
	double x,y,p,q;
	LPBYTE ib;
	// 画像確保
	ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	FillMemory(inBuf,sizeof(BYTE)*iSize,0);
	// 変換
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++){
		y = i/zy;     x = j/zx;
		m = (int)(y);  n = (int)(x);
		q = y-m;      p = x-n;
		if(ip_safe(n,m,iWidth,iHeight) && 
			ip_safe(n+1,m+1,iWidth,iHeight)){
			p1 = ib[XY2PI(n  ,m  ,iWidth)];
			p2 = ib[XY2PI(n+1,m  ,iWidth)];
			p3 = ib[XY2PI(n  ,m+1,iWidth)];
			p4 = ib[XY2PI(n+1,m+1,iWidth)];
			d = rint((1.0-q)*((1.0-p)*p1+p*p2) + q*((1.0-p)*p3+p*p4));
			if(d<0) d=0; if(d>iWHITE) d=iWHITE;
		} else d = 0;
		inBuf[XY2PI(j,i,xsize)] = d;
	}
	GlobalFree(ib);
}
//-----------------------------------------------------------------
// バイキュービック[3次補間]法による画像の拡大・縮小
//-----------------------------------------------------------------
void bcc_scale(LPBYTE inBuf,double zx,double zy)
{
	int i,j,k,m,n,d;
	int xsize=(zx>1.0)?iWidth : rint(iWidth*zx);
	double x,y,c[8],sum;
	LPBYTE ib;
	// 画像確保
	ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	FillMemory(inBuf,sizeof(BYTE)*iSize,0);
	// 変換
	for(i=0;i<iHeight;i++) for(j=0;j<iWidth;j++){
		y = i/zy;     x = j/zx;
		m = (int)(y); n = (int)(x);
		if(ip_safe(n-1,m-1,iWidth,iHeight) && 
			ip_safe(n+2,m+2,iWidth,iHeight)){
			// 16個の格子点を用いた補間(３次補間)
			c[0]=(m-1)-y; c[1]=m-y;
			c[2]=(m+1)-y; c[3]=(m+2)-y;
			c[4]=(n-1)-x; c[5]=n-x;
			c[6]=(n+1)-x; c[7]=(n+2)-x;
			// 補間関数
			for(k=0;k<8;k++){
				if(c[k]<0) c[k]=-c[k];
				if(0<=c[k] && c[k]<1)
					c[k]=1-2*c[k]*c[k]+c[k]*c[k]*c[k];
				else if(1<=c[k] && c[k]<2)
					c[k]=4-8*c[k]+5*c[k]*c[k]-c[k]*c[k]*c[k];
				else if(2<=c[k]) c[k]=0;
			}
			// 補間濃度
			sum = 0;
			sum += c[0]*c[4]*ib[XY2PI(n-1,m-1,iWidth)];
			sum += c[0]*c[5]*ib[XY2PI(n  ,m-1,iWidth)];
			sum += c[0]*c[6]*ib[XY2PI(n+1,m-1,iWidth)];
			sum += c[0]*c[7]*ib[XY2PI(n+2,m-1,iWidth)];
			sum += c[1]*c[4]*ib[XY2PI(n-1,m  ,iWidth)];
			sum += c[1]*c[5]*ib[XY2PI(n  ,m  ,iWidth)];
			sum += c[1]*c[6]*ib[XY2PI(n+1,m  ,iWidth)];
			sum += c[1]*c[7]*ib[XY2PI(n+2,m  ,iWidth)];
			sum += c[2]*c[4]*ib[XY2PI(n-1,m+1,iWidth)];
			sum += c[2]*c[5]*ib[XY2PI(n  ,m+1,iWidth)];
			sum += c[2]*c[6]*ib[XY2PI(n+1,m+1,iWidth)];
			sum += c[2]*c[7]*ib[XY2PI(n+2,m+1,iWidth)];
			sum += c[3]*c[4]*ib[XY2PI(n-1,m+2,iWidth)];
			sum += c[3]*c[5]*ib[XY2PI(n  ,m+2,iWidth)];
			sum += c[3]*c[6]*ib[XY2PI(n+1,m+2,iWidth)];
			sum += c[3]*c[7]*ib[XY2PI(n+2,m+2,iWidth)];
			if(sum<0) sum = 0; if(sum>iWHITE) sum = iWHITE;
			d = rint(sum);
		} else d=0;
		inBuf[XY2PI(j,i,xsize)] = d;
	}
}
//-----------------------------------------------------------------
// 平均画素法による画像サイズの縮小
// [兼ローパスフィルタ]
//-----------------------------------------------------------------
void sqm_scale(LPBYTE inBuf,double zx,double zy)
{
	int i,x,y,x1,x2,y1,y2;
	int lcm_x,lcm_y,p_count,meanp;
	int xsize,ysize,oxsize,oysize;
	LPBYTE oneline;
	// 画像サイズ[偶数に調節]
	xsize = oxsize = rint(iWidth*zx);
	ysize = oysize = rint(iHeight*zy);
	if(xsize%2) xsize--; if(ysize%2) ysize--;
	// 最小公倍数の計算
	lcm_x = lcm(iWidth ,xsize);
	lcm_y = lcm(iHeight,ysize);
	x1 = lcm_x / iWidth; y1 = lcm_y / iHeight;
	x2 = lcm_x / xsize;  y2 = lcm_y / ysize;
	// 描画サイズは元のまま
	if(xsize>iWidth ) xsize = oxsize = iWidth;
	if(ysize>iHeight) ysize = oysize = iHeight;
	//--------------------------------------------------
	// x 方向
	//--------------------------------------------------
	oneline=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iWidth*x1);
	for(y=0;y<iHeight;y++){
		for(x=0;x<iWidth;x++)
			for(i=0;i<x1;i++)
				oneline[x*x1+i] = inBuf[XY2PI(x,y,iWidth)];
		for(x=0;x<xsize;x++) {
			meanp = p_count = 0;
			for(i=0;i<x2;i++){
				if(x*x2+i < x1*iWidth){
					p_count ++;
					meanp += oneline[x*x2+i];
				}
			}
			if(p_count) meanp /= p_count;
			inBuf[XY2PI(x,y,oxsize)] = meanp;
		}
		for(x=xsize;x<iWidth;x++)
			inBuf[XY2PI(x,y,oxsize)] = 0;
	}
	GlobalFree(oneline);
	//--------------------------------------------------
	// y 方向
	//--------------------------------------------------
	oneline=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iHeight*y1);
	for(x=0;x<xsize;x++){
		for(y=0;y<iHeight;y++)
			for(i=0;i<y1;i++)
				oneline[y*y1+i] = inBuf[XY2PI(x,y,oxsize)];
		for(y=0;y<ysize;y++) {
			meanp = p_count = 0;
			for(i=0;i<y2;i++){
				if(y*y2+i < y1*iHeight){
					p_count ++;
					meanp += oneline[y*y2+i];
				}
			}
			if(p_count) meanp /= p_count;
			inBuf[XY2PI(x,y,oxsize)] = meanp;
		}
		for(y=ysize;y<iHeight;y++)
			inBuf[XY2PI(x,y,oxsize)] = 0;
	}
	GlobalFree(oneline);
}
void c_sqm_scale(LPBYTE inBuf,double zx,double zy)
{
	int i,x,y,x1,x2,y1,y2;
	int lcm_x,lcm_y,p_count,meanp[3];
	int bxsize,xsize,ysize,oxsize,oysize;
	LPBYTE oneline;
	// 画像サイズ[偶数に調節]
	xsize = oxsize = rint(iWidth*zx);
	ysize = oysize = rint(iHeight*zy);
	if(xsize%2) xsize--; if(ysize%2) ysize--;
	// 最小公倍数の計算
	lcm_x = lcm(iWidth ,xsize);
	lcm_y = lcm(iHeight,ysize);
	x1 = lcm_x / iWidth; y1 = lcm_y / iHeight;
	x2 = lcm_x / xsize;  y2 = lcm_y / ysize;
	// 描画サイズは元のまま
	if(xsize>iWidth ) xsize = oxsize = iWidth;
	if(ysize>iHeight) ysize = oysize = iHeight;
	// ビットマップ用サイズ
	if((oxsize*3)%4==0) bxsize = oxsize*3;
	else bxsize = oxsize*3+(4-(oxsize*3)%4);
	//--------------------------------------------------
	// x 方向
	//--------------------------------------------------
	oneline=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iLength*x1);
	for(y=0;y<iHeight;y++){
		for(x=0;x<iWidth;x++)
			for(i=0;i<x1;i++){
				oneline[(x*x1+i)*3+0] = inBuf[XY2PI(x*3,y,iLength)+0];
				oneline[(x*x1+i)*3+1] = inBuf[XY2PI(x*3,y,iLength)+1];
				oneline[(x*x1+i)*3+2] = inBuf[XY2PI(x*3,y,iLength)+2];
			}
		for(x=0;x<xsize;x++) {
			p_count = 0;
			for(i=0;i<3;i++) meanp[i] = 0;
			for(i=0;i<x2;i++){
				if((x*x2+i)*3 < x1*iLength){
					p_count ++;
					meanp[0] += oneline[(x*x2+i)*3+0];
					meanp[1] += oneline[(x*x2+i)*3+1];
					meanp[2] += oneline[(x*x2+i)*3+2];
				}
			}
			if(p_count) {
				meanp[0] /= p_count;
				meanp[1] /= p_count;
				meanp[2] /= p_count;
			}
			inBuf[XY2PI(x*3,y,bxsize)+0] = meanp[0];
			inBuf[XY2PI(x*3,y,bxsize)+1] = meanp[1];
			inBuf[XY2PI(x*3,y,bxsize)+2] = meanp[2];
		}
		for(x=xsize*3;x<bxsize;x++)
			inBuf[XY2PI(x,y,bxsize)] = 0;
	}
	GlobalFree(oneline);
	//--------------------------------------------------
	// y 方向
	//--------------------------------------------------
	oneline=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iHeight*y1*3);
	for(x=0;x<xsize;x++){
		for(y=0;y<iHeight;y++)
			for(i=0;i<y1;i++){
				oneline[(y*y1+i)*3+0] = inBuf[XY2PI(x*3,y,bxsize)+0];
				oneline[(y*y1+i)*3+1] = inBuf[XY2PI(x*3,y,bxsize)+1];
				oneline[(y*y1+i)*3+2] = inBuf[XY2PI(x*3,y,bxsize)+2];
			}
		for(y=0;y<ysize;y++) {
			for(i=0;i<3;i++) meanp[i] = 0;
			p_count = 0;
			for(i=0;i<y2;i++){
				if((y*y2+i)*3 < y1*iHeight*3){
					p_count ++;
					meanp[0] += oneline[(y*y2+i)*3+0];
					meanp[1] += oneline[(y*y2+i)*3+1];
					meanp[2] += oneline[(y*y2+i)*3+2];
				}
			}
			if(p_count){
				meanp[0] /= p_count;
				meanp[1] /= p_count;
				meanp[2] /= p_count;
			}
			inBuf[XY2PI(x*3,y,bxsize)+0] = meanp[0];
			inBuf[XY2PI(x*3,y,bxsize)+1] = meanp[1];
			inBuf[XY2PI(x*3,y,bxsize)+2] = meanp[2];
		}
		for(y=ysize;y<iHeight;y++) {
			inBuf[XY2PI(x*3,y,bxsize)+0] = 0;
			inBuf[XY2PI(x*3,y,bxsize)+1] = 0;
			inBuf[XY2PI(x*3,y,bxsize)+2] = 0;
		}
	}
	GlobalFree(oneline);
}
//-----------------------------------------------------------------
// OR 処理
//-----------------------------------------------------------------
void IPfunc_OR(LPBYTE iDst,LPBYTE iAdd)
{
	int x,y;
	for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++) {
		iDst[XY2PI(x,y,iWidth)] = (iDst[XY2PI(x,y,iWidth)]  |
								   iAdd[XY2PI(x,y,iWidth)]) ? iWHITE:iBLACK;
	}
}
//-----------------------------------------------------------------
// メディアンフィルタ
//-----------------------------------------------------------------
void IPfunc_median(LPBYTE inBuf)
{
	int i,x,y,px,py;
	BYTE c[9];
	LPBYTE ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
    CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	// フィルタリング開始
	for(y=1;y<iHeight-1;y++) for(x=1;x<iWidth-1;x++){
		for(i=0;i<8;i++){
			px = x+offset8[i].x;
			py = y+offset8[i].y;
			c[i] = ib[XY2PI(px,py,iWidth)];
		}
		inBuf[XY2PI(x,y,iWidth)] = median(c);
	}
	GlobalFree(ib);
}
BYTE median(BYTE c[9])
{
	int i,j;
	for(j=0;j<8;j++) for(i=0;i<8;i++){
		if(c[i+1]<c[i]) swapv(c[i],c[i+1]);
	}
	return c[4];
}
//-----------------------------------------------------------------
// 組織的ディザ法によるディザ画像の作成
//-----------------------------------------------------------------
void toDither(LPBYTE inBuf)
{
	static int DITHER_BSZ = 4;  // ブロックの画素数
	static int DITHER_NL  = 16; // 擬似階調数(=4^2)
    double width;               // 16段階画像の階調値の単位幅
    BYTE new_gray;              // 新しい階調(16階調での値)
    int x_block,y_block;        // 横・縦のブロック数
    int x,y,i,j,m,n,px,py;
    static int dither_matrix[4][4]={ // Bayer 型ディザ行列
        { 0, 8, 2,10},
        {12, 4,14, 6},
        { 3,11, 1, 9},
        {15, 7,13, 5}
    };
	//--------------------------------------------------
    // 16階調の画像の作成
	//--------------------------------------------------
    width = iWHITE/(double)DITHER_NL;
    for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++){
		new_gray=rint(inBuf[XY2PI(x,y,iWidth)]/width);
		if(new_gray > DITHER_NL-1)
			new_gray = DITHER_NL-1;
		inBuf[XY2PI(x,y,iWidth)] = new_gray;
    }
	//--------------------------------------------------
    // ディザ画像の作成
	//--------------------------------------------------
    x_block=iWidth  / DITHER_BSZ;  // 横のブロック数
    y_block=iHeight / DITHER_BSZ;  // 縦のブロック数
	if(iWidth  % DITHER_BSZ) x_block++;
	if(iHeight % DITHER_BSZ) y_block++;
    for(i=0;i<y_block;i++) for(j=0;j<x_block;j++){
		x=j*DITHER_BSZ;
		y=i*DITHER_BSZ;
		for(m=0;m<DITHER_BSZ;m++) for(n=0;n<DITHER_BSZ;n++){
			px=x+n; py=y+m;
			if(ip_safe(px,py,iWidth,iHeight)){
				if(inBuf[XY2PI(px,py,iWidth)] <= dither_matrix[m][n])
					inBuf[XY2PI(px,py,iWidth)]=0;
				else inBuf[XY2PI(px,py,iWidth)]=iWHITE;
			}
		}
	}
}
//-----------------------------------------------------------------
// 画素の階調分布に対して線形変換
//-----------------------------------------------------------------
void linear_transform(LPBYTE inBuf)
{
	int x,y,minv,maxv; BYTE gray;
    // 階調値の最小値，最大値を求める
    minv = 10000;  maxv = 0;
    for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++){
		if(inBuf[XY2PI(x,y,iWidth)]<minv) minv=inBuf[XY2PI(x,y,iWidth)];
		if(inBuf[XY2PI(x,y,iWidth)]>maxv) maxv=inBuf[XY2PI(x,y,iWidth)];
    }
    // 階調を線形変換して image2[y][x] に代入
    for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++){
		gray = rint(iWHITE*(inBuf[XY2PI(x,y,iWidth)]-minv)/
					(double)(maxv-minv));
		inBuf[XY2PI(x,y,iWidth)] = gray;
	}
}
//-----------------------------------------------------------------
// 表現階調数の変更
//    GetGray(),GetColor(): 等で計算したヒストグラムを使用
//-----------------------------------------------------------------
void change_graylevel(LPBYTE inBuf,int iLevel)
{
    int trhist[512],trans_table[512];
    int i,x,y,gray;
    int target_value; // 変換後の頻度の目標値
    double gray_step; // 表現階調間隔
    // ヒストグラム変換表の作成
    target_value = rint(iSize/(double)iLevel);
	FillMemory(trhist,sizeof(int)*512,0);
	FillMemory(trans_table,sizeof(int)*512,0);
    gray = 0;
    for(i=0;i<iMAX_GRAY;i++){
        if(abs(target_value-trhist[gray]) <
		   abs(target_value-(trhist[gray] + iGHist[i]))){
            gray ++;
            if(gray>=iLevel) gray=iLevel-1;
        }
        trans_table[i] = gray;
        trhist[gray]  += iGHist[i];
	}
	// 平滑化した画像の作成
    gray_step = (double)(iWHITE+1) / iLevel;
    for(y=0;y<iHeight;y++) for(x=0;x<iWidth;x++)
		inBuf[XY2PI(x,y,iWidth)] = (BYTE)rint(
			trans_table[inBuf[XY2PI(x,y,iWidth)]] * gray_step);
}
//-----------------------------------------------------------------
// 画像の空間フィルタリングを行う
//-----------------------------------------------------------------
void spacial_filtering(LPBYTE inBuf,int iFlt[9],double iMag)
{
	LPBYTE ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
    int i,x,y,px,py;
	double new_value,minv,maxv;
	// 画像の複製
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	// 最小・最大値の計算[線形変換用]
	minv = 10000.; maxv = 0.;
    for(y=1;y<iHeight-1;y++) for(x=1;x<iWidth-1;x++){
		new_value = 0.0;
		for(i=0;i<9;i++){
			px = x+offset9[i].x;
			py = y+offset9[i].y;
			new_value += iFlt[i] * ib[XY2PI(px,py,iWidth)];
		}
		new_value = new_value / iMag;
		if(new_value < minv) minv = new_value;
		if(new_value > maxv) maxv = new_value;
	}
    // フィルタリングを行う
    for(y=1;y<iHeight-1;y++) for(x=1;x<iWidth-1;x++){
		new_value = 0.0;
		for(i=0;i<9;i++){
			px = x+offset9[i].x;
			py = y+offset9[i].y;
			new_value += iFlt[i]*ib[XY2PI(px,py,iWidth)];
		}
		new_value = new_value / iMag;
		new_value = iWHITE / (maxv-minv) * (new_value-minv);
		inBuf[XY2PI(x,y,iWidth)] = (BYTE)rint(new_value);
	}
	GlobalFree(ib);
}
//-----------------------------------------------------------------
// 誤差拡散画像の作成
// 誤差分散 |  *  | 3/8 |    * : 注目点,
//          | 3/8 | 2/8 |    注）最右列, 最下行は誤差拡散しない
//-----------------------------------------------------------------
void make_error_difusion_image(LPBYTE inBuf)
{
    int x, y;
    double err1,err2,err3; // 階調誤差
    double* buf=(double*)GlobalAlloc(GPTR,sizeof(double)*iWidth*2);
    double modified_gray;
	// ２ライン分のバッファのクリア
	FillMemory(buf,sizeof(double)*iWidth*2,0);
	// 誤差拡散画像の作成
    for(y=0;y<iHeight;y++){
        for(x=0;x<iWidth;x++){
            // 近傍が外に出るときは単純２値化する
            if(x==iWidth-1 || y==iHeight-1){
                if(inBuf[XY2PI(x,y,iWidth)] <= iMAX_GRAY/2) 
                    inBuf[XY2PI(x,y,iWidth)]=iBLACK;
                else
					inBuf[XY2PI(x,y,iWidth)]=iWHITE;
            } else {
                modified_gray=inBuf[XY2PI(x,y,iWidth)]+buf[x];
                err1=modified_gray - 0.0;
                err2=modified_gray - iWHITE;
                if(err1 * err1<err2 * err2){
                    inBuf[XY2PI(x,y,iWidth)]=0;   
                    err3=err1;
                } else {
                    inBuf[XY2PI(x,y,iWidth)]=iWHITE; 
                    err3=err2;
                }
                // 近傍に応じて次の式を変える
                buf[x+1]       =buf[x+1]       + err3 * 3.0 / 8.0;
                buf[x+1+iWidth]=buf[x+1+iWidth]+ err3 * 2.0 / 8.0;
                buf[x  +iWidth]=buf[x  +iWidth]+ err3 * 3.0 / 8.0;
            }
        }
        for(x=0;x<iWidth;x++){
            buf[x]=buf[x+iWidth];
            buf[x+iWidth]=0.0;
        }
    }
	GlobalFree(buf);
}
//-----------------------------------------------------------------
// 画像のヒストグラムを作り画像として表示
//-----------------------------------------------------------------
void make_histogram_image(LPBYTE inBuf)
{
	int i,j,data,max_frequency;
    // ヒストグラムの作成
	make_histogram(inBuf);
    // 頻度の最大値を求める
    max_frequency=0;
    for(i=0;i<iMAX_GRAY;i++)
        if(iGHist[i] > max_frequency) 
            max_frequency = iGHist[i];
    // ヒストグラム画像の作成
    FillMemory(inBuf,sizeof(BYTE)*iSize,0);
    for(i=0;i<iMAX_GRAY;i++){
		// グラフ長さ
        data=rint(iWHITE*(double)iGHist[i]/(double)max_frequency);
        for(j=0;j<data;j++)
			inBuf[XY2PI(i,j,iWidth)] = iWHITE;
    }
}
//-----------------------------------------------------------------
// 画像を簡易3次元表示
//   plot_step:  線の描画間隔 [画素] 2〜10 程度
//   plot_ratio: 階調方向の倍率 0.0〜1.0(=255画素)
//-----------------------------------------------------------------
void draw_3D_graph(LPBYTE inBuf,int plot_step,double plot_ratio)
{
    int x,y,lx1,lx2,ly1,ly2;
	int x_size,y_size;
	int plot_height; // 階調方向の高さ=256.0 * plot_ratio
	double xmag,ymag;
	LPBYTE ib=(LPBYTE)GlobalAlloc(GPTR,sizeof(BYTE)*iSize);
	CopyMemory(ib,inBuf,sizeof(BYTE)*iSize);
	//-------------------------------------------------------
	// 描画パラメータ(plot_step,plot_ratio)の指定
	//-------------------------------------------------------
	x_size=iWidth  + iHeight / 2;  // グラフの想定横幅(画素)
	plot_height=rint(iMAX_GRAY*plot_ratio);
	y_size=iHeight + plot_height;  // グラフの想定縦幅(画素)
	xmag = (double)iWidth /(double)x_size;
	ymag = (double)iHeight/(double)y_size;
	//-------------------------------------------------------
	// 描画エリアの初期化
	//-------------------------------------------------------
	FillMemory(inBuf,sizeof(BYTE)*iSize,0);
	//-------------------------------------------------------
	// グラフの表示
	//-------------------------------------------------------
	for(y=0;y<iHeight;y++){
		if(y%plot_step == 0){
			// y が plot_step で割り切れるときだけ描く
			// 始めに曲線部分を黒で描いて隠面消去する
			for(x=0;x<iWidth;x++){
				lx1 = rint(xmag*(x+y/2));
				ly1 = rint(ymag*((y_size-1)-(plot_height+y)));
				lx2 = rint(xmag*(x+y/2));
				ly2 = rint(ymag*((y_size-1)-
								 (plot_height+y-
								  (ib[XY2PI(x,y,iWidth)]*plot_ratio))));
				draw_line(inBuf,lx1,ly1,lx2,ly2,0);
			}
			// 曲線部分を白で描く
			for(x=0;x<iWidth-1;x++){
				lx1 = rint(xmag*(x+y/2));
				ly1 = rint(ymag*((y_size-1)-
								 (plot_height+y-
								  (ib[XY2PI(x,y,iWidth)]*plot_ratio))));
				lx2 = rint(xmag*(x+y/2+1));
				ly2 = rint(ymag*((y_size-1)-
								 (plot_height+y-
								  (ib[XY2PI(x+1,y,iWidth)]*plot_ratio))));
				draw_line(inBuf,lx1,ly1,lx2,ly2,iWHITE);
			}
		}
	}
	GlobalFree(ib);
}