2018-11-07 22:57:28 Godsolve 阅读数 465
• ###### MATLAB图像处理

全面系统的学习MATLAB在图像处理中的应用

20789 人正在学习 去看看 魏伟

• 1.欧几里德距离
• 2.曼哈顿距离(City Block Distance)
Distance = |x2-x1|+|y2-y1|
• 3.象棋格距离(Chessboard Distance)
Distance = max(|x2-x1|,|y2-y1|)

##### 同时也欢迎各位关注我的微信公众号 南木的下午茶

// CVE6.2.cpp: 定义控制台应用程序的入口点。
//

#include "stdafx.h"
#include "opencv2/imgproc/imgproc.hpp"
#include "opencv2/highgui/highgui.hpp"
#include <stdio.h>

using namespace cv;

int voronoiType = -1;
int edgeThresh = 100;
int distType0 = CV_DIST_L1;

// The output and temporary images
Mat gray;

// threshold trackbar callback
static void onTrackbar(int, void*)
{
static const Scalar colors[] =
{
Scalar(0,0,0),
Scalar(255,0,0),
Scalar(255,128,0),
Scalar(255,255,0),
Scalar(0,255,0),
Scalar(0,128,255),
Scalar(0,255,255),
Scalar(0,0,255),
Scalar(255,0,255)
};

int distType = voronoiType >= 0 ? CV_DIST_L2 : distType0;

Mat edge = gray >= edgeThresh, dist, labels, dist8u;

if (voronoiType < 0)
else
distanceTransform(edge, dist, labels, distType, maskSize, voronoiType);

if (voronoiType < 0)
{
// begin "painting" the distance transform result
dist *= 5000;
pow(dist, 0.5, dist);

Mat dist32s, dist8u1, dist8u2;

dist.convertTo(dist32s, CV_32S, 1, 0.5);
dist32s &= Scalar::all(255);

dist32s.convertTo(dist8u1, CV_8U, 1, 0);
dist32s *= -1;

dist32s += Scalar::all(255);
dist32s.convertTo(dist8u2, CV_8U);

Mat planes[] = { dist8u1, dist8u2, dist8u2 };
merge(planes, 3, dist8u);
}
else
{
dist8u.create(labels.size(), CV_8UC3);
for (int i = 0; i < labels.rows; i++)
{
const int* ll = (const int*)labels.ptr(i);
const float* dd = (const float*)dist.ptr(i);
uchar* d = (uchar*)dist8u.ptr(i);
for (int j = 0; j < labels.cols; j++)
{
int idx = ll[j] == 0 || dd[j] == 0 ? 0 : (ll[j] - 1) % 8 + 1;
float scale = 1.f / (1 + dd[j] * dd[j] * 0.0004f);
int b = cvRound(colors[idx][0] * scale);
int g = cvRound(colors[idx][1] * scale);
int r = cvRound(colors[idx][2] * scale);
d[j * 3] = (uchar)b;
d[j * 3 + 1] = (uchar)g;
d[j * 3 + 2] = (uchar)r;
}
}
}

imshow("Distance Map", dist8u);
}

static void help()
{
printf("\nProgram to demonstrate the use of the distance transform function between edge images.\n"
"Usage:\n"
"./distrans [image_name -- default image is stuff.jpg]\n"
"\nHot keys: \n"
"\tESC - quit the program\n"
"\tC - use C/Inf metric\n"
"\tL1 - use L1 metric\n"
"\tL2 - use L2 metric\n"
"\t0 - use precise distance transform\n"
"\tv - switch to Voronoi diagram mode\n"
"\tp - switch to pixel-based Voronoi diagram mode\n"
"\tSPACE - loop through all the modes\n\n");
}

const char* keys =
{
"{1| |stuff.jpg|input image file}"
};

int main(int argc, const char** argv)
{
help();
CommandLineParser parser(argc, argv, keys);
//string filename = parser.get<string>("1");
imshow("原图", img);
imshow("灰度图", gray);
resize(gray, gray, Size(), 0.25, 0.25, 1);
//if (gray.empty())
//{
//	printf("Cannot read image file: %s\n", filename.c_str());
//	help();
//	return -1;
//}

namedWindow("Distance Map", CV_WINDOW_NORMAL);
createTrackbar("Brightness Threshold", "Distance Map", &edgeThresh, 255, onTrackbar, 0);

for (;;)
{
// Call to update the view
onTrackbar(0, 0);

int c = waitKey(0) & 255;

if (c == 27)
break;

if (c == 'c' || c == 'C' || c == '1' || c == '2' ||
c == '3' || c == '5' || c == '0')
voronoiType = -1;

if (c == 'c' || c == 'C')
distType0 = CV_DIST_C;
else if (c == '1')
distType0 = CV_DIST_L1;
else if (c == '2')
distType0 = CV_DIST_L2;
else if (c == '3')
else if (c == '5')
else if (c == '0')
else if (c == 'v')
voronoiType = 0;
else if (c == 'p')
voronoiType = 1;
else if (c == ' ')
{
if (voronoiType == 0)
voronoiType = 1;
else if (voronoiType == 1)
{
voronoiType = -1;
distType0 = CV_DIST_C;
}
else if (distType0 == CV_DIST_C)
distType0 = CV_DIST_L1;
else if (distType0 == CV_DIST_L1)
distType0 = CV_DIST_L2;
voronoiType = 0;
}
}

return 0;
}



2019-10-29 19:16:14 qq_33446100 阅读数 76
• ###### MATLAB图像处理

全面系统的学习MATLAB在图像处理中的应用

20789 人正在学习 去看看 魏伟

# 说在前面

• 编译环境：vs2017
• opencv版本：4.0.1
• 参考：《图像处理、分析与机器视觉(第4版)》

# 概念

• ## 数字图像

对于一个图像函数$f(x,y)$，若其定义域以及值域都是离散的，那么称之为数字的；
通常，通过采样将连续图像进行数字化。

一个连续图像在采样点处被数字化；采样点之间构成的关系为栅格。
• ## 距离

• 定义
满足以下条件的函数：
$D(\bold p,\bold q)\geq 0,当且仅当\bold p=\bold q时,D(\bold p,\bold q)=0\\ D(\bold p,\bold q)=D(\bold q,\bold p)\\ D(\bold p,\bold r)\leq D(\bold p,\bold q)+D(\bold q,\bold r)$
• 欧式距离
$D_E[(i,j),(h,k)]=\sqrt{(i-h)^2+(j-k)^2}$
• 城市街区距离
只允许横向以及纵向的移动

$D_4[(i,j),(h,k)]=|i-h|+|j-k|$
• 棋盘距离
允许横向、纵向以及对角线上的移动

$D_8[(i,j),(h,k)]=max\{|i-h|+|j-k|\}$
• ## 距离变换算法

• 按照一种距离度量D，D是D4或者D8，对大小为MxN的图像的一个子集S进行距离变换，建立一个MxN的数组F并作初始化：子集S中的元素置为0，其他位置为无穷。
• 按行遍历图像，从上到下，从左到右。对于上方和左边的邻接像素，设
$F(\bold p)=min[F(\bold p), D(p,q)+F(q)]$
• 按行遍历图像，从下到上，从右到左。对于下方和右边的邻接像素，设
$F(\bold p)=min[F(\bold p), D(p,q)+F(q)]$
• 数组F中得到的是子集S的斜切

# Code

• ## c++实现

并未处理图像，仅用于展示算法；
时间复杂度：$O(n^2)$
#include <iostream>
#include <vector>
#include <algorithm>
using namespace std;
/*
@{func} 判断i,j是否在范围内
*/
bool InArea(int i, int j, int rows, int cols)
{
if (i<0 || i>=rows)
return false;
if (j<0 || j>=cols)
return false;

return true;
}
/*
@{param: i j} 点位置
@{param: rows cols} 图像大小
@{param: f} 目标数组
@{func} 处理上以及左边的近邻像素
*/
void D4AL(int i,int j, int rows, int cols, vector<vector<int>> &f)
{
//上
if (InArea(i - 1, j, rows, cols))
f[i][j] = min(f[i][j], 1 + f[i - 1][j]);
//左上
if (InArea(i - 1, j - 1, rows, cols))
f[i][j] = min(f[i][j], 2 + f[i - 1][j - 1]);
//左
if (InArea(i, j - 1, rows, cols))
f[i][j] = min(f[i][j], 1 + f[i][j - 1]);
//左下
if (InArea(i + 1, j - 1, rows, cols))
f[i][j] = min(f[i][j], 2 + f[i + 1][j - 1]);
}
/*
@{param: i j} 点位置
@{param: rows cols} 图像大小
@{param: f} 目标数组
@{func} 处理下以及右边的近邻像素
*/
void D4BR(int i, int j, int rows, int cols, vector<vector<int>> &f)
{
//下
if (InArea(i + 1, j, rows, cols))
f[i][j] = min(f[i][j], 1 + f[i + 1][j]);
//右下
if (InArea(i + 1, j + 1, rows, cols))
f[i][j] = min(f[i][j], 2 + f[i + 1][j + 1]);
//右
if (InArea(i, j + 1, rows, cols))
f[i][j] = min(f[i][j], 1 + f[i][j + 1]);
//右上
if (InArea(i - 1, j + 1, rows, cols))
f[i][j] = min(f[i][j], 2 + f[i - 1][j + 1]);
}

/*
@{param:src} 源图像
@{param:f}   目标数组
*/
void DistanceTransformD4(vector<vector<int>> &src, vector<vector<int>> &f)
{
int cols = src[0].size();
int rows = src.size();

//初始化
for (int i = 0; i < rows; ++i)
for (int j = 0; j < cols; ++j)
if (src[i][j] == 1)
f[i][j] = 0;
else
f[i][j] = INT_MAX - 2;//简单的防止溢出
//按行遍历图像，从上到下，从左到右
for (int i = 0; i < rows; ++i)
for (int j = 0; j < cols; ++j)
D4AL(i, j, rows, cols, f);

//按行遍历图像，从下到上，从右到左
for (int i = rows - 1; i >= 0; --i)
for (int j = cols - 1; j >= 0; --j)
D4BR(i, j, rows, cols, f);
}

int main()
{
vector<vector<int>> src = {
{0,1,1,0,0,0,1,0 },
{0,1,0,0,0,0,0,1 },
{0,1,0,0,0,0,0,0 },
{0,1,0,0,0,0,0,0 }
};
int rows = src.size();
int cols = src[0].size();
vector<vector<int>> f(rows, vector<int>(cols, 0));

cout << "SRC:" << endl;
for (int i = 0; i < rows; ++i)
{
for (int j = 0; j < cols; ++j)
cout << src[i][j] << " ";
cout << endl;
}

DistanceTransformD4(src, f);

cout << "\nResult:" << endl;
for (int i = 0; i < rows; ++i)
{
for (int j = 0; j < cols; ++j)
cout << f[i][j] << " ";
cout << endl;
}

return 0;
}

• ## opencv

void cv::distanceTransform	(
InputArray 	src,
OutputArray 	dst,
OutputArray 	labels,
int 	distanceType,
int 	labelType = DIST_LABEL_CCOMP
)

栗子?

# 应用

• 移动机器人领域的路径规划以及障碍躲避
• 图像中寻找最近特征，骨架抽取
• 待补充

2019-04-11 02:10:11 qq_35306281 阅读数 290
• ###### MATLAB图像处理

全面系统的学习MATLAB在图像处理中的应用

20789 人正在学习 去看看 魏伟

@图像处理_距离变换_c++代码实现

1. 输入：二值图像。
2. 从图像左上角第二行开始，从左向右、从上到下移动窗口扫描每个像素，类似滤波过程，检测在中心像素P的周围四个像素与中心像素的距离，若中心像素为0，则跳过。保存最小距离与位置作为结果。
3. 从右下角倒数第二行开始，从右向左，从下到上扫描图像。其它同2。
4. 输出图像。

/*created at 2019/4/11*/
#include <iostream>
#include <algorithm>
#include <opencv2/opencv.hpp>

void distance_transform_3x3(cv::Mat& src)
{
int rows = src.rows;
int cols = src.cols;
float sum[5];
/*第一次扫描*/
for(int r = 1; r < rows - 1; ++r)
{
for(int c = 1; c < cols - 1; ++c)
{
if(src.at<uchar>(r,c))
{
sum[0] = 1.4142 + src.at<uchar>(r-1,c-1);
sum[1] = 1      + src.at<uchar>(r-1,  c);
sum[2] = 1.4142 + src.at<uchar>(r-1,c+1);
sum[3] = 1	    + src.at<uchar>(r,  c-1);
sum[4] =          src.at<uchar>(r,    c);
std::sort(sum, sum+5);
src.at<uchar>(r,c) = sum[0];
}
}
}
/*第二次扫描*/
for(int r = rows - 1; r > 0; --r)
{
for(int c = cols - 1; c > 0; --c)
{
if(src.at<uchar>(r,c))
{
sum[0] =          src.at<uchar>(r,    c);
sum[1] = 1	    + src.at<uchar>(r,  c+1);
sum[2] = 1.4142 + src.at<uchar>(r+1,c-1);
sum[3] = 1      + src.at<uchar>(r+1,  c);
sum[4] = 1.4142	+ src.at<uchar>(r+1,c+1);
std::sort(sum, sum+5);
src.at<uchar>(r,c) = sum[0];
}
}
}
}

int main(int argc, char* argv[])
{
distance_transform_3x3(src );
cv::imshow("out.jpg", src );
cv::waitKey(0);
return 0;
}


2017-04-02 15:47:39 qq_34784753 阅读数 6576
• ###### MATLAB图像处理

全面系统的学习MATLAB在图像处理中的应用

20789 人正在学习 去看看 魏伟

1.将输入图片转换为二值图像，前景设置为1，背景设置为0；

//距离变换的扫描实现

#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;

//计算欧氏距离的函数
float calcEuclideanDiatance(int x1, int y1, int x2, int y2)
{
return sqrt(float((x2 - x1)*(x2 - x1) + (y2 - y1)*(y2 - y1)));
}

//计算棋盘距离的函数
int calcChessboardDistance(int x1, int y1, int x2, int y2)
{
return cv::max(abs(x1 - x2), (y1 - y2));
}

//计算麦哈顿距离(街区距离)
int calcBlockDistance(int x1, int y1, int x2, int y2)
{
return abs(x1 - x2) + abs(y1 - y2);
}

//距离变换函数的实现
void distanceTrans(Mat &srcImage, Mat &dstImage)
{
CV_Assert(srcImage.data != nullptr);
//CV_Assert（）若括号中的表达式值为false，则返回一个错误信息。
//定义灰度图像的二值图像
Mat grayImage, binaryImage;
//将原图像转换为灰度图像
cvtColor(srcImage, grayImage, CV_BGR2GRAY);
//将灰度图像转换为二值图像
threshold(grayImage, binaryImage, 100, 255, THRESH_BINARY);
imshow("二值化图像", binaryImage);

int rows = binaryImage.rows;
int cols = binaryImage.cols;
uchar *pDataOne;
uchar *pDataTwo;
float disPara = 0;
float fDisMIn = 0;

//第一遍遍历图像，使用左模板更新像素值
for (int i = 1; i < rows - 1; i++)
{
//图像的行指针的获取
pDataOne = binaryImage.ptr<uchar>(i);
for (int j = 1; j < cols; j++)
{
//分别计算左模板掩码的相关距离
//PL  PL
//PL  P
//PL
pDataTwo = binaryImage.ptr<uchar>(i - 1);
disPara = calcEuclideanDiatance(i, j, i - 1, j - 1);
fDisMIn = cv::min((float)pDataOne[j], disPara + pDataTwo[j - 1]);
disPara = calcEuclideanDiatance(i, j, i - 1, j);
fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j]);
pDataTwo = binaryImage.ptr<uchar>(i);
disPara = calcEuclideanDiatance(i, j, i, j - 1);
fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j-1]);
pDataTwo = binaryImage.ptr<uchar>(i+1);
disPara = calcEuclideanDiatance(i, j, i+1,j-1);
fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j - 1]);
pDataOne[j] = (uchar)cvRound(fDisMIn);
}
}

//第二遍遍历图像，使用右模板更新像素值
for (int i = rows - 2; i > 0; i--)
{
pDataOne = binaryImage.ptr<uchar>(i);
for (int j = cols - 1; j >= 0; j--)
{
//分别计算右模板掩码的相关距离
//pR  pR
//pR  p
//pR
pDataTwo = binaryImage.ptr<uchar>(i + 1);
disPara = calcEuclideanDiatance(i, j, i + 1, j);
fDisMIn = cv::min((float)pDataOne[j], disPara + pDataTwo[j]);
disPara = calcEuclideanDiatance(i, j, i + 1, j + 1);

fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j+1]);
pDataTwo = binaryImage.ptr<uchar>(i);
disPara = calcEuclideanDiatance(i, j, i, j +1);
fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j + 1]);
pDataTwo = binaryImage.ptr<uchar>(i - 1);
disPara = calcEuclideanDiatance(i, j, i - 1, j + 1);
fDisMIn = cv::min(fDisMIn, disPara + pDataTwo[j + 1]);
pDataOne[j] = (uchar)cvRound(fDisMIn);
}
}
dstImage = binaryImage.clone();
}

//主函数
int main()
{
if (!srcImage.data)
{
cout << "读入图片错误！" << endl;
system("pause");
return -1;
}
Mat dstImage;
distanceTrans(srcImage, dstImage);
imshow("距离变换图像", dstImage);
waitKey();
return 0;
}

CV_EXPORTS_AS(distanceTransformWithLabels) void distanceTransform( InputArray src,
OutputArray dst, OutputArray labels, int distanceType, int maskSize,int labelType=DIST_LABEL_CCOMP );
//! computes the distance transform map
CV_EXPORTS_W void distanceTransform( InputArray src, OutputArray dst, int distanceType, int maskSize );


 DIST_USER   = ⑴,  //!< User defineddistance
DIST_L1      = 1,   //!< distance = |x1-x2| + |y1-y2|
DIST_L2      = 2,   //!< the simple euclidean distance
DIST_C       = 3,   //!< distance = max(|x1-x2|,|y1-y2|)
DIST_L12     = 4,   //!< L1-L2 metric: distance =2(sqrt(1+x*x/2) - 1))
DIST_FAIR    = 5,   //!< distance = c^2(|x|/c-log(1+|x|/c)),c = 1.3998
DIST_WELSCH  = 6,   //!< distance = c^2/2(1-exp(-(x/c)^2)), c= 2.9846
DIST_HUBER   = 7    //!< distance = |x|<c ? x^2/2 :c(|x|-c/2), c=1.345


//使用OpenCV中的相应的函数实现距离变换

#include <iostream>
#include <opencv2\core\core.hpp>
#include <opencv2\highgui\highgui.hpp>
#include <opencv2\imgproc\imgproc.hpp>

using namespace cv;
using namespace std;

int main()
{
if (!srcImage.data)
{
cout << "读入图片错误！" << endl;
system("pause");
return -1;
}
//转换为灰度图像
Mat grayImage;
cvtColor(srcImage, grayImage, CV_BGR2GRAY);
//转换为二值图像
Mat binaryImage;
threshold(grayImage, binaryImage, 100, 255, THRESH_BINARY);
//距离变换
Mat dstImage;
//使用欧式距离进行计算
//归一化矩阵
normalize(dstImage, dstImage, 0, 1, NORM_MINMAX);
imshow("二值化图像", binaryImage);
imshow("距离变换后的图像", dstImage);
waitKey();
return 0;
}






2018-05-11 17:07:53 u013921430 阅读数 1220
• ###### MATLAB图像处理

全面系统的学习MATLAB在图像处理中的应用

20789 人正在学习 去看看 魏伟

## 前言

距离变换（distance transform，DT）在图像处理、计算机视觉等领域有非常多的运用，例如骨架化、目标细化、黏连目标分离等。在Matlab和OpenCV中都有进行距离变换的函数。Matlab中有bwdist函数用于二值化图像的距离变换， OpenCV中有cvDistTransform函数用于进行二值化图像的距离变换。

除此之外，Matlab中还提供了graydist函数用于灰度加权距离变换，但是graydist居然要输入一个种子点，这让我觉得很神奇，因为在我的印象中灰度加权距离变换是不需要种子点的！自己测试运行后就发现，这是扯淡呢！算出来的结果是图像中各个点到种子点的距离加权变换。所以，我必须来讲讲我所理解的不需要种子点的灰度加权距离变换。

## 灰度加权距离变换

一般的距离变换需要首先对图像进行二值化，这样的过程会损失掉图像中某些有用的信息，因为单纯的距离变换获得的是前景点到背景的最短距离，损失了灰度信息所代表的物理意义。因此在距离变换的基础上加上图像灰度信息，就有了灰度加权距离变换——GWDT（gray-weighted image distance transform），也可以称之为灰度尺度距离变换——GSDT（gray-scale image distance transform）。

理解了距离变换，灰度加权距离变换就很好理解了，就是在距离变换的基础上加上了灰度信息。那么如何操作呢？我采用类似于区域增长的方法，在快速行进法的框架下来进行实现它，基本步骤如下；

1.    设定阈值分离前景和背景，将背景点的状态设置为ALIVE（其实所有的背景点都是seed point），将前景设为FAR，前景点的值初始化为无限大；

2.    初始化边缘点，在FAR点中选取8邻域（对二维图像进行处理）中有背景点的点，这些点就是边缘点，将边缘点的值赋予原图中的值，更改其状态为TRAIL，并将其放入一个队列Q中；

3.    在队列Q中不放回地取出点，假设取出的点记为y，计算y周围的8个点x，根据以下公式计算x，如果x的值被更新了，判断其是否已经在队列中，如果不在，那么就将其放入队列中。（在我实现的时候，考虑到效率的因素，并没有判断其是否已经在队列中，而是将其直接放入队列中，这个跟queue类和C++函数传参有关系

4.    循环第三步，直到队列Q为空。

从上面的步骤中可以看出，实现过程与快速行进算法还是存在一定的差异的，如果按照快速行进算法的思想，那么对于点y，只会计算周围8临域中属性是FAR的点的值，不会对ALIVE的点进行更新。但是上述的过程，会对ALIVE点进行更新。我这样做是为了尽量保证每个点在变换过后值尽量小，大家可以自行修改代码来决定要不要对ALIVE点进行更新。

如下图所示，假设图中蓝色的点为初始化后的边缘点，他的状态为TRAIL，点A的状态为FAR；依次遍历S1-S4周围的8个点，很显然A是这四个点的共同邻域，A取四个点对其进行变换后的最小值；

图 1 迭代过程示意图标题标题

## 处理结果

图 2 原始图像

图 3  处理后的图像

### 结果分析

原图中我加入了一些噪声，在GWDT后，前景点特别是远离背景的前景点的灰度值相对被增强了，这样背景点就变相的被抑制了。此外，通过设定阈值可以分离区域内的黏连的大米，从而自动进行记数。等等等等，GWDT的应用还有很多。总的来说，GWDT是在传统的距离变换的基础上保留和利用了灰度信息，更加具有物理意义和价值。

## 代码

为了方便大家学习，为大家提供源码；

//-------------------------
//潘正宇 ，2018.05.10
//灰度尺度距离变换，GSDT/GWDT
//-------------------------

#include <iostream>
#include <vector>
#include <opencv2/opencv.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <queue>
#include <string>

using namespace std;
using namespace cv;

#define INF 3.4e+38                        // 无穷大的点
enum{ ALIVE = -1, TRIAL = 0, FAR = 1 };    //设定三种状态；

bool find(queue<int> que,int value)
{
int size = que.size();
for (int i = 0; i <size; ++i)
{
int a = que.front();
if (a == value)
{
return true;
}
que.pop();
}
return false;
}

bool GSDT(string &picPath)
{
Mat Scr = imread(picPath, 0);    //以灰度图的模式读入图像

imshow("原始图", Scr);

Mat gsdtImg;                     //创建GSDT后的图
gsdtImg.create(Scr.size(), CV_32FC1);

int row = Scr.rows;
int col = Scr.cols;

Mat  mat_mean, mat_stddev;
meanStdDev(Scr, mat_mean, mat_stddev);
double mean;
mean = mat_mean.at<double>(0, 0);   //获取图像均值

int *state= new int[col*row];

for (int i = 0; i < row; i++)
{
for (int j = 0; j < col; j++)
{
if (Scr.at<uchar>(i, j)<(mean))
{
gsdtImg.at<float>(i, j) = Scr.at<uchar>(i, j);       //背景点是ALIVE的
state[i*col + j] = ALIVE;
}
else
{
gsdtImg.at<float>(i, j) = INF;
state[i*col + j] = FAR;
}

}
}

queue<int> TrailQue;                     //定义队列用于存放TRAIL的点
for (int i = 1; i < row-1; i++)
{
for (int j = 1; j < col-1; j++)
{
if (state[i*col + j] == FAR)     //如果这个点是前景点，那么搜寻这个点周围是否有背景点，如果是那么它就是边缘点；
{
for (int o = -1; o <= 1; o++)
{
for (int p = -1; p <= 1; p++)
{
if (state[(i + o)*col + j + p] == ALIVE)        //找到所有的边缘点，并且将其放入队列TrailQue；
{
state[i*col + j] = TRIAL;
gsdtImg.at<float>(i, j) = Scr.at<uchar>(i, j);
TrailQue.push(i*col + j);
break;
}
}
}
}

}
}

while (!TrailQue.empty())
{
int P_row = TrailQue.front() / col;    ///获取TrailQue中点的坐标
int P_col = TrailQue.front() % col;

if (P_row < 1){ P_row = 1; }if (P_row>row - 2){ P_row = row - 2; }
if (P_col < 1){ P_col = 1; }if (P_col>col - 2){ P_col = col - 2; }

for (int o = -1; o <= 1; o++)
{
for (int p = -1; p <= 1; p++)
{
int N = (P_row + o)*col + P_col + p;

double len = sqrt(o*o + p*p);
float gs = gsdtImg.at<float>(P_row, P_col) + Scr.at<uchar>(P_row + o, P_col + p)*len;

//---比较该点现有的GSDT值与P点给它的值；
if (gsdtImg.at<float>(P_row + o, P_col + p) > gs)
{
state[N] = TRIAL;
gsdtImg.at<float>(P_row + o, P_col + p) = gs;

TrailQue.push(N);

/*if (!find(TrailQue, N))
{
TrailQue.push(N);
}*/

}

}
}

TrailQue.pop();                        //从TrailQue中删除这个点
}

//寻找最大值；
float Max = 0;
for (int i = 0; i < row; i++)
{
for (int j = 0; j < col; j++)
{
if (gsdtImg.at<float>(i, j)>(INF / 100))
{
gsdtImg.at<float>(i, j) = 0;
}
if (gsdtImg.at<float>(i, j)>Max)//&&gsdtImg.at<float>(i, j)<(INF/10)
{
Max = gsdtImg.at<float>(i, j);
}

}
}

//将图像像素区间压缩到0-255
Mat Dst;
Dst = Scr.clone();
for (int i = 0; i < row; i++)
{
for (int j = 0; j < col; j++)
{
Dst.at<uchar>(i, j) = 255 * gsdtImg.at<float>(i, j) / (Max+1);
}
}

imwrite("temp.bmp", Dst);

//归一化
gsdtImg = gsdtImg / (Max+1);

cout << gsdtImg.at<float>(20, 206);
imshow("GSDT",gsdtImg);
waitKey(0);

return true;
}

void main()
{
string picPath = "G:\\博客\\图像处理\\灰度尺度距离变换\\rice.jpg";
GSDT(picPath);

system("pause");
return;
}