边缘检测

float fSqrtByCarmack(float number)
{
	int i;
	float x2, y;
	const float threehalfs = 1.5F;

	x2 = number * 0.5F;
	y = number;
	i = *(int *)&y;            
	i = 0x5f375a86 - (i >> 1); 
	y = *(float *)&i;
	y = y * (threehalfs - (x2 * y * y));   
	
	return 1.0f / y;
}

void sqrt_cmp()
{
	float a;
	time_t start;
	int cishu = 1000000;

	cin >> a >> cishu;
	start = clock();
	for (int i = 0; i < cishu; i++)
	{
		fSqrtByCarmack(a);
	}
	cout << "fSqrtByCarmack: " << clock() - start << endl;
	start = clock();
	for (int i = 0; i < cishu; i++)
	{
		sqrtf(a);
	}
    start = clock();
	cout << "sqrtf: " << clock() - start << endl;	
	
/*
平台配置x86 debug
结果：
fSqrtByCarmack: 33
sqrtf: 27
结论：
自带的更快，原因可能跟硬件有关
*/

}

main.cpp调用homework.h中特定作业的题目的函数
homework.cpp调用bmpConverter.h的类来完成任务
bmpConverter.cpp调用BMPCore.h的方法来完成相应功能
    优点：结构更加明确，相较之前，bug变少了（可能是最近写代码状态不错），调试时间明显缩减
    缺点：修改参数会变得非常繁琐，因为文件变多了

// main.cpp
int main()
{
	//hw4_grad1_grad2_shenjun();
	//hw4_num_dect();
	//test();
	//hw4_grad1_grad2_shenjun();
	return 0;
}
// homework.cpp


void hw4_grad1_grad2_shenjun()
{
	bmpConverter bmpCvt("./pic/H0401Gry.bmp");  
	bmpCvt.AEadgeDectGrad_SJ2d1c(0.5);
	bmpCvt.Img2Bmp("./pic/H0401js_gd.bmp");
	bmpConverter bmpCvt2("./pic/H0401Gry.bmp");  
	bmpCvt2.AEadgeDectLaplacain_Gd2d1c(0.5);
	bmpCvt2.Img2Bmp("./pic/H0401la_gd.bmp");
}

void test()
{
	bmpConverter bmpCvt("./pic/treeGry.bmp");  
	bmpCvt.test();
	//bmpCvt.Img2Bmp("./pic/tree.bmp");
}

void hw4_num_dect()
{
	// 缩小
	// 得到边缘强度
	// 积分图
	// 得到最大强度区域
	// 放大
	bmpConverter bmpCvt("./pic/H0402Gry.bmp");  
	bmpCvt.num_dect();
	bmpCvt.Img2Bmp("./pic/H0402_det_sb.bmp");
	bmpConverter bmpCvt1("./pic/H0403Gry.bmp");  
	bmpCvt1.num_dect();
	bmpCvt1.Img2Bmp("./pic/H0403_det_sb.bmp");
	bmpConverter bmpCvt2("./pic/H0404Gry.bmp");  
	bmpCvt2.num_dect();
	bmpCvt2.Img2Bmp("./pic/H0404_det_sb.bmp");

}

void hw4_canny()
{
	bmpConverter bmpCvt("./pic/H0401Gry.bmp");  
	bmpCvt.AEadgeCanny2d1c(0.4, 6, 125, 75);
	bmpCvt.Img2Bmp("./pic/H0401_canny.bmp");
}

//bmpConverter.cpp

void getEdg(BYTE * pSrc, int width, int height, BYTE threshold)
{
	BYTE * pCur, *pEnd;
	for (pCur = pSrc, pEnd = pSrc + width * height; pCur < pEnd; pCur++)
	{
		*pCur = *pCur > threshold ? 255 : 0;
	}
}
// 边缘检测，一阶梯度+沈俊算子
Img bmpConverter::AEadgeDectGrad_SJ2d1c(double a0, bool inplace)
{
	if (!pImg)return Img();
	BYTE *pDstImg, *pTmpImg;
	BYTE threshold;
	int i;

	threshold = 10;
	pDstImg = new BYTE[width * height];
	pTmpImg = new BYTE[width * height];
	a0 = (max)(0.01, (min)(a0, 0.99));

	AEDOpShenJun4_8bit(pImg, pTmpImg, width, height, a0, pDstImg);
	AEDOpGrad1_8bit(pImg, width, height, pTmpImg);
	getEdg(pTmpImg, width, height, threshold);

	for (i = 0; i < width * height; i++)
	{
		pDstImg[i] = pDstImg[i] & pTmpImg[i];
	}

	delete pTmpImg;
	if (inplace)
	{
		delete pImg;
		pImg = pDstImg;
		return Img();
	}
	else
	{
		return Img(pDstImg, width, height);
	}

	return Img();
}

Img bmpConverter::AEadgeDectLaplacain_Gd2d1c(double a0, BYTE threshold, bool inplace)
{
	if (!pImg)return Img();
	BYTE *pDstImg, *pTmpImg;
	int i;

	pDstImg = new BYTE[width * height];
	pTmpImg = new BYTE[width * height];
	a0 = (min)(0.01, (max)(a0, 0.99));

	AEDOpShenJun4_8bit(pImg, pTmpImg, width, height, a0, pDstImg);
	AEDOpGrad1_8bit(pImg, width, height, pTmpImg);
	getEdg(pTmpImg, width, height, threshold);

	for (i = 0; i < width * height; i++)
	{
		pDstImg[i] = pDstImg[i] & pTmpImg[i];
	}

	delete pTmpImg;
	if (inplace)
	{
		delete pImg;
		pImg = pDstImg;
		return Img();
	}
	else
	{
		return Img(pDstImg, width, height);
	}

	return Img();
}


Img bmpConverter::AEadgeDectsobel_ShenJun2d1c(double a0, BYTE threshold, bool inplace)
{
	if (!pImg)return Img();
	BYTE *pDstImg, *pTmpImg;
	int i;

	pDstImg = new BYTE[width * height];
	pTmpImg = new BYTE[width * height];
	a0 = (max)(0.01, (min)(a0, 0.99));

	AEDOpShenJun4_8bit(pImg, pTmpImg, width, height, a0, pDstImg);
	AEDOpSobel_8bit(pImg, width, height, pTmpImg);

	getEdg(pTmpImg, width, height, threshold);
	for (i = 0; i < width * height; i++)
	{
		pDstImg[i] = pDstImg[i] & pTmpImg[i];
	}


	delete []pTmpImg;
	if (inplace)
	{
		delete pImg;
		pImg = pDstImg;
		return Img();
	}
	else
	{
		return Img(pDstImg, width, height);
	}
}
// 缩小到1/16，使用均值，已知宽高
void bmpConverter::shrink16(BYTE *pSrcImg, int width, int height, BYTE *pDstImg)
{
	int sum;
	BYTE *pCur, *pDes;
	int x, y;

	for (y = 0, pCur = pSrcImg, pDes = pDstImg; y < height; y += 4)
	{
		for (x = 0; x < width; x += 4)
		{
			// 求和
			//sum = 0, y = 0;
			//for (y = 0; y < 4 * width; y += width)
			//{
			//	sum += pCur[y];
			//	sum += pCur[y + 1];
			//	sum += pCur[y + 2];
			//	sum += pCur[y + 3];
			//}
			//*pDes++ = sum / 16;
			*pDes++ = *pCur;
			pCur += 4;
		}
		pCur += 3 * width;
	}

}

void getMaxAreaCalGraph(int *pSum, int m, int n, int width, int height, int &mx, int &my)
{
	int *pY1, *pY2, *pY3;
	int halfx, halfy, w;
	int y, x1, x2;
	int curSum, befSum, nxtSum, maxSum, c, tempSum;
	int rowSum[1024];

	// 初始化
	m = m | 1; 
	n = n | 1; 
	halfx = m / 2; 
	halfy = n / 2; 
	w = width + 1;
	c = (1 << 23) / (m * n); 
	mx = 0, my = 0;
	maxSum = -1e9;
	for (x1 = width - m; x1 <= width + m; x1++)rowSum[x1] = 0;

	// 最大值区域求取
	//for (y = 2 * n, pY3 = pSum, pY1 = pY3 + w * n, pY2 = pY1 + w * n;	y < height;	y++, pY1 += w, pY2 += w)
	for (y = 2 * n, pY3 = pSum, pY1 = pY3 + w * n, pY2 = pY1 + w * n;	y < height;	y++, pY1 += w, pY2 += w, pY3 += w)
	{
		befSum = 0, curSum = 0;
		for (x1 = 0, x2 = m; x2 <= width; x1++, x2++)
		{
			rowSum[x1] = pY1[x1] - pY1[x2] - pY2[x1] + pY2[x2];
		}
		for (x1 = m; x1 <= width; x1++)
		{
			//tempSum = rowSum[x1] - rowSum[x1 - m] - rowSum[x1 + m];
			tempSum = rowSum[x1] - rowSum[x1 - m] - rowSum[x1 + m] -   (pY3[x1] - pY3[x1 + w] - pY1[x1] + pY1[x1 + w]);
			if (tempSum > maxSum)
			{
				maxSum = tempSum;
				mx = x1, my = y;

			}
		}
	}
	return;
}

void drawSquare(BYTE *pSrcImg, int width, int height, int mx, int my, int w, int h)
{
	BYTE *pCur, *pY, *pX;

	pCur = pSrcImg;
	// 上横线
	for (pCur = pSrcImg + width * my + mx, pX = pCur + w; pCur < pX; ) *pCur++ = 255;
	//下横线
	for (pCur = pSrcImg + width * (my + h) + mx, pX = pCur + w; pCur < pX; ) *pCur++ = 255;
	// 左竖线
	for (pCur = pSrcImg + width * my + mx, pY = pCur + h * width; pCur < pY; pCur += width) *pCur = 255;
	// 右竖线
	for (pCur = pSrcImg + width * (my)+ mx + w, pY = pCur + h * width; pCur <= pY; pCur += width) *pCur = 255;
}
void bmpConverter::num_dect()
{
// 缩小
// 得到边缘强度
// 积分图
// 得到最大强度区域
// 放大

	BYTE *pDstImg, *pTmpImg, *pSrcImg, *pOrgImg;
	int *pCalGraph;
	int winSize, mx, my;
	BYTE threshold;
	double a0;

	// 初始化
	a0 = 0.1;
	threshold = 80;
	winSize = 25;
	pCalGraph = new int[(width / 4 + (winSize | 1)) * (height / 4 + (winSize | 1))];
	pSrcImg = new BYTE[width * height / 16];
	//BYTE * pSbImg = new BYTE[width * height / 16];
	pOrgImg = pImg;

	// 缩小图片
	shrink16(pImg, width, height, pSrcImg);
	width /= 4, height /= 4;
	pImg = pSrcImg;
	// 得到边缘数目
	AEadgeDectsobel_ShenJun2d1c(a0, threshold);
	//AEadgeCanny2d1c(1, 5, 60, 80);
	//AEDOpSobel_8bit(pImg, width, height, pSbImg);
	//delete[] pSrcImg;
	//pImg = pSbImg;
	Img2Bmp("./pic/temp.bmp");
	// 获得边缘数目的积分图
	AgetCalGraph(pImg, width, height, 1, 1, pCalGraph);
	// 原图像100x100的区域
	// 4x4缩小之后，25x25
	// 寻找计算图中值最大区域
	getMaxAreaCalGraph(pCalGraph, winSize, winSize, width, height, mx, my);
	//drawSquare(pImg, width, height, mx, my - winSize, winSize, winSize);
	mx *= 4, my = my * 4, width = width * 4, height = height * 4, winSize = winSize * 4;
	pSrcImg = pImg;
	pImg = pOrgImg;
	// 画出框框
	drawSquare(pImg, width, height, mx, my - winSize, winSize, winSize);
	delete[] pCalGraph;
	delete[] pSrcImg;
}

// canny算子边缘检测
Img bmpConverter::AEadgeCanny2d1c(double sigma, int Fsize, int uBound, int lBound, bool inplace)
{
	BYTE *pDstImg, *pTmpImg;
	int *pFilter, *dx, *dy;

	// 初始化
	sigma = 0.4;
	Fsize = 5;
	uBound = 125;
	lBound = 75;
	pDstImg = new BYTE[width * height];
	dx = new int[width * height];
	dy = new int[width * height];
	pFilter = new int[width * height];
	pTmpImg = new BYTE[width * height];

	AEDOpCanny_8bit(
		pImg,
		pTmpImg,
		dx, dy,
		width, height,
		pFilter, sigma, Fsize,
		uBound, lBound,
		pDstImg
	);

	// 结束
	delete[] pTmpImg;
	delete[] pFilter;
	delete[] dx;
	delete[] dy;
	if (inplace)
	{
		delete[] pImg;
		pImg = pDstImg;
		return Img();
	}
	else
	{
		return Img(pDstImg, width, height);
	}
}
void bmpConverter::test()
{

// 缩小
// 得到边缘强度
// 积分图
// 得到最大强度区域
// 放大

	BYTE *pDstImg, *pTmpImg, *pSrcImg, *pOrgImg, *pTansImg;
	int *dx, *dy;
	int Fsize, uBound, lBound;
	BYTE threshold;
	double sigma;

	pDstImg = new BYTE[height * width];
	double pFilter[]{
		1 , 1 , 1 ,
		1 , -8 , 1 ,
		1 , 1 , 1 , };
	for (int i = 0; i < 9; i++)pFilter[i] *= -1;
	AConv2d3_8bit(pImg, width, height, pFilter, pDstImg);
	pImg = pDstImg;
	Img2Bmp("./pic/tree6.bmp");

}

// 邻域运算，边缘检测梯度算子4邻域8bit图像
void AEDOpGrad1_8bit(BYTE *pSrc, int width, int height, BYTE *pDst)
{
	BYTE *pCur, *pDes;
	int dx, dy;
	int x, y, i, m;
	BYTE LUT[1024];  // 局部查表

	for (i = 0; i < 1024; i++)LUT[i] = (BYTE)sqrt(double(i));

	for (y = 0, pCur = pSrc, pDes = pDst; y < height - 1; y++)
	{
		for (x = 0; x < width - 1; x++, pCur++)
		{
			dx = *pCur - pCur[1];
			dy = *pCur - pCur[width];
			m = dx * dx + dy * dy;
			*(pDes++) = m < 1024 ? LUT[m] : (min)(255, (int)(sqrt((double) m)));
		}
		*(pDes++) = 0;  // 尾列不做,边缘强度赋0
		pCur++;
	}
	memset(pDes, 0, width);  // 尾行不做,边缘强度赋0
	return;
}
// 邻域运算，边缘检测索贝尔算子4邻域8bit图像
void AEDOpSobel_8bit(BYTE *pGryImg, int width, int height, BYTE *pSbImg)
{
	BYTE *pGry, *pSb;
	int dx, dy;
	int x, y;

	memset(pSbImg, 0, width); //首行不做,边缘强度赋0
	for (y = 1, pGry = pGryImg + width, pSb = pSbImg + width; y < height - 1; y++)
	{
		*(pSb++) = 0;  //首列不做,边缘强度赋0
		pGry++;
		for (x = 1; x < width - 1; x++, pGry++)
		{
			//求dx
			dx = *(pGry - 1 - width) + (*(pGry - 1) * 2) + *(pGry - 1 + width);
			dx -= *(pGry + 1 - width) + (*(pGry + 1) * 2) + *(pGry + 1 + width);
			//求dy
			dy = *(pGry - width - 1) + (*(pGry - width) * 2) + *(pGry - width + 1);
			dy -= *(pGry + width - 1) + (*(pGry + width) * 2) + *(pGry + width + 1);
			//结果
			*(pSb++) = min(255, abs(dx) + abs(dy));
		}
		*(pSb++) = 0; //尾列不做,边缘强度赋0
		pGry++;
	}
	memset(pSb, 0, width); //尾行不做,边缘强度赋0
	return;
}


// 邻域运算，边缘检测Laplacian算子4邻域8bit图像
void AEDOpLaplacian4_8bit(BYTE *pSrc, BYTE *pTmp, int width, int height, BYTE *pDst)
{
	BYTE *pCur, *pDes;
	int grad;
	int x, y;

	memset(pTmp, 0, width);  // 首行不做,边缘强度赋0
	for (y = 1, pCur = pSrc + width, pDes = pTmp + width; y < height - 1; y++)
	{
		*(pDes++) = 0;  // 首列不做,边缘强度赋0
		pCur++;
		for (x = 1; x < width - 1; x++, pCur++)
		{
			// grad更新
			grad = -((int)*pCur) * 4;
			grad += pCur[-width];
			grad += pCur[-1];
			grad += pCur[+1];
			grad += pCur[+width];
			*(pDes++) = grad > 0;
		}
		*(pDes++) = 0; //尾列不做,边缘强度赋0
		pCur++;
	}
	memset(pDes, 0, width); //尾行不做,边缘强度赋0

	pDes = pDst + width, pCur = pTmp + width;
	memset(pDst, 0, width);
	for (y = 1; y < height - 1; y++)
	{
		*pDes++ = 0, pCur++;
		for (x = 1; x < width - 1; x++, pCur++)
		{
			if (*pCur && (!pCur[-1] || !pCur[1] || !pCur[-width] || !pCur[width]))*pDes++ = 255;
			else *pDes++ = 0;
		}
		*pDes++ = 0, pCur++;
	}
	memset(pDes, 0, width);
	return;
}


// 邻域运算，边缘检测沈俊算子4邻域8bit图像
void AEDOpShenJun4_8bit(
	BYTE *pSrcImg, //原始灰度图像 
	BYTE *pTmpImg, //辅助图像
	int width, int height,
	double a0,
	BYTE *pDstImg
)
{
	BYTE * pCur, * pDes, * pEnd;
	int LUT[511], *A;
	int x, y, i, pre;
	// 初始化
	A = LUT + 255;
	for (A[0] = 0, i = 1; i < 256; i++)  // 向上取整
	{
		A[i] = a0 * i + 0.5;	
		A[-i] = -a0 * i - 0.5;  
	}
	// 差分
	// 横向
	for (y = 0, pCur = pSrcImg, pDes = pTmpImg; y < height; y++)
	{
		pre = *pDes++ = *pCur++;
		// 从左到右
		for (x = 1; x < width; x++)
		{
			pre = *pDes++ = pre + A[*pCur++ - pre];
		}
		// 从右到左
		pDes--;
		for (x = 1; x < width; x++)
		{
			pDes--;
			pre = *pDes = pre + A[*pDes - pre];
		}
		pDes += width;
	}
	// 纵向
	for (x = 0, pDes=pTmpImg; x < width; x++)
	{
		// 从上到下
		pDes = pTmpImg + x, pre = *pDes;
		for (y = 1; y < height; y++)
		{
			pDes += width;
			pre = *pDes = pre + A[*pDes - pre];
		}
		// 从下到上
		for (y = 1; y < height; y++)
		{
			pDes -= width;
			pre = *pDes = pre + A[*pDes - pre];
		}
	}


	// 求二阶梯度的符号
	for (pCur = pTmpImg, pEnd = pTmpImg + width * height, pDes = pSrcImg; pCur<pEnd; pDes++)
	{
		*(pCur++) = (*pCur <= *pDes);
	}

	// 求边缘
	pDes = pDstImg + width, pCur = pTmpImg + width;
	memset(pDstImg, 0, width);
	for (y = 1; y < height - 1; y++)
	{
		*pDes++ = 0, pCur++;
		for (x = 1; x < width - 1; x++, pCur++)
		{
			if (!*pCur && (pCur[-1] ||pCur[1] || pCur[-width] || pCur[width]))*pDes++ = 255;
			else *pDes++ = 0;
		}
		*pDes++ = 0, pCur++;
	}
#ifdef DEBUG
	//printf("求边缘ok\n");
#endif // DEBUG
	memset(pDes, 0, width);
#ifdef DEBUG
	//printf("return ok\n");
#endif // DEBUG

	return;
}

void AGuassFilter1d_8bit(BYTE * pSrc, int width, int height, int * pFilter, int m, BYTE * pDes)
{
	int offset = 17;
	int halfm, y, x, f, sum = 0;
	BYTE  *pCur, *pRes, *pY;

	halfm = m / 2;

	for (y = 0, pY = pSrc, pRes = pDes; y < height; y++, pY += width)
	{
		// 处理左侧
		for (x = 0, sum = 0; x < halfm; x++, sum = 0)
		{
			for (f = halfm - x, pCur = pY; f < m; f++)
			{
				sum += pFilter[f] * *pCur++;
			}
			sum >>= offset;
			*pRes++ = sum;
		}
		// 处理中间
		for (x = halfm, sum = 0; x < width - halfm; x++, sum = 0)
		{
			for (f = 0, pCur = pY + x - halfm; f < m; f++)
			{
				sum += pFilter[f] * *pCur++;
			}
			sum >>= offset;
			*pRes++ = sum;
		}
		// 处理右侧
		for (x = width - halfm, sum = 0; x < width; x++, sum = 0)
		{
			for (f = 0, pCur = pY + x - halfm; pCur < pY + width; f++)
			{
				sum += pFilter[f] * *pCur++;
			}
			sum >>= offset;
			*pRes++ = sum;
		}
	}
}

void Transpose(BYTE * pSrc, int w, int h, BYTE * pDes)
{
	BYTE * pX, *pCur;
	int x, y;

	pCur = pDes;

	for (x = 0; x < w; x++)
	{
		for (y = 0, pX = pSrc + x; y < h; y++, pX += w)
		{
			*pCur++ = *pX;
		}
	}

	return;
}

void getGuassFilter_8bit(double dev, int m, int * pFilter, int offset)
{
	int halfm;
	double sum;
	double dF[1024];

	halfm = m / 2;
	sum = dF[halfm] = 1;

	for (int i = 1; i <= halfm; i++)
	{
		dF[halfm + i] = dF[halfm - i] = exp(-i * i / 2.0 / dev / dev);
		sum += dF[halfm + i] + dF[halfm + i];
	}
	for (int i = 0; i <= halfm; i++)
	{
		pFilter[halfm + i] = pFilter[halfm - i] = (1 << offset) * dF[halfm + i] / sum;
	}
}

void AEDOpSobe4Canny(BYTE *pGryImg, int width, int height, int *pDX, int *pDY, BYTE *pSbImg)
{
	BYTE *pGry, *pSb;
	int *pDx, *pDy;
	int dx, dy;
	int x, y;

	pDx = pDX;
	pDy = pDY;
	memset(pSbImg, 0, width); //首行不做,边缘强度赋0
	for (y = 1, pGry = pGryImg + width, pSb = pSbImg + width; y < height - 1; y++)
	{
		*(pSb++) = 0;  //首列不做,边缘强度赋0
		*(pDx++) = 0;
		*(pDy++) = 0;
		pGry++;
		for (x = 1; x < width - 1; x++, pGry++)
		{
			//求dx
			dx = *(pGry - 1 - width) + (*(pGry - 1) * 2) + *(pGry - 1 + width);
			dx -= *(pGry + 1 - width) + (*(pGry + 1) * 2) + *(pGry + 1 + width);
			//求dy
			dy = *(pGry - width - 1) + (*(pGry - width) * 2) + *(pGry - width + 1);
			dy -= *(pGry + width - 1) + (*(pGry + width) * 2) + *(pGry + width + 1);
			//结果
			*(pSb++) = min(255, abs(dx) + abs(dy));
			*(pDx++) = dx;
			*(pDy++) = dy;

		}
		*(pSb++) = 0; //尾列不做,边缘强度赋0
		*(pDx++) = 0;
		*(pDy++) = 0;
		pGry++;
	}
	memset(pSb, 0, width); //尾行不做,边缘强度赋0
	return;
}

// 对梯度幅值进行非极大值抑制
void GradMagThresholding(BYTE *pSrcImg, int width, int height, int *pDx, int *pDy)
{
	BYTE *pCur;
	int *pDX, *pDY;
	BYTE Cur;
	int DX, DY;
	int x, y;
	BYTE btan, binvtan;
	// 不做边缘
	for (y = 1, pCur = pSrcImg + width, pDX = pDx + width, pDY = pDy + width; y < height - 1; y++)
	{
		pCur++, pDX++, pDY++; 
		for (x = 1; x < width - 1; x++, pCur++, pDX++, pDY++)
		{
			Cur = *pCur, DX = *pDX, DY = *pDY;
			// 取周围的8个点作比较
			if (!Cur)continue; // 模为0，则继续
			if (!DX)  // 纵向梯度
			{
				if (Cur < pCur[-width] || Cur < pCur[width])*pCur = 0;
			}
			else if (!DY) // 纵向梯度
			{
				if (Cur < pCur[-1] || Cur < pCur[1])*pCur = 0;
			}
			else 
			{
				btan = DY / DX;
				if (!btan)  // abs(tan)小于1
				{
					binvtan = DX / DY;
					if (binvtan > 2 || binvtan < -2)
					{
						if (Cur < pCur[-1] || Cur < pCur[1])*pCur = 0;
					}
					else if (binvtan > 0)
					{
						if (Cur < pCur[-width - 1] || Cur < pCur[width + 1])*pCur = 0;
					}
					else
					{
						if (Cur < pCur[-width + 1] || Cur < pCur[width - 1])*pCur = 0;
					}
				}
				else
				{
					if (btan >= 2 || btan < -2)
					{
						if (Cur < pCur[-width] || Cur < pCur[width])*pCur = 0;
					}
					else if (btan > 0)
					{
						if (Cur < pCur[-width - 1] || Cur < pCur[width + 1])*pCur = 0;
					}
					else
					{
						if (Cur < pCur[-width + 1] || Cur < pCur[width - 1])*pCur = 0;
					}
				}
			}
		}
		pCur++, pDX++, pDY++;
	}
}


void DoubThre_EdgTrack(BYTE *pSrcImg, int width, int height, int uBound, int lBound)
{
	BYTE *pCur;
	int x, y;

	for (y = 1, pCur = pSrcImg + width; y < height - 1; y++)  // 不处理边缘
	{
		pCur++;
		for (x = 1; x < width - 1; x++)
		{
			BYTE &Cur = *pCur++;
			if (Cur > uBound)Cur = 255;  // 该点是强连接点
			else if (Cur < lBound)Cur = 0;  // 该点无连接
			else  // 该点时弱连接点
			{
				// 如果周围有强连接点，则弱->强
				if (pCur[-width - 1] == 255 || pCur[-width] == 255 || pCur[-width + 1] == 255 ||
					pCur[-1] == 255 || pCur[1] > uBound ||
					pCur[width - 1] > uBound || pCur[width] > uBound || pCur[width + 1] > uBound)
				{
					Cur = 255;
				}
				else Cur = 0;  // 否则弱->无
			}
		}
		pCur++;
	}
}

void AEDOpCanny_8bit(
	BYTE *pSrcImg,							// 原图像
	BYTE * pTmpImg,		
	int *dx, int *dy,						// sobel算子处理得到的梯度
	int width, int height,					// 图像宽高
	int *pFilter, double sigma, int Fsize,	// 高斯滤波器的参数
	int uBound, int lBound,					// 用双阈值算法检测和连接边缘
	BYTE * pDstImg							// 目标图像
)
{
	// 高斯平滑
	// sobel算子处理
	// 对梯度幅值进行非极大值抑制
	// 用双阈值算法检测和连接边缘


	// 初始化
	Fsize |= 1;
	getGuassFilter_8bit(sigma, Fsize, pFilter);

	// 施展高斯平滑
	AGuassFilter1d_8bit(pSrcImg, width, height, pFilter, Fsize, pDstImg);
	Transpose(pDstImg, width, height, pTmpImg);
	AGuassFilter1d_8bit(pTmpImg, height, width, pFilter, Fsize, pDstImg);
	Transpose(pDstImg, height, width, pTmpImg);
	
	// sobel算子处理
	AEDOpSobe4Canny(pTmpImg, width, height, dx, dy, pDstImg);

	// 对梯度幅值进行非极大值抑制
	GradMagThresholding(pDstImg, width, height, dx, dy);

	// 用双阈值算法检测和连接边缘
	DoubThre_EdgTrack(pDstImg, width, height, uBound, lBound);
}

void AConv2d3_8bit(BYTE * pSrcImg, int width, int height, double * pFilter, BYTE * pDstImg)
{
	BYTE *pCur, *pDes;
	int x, y, i, j;
	double *pF, sum;

	pCur = pSrcImg;
	pDes = pDstImg;
	pF = pFilter;

	for (x = 0; x < width; x++)*pDes++ = *pCur++;  // 边缘处理
	for (y = 1; y < height - 1; y++)
	{
		*pDes++ = *pCur++;
		for (x = 1; x < width - 1; x++, pCur++)
		{
			sum = pCur[-1 - width] * pF[0] + pCur[-width] * pF[1] + pCur[1 - width] * pF[2] +
				pCur[-1] * pF[3] + pCur[0] * pF[4] + pCur[1] * pF[5] +
				pCur[-1 + width] * pF[6] + pCur[width] * pF[7] + pCur[width + 1] * pF[8];
			//cout << pCur[-1 - width] * pF[0] << endl;
			//cout << sum << endl;
			*pDes++ = max(0, min(255, (sum)));
		}
		*pDes++ = *pCur++;

	}
	for (x = 0; x < width; x++)*pDes++ = *pCur++;  // 边缘处理

}