跳轉到內容

分形/複平面迭代/等勢

來自華夏公益教科書,開放的世界,開放的書籍

邊界追蹤,也稱為輪廓追蹤,指的是二進位制數字區域的邊界追蹤可以看作是一種分割技術,它識別數字區域的邊界畫素。邊界追蹤是分析該區域的重要第一步。邊界是一個拓撲概念。然而,數字影像不是拓撲空間。因此,不可能在數學上精確地定義數字影像中邊界的概念。關於追蹤數字影像 I 子集 S 邊界的多數出版物都描述了演算法,這些演算法找到屬於 S 的一組畫素,這些畫素在其直接鄰域內具有屬於 S 和其補集 I - S 的畫素。根據這個定義,子集 S 的邊界不同於補集 I - S 的邊界,這是一個拓撲悖論。為了正確地定義邊界,有必要引入一個與給定數字影像相對應的拓撲空間。這樣的空間可以是二維抽象細胞復形。它包含三個維度的細胞:與數字影像畫素相對應的二維細胞,表示位於兩個相鄰畫素之間的短線的單維細胞或“裂縫”,以及對應於畫素角點的零維細胞或“點”。然後,子集 S 的邊界是裂縫和點的序列,而這些裂縫和點的鄰域與子集 S 及其補集 I - S 相交。以這種方式定義的邊界完全對應於拓撲定義,也對應於我們對邊界的直觀想象,因為 S 的邊界既不包含 S 的元素,也不包含其補集的元素。它應該只包含位於 S 和補集之間的元素。這正是復形的裂縫和點。這種邊界追蹤方法在弗拉基米爾·A·科瓦列夫斯基[1] 的書和網站上有所描述。 [2]


  • 邊界追蹤
  • 函式 f 的形狀 = 臨界水平曲線(包含 f 的臨界點的水平曲線)的形狀 [3]

演算法

[編輯 | 編輯原始碼]

用於邊界追蹤的演算法:[4]

  • 方形追蹤演算法 [5]
  • Moore 鄰域追蹤演算法
  • 徑向掃描 [6]
  • Theo Pavlidis 演算法 [7]
  • 可以使用向量代數追蹤邊界的通用方法可以在以下位置找到。 [8]
  • 邊界追蹤的擴充套件,用於將追蹤的邊界分割成開放和封閉子部分,在以下位置有所描述。 [9]
  • 水平集邊界的邊緣檢測:水平集方法 (LSM) + 邊緣檢測(例如 Sobel 濾波器)
  • 追蹤等值曲線(等值線)
  • 吸引子周圍的封閉曲線(圓)的反向迭代


  • 函式 f(x,y)
    • 數學方程(用於符號計算)
    • 程式子程式(用於數值計算)
  • 平面
    • 浮點數矩陣(標量場) [10] [11] (用於數值計算)



Moore 鄰域追蹤演算法

[編輯 | 編輯原始碼]
Moore 鄰域由九個細胞組成:一箇中心細胞和圍繞它的八個細胞。

在元胞自動機中,Moore 鄰域定義在二維方形格子上,由一箇中心細胞和圍繞它的八個細胞組成。


這個鄰域以元胞自動機理論的先驅愛德華·F·摩爾命名。

重要性

[編輯 | 編輯原始碼]

它是兩種最常用的鄰域型別之一,另一種是馮·諾依曼鄰域,它不包括角點細胞。例如,著名的康威生命遊戲就使用了 Moore 鄰域。它類似於計算機圖形學中 8 連通畫素的概念。

一個細胞的 Moore 鄰域是指該細胞本身以及與它切比雪夫距離為 1 的細胞。

這個概念可以擴充套件到更高維度,例如,為三維元胞自動機形成一個 26 細胞立方鄰域,就像 3D 生命 所用的一樣。在維度 d 中,其中 ,鄰域的大小為 3d - 1。

在二維中,給定範圍 r擴充套件 Moore 鄰域中細胞的數量為 (2r + 1)2

演算法

[編輯 | 編輯原始碼]

Moore 鄰域公式背後的想法是找到給定圖的輪廓。這個想法對 18 世紀的大多數分析師來說是一個巨大的挑戰,因此從Moore 圖 推匯出了一種演算法,後來被稱為 Moore 鄰域演算法。

Moore 鄰域追蹤演算法的虛擬碼如下:

Input: A square tessellation, T, containing a connected component P of black cells.
Output: A sequence B (b1, b2, ..., bk) of boundary pixels i.e. the contour.
Define M(a) to be the Moore neighborhood of pixel a.
Let p denote the current boundary pixel.
Let c denote the current pixel under consideration i.e. c is in M(p).
Let b denote the backtrack of c (i.e. neighbor pixel of p that was previously tested)
 
Begin
  Set B to be empty.
  From bottom to top and left to right scan the cells of T until a black pixel, s, of P is found.
  Insert s in B.
  Set the current boundary point p to s i.e. p=s
  Let b = the pixel from which s was entered during the image scan.
  Set c to be the next clockwise pixel (from b) in M(p).
  While c not equal to s do
    If c is black
      insert c in B
      Let b = p
      Let p = c
      (backtrack: move the current pixel c to the pixel from which p was entered)
      Let c = next clockwise pixel (from b) in M(p).
    else
      (advance the current pixel c to the next clockwise pixel in M(p) and update backtrack)
      Let b = c
      Let c = next clockwise pixel (from b) in M(p).
    end If
  end While
End

終止條件

[編輯 | 編輯原始碼]

最初終止條件是在第二次訪問起始畫素後停止。 這限制了演算法將完全遍歷的輪廓集。 Jacob Eliosoff 提出的改進停止條件是在您最初進入相同方向的第二個時間進入起始畫素後停止。

David Rand 跟蹤

[編輯 | 編輯原始碼]

David Rand 描述[12]

The plotting algorithm can be summarized thus:

a) Given a matrix of floating point values which are the values of a function z = f(x,y) given at the nodes of a grid of x and y values (the grid values are assumed equally spaced, although the horizontal and vertical spacing may differ), the program determines the minimum and maximum values of z and then computes a number of contour values (in this implementation, 10 values) by linear or logarithmic interpolation between the extrema.

b) The program "walks" about the grid of points looking for any segment (i.e. a line joining two adjacent nodes in the grid) which must be crossed by one of the contours because some contour value lies between the values of z at the nodes.

c) Having found such a segment, it finds the intersection point of the contour and the segment by linear interpolation between the nodes. It also stores the information that the current contour value has been located on the current segment, so that this operation will not be repeated.

d) The program then attempts to locate a neighbouring segment having a similar property - that is, crossed by the same contour. If it finds one, it determines the intersection point as in step c) and then draws a straight line segment joining the previous intersection point with the current one. This step is repeated until no such neighbour can be found, taking care to exclude any segment which has already been dealt with.

e) Steps b), c) and d) are repeated until no segment can be found whose intersection with any contour value has not already been processed.



方形跟蹤演算法

[編輯 | 編輯原始碼]

方形跟蹤演算法簡單而有效。 它的行為完全基於您是在黑色還是白色單元格上(假設白色單元格是形狀的一部分)。 首先,從左上角掃描到右側,並逐行掃描。 進入第一個白色單元格後,演算法的核心就開始了。 它主要包含兩個規則

  • 如果您在白色單元格中,則向左移動。
  • 如果您在黑色單元格中,則向右移動。

請記住,您如何進入當前單元格很重要,因此可以定義左右。

public void GetBoundary(byte[,] image)
{
    for (int j = 0; j < image.GetLength(1); j++)
        for (int i = 0; i < image.GetLength(0); i++)
            if (image[i, j] == 255)               // Found first white pixel
                SquareTrace(new Point(i, j));
}

public void SquareTrace(Point start)
{
    HashSet<Point> boundaryPoints = new HashSet<Point>();  // Use a HashSet to prevent double occurrences
    // We found at least one pixel
    boundaryPoints.Add(start);

    // The first pixel you encounter is a white one by definition, so we go left. 
    // Our initial direction was going from left to right, hence (1, 0)
    Point nextStep = GoLeft(new Point(1, 0));               
    Point next = start + nextStep;
    while (next != start)
    {
        // We found a black cell, so we go right and don't add this cell to our HashSet
        if (image[next.x, next.y] == 0)
        {
            next = next - nextStep;
            nextStep = GoRight(nextStep);
            next = next + nextStep;
        }
        // Alternatively we found a white cell, we do add this to our HashSet
        else
        {
            boundaryPoints.Add(next);
            nextStep = GoLeft(nextStep);
            next = next + nextStep;
        }
    }
}

private Point GoLeft(Point p) => new Point(p.y, -p.x);
private Point GoRight(Point p) => new Point(-p.y, p.x);

測試值

[編輯 | 編輯原始碼]
  • 在引數平面上
    • c = -0.752174630125934 +0.057670473174377*i。 勢能 = 2.108396223015232e-16
    • c = -0.750087705298239 +0.014275963195317 *i,(Mandel 不能繪製這樣的等勢線)

程式碼

[編輯 | 編輯原始碼]

來自 Maxima CAS 的 Lisp 程式碼

函式 ploteq

 ploteq (exp, ...options...)

繪製 exp 的等勢曲線,exp 應該是依賴於兩個變數的表示式。 曲線是透過積分定義正交軌跡的微分方程得到的,該正交軌跡是透過給定表示式的梯度得到的自治系統的解。 圖表也可以顯示該梯度系統的積分曲線(選項 fieldlines)。

即使從控制檯中的 Maxima 會話執行,該程式也需要 Xmaxima,因為繪圖將由 Xmaxima 中的 Tk 指令碼建立。 預設情況下,繪圖區域將為空,直到使用者單擊一個點(或在設定選單中或透過 trajectory_at 選項提供其座標)。

plotdf 接受的大多數選項也可以用於 ploteq,並且繪圖介面與 plotdf 中描述的相同。

示例

  load("plotdf.lisp")$
  V: 900/((x+1)^2+y^2)^(1/2)-900/((x-1)^2+y^2)^(1/2)$
  ploteq(V,[x,-2,2],[y,-2,2],[fieldlines,"blue"])$


單擊一個點將繪製經過該點的等勢曲線(紅色)和正交軌跡(藍色)。


原始碼檔案

  • share/dynamics/dynamics.mac
  • share/diffequations/drawdf.mac
  • share/dynamics/plotdf.lisp


;; plotdf.mac - Adds a function plotdf() to Maxima, which draws a Direction
;;              Field for an ordinary 1st order differential equation,
;;              or for a system of two autonomous 1st order equations.
;;   
;; Copyright (C) 2004, 2008, 2011 Jaime E. Villate <villate@fe.up.pt>
;;


  
;; plot equipotential curves for a scalar field f(x,y)
(defun $ploteq (fun &rest options)
  
  (let (file cmd mfun (opts " ") (s1 '$x) (s2 '$y))
    (setf mfun `((mtimes) -1 ,fun))
    ;; if the variables are not x and y, their names must be given right
    ;; after the expression for the ode's
    (when
      (and (listp (car options)) (= (length (car options)) 3)
            (or (symbolp (cadar options)) ($subvarp (cadar options)))
            (or (symbolp (caddar options)) ($subvarp (caddar options))))
      (setq s1 (cadar options))
      (setq s2 (caddar options))
      (setq options (cdr options)))
    (defun subxy (expr)
      (if (listp expr)
          (mapcar #'subxy expr)
          (cond ((eq expr s1) '$x) ((eq expr s2) '$y) (t expr))))
    (setf mfun (mapcar #'subxy mfun))
;; the next two lines should take into account parameters given in the options
;;    (if (delete '$y (delete '$x (rest (mfuncall '$listofvars ode))))
;;        (merror "The equation(s) can depend only on 2 variable which must be specified!"))
    (setq cmd (concatenate 'string " -dxdt \""
			   (expr_to_str (mfuncall '$diff mfun '$x))
			   "\" -dydt \""
			   (expr_to_str (mfuncall '$diff mfun '$y)) 
			   "\" "))
    
    ;; parse options and copy them to string opts
    (cond (options
           (dolist (v options) 
             (setq opts (concatenate 'string opts " "
                                  (plotdf-option-to-tcl v s1 s2))))))

    (unless (search "-vectors " opts)
      (setq opts (concatenate 'string opts " -vectors {}")))
    (unless (search "-fieldlines " opts)
      (setq opts (concatenate 'string opts " -fieldlines {}")))
    (unless (search "-curves " opts)
      (setq opts (concatenate 'string opts " -curves {red}")))
    (unless (search "-windowtitle " opts)
      (setq opts (concatenate 'string opts " -windowtitle {Ploteq}")))
    (unless (search "-xaxislabel " opts)
      (setq opts (concatenate 'string opts " -xaxislabel " (ensure-string s1))))
    (unless (search "-yaxislabel " opts)
      (setq opts (concatenate 'string opts " -yaxislabel " (ensure-string s2))))
							      
    (setq file (plot-temp-file (format nil "maxout~d.xmaxima" (getpid))))
    (show-open-plot
     (with-output-to-string
         (st)
       (cond ($show_openplot (format st "plotdf ~a ~a~%" cmd opts))
             (t (format st "{plotdf ~a ~a}" cmd opts))))
     file)
    file))

參考資料

[編輯 | 編輯原始碼]
  1. Kovalevsky,V.,影像處理與蜂窩拓撲,施普林格 2021,ISBN 978-981-16-5771-9
  2. http://www.kovalevsky.de,講座“跟蹤二維影像中的邊界”
  3. mathoverflow.net 問題:有理函式的形狀是什麼 ?
  4. 輪廓跟蹤演算法
  5. Abeer George Ghuneim:方形跟蹤演算法
  6. Abeer George Ghuneim:徑向掃描演算法
  7. Abeer George Ghuneim:Theo Pavlidis 演算法
  8. 基於向量代數的二值影像中物體外部和內部邊界的跟蹤,工程科學進展雜誌第 3 卷第 1 期,2010 年 1 月至 6 月,第 57-70 頁 [1]
  9. 基於圖論的跟蹤邊界到開放和封閉子部分的分割,計算機視覺與影像理解,第 115 卷第 11 期,2011 年 11 月,第 1552-1558 頁 [2]
  10. David Rand 編寫的 Java 輪廓繪圖
  11. stackoverflow 問題:實現輪廓繪圖的方法
  12. ContourPlottinginJava 由 D Rand 開發
  13. 由 Josef Boehm 和 Leon Magiera 解決的...
華夏公益教科書