【主动轮廓模型（一）】《Snakes: Active Contour Models》算法原理与OpenCV实现

傷城~ 2022-12-11 11:28 714阅读 0赞

文章目录

  • 1 概述
  • 2 算法原理
    • 2.1 内部能量 E i n t E_{int} Eint​
    • 2.2 图像能量 E i m a g e E_{image} Eimage​
      • （1）线函数 E l i n e E_{line} Eline​
      • （2）边函数 E e d g e E_{edge} Eedge​
      • （3）末端函数 E t e r m E_{term} Eterm​
    • 2.3 外部能量 E c o n E_{con} Econ​
  • 3 模型求解
  • 4 算法实现（OpenCV3）

1 概述

主动轮廓模型（也称Active Contour Model、Snake）是Kass等人在1988年提出的，该算法将图像分割问题转换为求解能量泛函最小值的问题。主要思路是通过构造能量泛函，经过算法迭代，轮廓曲线由初始位置逐渐向使能量函数最小（或局部极小）的图像边缘逼近，最终分割出目标。

2 算法原理

首先需要人为地在图像上给出初始轮廓曲线，确切的说是一组用于控制曲线形状的控制点： v ( s ) = [ x ( s ) , y ( s ) ] s ∈ [ 0 , 1 ] v(s)=[x(s),y(s)] s\in[0,1] v(s)=[x(s),y(s)]s∈[0,1]，这些点收尾相连构成一个封闭的轮廓线。其中 x ( s ) x(s) x(s)和 y ( s ) y(s) y(s)分别表示每个控制点在图像中的坐标位置， s s s是以傅立叶变换形式描述边界的自变量，也可以理解为弧长。则Snake曲线的能量函数表示为：

E s n a k e ∗ = ∫ 0 1 E s n a k e ( v ( s ) ) d s = ∫ 0 1 E i n t ( v ( s ) ) + E i m a g e ( v ( s ) ) + E c o n ( v ( s ) ) d s = ∫ 0 1 E i n t ( v ( s ) ) + E e x t ( v ( s ) ) d s (1) \begin{aligned} E_{snake}^* &= \int_0^1 E_{snake}(v(s))ds \\ &= \int_0^1 E_{int}(v(s))+E_{image}(v(s))+E_{con}(v(s))ds \\ &= \int_0^1 E_{int}(v(s))+E_{ext}(v(s))ds \end{aligned} \tag{1} Esnake∗​​=∫01​Esnake​(v(s))ds=∫01​Eint​(v(s))+Eimage​(v(s))+Econ​(v(s))ds=∫01​Eint​(v(s))+Eext​(v(s))ds​(1)

其中， E i n t E_{int} Eint​为内部能量， E i m a g e E_{image} Eimage​为图像能量， E c o n E_{con} Econ​为外部约束能量。其中图像能量和外部约束能量统称为外部能量，即： E e x t = E i m a g e + E c o n E_{ext}=E_{image}+E_{con} Eext​=Eimage​+Econ​。

2.1 内部能量 E i n t E_{int} Eint​

内部能量 E i n t E_{int} Eint​由保证曲线连续性的一阶导和保证曲线平滑的二阶导组成，表示为：

E i n t = 1 2 ( α ( s ) ∣ v s ( s ) ∣ 2 + β ( s ) ∣ v s s ( s ) ∣ 2 ) (2) E_{int}=\frac{1}{2} (\alpha(s)|v_s(s)|^2+\beta(s)|v_{ss}(s)|^2) \tag{2} Eint​=21​(α(s)∣vs​(s)∣2+β(s)∣vss​(s)∣2)(2)

通过调整权值 α ( s ) \alpha(s) α(s)和 β ( s ) \beta(s) β(s)可以控制曲线的形状。例如将 β ( s ) \beta(s) β(s)置为0可以让曲线最终出现拐角，即曲线二阶不连续。

2.2 图像能量 E i m a g e E_{image} Eimage​

在低层次的计算机视觉中，需要能够将轮廓吸引到特定图像特征的能量函数。原始的Snake轮廓模型中提出了3种不同的能量函数，分别将snake轮廓吸引到线、边和末端。完整的图像能量 E i m a g e E_{image} Eimage​可以表示为这3个能量函数的权值组合，通过调整这3个权值，可以形成不同的轮廓形状。

E i m a g e = ω l i n e E l i n e + ω e d g e E e d g e + ω t e r m E t e r m (3) E_{image}=\omega_{line}E_{line}+\omega_{edge}E_{edge}+\omega_{term}E_{term} \tag{3} Eimage​=ωline​Eline​+ωedge​Eedge​+ωterm​Eterm​(3)

（1）线函数 E l i n e E_{line} Eline​

最简单直接且有用的图像能量函数是图像本身，即图像本身的灰度。令：

E l i n e = I ( x , y ) (4) E_{line}=I(x,y) \tag{4} Eline​=I(x,y)(4)

控制 ω l i n e \omega_{line} ωline​的正负号可以控制轮廓被吸引到较暗的线或是较亮的线，也就是使轮廓试图靠近轮廓的最暗或最亮处。

然而，如果snake轮廓的一部分到达了一个低能量的图像特征位置，这个样条项将推动snake轮廓临近的部分朝着这个特征可能的延续方向移动，这会使其在一个最优的局部最小位置引入一个较大的能量。一种解决方案是允许snake轮廓和模糊能量函数平衡，然后慢慢降低模糊程度。

Marr和Hildreth在他们的论文里证明了“灰度的突然变化会在一阶导数中引起波峰或波谷，或在二阶导数中等效地引起零交点”。为了显示图像尺度空间连续性和Marr-Hildreth边缘检测理论的关系，snake模型中采用了模糊的边能量函数：

E l i n e = − ( G σ ∗ ∇ 2 I ) 2 (6) E_{line}=-(G_{\sigma}*\nabla^2I)^2 \tag{6} Eline​=−(Gσ​∗∇2I)2(6)

其中 G σ G_{\sigma} Gσ​是标准差为 σ \sigma σ的高斯函数，该函数的最小值位于在Marr-Hildreth理论中被定义的 G σ ∗ ∇ 2 I G_{\sigma}*\nabla^2I Gσ​∗∇2I的零交点处。在能量函数中加入这一项意味着snake轮廓在被吸引到零交点的同时，仍然受到它自己的平滑限制。

（2）边函数 E e d g e E_{edge} Eedge​

在图像上找边缘可以通过梯度来实现，令：

E e d g e = − ∣ ∇ I ( x , y ) ∣ 2 (5) E_{edge}=-|\nabla I(x,y)|^2 \tag{5} Eedge​=−∣∇I(x,y)∣2(5)

在某个点上，梯度越大上式的能量越小，则snake轮廓将被吸引到梯度较大的区域。

（3）末端函数 E t e r m E_{term} Eterm​

为了找到轮廓的终止位置，作者将平滑过图像中等高线的曲率加入到能量函数中。令 C ( x , y ) = G σ ( x , y ) ∗ I ( x , y ) C(x,y)=G_{\sigma}(x,y)*I(x,y) C(x,y)=Gσ​(x,y)∗I(x,y)是模糊后的图像， θ = t a n − 1 ( C y C x ) \theta=tan^{-1}(\frac{C_y}{C_x}) θ=tan−1(Cx​Cy​​)是梯度角， n = ( c o s θ , s i n θ ) n=(cos\theta,sin\theta) n=(cosθ,sinθ)、 n ⊥ = ( − s i n θ , c o s θ ) n_{\perp}=(-sin\theta,cos\theta) n⊥​=(−sinθ,cosθ)是沿着和垂直梯度方向的单位向量。则 C ( x , y ) C(x,y) C(x,y)中等高线的曲率可以表示为：

E t e r m = ∂ θ ∂ n ⊥ = ∂ C 2 / ∂ n ⊥ 2 ∂ C / ∂ n = C y y C x 2 − 2 C x y C x C y + C x x C y 2 ( C x 2 + C y 2 ) 3 / 2 (6) \begin{aligned} E_{term}&=\frac{\partial \theta}{\partial n_{\perp}} \\ &=\frac{\partial C^2 / \partial n_{\perp}^2}{\partial C / \partial n} \\ &=\frac{C_{yy}C_x^2-2C_{xy}C_xC_y+C_{xx}C_y^2}{(C_x^2+C_y^2)^{3/2}} \end{aligned} \tag{6} Eterm​​=∂n⊥​∂θ​=∂C/∂n∂C2/∂n⊥2​​=(Cx2​+Cy2​)3/2Cyy​Cx2​−2Cxy​Cx​Cy​+Cxx​Cy2​​​(6)

通过组合 E e d g e E_{edge} Eedge​和 E t e r m E_{term} Eterm​，我们能够创建一个被边和末端吸引的snake曲线。

2.3 外部能量 E c o n E_{con} Econ​

外部能量 E c o n E_{con} Econ​来自于外部的约束力。在论文原文中，作者给了一个外部约束的例子， 包括外部的固定点、连接两条snake曲线的锚点或鼠标拖动的点。例如，为了在 x 1 x_1 x1​和 x 2 x_2 x2​点之间创建一个连接的弹性力，就可以将 − k ( x 1 − x 2 ) 2 -k(x_1-x_2)^2 −k(x1​−x2​)2添加到外部能量 E c o n E_{con} Econ​中。

在这里插入图片描述

3 模型求解

由欧拉方程，求解能量 E s n a k e ∗ E_{snake}^* Esnake∗​的最小值（局部极小值），公式（1）的导数必须满足：

α v s s ( s ) − β v s s s s ( s ) − ∇ E i m a g e ( v ( s ) ) − ∇ E c o n ( v ( s ) ) = 0 (7) \alpha v_{ss}(s) - \beta v_{ssss}(s) - \nabla E_{image}(v(s)) - \nabla E_{con}(v(s)) = 0 \tag{7} αvss​(s)−βvssss​(s)−∇Eimage​(v(s))−∇Econ​(v(s))=0(7)

由 v ( s ) = [ x ( s ) , y ( s ) ] v(s)=[x(s),y(s)] v(s)=[x(s),y(s)]，将上式改写成x和y两个方向有：

{ α x s s ( s ) − β x s s s s ( s ) − ∂ E e x t ∂ x = 0 α y s s ( s ) − β y s s s s ( s ) − ∂ E e x t ∂ y = 0 (8) \begin{cases} \alpha x_{ss}(s) - \beta x_{ssss}(s) - \frac{\partial{E}_{ext}}{\partial x} = 0 \\ \alpha y_{ss}(s) - \beta y_{ssss}(s) - \frac{\partial{E}_{ext}}{\partial y} = 0 \end{cases} \tag{8} { αxss​(s)−βxssss​(s)−∂x∂Eext​​=0αyss​(s)−βyssss​(s)−∂y∂Eext​​=0​(8)

利用顺序连接的锚点的坐标差分来近似导数：

{ x s s ( s ) = x ( s + 1 ) + x ( s − 1 ) − 2 x ( s ) x s s s s ( s ) = ( x ( s + 2 ) + x ( s ) − 2 x ( s + 1 ) ) + ( x ( s ) + x ( s − 2 ) − 2 x ( s − 1 ) ) − 2 ( x ( s + 1 ) + x ( s − 1 ) − 2 x ( s ) ) (9) \begin{cases} x_{ss}(s)=x(s+1)+x(s-1)-2x(s) \\ \begin{aligned} x_{ssss}(s)=&(x(s+2)+x(s)-2x(s+1)) \\ &+(x(s)+x(s-2)-2x(s-1)) \\ &-2(x(s+1)+x(s-1)-2x(s)) \end{aligned} \end{cases} \tag{9} ⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧​xss​(s)=x(s+1)+x(s−1)−2x(s)xssss​(s)=​(x(s+2)+x(s)−2x(s+1))+(x(s)+x(s−2)−2x(s−1))−2(x(s+1)+x(s−1)−2x(s))​​(9)

将（9）带入（8），同时令 f x = ∂ E e x t ∂ x , f y = ∂ E e x t ∂ y f_x=\frac{\partial E_{ext}}{\partial x}, f_y=\frac{\partial E_{ext}}{\partial y} fx​=∂x∂Eext​​,fy​=∂y∂Eext​​有：

{ β x ( s − 2 ) − ( α + 4 β ) x ( s − 1 ) + ( 2 α + 6 β ) x ( s ) − ( α + 4 β ) x ( s + 1 ) + β x ( s + 2 ) + f x = 0 β y ( s − 2 ) − ( α + 4 β ) y ( s − 1 ) + ( 2 α + 6 β ) y ( s ) − ( α + 4 β ) y ( s + 1 ) + β y ( s + 2 ) + f y = 0 (10) \begin{cases} \beta x(s-2)-(\alpha+4\beta)x(s-1)+(2\alpha+6\beta)x(s)-(\alpha+4\beta)x(s+1)+\beta x(s+2)+f_x=0 \\ \beta y(s-2)-(\alpha+4\beta)y(s-1)+(2\alpha+6\beta)y(s)-(\alpha+4\beta)y(s+1)+\beta y(s+2)+f_y=0 \end{cases} \tag{10} { βx(s−2)−(α+4β)x(s−1)+(2α+6β)x(s)−(α+4β)x(s+1)+βx(s+2)+fx​=0βy(s−2)−(α+4β)y(s−1)+(2α+6β)y(s)−(α+4β)y(s+1)+βy(s+2)+fy​=0​(10)

令：

{ a = 2 α + 6 β b = − ( α + 4 β ) c = β (11) \begin{cases} a=2\alpha + 6\beta \\ b=-(\alpha+4\beta) \\ c=\beta \end{cases} \tag{11} ⎩⎪⎨⎪⎧​a=2α+6βb=−(α+4β)c=β​(11)

则将（10）写成矩阵形式为：

{ A x + f x ( x , y ) = 0 A y + f y ( x , y ) = 0 (12) \begin{cases} Ax+f_x(x,y)=0 \\ Ay+f_y(x,y)=0 \end{cases} \tag{12} { Ax+fx​(x,y)=0Ay+fy​(x,y)=0​(12)

其中A为五对角带状矩阵：

A = [ a b c ⋯ c b b a b c ⋯ c c b a b c ⋯ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋯ c b a b c c ⋯ c b a b b c ⋯ c b a ] , x = [ x 1 x 2 x 3 ⋮ x n − 1 x n x 1 ] , f x = [ f x 1 f x 2 f x 3 ⋮ f x n − 1 f x n f x 1 ] (13) A= \begin{bmatrix} a & b & c & \cdots & c & b \\ b & a & b & c & \cdots & c \\ c & b & a & b & c & \cdots \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \cdots & c & b & a & b & c \\ c & \cdots & c & b & a & b \\ b & c & \cdots & c & b & a \end{bmatrix}, x= \begin{bmatrix} x_1 \\ x_2 \\ x_3 \\ \vdots \\ x_{n-1} \\ x_n \\ x_1 \end{bmatrix}, f_x= \begin{bmatrix} f_{x_1} \\ f_{x_2} \\ f_{x_3} \\ \vdots \\ f_{x_{n-1}} \\ f_{x_n} \\ f_{x_1} \end{bmatrix} \tag{13} A=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡​abc⋮⋯cb​bab⋮c⋯c​cba⋮bc⋯​⋯cb⋮abc​c⋯c⋮bab​bc⋯⋮cba​⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤​,x=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡​x1​x2​x3​⋮xn−1​xn​x1​​⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤​,fx​=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡​fx1​​fx2​​fx3​​⋮fxn−1​​fxn​​fx1​​​⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤​(13)

利用梯度下降法，在第t次迭代有：

{ A x t + f x ( x t − 1 , y t − 1 ) = − γ ( x t − x t − 1 ) A y t + f y ( y t − 1 , y t − 1 ) = − γ ( y t − y t − 1 ) (14) \begin{cases} Ax_t+f_x(x_{t-1},y_{t-1})=-\gamma(x_t-x_{t-1}) \\ Ay_t+f_y(y_{t-1},y_{t-1})=-\gamma(y_t-y_{t-1}) \end{cases} \tag{14} { Axt​+fx​(xt−1​,yt−1​)=−γ(xt​−xt−1​)Ayt​+fy​(yt−1​,yt−1​)=−γ(yt​−yt−1​)​(14)

其中 γ \gamma γ是迭代步长，公式（14）的解为：

{ x t = ( A + γ I ) − 1 ( x t − 1 − f x ( x t − 1 , y t − 1 ) ) y t = ( A + γ I ) − 1 ( y t − 1 − f y ( x t − 1 , y t − 1 ) ) (15) \begin{cases} x_t=(A+\gamma I)^{-1}(x_{t-1}-f_x(x_{t-1},y_{t-1})) \\ y_t=(A+\gamma I)^{-1}(y_{t-1}-f_y(x_{t-1},y_{t-1})) \end{cases} \tag{15} { xt​=(A+γI)−1(xt−1​−fx​(xt−1​,yt−1​))yt​=(A+γI)−1(yt−1​−fy​(xt−1​,yt−1​))​(15)

由于 A A A为五对角带状矩阵，因此 A + γ I A+\gamma I A+γI也是一个五对角带状矩阵，原文中作者用LU分解来求其逆矩阵。

4 算法实现（OpenCV3）

网上关于snake算法的实现多是matlab代码。OpenCV2中可用cvSnakeImage来实现，但这个函数在OpenCV3中已经被删除。本文将在OpenCV3.14中实现snake算法代码。

  1. cv::Mat Interate(
  2.     cv::Mat image,
  3.     cv::Mat xs,
  4.     cv::Mat ys,
  5.     double alpha,
  6.     double beta,
  7.     double gamma,
  8.     double kappa,
  9.     double wl,
  10.     double we,
  11.     double wt,
  12.     int iterations
  13. )
  14. { 
  15.     // 相关参数
  16.     int N = iterations;
  17.     cv::Mat smth = image.clone();
  19.     // 图像大小
  20.     qDebug() << "Calculating size of image";
  21.     cv::Size size = image.size();
  22.     int row = size.height;
  23.     int col = size.width;
  25.     // 计算外部力（图像力）
  26.     qDebug() << "Computing external forces";
  27.     cv::Mat E_line = smth.clone(); // E_line is simply the image intensities
  29.     cv::Mat gradx, grady;
  30.     cv::Sobel(smth, gradx, smth.depth(), 1, 0, 1, 1, 0, cv::BORDER_CONSTANT);
  31.     cv::Sobel(smth, grady, smth.depth(), 0, 1, 1, 1, 0, cv::BORDER_CONSTANT);
  33.     qDebug() << "Computing gradx and grady";
  34.     cv::Mat E_edge(row, col, CV_32FC1);
  35.     for (int i = 0; i < gradx.rows; i++)
  36.     { 
  37.         for (int j = 0; j < gradx.cols; j++)
  38.         { 
  39.             float v_gradx = gradx.at<float>(i, j);
  40.             float v_grady = grady.at<float>(i, j);
  42.             E_edge.at<float>(i, j) = -1 * std::sqrt(v_gradx * v_gradx + v_grady * v_grady); // E_edge is measured by gradient in the image
  43.         }
  44.     }
  46.     // 导数mask
  47.     qDebug() << "masks for taking various derivatives";
  48.     cv::Mat m1 = (cv::Mat_<float>(1, 2) << -1, 1);
  49.     cv::Mat m2 = (cv::Mat_<float>(2, 1) << -1, 1);
  50.     cv::Mat m3 = (cv::Mat_<float>(1, 3) << 1, -2, 1);
  51.     cv::Mat m4 = (cv::Mat_<float>(3, 1) << 1, -2, 1);
  52.     cv::Mat m5 = (cv::Mat_<float>(2, 2) << 1, -1, -1, 1);
  54.     cv::Mat cx, cy, cxx, cyy, cxy;
  55.     filter2D(smth, cx, -1, m1);
  56.     filter2D(smth, cy, -1, m2);
  57.     filter2D(smth, cxx, -1, m3);
  58.     filter2D(smth, cyy, -1, m4);
  59.     filter2D(smth, cxy, -1, m5);
  61.     // 计算 E_term
  62.     cv::Mat E_term(row, col, CV_32FC1);
  63.     for (int i = 0; i < row; i++)
  64.     { 
  65.         for (int j = 0; j < col; j++)
  66.         { 
  67.             int v_cx = cx.at<float>(i, j);
  68.             int v_cy = cy.at<float>(i, j);
  69.             int v_cxx = cxx.at<float>(i, j);
  70.             int v_cyy = cyy.at<float>(i, j);
  71.             int v_cxy = cxy.at<float>(i, j);
  73.             E_term.at<float>(i, j) = (v_cyy*v_cx*v_cx - 2 * v_cxy*v_cx*v_cy + v_cxx * v_cy*v_cy) / (std::pow((1 + v_cx * v_cx + v_cy * v_cy), 1.5));
  74.         }
  75.     }
  77.     // 计算E_ext
  78.     cv::Mat E_ext = (wl*E_line + we * E_edge - wt * E_term);
  80.     // 计算梯度
  81.     cv::Mat fx, fy;
  82.     cv::Sobel(E_ext, fx, E_ext.depth(), 1, 0, 1, 0.5, 0, cv::BORDER_CONSTANT);
  83.     cv::Sobel(E_ext, fy, E_ext.depth(), 0, 1, 1, 0.5, 0, cv::BORDER_CONSTANT);
  85.     cv::transpose(xs, xs);
  86.     cv::transpose(ys, ys);
  88.     int m = xs.rows;
  89.     int n = 1;
  91.     int mm = fx.cols;
  92.     int nn = fx.rows;
  94.     // 计算五对角状矩阵，b(i)表示vi系数(i = i - 2 到 i + 2)
  95.     double b[5];
  96.     b[0] = beta;
  97.     b[1] = -(alpha + 4 * beta);
  98.     b[2] = 2 * alpha + 6 * beta;
  99.     b[3] = b[1];
  100.     b[4] = b[0];
  102.     cv::Mat A = cv::Mat::eye(m, m, CV_32FC1);
  103.     cv::Mat eyeMat0 = cv::Mat::eye(m, m, CV_32FC1);
  104.     circRowShift(eyeMat0, 2);
  105.     eyeMat0.convertTo(eyeMat0, CV_32FC1);
  106.     A = b[0] * eyeMat0;
  108.     cv::Mat eyeMat1 = cv::Mat::eye(m, m, CV_32FC1);
  109.     circRowShift(eyeMat1, 1);
  110.     eyeMat1.convertTo(eyeMat1, CV_32FC1);
  111.     A = A + b[1] * eyeMat1;
  113.     cv::Mat eyeMat2 = cv::Mat::eye(m, m, CV_32FC1);
  114.     circRowShift(eyeMat2, 0);
  115.     eyeMat2.convertTo(eyeMat2, CV_32FC1);
  116.     A = A + b[2] * eyeMat2;
  118.     cv::Mat eyeMat3 = cv::Mat::eye(m, m, CV_32FC1);
  119.     circRowShift(eyeMat3, -1);
  120.     eyeMat3.convertTo(eyeMat3, CV_32FC1);
  121.     A = A + b[3] * eyeMat3;
  123.     cv::Mat eyeMat4 = cv::Mat::eye(m, m, CV_32FC1);
  124.     circRowShift(eyeMat4, -2);
  125.     eyeMat4.convertTo(eyeMat4, CV_32FC1);
  126.     A = A + b[4] * eyeMat4;
  128.     // 计算矩阵的逆
  129.     cv::Mat Ainv(A.size(), CV_32FC1);
  130.     A = A + gamma * cv::Mat::eye(m, m, CV_32FC1);
  131.     cv::invert(A, Ainv); // Computing Ainv 
  133.     cv::Mat srcImg = cv::imread("D:/endo_image.jpg");
  134.     // 迭代更新曲线
  135.     for (int i = 0; i < N; i++)
  136.     { 
  137.         cv::Mat intFx(fx.size(), CV_32FC1);
  138.         cv::Mat intFy(fy.size(), CV_32FC1);
  140.         cv::remap(fx, intFx, xs, ys, cv::INTER_LINEAR, cv::BORDER_CONSTANT);
  141.         cv::remap(fy, intFy, xs, ys, cv::INTER_LINEAR, cv::BORDER_CONSTANT);
  143.         cv::Mat ssx(xs.size(), CV_32FC1);
  144.         cv::Mat ssy(ys.size(), CV_32FC1);
  145.         for (int k = 0; k < xs.rows; k++)
  146.         { 
  147.             for (int l = 0; l < xs.cols; l++)
  148.             { 
  149.                 ssx.at<float>(k, l) = gamma * xs.at<float>(k, l) - kappa * intFx.at<float>(k, l);
  150.                 ssy.at<float>(k, l) = gamma * ys.at<float>(k, l) - kappa * intFy.at<float>(k, l);
  151.             }
  152.         }
  154.         // 更新曲线位置
  155.         xs = Ainv * ssx;
  156.         ys = Ainv * ssy;
  158.         cv::Mat resultImg = srcImg.clone();
  159.         for (int j = 0; j < xs.rows; j++)
  160.         { 
  161.             cv::Point center = cv::Point(xs.at<float>(j, 0), ys.at<float>(j, 0));
  162.             cv::circle(resultImg, center, 4, cv::Scalar(0, 255, 0), 2);
  163.         }
  165.         // 显示
  166.         cv::imshow("result", resultImg);
  167.         cv::waitKey(30);
  168.     }
  170.     return image;
  171. }

运行结果：

在这里插入图片描述

发表评论

表情：
评论列表 （有 0 条评论，714人围观）

还没有评论，来说两句吧...

相关阅读

    相关 Snake 算法

    分享一下我老师大神的人工智能教程！零基础，通俗易懂！[http://blog.csdn.net/jiangjunshow][http_blog.csdn.net_jiangju

    梦里梦外;梦里梦外;/ 0 赞/ 195 阅读