基于八叉樹的網格生成算法剖析

前言

　　對于網格生成這個主題，之前的網格生成系列的三篇部落格文章分别介紹了MC算法，SMC算法以及Cuberille算法三種方法。同時還有一篇介紹網格生成與種子點生長算法高效結合的算法。本篇文章繼續這一主題，介紹采用八叉樹結構來進行網格生成的算法，這個算法并不是獨立于之前介紹的算法，而是基于了SMC算法，同時也采納了一些MC算法範疇内的思想。換句話說，就是使用八叉樹對上述網格生成進行改良，進而進一步減少網格的規模。

研究動機與SMC算法的網格規模

　　在SMC算法那一篇中已經指出，SMC算法是相對于MC算法的簡化改良，它的一大特點就是生成網格的規模相對于MC算法大大減少。但是SMC算法也畢竟是從體元中産生三角片，三角片的大小是不會超過一個體元的範圍的。而三維圖像往往會有很大的尺寸這就意味着體元的數量會非常大。因而即便不采用MC算法而采用SMC算法，仍然會産生大量密集的小三角片。比如下圖是Engine資料的三維圖像模型，可以看出三角網格是由很多小三角片構成的。

基于八叉樹的網格生成算法剖析

　　進一步觀察這樣的模型，我們發現：模型中實際上存在大量的共面三角片。了解網格型的話就會意識到這不是一個好的網格模型，因為網格模型通過小三角形片去拟合模型表面的曲面。而像上圖那樣的平面，完全可以使用更少的三角片去表示。比如一個具有下圖左方結構的Mesh用8三角片個表示一個正方形，完全等價于右圖用2個三角片表示的正方形。

基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
8個三角片組成一個正方形	2個也能組個一樣的

　　是以這就引入了進一步減少網格規模的主題，這也是本文使用八叉樹來簡化網格生成的研究動機。無論是SMC算法還是MC算法，其生成的網格都有存在大量共面三角形的可能，而這些數量衆多的小三角形也有合并成較少的大三角形的可能。尤其是SMC算法，由于其小三角形的法向和形狀具有良好的可枚舉性，因而在局部非常容易出現能夠合并的小三角形。

　　為此還必須指出，本文即将讨論的方法與針對Mesh進行網格削減的一系列算法（如之前的文章提到的頂點聚簇）并不是同一類算法。後者旨在脫離圖像的概念，僅從幾何意義上去對Mesh進行削減。而本文的算法是在網格生成之前，通過八叉樹結構進行網格生成的優化，進而直接在Mesh被構造時就減少網格的規模。是以其本質上還是算作Mesh Generation算法而不是Mesh Processing算法。

空間劃分樹的分類-均等樹和BON樹

　　熟悉計算機的都知道樹這個資料結構，尤其是二叉樹。在計算機科學中，樹形結構經常被用來組織和檢索資料。當然樹的用處不止是這些。在計算機圖形學中，經常使用樹來對空間進行劃分，也就是把空間組織成一個層次結構。這樣樹的每一個節點就被賦予了空間上的意義，可以用來代表一個區域，而樹的父子關系也正好能表示區域之間的包含關系。一般的，一維空間可以使用二叉樹來劃分，而平面和空間區域分别使用四叉樹和八叉樹來進行劃分。

　　對于一維的資料，比如一個下面的範圍[0-12]，可以使用二叉樹劃分。最終使得範圍内的每一個子元素都在二叉樹的節點上。這樣的劃分方式也比較常見，簡單的對範圍空間除以2即可。是以擴充到平面和空間就是對各個軸坐标範圍不停的除以2，這種方式與快速排序中遞歸樹的劃分方式有一點相似（注意不是完全一樣的），我們把這種方式劃分的樹叫做均等樹。比如下圖展示了對0-12範圍的資料進行均等劃分到底所建立的樹。均等樹對于[A,B]這樣的範圍，首先找到(A+B)/2，然後就将區域分成[A,(A+B)/2]，[(A+B)/2,B]。之後再繼續按相同的方式劃分到不能再分為止。

基于八叉樹的網格生成算法剖析

均等樹的劃分方式

　　BON樹是一種不大同于上述均等樹劃分方式的樹。從上圖可以看出均等樹的算法在除2的時候會遇到奇數偶數的問題。實際上在遇到奇數個元素的範圍内，均等樹所劃分的子樹實際上是不均等的，按上文說的辦法，遇到奇數個元素，就會把中位數分在靠前的範圍中，而靠後的範圍實際上會少一個元素。這樣劃分到最後會出現單枝葉子節點的情況。而BON樹為了保證樹形，會在一開始就把最初的範圍變成2的幂，進而能讓劃分的範圍一直處于偶數個元素。例如下圖展示使用BON樹的分法來劃分0-12的範圍，首先需要找到包含這個範圍的最小的2的幂，即16，然後再對這個範圍進行劃分，超出範圍的節點不再建立。

基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
BON樹實際能表示到最近的2的幂	BON樹超出部分不建立分支

　　BON樹的劃分邏輯相比于均等樹，并不更加複雜，而且相比于均等樹存在一個潛在的好處。例如上圖中對0-12範圍所建立的均等樹，可以從節點的編号上就能夠獲得節點的位置資訊。例如9這個節點，9的二進制位為1001，正好對應這這個9的節點被BON樹所分的層次資訊。從圖上0對應左子樹，1對應右子樹，這樣從根節點到1001這個數就正好是右→左→左→右。這樣同時提供了樹節點插入的一個思路：在表示0-M範圍的BON樹中，若要插入A，則可以跟據A的二進制位和樹的層數來确定這個A在樹中的位置。同時，BON樹的非葉子節點都能指代一個範圍。例如9上面的父節點，就能代表100X（X為0或1）這個由兩個數組成的範圍，而這個父節點的父親又能指代10XX的範圍。

基于八叉樹的網格生成算法剖析

　　以上是對兩種具有不同空間劃分方式的樹的介紹，本文所要繼續介紹的八叉樹劃分方法，是基于BON樹的。同時在下文中會詳細說明這樣的樹結構究竟是如何與具體的網格生成算法結合的。

樹結構與資訊的概括

　　在介紹具體算法之前，還需要繼續對空間劃分樹的一大特點進行介紹，而對這個特點的利用也是本文算法的基本思想。這個特點就是樹結構具有對局部資訊的概括能力。

　　了解算法的人都比較清楚，像樹這種具有層次的資料結構，很好的實作了對資訊的概括。比如二叉排序樹組織資料，相比與使用數組來說，樹的每一個節點都概括了自己子樹的資訊，就是：左子樹下面的都比自己小，右子樹下面的都比自己大，這就是一種概括。那麼當有資料來檢索的時候，與一個節點的比較，就能判斷出數組通路一個元素更多的大小關系。這就是層次結構之是以能實作高效檢索的原因。而數組的每一個位置不像樹那樣具有概括資訊的能力，通路每個位置都隻能知道這一個位置的元素的資訊。

　　而對于本文說的空間劃分樹，可以通過下面的例子來說明其是如何概括資訊的。例如下圖中展示了一個編号為0-12的範圍，這個範圍分布有紅綠兩種顔色。也就是有的編号處是紅色，有的是綠色。

基于八叉樹的網格生成算法剖析

　　使用空間劃分樹來對這個範圍組織一個BON樹結果如下：

基于八叉樹的網格生成算法剖析

　　假如我們想用最少的節點來表示這樣一個紅綠數組，那麼就可以對這個樹執行一個shrink操作。Shrink代表收縮是expand的反義詞，類似于我們在treeview菜單上HEAD的那個節點上點選了”-”的動作一樣。

基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
Expand後的樹	Shrink後的樹

　　這樣每一個節點都會概括自己的子節點，當發現他們的顔色一緻的時候，就會把這些資訊概括到自己身上。下圖展示了這種自底向上概括的過程，最後形成的樹充分概括了這段紅綠數組的資訊。

基于八叉樹的網格生成算法剖析	原始BON樹
基于八叉樹的網格生成算法剖析	shrink倒數第一層
基于八叉樹的網格生成算法剖析	shrink倒數第二層
基于八叉樹的網格生成算法剖析	shrink倒數第三層，不能再shrink完畢

　　當對最後shrink完畢的樹檢索到最右的節點時，在這個節點上能夠得知他所表示的範圍即8-12都是綠色。

　　了解了樹的自底向上的shrink操作，就會在大方向上清楚了本文的算法将如何減少三角形的面片數。如下圖所示，如果左圖中這樣八個體元中的三角形都被抽取出來，這樣形成的一個大三角形平面是由4個三角形組成的，而其實如果能夠将這些體元合并後再進行抽取。那麼抽取出的三角形就隻一個。這樣實作了用更少的三角形表示相同的平面區域。這也就是本文所述的算法所要達到的直接目的。

基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
相鄰的八個體元中的三角片	體元合并在一起後再抽取，隻有一個三角片了

基于八叉樹的網格生成

　　在介紹MC算法和SMC算法的時候我們已經知道，三角形片都是從邊界體元中抽取出來的。三維圖像中一般都有大量的實體元和空體元，這些體元對三角形抽取是無用的，隻有邊界體元是我們需要的。MC和SMC算法執行過程中，對三維圖像的三個軸的三重循環掃描過程實際上就是一個搜尋邊界體元的過程，找出所有體元配置不為0和255的體元，就是邊界體元，每找到這樣的體元，就抽取其中的三角形。在種子點生長算法與網格生成算法結合的那文章裡面，也是這樣的方式，無非就是尋找邊界體元的方式變成基于種子點生長的方式而不是全部掃描。是以我們可以把使用MC、SMC算法生成網格的步驟總結成三步：

尋找邊界體元；
抽取體元三角片；
組合三角片成Mesh；

　　而本文的基于八叉樹的網格算法，通過樹的shrink來減少需要抽取三角片的體元數量，這樣在上面步驟的基礎上有如下的改變：

尋找邊界體元；
将邊界體元插入八叉樹；
對八叉樹執行shrink；
抽取shrink後的各體元三角片；
組合三角片成Mesh；

　　從上文的步驟可以想到，shrink之前樹中存着都是同一層的等大小的體元，而shrink之後，就一部分節點被父節點吸收概括成一個超體元，那麼這棵樹裡的體元就會由各種大小不同的體元構成，既有邊長為1的機關體元，也有邊長為2的幂如2、4、8等的超體元。

　　上文的步驟是粗略的一個步驟，其中最為關鍵的是中間三步，每一步都還有具體需要關心的問題：第二步将邊界體元插入八叉樹，是如何的插法；第三步對樹進行shrink，那麼符合什麼樣條件的節點可以shrink，shrink到什麼程度為止；第四步中提取體元的三角片，如果是機關體元，提取的方式和SMC算法一樣，但如果是超體元将如何提取三角片。解決好了上述問題，就完整實作了本文所介紹的算法。

八叉樹的建立以及邊界體元的插入

　　首先需要解決建樹以及插入節點的問題。本算法的實作使用的是BON八叉樹，空間樹的建立必然關聯着一個空間範圍，邊界體元集合是有範圍的，而且至多也不會超出三維圖像的範圍。為了建樹，首先定義樹的節點類型OctreeNode，聲明其為模版類型友善于攜帶不同的參數，更具一般性。

基于八叉樹的網格生成算法剖析

public class OctreeNode<T>
{
    public OctreeNode<T>[] Children;//孩子指針,數組大小為8
    public OctreeNode<T> Parent;//父節點指針
    public T Parms;//攜帶的參數
    public int XMin;//所代表範圍的X軸下界
    public int YMin;//所代表範圍的Y軸下界
    public int ZMin;//所代表範圍的Z軸下界
    public int XMax;//所代表範圍的X軸下界
    public int YMax;//所代表範圍的Y軸下界
    public int ZMax;//所代表範圍的Z軸下界
    public int IndexInParent;//自己在父節點孩子數組中的索引
    public int LayerIndex;//自己所在的層索引
    public bool IsLeaf()
    {
        return (XMin == XMax)&&(YMin==YMax)&&(ZMin==ZMax);
    }//傳回是否是葉子節點
    public OctreeNode()
    {
    }
    public override string ToString()
    {
        if (IsLeaf())
        {
            return string.Format("[{0},{1}][{2},{3}][{4},{5}] {6} Leaf", XMin, XMax, YMin, YMax, ZMin, ZMax,Parms);
        }
        if (Parms == null)
            return string.Format("[{0},{1}][{2},{3}][{4},{5}] {6}", XMin, XMax, YMin, YMax, ZMin, ZMax, "not simple"); ;
        return string.Format("[{0},{1}][{2},{3}][{4},{5}] {6}", XMin, XMax, YMin, YMax, ZMin, ZMax,Parms);
    }
}//BON八叉樹節點

基于八叉樹的網格生成算法剖析

　　上述定義中為了友善附加了很多資訊如節點的各軸範圍，節點的層次以及在父節點的索引等，這些資訊有的是可以動态擷取的，開辟空間記錄下來是用空間換時間，本文為友善實作，附加了比較多的這些輔助資訊。

　　那麼下一步就是定義一顆BON樹RegionOctree：

基于八叉樹的網格生成算法剖析

public class RegionOctree<T>
{
    private static int GetMax2Power(int xmax,int ymax,int zmax,ref int log)
    {
        int max = xmax;
        if (ymax > max)
            max = ymax;
        if (zmax > max)
            max = zmax;
        if ((max & (max - 1)) == 0)
        {
            double L = Math.Log(max, 2);
            log = (int)L + 1;
            return max;
        }
        else
        {
            double L = Math.Log(max, 2);
            log = (int)L + 2;
            return (int)Math.Pow(2, log - 1);
        }
    }
    private int Width;//樹所關聯空間範圍的X上界
    private int Height;//樹所關聯空間範圍的Y上界
    private int Depth;//樹所關聯空間範圍的Z上界
    public OctreeNode<T> Root;//樹根節點
    public  int NodeCount;//所有節點總數
    public  int LeafCount;//葉子節點
    private int Scale;//2的幂包圍盒邊長
    private int LayerNum;//層次數
    private OctreeNode<T>[] NodeLayers;//指代一條由根通往葉子的路徑

    public RegionOctree(int width,int height,int depth)//使用範圍構造BON樹
    {
        this.Width = width;
        this.Height = height;
        this.Depth = depth;
        Scale = GetMax2Power(Width,Height,Depth,ref LayerNum);
        NodeCount = 0;
        Root = new OctreeNode<T>();
        Root.XMin = 0;
        Root.XMax = Scale-1;
        Root.YMin = 0;
        Root.YMax = Scale-1;
        Root.ZMin = 0;
        Root.ZMax = Scale-1;
        Root.Parent = null;
        Root.IndexInParent = -1;
        Root.LayerIndex = LayerNum - 1;
        Root.Children = new OctreeNode<T>[8];
        NodeLayers = new OctreeNode<T>[LayerNum];
        NodeLayers[0] = Root;
    }

    public OctreeNode<T> CreateToLeafNode(int x,int y,int z)
    {
        LeafCount++;
        for (int i = 1; i <= LayerNum - 1; i++)
        {
            int index = GetIndexOn(x, y, z, LayerNum - i-1);
            if (NodeLayers[i - 1].Children[index] == null)
            {
                OctreeNode<T> node = new OctreeNode<T>();
                NodeCount++;
                node.Parent = NodeLayers[i - 1];
                node.IndexInParent = index;
                node.Children = new OctreeNode<T>[8];
                node.LayerIndex = NodeLayers[i - 1].LayerIndex - 1;
                InitRangeByParentAndIndex(node, NodeLayers[i - 1], index);
                NodeLayers[i - 1].Children[index] = node;
            }
            NodeLayers[i]=NodeLayers[i-1].Children[index];
        }
        return NodeLayers[NodeLayers.Length - 1];
    }//将關聯着坐标（x，y，z）處元素一路插入到底層為葉子節點

    private int GetIndexOn(int x, int y, int z, int bitindex)
    {
        int ret = 0;
        if ((x & (1 << bitindex)) != 0)
        {
            ret |= 1;
        }
        if ((y & (1 << bitindex)) != 0)
        {
            ret |= 2;
        }
        if ((z & (1 << bitindex)) != 0)
        {
            ret |= 4;
        }
        return ret;
    }

    private void InitRangeByParentAndIndex(OctreeNode<T> node,OctreeNode<T> pnode, int index)
    {
        int deltaX = (pnode.XMax - pnode.XMin + 1) / 2;
        int deltaY = (pnode.YMax - pnode.YMin + 1) / 2;
        int deltaZ = (pnode.ZMax -  pnode.ZMin + 1) / 2;
        if ((index & 1) == 0)
        {
            node.XMin = pnode.XMin;
            node.XMax = pnode.XMin + deltaX - 1;
        }
        else
        {
            node.XMin = pnode.XMin + deltaX;
            node.XMax = pnode.XMax;
        }
        if ((index & 2) == 0)
        {
            node.YMin = pnode.YMin;
            node.YMax = pnode.YMin + deltaY - 1;
        }
        else
        {
            node.YMin = pnode.YMin + deltaY;
            node.YMax = pnode.YMax;
        }
        if ((index & 4) == 0)
        {
            node.ZMin = pnode.ZMin;
            node.ZMax = pnode.ZMin + deltaZ - 1;
        }
        else
        {
            node.ZMin = pnode.ZMin + deltaZ;
            node.ZMax = pnode.ZMax;
        }
    }//使用父節點的資訊初始化子節點的範圍
}

基于八叉樹的網格生成算法剖析

　　從上述代碼可以看出，對三維圖像範圍(0-width,0-height,0-depth)關聯一顆BON樹，那麼這顆BON樹就首先需要明确自己的層數，也就是需要找到最小的2的幂的包圍盒。GetMax2Power函數實作了這一功能，例如一個具有範圍(200,300,400)的空間，尋找到的最小2的幂包圍盒就是(512,512,512)。

　　在了解插入方法之前還需要對子節點的順序進行固定的編号，在這裡采用的方式是按0-7的二進制位000-111進行配置設定：從右數第一位1表示X軸正方向；0表示X軸負方向。從右數第二位1表示Y軸正方向；0表示Y軸負方向。從右數第三位1表示Z軸正方向；0表示Z軸負方向。是以8個子體元的編号位置對應關系就如下表所示。

基于八叉樹的網格生成算法剖析

X範圍：（XMIN~XMAX）

Y範圍：（YMIN~YMAX）

Z範圍：（ZMIN~ZMAX）

設各方向範圍中點分别為：XMID,YMID,ZMID

孩子編号	對應範圍
（XMIN-XMID，YMIN-YMID，ZMIN-ZMID）
1	（XMID-XMAX，YMIN-YMID，ZMIN-ZMID）
2	（XMIN-XMID，YMID-YMAX，ZMIN-ZMID）
3	（XMID-XMAX，YMID-YMAX，ZMIN-ZMID）
4	（XMIN-XMID，YMIN-YMID，ZMID-ZMAX）
5	（XMID-XMAX，YMIN-YMID，ZMID-ZMAX）
6	（XMIN-XMID，YMID-YMAX，ZMID-ZMAX）
7	（XMID-XMAX，YMID-YMAX，ZMID-ZMAX）

圖示

父節點資訊

孩子資訊

　　為了了解與插入相關的三個函數CreateToLeafNode、GetIndexOn與InitRangeByParentAndIndex。下面使用一組資料來舉例說明一個具有坐标(80,85,22)的元素是如何插入一個範圍為(200,200,210)的BON樹的。首先根據這三個軸坐标的最大值210，找到對應比它大的最小的2的幂為256，是以這個BON樹所能表示的的最大範圍是(256，256，256)而層數（或者說深度）為8（256為2的8次方）。

　　确定了BON樹的層數為8後，将這(80,85,22)的三個坐标的二進制位像下圖一樣填在3個長度為8（等于層數）的表裡，我們就可以從這個表中得出插入這個節點的方法：

基于八叉樹的網格生成算法剖析

　　圖中的層内索引是每一豎列的二進制位組合起來的結果，也是X行的值*4與Y行的值*2與Z行的值*1相加的結果。這個層内索引表達了這個節點在每一層應當被插入到哪一顆子樹。從上圖的結果看，層内索引分别為第0、第3、第0、第7、第0、第6、第2。這就意味着将(80,85,22)插入這課BON樹的步驟為：

将[80,85,22]插入0層第0個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入1層第3個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入2層第0個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入3層第7個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入4層第0個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入5層第6個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入6層第4個孩子節點，若該孩子節點為NULL則建立之。
将[80,85,22]插入7層第2個孩子節點，若該孩子節點為NULL則建立之。

　　CreateToLeafNode函數主邏輯便是這樣一個過程，這樣的過程結束後，私有成員NodeLayers内的資訊包含了所有經過的節點。GetIndexOn函樹通過檢查位來找到對應層的層内索引。InitRangeByParentAndIndex函數則負責初始化每一個節點的建立資訊。例如一個父節點所表示的範圍為[X0-X1,Y0-Y1,Z0-Z1]，那麼他的子節點所表示的範圍就正好将每個軸的範圍平均劈成兩半。對應的資訊可以參見上文的表中孩子資訊中的範圍。

　　這樣八叉樹中如何插入邊界體元資訊的問題就解決了，體元集合的插入就相當于是三維坐标的依次插入，這樣插入之後每一個邊界體元都處在八叉樹的葉子節點上。當然這個八叉樹隻是可能的最大分支數為8。一般來說，會有很多節點不滿8個子節點。是以這樣的樹就叫做Branch On Need樹，簡稱為BON樹。

Shrink操作的合并原則

　　在之前介紹算法目标的時候就已經提到，本算法的目的是減少共面三角形數，也就是企圖合并那些在一個平面上的小三角片，将其合成大三角片。如下圖所示的4個例子情況，都是可以合并的情況：

合并前	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
合并後	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析

　　從上圖中可見，假如一個節點的子節點裡面的三角片全部是共平面的，這樣的節點就可以shrink。再進一步考慮，有的體元配置中隻有一個平面，而有的體元卻包含了了多個不共面的三角片，如下圖所示：

基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
不能shrink的例子之一	不能shrink的例子之二	能shrink的例子之一	能shrink的例子之二

　　那麼後者就絕對沒有被shrink的可能。而前者的體元若作為某節點的孩子，假如他的兄弟節點和他一樣具有共面的三角片，那麼這些兄弟就可以與之合并成一個大的超體元，那麼現在的問題就是如何判斷三角形的共面情況。考慮到SMC算法的三角形片形狀和方向是可以枚舉的那幾種形狀，我們可以把這些資訊結合體元配置做一個映射，進而快速判斷三角形片是否是共面的。使用空間解析幾何的知識我們可以知道：空間中的所有平面都可以使用方程Ax+By+Cz=D來表示。其中（A,B,C）的組合表征平面的方向，而D值表征位移，确定了這四個系數，平面就被唯一确定。

　　在SMC算法介紹中，我們已經知道三角形片的頂點都是體元的體素。下面以單個體元為例子來說明三角片的平面方程。下圖是筆者使用WPF技術制作的訓示SMC體元配置和三角片關系的軟體截圖，可以通過http://files.cnblogs.com/chnhideyoshi/Debug.rar下載下傳。這個軟體通過拖動滑塊來觀看不同的體元配置的三角片情況。

共面配置之一	基于八叉樹的網格生成算法剖析
共面配置之二	基于八叉樹的網格生成算法剖析
共面配置之三	基于八叉樹的網格生成算法剖析
非共面配置	基于八叉樹的網格生成算法剖析

　　我們可以看出，凡是具有多個不共面三角形的情況，都被标注為“非共面配置”；而隻有一個正三角形或者是由兩個三角形組合成長方形這樣的一類的體元配置被标注為“共面配置”。具有共面配置的體元就存在着與兄弟節點合并的可能。同時還可以發現，這些三角形的方程的形式是有限的，ABC的值隻在0、1、-1三個數中變化。D值的變化也是有限的，但由于這隻是基準點坐标為（0,0,0）的體元，一旦将這個體元坐标一般化為（Dx，Dy，Dz），那麼這個D值就會跟随這個坐标進行變化，而ABC是不變的。

　　所謂共面配置隻是描述了體元内三角形的形式，或者說三角形平面方向，即方程中的(ABC)三元組，而三角形所在平面的具體位置還依賴于D，D不同即使ABC相同，這樣的三角形也隻能是平行的三角形，而不是共面的。是以能夠合并的體元，其中的三角形不僅要具有相同的法向，也要有相同的D值。也就是說，八叉樹的節點可以記錄下所代表體元中三角形片的法向和D值，這兩個資訊關系到體元能否與兄弟節點合并。

　　是以綜合以上的陳述，總結一下對于八叉樹的一個節點P，他的孩子C0-C7能夠shrink到P必須遵循以下的原則：

C0-C7要麼為NULL（即非邊界體元區域），要麼為一種共面配置。
C0-C7不為NULL的對應體元的共面配置必須均為同一種
C0-C7的不為NULL的節點均有相同的D值
C0-C7的後代節點若不為NULL則需均不含合并失敗的節點，也就是一個節點若shrink失敗，其父節點不再shrink。

　　由于使用三元組來表達三角片方向這個資訊會有諸多備援，考慮到法向量是有限可枚舉的，是以可以為這些法向量進行編号，設定查找數組，以編号來代替三元組來進行存儲。設定的查找數組如下：

基于八叉樹的網格生成算法剖析

public struct Int32Quad
{
    public int A;
    public int B;
    public int C;
    public int D;
    public Int32Quad(int a, int b, int c, int d)
    {
        this.A = a;
        this.B = b;
        this.C = c;
        this.D = d;
    }
    public override string ToString()
    {
        return string.Format("{0}x+{1}y+{2}z={3}", A, B, C, D);
    }
}//代表一個平面方程對象
public class OctreeTable
{
    public static byte[] ConfigToNormalTypeId = new byte[256]
    {
        13,0,1,2,3,13,4,13,5,6,13,13,7,13,13,8,3,9,13,13,13,13,13,13,13,13,13,0,13,13,13,13,5,13,10,13,13,13
        ,13,1,13,13,13,13,13,13,13,13,7,13,13,11,13,13,13,13,13,13,13,13,13,13,13,2,0,13,13,13,9,13,13,13,13
        ,13,13,13,13,13,3,13,13,13,13,13,13,13,13,13,13,13,
        13,13,13,13,13,13,6,13,13,13,13,13,12,13,13,13,13,13,13,13,13,4,13,13,5,13,13,13,13,10,13,13,13,13
        ,13,13,13,1,1,13,13,13,13,13,13,13,10,13,13,13,13,5,13,
        13,4,13,13,13,13,13,13,13,13,12,13,13,13,13,13,6,13,13,13,13,13,13,
        13,13,13,13,13,13,13,13,13,13,13,3,13,13,13,13,13,13,13,13,13,9,13,13,13,0,2,13,13,13,
        13,13,13,13,13,13,13,13,11,13,13,7,13,13,13,13,13,13,13,13,1,13,13,13,13,10,13,5,13,13,
        13,13,0,13,13,13,13,13,13,13,13,13,9,3,8,13,13,7,13,13,6,5,13,4,13,3,2,1,0,13
    };//體元配置對應的法向量索引
    public static Int16Triple[] NormalTypeIdToNormal = new Int16Triple[13]
    {
        new Int16Triple(1,-1,-1),
        new Int16Triple(1,-1,1),
        new Int16Triple(1,-1,0),
        new Int16Triple(1,1,1),
        new Int16Triple(1,0,1),
        new Int16Triple(1,1,-1),
        new Int16Triple(1,0,-1),
        new Int16Triple(1,1,0),
        new Int16Triple(1,0,0),
        new Int16Triple(0,1,1),
        new Int16Triple(0,1,-1),
        new Int16Triple(0,1,0),
        new Int16Triple(0,0,1)
    };//枚舉的法向量集合

    public static byte[] ConfigToEqType = new byte[256]
    {
        55,0,1,2,3,55,4,55,5,6,55,55,7,55,
        55,8,9,10,55,55,55,55,55,55,55,
        55,55,11,55,55,55,55,12,55,13,55,
        55,55,55,14,55,55,55,55,55,55,55
        ,55,15,55,55,16,55,55,55,55,55,55
        ,55,55,55,55,55,17,18,55,55,55,19,
        55,55,55,55,55,55,55,55,55,20,55,
        55,55,55,55,55,55,55,55,55,55,55,
        55,55,55,55,0,21,55,55,55,55,55,22,
        55,55,55,55,55,55,55,55,23,55,55,24,55,
        55,55,55,25,55,55,55,0,55,0,0,26,27,55,55,
        55,55,55,55,55,28,55,55,55,55,29,55,55,
        30,55,55,55,55,55,55,55,55,31,55,55,55,
        55,55,32,55,55,55,55,55,55,55,55,55,55,
        55,55,55,55,55,0,55,33,55,55,55,55,
        55,0,55,55,55,34,55,0,0,35,36,55,55,55,
        55,55,55,55,55,55,55,55,37,55,55,38,55,
        55,55,55,55,55,55,0,39,55,55,0,55,40,0,41,55,
        55,55,55,42,55,55,0,55,55,55,0,55,0,43,44,45,55
        ,55,46,55,0,47,48,55,49,0,50,51,52,53,55
    };//體元配置對應的方程索引，55表示非共面配置
    public static Int32Quad[] EqTypeToEqQuad = new Int32Quad[54]
    {
        new Int32Quad(1,-1,-1,-1),
        new Int32Quad(1,-1,1,0),
        new Int32Quad(1,-1,0,0),
        new Int32Quad(1,1,1,1),
        new Int32Quad(1,0,1,1),
        new Int32Quad(1,1,-1,0),
        new Int32Quad(1,0,-1,0),
        new Int32Quad(1,1,0,1),
        new Int32Quad(1,0,0,1),
        new Int32Quad(1,1,1,2),
        new Int32Quad(0,1,1,1),
        new Int32Quad(1,-1,-1,0),
        new Int32Quad(1,1,-1,1),
        new Int32Quad(0,1,-1,0),
        new Int32Quad(1,-1,1,1),
        new Int32Quad(1,1,0,1),
        new Int32Quad(0,1,0,0),
        new Int32Quad(1,-1,0,1),
        new Int32Quad(1,-1,-1,0),
        new Int32Quad(0,1,1,1),
        new Int32Quad(1,1,1,2),
        new Int32Quad(1,0,-1,0),
        new Int32Quad(0,0,1,1),
        new Int32Quad(1,0,1,2),
        new Int32Quad(1,1,-1,0),
        new Int32Quad(0,1,-1,-1),
        new Int32Quad(1,-1,1,2),
        new Int32Quad(1,-1,1,1),
        new Int32Quad(0,1,-1,0),
        new Int32Quad(1,1,-1,1),
        new Int32Quad(1,0,1,1),
        new Int32Quad(0,0,1,0),
        new Int32Quad(1,0,-1,1),
        new Int32Quad(1,1,1,1),
        new Int32Quad(0,1,1,0),
        new Int32Quad(1,-1,-1,1),
        new Int32Quad(1,-1,0,0),
        new Int32Quad(0,1,0,1),
        new Int32Quad(1,1,0,2),
        new Int32Quad(1,-1,1,0),
        new Int32Quad(0,1,-1,1),
        new Int32Quad(1,1,-1,2),
        new Int32Quad(1,-1,-1,-1),
        new Int32Quad(0,1,1,2),
        new Int32Quad(1,1,1,3),
        new Int32Quad(1,0,0,0),
        new Int32Quad(1,1,0,0),
        new Int32Quad(1,0,-1,-1),
        new Int32Quad(1,1,-1,-1),
        new Int32Quad(1,0,1,0),
        new Int32Quad(1,1,1,0),
        new Int32Quad(1,-1,0,-1),
        new Int32Quad(1,-1,1,-1),
        new Int32Quad(1,-1,-1,-2),
    };//體元配置對應的三角片方程集合
    public const byte NormalNotSimple = 13;
}

基于八叉樹的網格生成算法剖析

　　這樣我們就可以為八叉樹節點Node的T類型成員指定一個類型NodeParms，這個類型的變量攜帶體元的共面資訊NormalTypeId和D，為友善起見還加上了體元配置Config，一共占據2個位元組加一個整型的空間。以便在shrink執行時提取該資訊進行判斷。

基于八叉樹的網格生成算法剖析

public class NodeParms
{
    public byte Config;
    public byte NormalTypeId;
    public int D;
    public override string ToString()
    {
        if (NormalTypeId < OctreeTable.NormalTypeIdToNormal.Length)
            return "\"" + NormalTypeId + "," + D + "\"";
        else
            return "not simple";
    }
}

基于八叉樹的網格生成算法剖析

　　這裡附上完整的OctreeSurfaceGenerator的代碼如下：

基于八叉樹的網格生成算法剖析

public static class OctreeSurfaceGenerator
{
    #region consts
    public static byte VULF = 1 << 0;
    public static byte VULB = 1 << 1;
    public static byte VLLB = 1 << 2;
    public static byte VLLF = 1 << 3;
    public static byte VURF = 1 << 4;
    public static byte VURB = 1 << 5;
    public static byte VLRB = 1 << 6;
    public static byte VLRF = 1 << 7;
    //以上為體素為實點的位标記
    public static Int16Triple[] PointIndexToPointDelta = new Int16Triple[8]
    {
        new Int16Triple(0, 1, 1 ),
        new Int16Triple(0, 1, 0 ),
        new Int16Triple(0, 0, 0 ),
        new Int16Triple(0, 0, 1 ),
        new Int16Triple(1, 1, 1 ),
        new Int16Triple(1, 1, 0 ),
        new Int16Triple(1, 0, 0 ),
        new Int16Triple(1, 0, 1 )
    };//體元内每個體素相對基準體素坐标的偏移
    public static byte[] PointIndexToFlag = new byte[8]
    {
        VULF,
        VULB,
        VLLB,
        VLLF,
        VURF,
        VURB,
        VLRB,
        VLRF
    };//每個體素對應的位标記
    #endregion
    public static Mesh GenerateSurface(BitMap3d bmp)
    {
        int width = bmp.width;
        int height = bmp.height;
        int depth = bmp.depth;
        Int16Triple[] tempArray = new Int16Triple[8];
        #region CreateTree
        RegionOctree<NodeParms> otree = new RegionOctree<NodeParms>(width, height, depth);
        Queue<OctreeNode<NodeParms>> nodequeue = new Queue<OctreeNode<NodeParms>>();
        for (int k = 0; k < depth-1; k++)
        {
            for (int j = 0; j < height-1; j++)
            {
                for (int i = 0; i < width-1; i++)
                {
                    byte value = 0;
                    for (int pi = 0; pi < 8; pi++)
                    {
                        tempArray[pi].X = i + PointIndexToPointDelta[pi].X;
                        tempArray[pi].Y = j + PointIndexToPointDelta[pi].Y;
                        tempArray[pi].Z = k + PointIndexToPointDelta[pi].Z;
                        if (InRange(bmp, tempArray[pi].X, tempArray[pi].Y, tempArray[pi].Z)&&
                            IsWhite(bmp, tempArray[pi].X, tempArray[pi].Y, tempArray[pi].Z))
                        {
                            value |= PointIndexToFlag[pi];
                        }
                    }
                    if (value != 0 && value != 255)
                    {
                        OctreeNode<NodeParms> leafnode = otree.CreateToLeafNode(i, j, k);
                        leafnode.Parms = new NodeParms();
                        leafnode.Parms.Config = value;
                        leafnode.Parms.NormalTypeId = OctreeTable.ConfigToNormalTypeId[value];
                        leafnode.Parms.D = CaculateDFromNormalAndCoord(i, j, k, value);
                        nodequeue.Enqueue(leafnode.Parent);
                    }
                }
            }
        }
        #endregion
        #region Shrink
        while (nodequeue.Count != 0)
        {
            OctreeNode<NodeParms> node = nodequeue.Dequeue();
            byte normalType = OctreeTable.NormalNotSimple;
            int D = int.MinValue;
            if (CanMergeNode(node, ref normalType, ref D))
            {
                node.Parms = new NodeParms();
                node.Parms.NormalTypeId = normalType;
                node.Parms.D = D;
                nodequeue.Enqueue(node.Parent);
            }
        }
        #endregion
        #region ExtractTriangles
        MeshBuilder_IntegerVertex mb = new MeshBuilder_IntegerVertex(width+1, height+1, depth+1);
        nodequeue.Enqueue(otree.Root);
        while (nodequeue.Count != 0)
        {
            OctreeNode<NodeParms> node = nodequeue.Dequeue();
            if (node.Parms == null)
            {
                for (int i = 0; i < 8; i++)
                {
                    if (node.Children[i] != null)
                        nodequeue.Enqueue(node.Children[i]);
                }
            }
            else
            {
                if (node.Parms.NormalTypeId != OctreeTable.NormalNotSimple)
                {
                    if (node.IsLeaf())
                    {
                        GenerateFaceLeaf(node, mb, ref tempArray, bmp);
                    }
                    else
                    {
                        GenerateFace(node, mb, ref tempArray, bmp);
                    }
                }
                else
                {
                    if (node.IsLeaf())
                    {
                        GenerateFaceLeaf(node, mb, ref tempArray, bmp);
                    }
                    else
                    {
                        for (int i = 0; i < 8; i++)
                        {
                            if (node.Children[i] != null)
                                nodequeue.Enqueue(node.Children[i]);
                        }
                    }
                }
            }
        }//采用層次周遊尋找需要抽取三角片的節點
        #endregion
        return mb.GetMesh();
    }

    private static bool InRange(BitMap3d bmp, int x, int y, int z)
    {
        return x > 0 && x < bmp.width && y > 0 && y < bmp.height && z > 0 && z < bmp.depth;
    }
    private static bool IsWhite(BitMap3d bmp, int x, int y, int z)
    {
        return bmp.GetPixel(x, y, z) == BitMap3d.WHITE;
    }

    private static void GenerateFace(OctreeNode<NodeParms> node, MeshBuilder_IntegerVertex mb, ref Int16Triple[] tempArray, BitMap3d bmp)
    {
        InitVoxelPositionForNodeRange(node.XMin, node.XMax, node.YMin, node.YMax, node.ZMin, node.ZMax, ref tempArray);
        //需要找到該節點的端點位置
        byte cubeConfig = 0;
        for (int pi = 0; pi < 8; pi++)
        {
            if (InRange(bmp, tempArray[pi].X, tempArray[pi].Y, tempArray[pi].Z)
                && IsWhite(bmp,tempArray[pi].X,tempArray[pi].Y ,tempArray[pi].Z))
            {
                cubeConfig |= PointIndexToFlag[pi];
            }
        }
        int index = 0;
        while (MCTable.TriTable[cubeConfig, index] != -1)
        {
            int ei1 = MCTable.TriTable[cubeConfig, index];
            int ei2 = MCTable.TriTable[cubeConfig, index + 1];
            int ei3 = MCTable.TriTable[cubeConfig, index + 2];

            Int16Triple p1 = GetIntersetedPointAtEdge(node, ei1, OctreeTable.NormalTypeIdToNormal[node.Parms.NormalTypeId], node.Parms.D);
            Int16Triple p2 = GetIntersetedPointAtEdge(node, ei2, OctreeTable.NormalTypeIdToNormal[node.Parms.NormalTypeId], node.Parms.D);
            Int16Triple p3 = GetIntersetedPointAtEdge(node, ei3, OctreeTable.NormalTypeIdToNormal[node.Parms.NormalTypeId], node.Parms.D);

            mb.AddTriangle(p1, p2, p3);
            index += 3;
        }
    }//對非葉子節點的超體元的抽取需要參考MCTable求取被截斷邊的資訊

    private static void GenerateFaceLeaf(OctreeNode<NodeParms> node, MeshBuilder_IntegerVertex mb, ref Int16Triple[] tempArray, BitMap3d bmp)
    {
        for (int k = 0; k < 8; k++)
        {
            tempArray[k].X = node.XMin + PointIndexToPointDelta[k].X;
            tempArray[k].Y = node.YMin + PointIndexToPointDelta[k].Y;
            tempArray[k].Z = node.ZMin + PointIndexToPointDelta[k].Z;
        }
        byte value = node.Parms.Config;
        int index = 0;
        while (SMCTable.TableFat[value, index] != -1)
        {
            Int16Triple t0 = tempArray[SMCTable.TableFat[value, index]];
            Int16Triple t1 = tempArray[SMCTable.TableFat[value, index + 1]];
            Int16Triple t2 = tempArray[SMCTable.TableFat[value, index + 2]];
            mb.AddTriangle(t0, t1,t2);
            index += 3;
        }
    }//對葉子節點的機關體元的抽取和SMC算法中的抽取一緻

    private static bool CanMergeNode(OctreeNode<NodeParms> node, ref byte normalType, ref int D)
    {
        for (int i = 0; i < 8; i++)
        {
            if (node.Children[i] != null)
            {
                if (node.Children[i].Parms == null)
                {
                    return false;//說明其下存在不能合并的節點 合并失敗
                }
                else
                {
                    if (node.Children[i].Parms.NormalTypeId != OctreeTable.NormalNotSimple)
                    {
                        normalType = node.Children[i].Parms.NormalTypeId;
                        D = node.Children[i].Parms.D;//記錄其中的共面配置資訊
                    }
                    else
                    {
                        return false;//遇到了非共面配置 合并失敗
                    }
                }
            }
        }
        for (int i = 0; i < 8; i++)
        {
            if (node.Children[i] != null)
            {
                if (node.Children[i].Parms.NormalTypeId != normalType || node.Children[i].Parms.D != D)
                {
                    return false;//體元配置均為共面類型但不是同一種的話也不中
                }
            }
        }
        return true;
    }

    private static Int16Triple GetIntersetedPointAtEdge(OctreeNode<NodeParms> node, int edgeIndex, Int16Triple normal, int d)
    {
        int x = 0, y = 0, z = 0;
        switch (edgeIndex)
        {
            case 0: { x = node.XMin; y = node.YMax + 1; return new Int16Triple(x, y, (d - normal.X * x - normal.Y * y) / normal.Z); }
            case 2: { x = node.XMin; y = node.YMin; return new Int16Triple(x, y, (d - normal.X * x - normal.Y * y) / normal.Z); }
            case 4: { x = node.XMax + 1; y = node.YMax + 1; return new Int16Triple(x, y, (d - normal.X * x - normal.Y * y) / normal.Z); }
            case 6: { x = node.XMax + 1; y = node.YMin; return new Int16Triple(x, y, (d - normal.X * x - normal.Y * y) / normal.Z); }

            case 8: { y = node.YMax + 1; z = node.ZMax + 1; return new Int16Triple((d - normal.Y * y - normal.Z * z) / normal.X, y, z); }
            case 9: { y = node.YMax + 1; z = node.ZMin; return new Int16Triple((d - normal.Y * y - normal.Z * z) / normal.X, y, z); }
            case 10: { y = node.YMin; z = node.ZMin; return new Int16Triple((d - normal.Y * y - normal.Z * z) / normal.X, y, z); }
            case 11: { y = node.YMin; z = node.ZMax + 1; return new Int16Triple((d - normal.Y * y - normal.Z * z) / normal.X, y, z); }

            case 1: { x = node.XMin; z = node.ZMin; return new Int16Triple(x, (d - normal.X * x - normal.Z * z) / normal.Y, z); }
            case 3: { x = node.XMin; z = node.ZMax + 1; return new Int16Triple(x, (d - normal.X * x - normal.Z * z) / normal.Y, z); }
            case 5: { x = node.XMax + 1; z = node.ZMin; return new Int16Triple(x, (d - normal.X * x - normal.Z * z) / normal.Y, z); }
            case 7: { x = node.XMax + 1; z = node.ZMax + 1; return new Int16Triple(x, (d - normal.X * x - normal.Z * z) / normal.Y, z); }

            default: throw new Exception();
        }
    }

    private static int CaculateDFromNormalAndCoord(int cx, int cy, int cz, byte config)
    {
        byte index = OctreeTable.ConfigToEqType[config];
        if (index >= OctreeTable.EqTypeToEqQuad.Length)
            return int.MinValue;
        Int32Quad eq = OctreeTable.EqTypeToEqQuad[index];
        return eq.D + eq.A * cx + eq.B * cy + eq.C * cz;
    }

    private static void InitVoxelPositionForNodeRange(int xmin, int xmax, int ymin, int ymax, int zmin, int zmax, ref Int16Triple[] temp)
    {
        temp[0].X = xmin;  //(0, 1, 1, VULF);
        temp[0].Y = ymax + 1;
        temp[0].Z = zmax + 1;

        temp[1].X = xmin;//(0, 1, 0, VULB);
        temp[1].Y = ymax + 1;
        temp[1].Z = zmin;

        temp[2].X = xmin;//(0, 0, 0, VLLB);
        temp[2].Y = ymin;
        temp[2].Z = zmin;

        temp[3].X = xmin;//(0, 0, 1, VLLF);
        temp[3].Y = ymin;
        temp[3].Z = zmax + 1;

        temp[4].X = xmax + 1;//(1, 1, 1, VURF);
        temp[4].Y = ymax + 1;
        temp[4].Z = zmax + 1;

        temp[5].X = xmax + 1; //(1, 1, 0, VURB);
        temp[5].Y = ymax + 1;
        temp[5].Z = zmin;

        temp[6].X = xmax + 1;//(1, 0, 0, VLRB);
        temp[6].Y = ymin;
        temp[6].Z = zmin;

        temp[7].X = xmax + 1;  //(1, 0, 1, VLRF);
        temp[7].Y = ymin;
        temp[7].Z = zmax + 1;
    }

    public static void Test1()
    {
        BitMap3d image = BitMap3d.CreateSampleTedVolume(100);
        Mesh m=GenerateSurface(image);
        PlyManager.Output(m, "D://VTKproj//Tree_Test6_T1.ply");
    }

    public static void Test2()
    {
        BitMap3d image = BitMap3d.CreateSampleEngineVolume();
        Mesh m = GenerateSurface(image);
        PlyManager.Output(m, "D://VTKproj//Tree_Test6_T2.ply");
    }

    public static void Test()
    {
        Test2();
    }
}

基于八叉樹的網格生成算法剖析

　　可以看出在主函數的函數體中的三大步邏輯：

建立樹
Shrink
三角形提取

　　在函數CanMergeNode中充分展現出了本段所述的合并原則。

體元三角形的抽取

　　雖然已經貼出了完整的代碼，但是最後還剩下一個需要解答的問題，就是如何抽取體元三角形。首先在SMC算法中我們知道如果是機關體元，直接根據體元配置查找三角形表，這樣就可以找到三角形的頂點和連接配接方式。三角形的頂點一定是體元的8個體素點之一。而在本文介紹的算法中，涉及到了對非機關體元、也就是超體元的三角片提取。不難想到，八叉樹的葉子節點存儲的都是機關體元的坐标，在shrink之後若其沒有被上一級吸收則仍然為機關體元，那麼這樣的節點抽取三角形仍然采用SMC算法的辦法；但是對于超體元中三角片的提取，則有必要描述一下這之中的技巧。

　　對于超體元首先必須明确這樣一個事實：雖然超體元中的三角形片的頂點也一定是體素點，但這些頂點不見得是超體元的8個端點，而是可能是其下層節點體元的端點。如下圖所示的節點為葉子節點的父親，包含了8個機關體元。可以看出，合并成一個超體元後，三角片穿過這個超體元的邊，而不一定經過這個超體元的端點。

基于八叉樹的網格生成算法剖析

　　是以這裡在提取三角形的時候，就不再是隻考慮使用SMC三角形表，因為SMC三角形表隻能提取機關體元的三角形。由于MC三角形表涉及到了對三角形穿過體元邊的資訊，是以這裡對超體元提取三角形還得用上MC的三角形表。

　　在提取三角形之前還需要考慮超體元配置的求法，由于超體元本質上表示了一個立方體範圍的體元坐标，是以其中包含了多個體元的體素，我們需要考察超體元端點的體素來決定這個超體元的配置。下面的表通過對比來顯示超體元和機關體元求取體元配置的不同

體素編号	體素名稱	體素位置
VULF	（x，y+1，z+1）
1	VULB	（x，y+1，z）
2	VLLB	（x，y，z）
3	VLLF	（x，y，z+1）
4	VURF	（x+1，y+1，z+1）
5	VURB	（x+1，y+1，z）
6	VLRB	（x+1，y，z）
7	VLRF	（x+1，y，z+1）

體素編号	體素名稱	體素位置
VULF	（xmin，ymax+1，zmax+1）
1	VULB	（xmin，ymax+1，zmin）
2	VLLB	（xmin，ymin，zmin）
3	VLLF	（xmin，ymin，zmax+1）
4	VURF	（xmax+1，ymax+1，zmax+1）
5	VURB	（xmax+1，ymax+1，zmin）
6	VLRB	（xmax+1，ymin，zmin）
7	VLRF	（xmax+1，ymin，zmax+1）

機關體元，坐标為（x,y,z）

超體元，指代範圍（xmin-xmax,ymin-ymax,zmin-zmax）

　　在求取超體元的插值點的時候，根據邊的方向和三角形的方程，可以枚舉出每條邊上插值點如下表所示：

基于八叉樹的網格生成算法剖析

邊代号	邊兩端頂點索引	邊方程	與三角片方程AX+BY+CZ=D的交點位置
e0	（0，1）	x=xmin,y=ymax+1	（xmin,ymax+1,(D-Axmin-B(ymax+1))/C）
e1	（1，2）	x=xmin,z=zmin	（xmin,(D-Axmin-Czmin)/B,zmin）
e2	（2，3）	x=xmin,y=ymin	（xmin,ymin,(D-Axmin-Bymin)/C）
e3	（3，0）	x=xmin,z=zmax+1	（xmin,(D-Axmin-C(zmax+1))/B,zmax+1）
e4	（4，5）	x=xmax+1,y=ymax+1	（xmax+1,ymax+1,(D-A(xmax+1)-B(ymax+1))/C）
e5	（5，6）	x=xmax+1,z=zmin	（xmax+1,(D-A(xmax+1)-Czmin)/B,zmin）
e6	（6，7）	x=xmax+1,y=ymin	（xmax+1,ymin,(D-A(xmax+1)-Bymin)/C）
e7	（7，4）	xmax+1,z=zmax+1	（xmax+1,(D-A(xmax+1)-C(zmax+1))/B,zmax+1）
e8	（0，4）	y=ymax+1,z=zmax+1	（(D-B(ymax+1)-C(zmax+1))/A,ymax+1,zmax+1）
e9	（1，5）	y=ymax+1,z=zmin	（(D-B(ymax+1)-Czmin)/A,ymax+1,zmin）
e10	（2，6）	y=ymin,z=zmin	（(D-Bymin-Czmin)/A,ymin,zmin）
e11	（3，7）	y=ymin,z=zmax+1	（(D-Bymin-C(zmax+1))/A,ymin,zmax+1）

圖示超體元，指代範圍（xmin-xmax,ymin-ymax,zmin-zmax）

邊資訊

　　這樣就解釋了第三個問題，同時也說明了RegionOctree類的其他函數的含義。

算法實作與測試

　　除了上文貼出的代碼之外，還涉及到哈希表類、Mesh類等在之前的博文中都有涉及到，這裡不再重複貼出。

　　算法結果的測試采用兩組資料，一組是TED三維圖像，其中的内容部分表現為一個三棱錐。還有一組是Engine資料。算法測試時會與SMC算法的結果進行對比。

　　首先是算法的局部結果對比：

TED資料	SMC結果預覽	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
八叉樹算法結果預覽	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析

Engine資料	SMC結果預覽	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析
八叉樹算法結果預覽	基于八叉樹的網格生成算法剖析	基于八叉樹的網格生成算法剖析

　　算法結果參數對比如下表所示：

SMC算法頂點數	八叉樹算法頂點數	SMC算法面片數	八叉樹算法面片數
TED資料	16471	3123	32938	5147
Engine資料	216147	144117	432370	259479

　　可見八叉樹算法針對有大片平面區域的模型，能夠一定程度上減少網格的規模。算法工程的下載下傳位址為：https://github.com/chnhideyoshi/SeededGrow2d

基于八叉樹的網格生成算法剖析

繼續閱讀

GitHub開源！《深度學習面試指南》

（三十六：2021.01.13）Pre-MICCAI 2019學習（一）《從植入粒子後的CT圖像自動定位植入的粒子》講在前面摘要論文内容

Jetson TX2 opencv3.4的編譯安裝

二維圖像三維化錄影機标定視覺計算理論（Computational Theory of Vision）

Photoshop最常見的問答集

Visual C++多媒體設計及圖形、圖像處理

DLP LightCrafter4500投影圖案步驟整理(四)------常見投影失敗問題：Pattern Number OOR

3d相機機器視覺檢測技術

使用自相關函數的紋理圖像周期計算

2014ECCV視覺跟蹤錄用文章

opencv視覺跟蹤——消除背景模組化

JPEG壓縮技術JPEG壓縮的步驟

JPEG壓縮技術的原理

Python3 色情圖檔識别 Python 3 色情圖檔識别

BMP檔案結構及圖像每行位元組計算方法

如何用Opencv求圖像的灰階投影曲線