您的位置：首页 > 其它

种子点生长算法下——三维种子点生长

2016-05-17 01:00 246 查看
上一篇文章“二维种子点生长研究”介绍了二维种子点生长，对于平面图像，二维种子点生长其实很容易理解，其典型应用是“油漆桶”。例如在一个具有多片联通区域的二维图上，我们能用油漆桶寻找所需要的单片联通区域。



　　至于三维种子点生长，其目的与二维种子点生长一样，都是要找出图像中的联通区域。三维图像由于其多了一维，所以三维的联通区域就不太好用二维图片来查看了。关于三位图片的介绍可以参考第一篇“图像数据的组织方式”。不过借助一些三维可视化软件，如ParaView，我们也可以观察三维图像的内容。

　　下面的图片是显示在ParaView下打开一张三维图片Lobster.raw的结果,我们通过调节不同体素值的透明度，就能凸显出三维图像中所需要关心的内容。
图像预览


像素透明度曲线


　　从上图看，由于三维图像本身是个长方体，各个像素有自己的颜色，所以只能看到外表面的样子。顺便提一下，在三维图像中，我们一般不使用“像素”这个词，而是使用“体素”来形容组成图像的这些单元。我们可以使用ParaView体绘制属性中调节透明度的选项，把一部分值较小的体素设为透明，就能看到这个三维图像中，值较大的体素就像一个个小砖头一样堆成了一个龙虾的形状。也就是说，如果我们在这些体素中找一个种子点，利用三维图像种子点生长算法，就能找到所有组成这个龙虾的体素。

　　上一篇文章已经介绍了二维种子点生长算法的相关知识，介绍了三种种子点生长算法—泛洪法、扫描线法和区段法，并在最后分析了他们的性能。本文的核心部分是三维种子点点生长的讨论，也就是把上一篇文章的算法扩展到三维，那么本文也按照相同的结构来说明这三种算法是如何扩展到三维的。

泛洪法
扫描线法
区段法
算法分析与实验

一、泛洪法

　　泛洪法的基本思想上篇文章详细介绍过，在三维图像中，需要扩展的主要是邻域的范围。邻域从二维扩展到三维后，就不再是四邻域和八邻域而是六邻域和二十六邻域。示意图如下：



　　其中绿色的点表示为当前点（蓝色点）的6邻域，绿色+红色的点为这个像素的26邻域。

　　对于三维点P，P的6邻域可以表示为：

　　　　S=Adj6(P)=｛P0(P.X-1,P.Y,P.Z)，P1(P.X+1,P.Y,P.Z)，P2(P.X,P.Y-1,P.Z)，P3(P.X,P.Y+1P.Z)，P4(P.X,P.Y,P.Z-1)，P5(P.X,P.Y,P.Z+1)｝。

　　P的26邻域可以表示为：

　　　　S=Adj26(P)=｛P0(P.X-1,P.Y,P.Z)，P1(P.X+1,P.Y,P.Z)，P2(P.X,P.Y-1,P.Z)，P3(P.X,P.Y+1,P.Z)，

　　　　　　　　　　　　P4(P.X-1,P.Y-1,P.Z)，P5(P.X+1,P.Y+1,P.Z)，P6(P.X+1,P.Y-1,P.Z)，P7(P.X-1,P.Y+1,P.Z)，

　　　　　　　　　　　　P8(P.X-1,P.Y-1,P.Z-1)，P9(P.X+1,P.Y+1,P.Z-1)，P10(P.X+1,P.Y-1,P.Z-1)，P11(P.X-1,P.Y+1,P.Z-1)

　　　　　　　　　　　　P12(P.X,P.Y-1,P.Z-1)，P13(P.X,P.Y+1,P.Z-1)，P14(P.X+1,P.Y,P.Z-1)，P15(P.X-1,P.Y,P.Z-1)，

　　　　　　　　　　　　P16(P.X,P.Y,P.Z-1)，P17(P.X,P.Y,P.Z+1)，P18(P.X+1,P.Y-1,P.Z+1)，P19(P.X-1,P.Y+1,P.Z+1)，

　　　　　　　　　　　　P20(P.X-1,P.Y-1,P.Z+1)，P21(P.X+1,P.Y+1,P.Z+1)，P22(P.X,P.Y-1,P.Z+1)，P23(P.X,P.Y+1,P.Z+1)

　　　　　　　　　　　　P24(P.X-1,P.Y,P.Z+1)，P25(P.X+1,P.Y,P.Z+1)｝。

　　那么实际上，三维泛洪法采用和二维泛洪法一样的逻辑。

FloodFill(seed,bitmap,includePredicate,process)
　　set all postions of flagMap to false
　　set container Q to empty.
　　push seed into Q
　　set seed position in flagMap to true
　　process(seed)
　　while Q is not empty
　　    pop a point P from Q
　　    foreach point T in Adj(P)
　　        if  includePredicate(T) is true
　　            set position of T in flagMap to true
　　            push T into Q
　　            process(T)
　　 end


　　注意在三维的情况下，相应的点类、图像类、位图标记表类等基础数据结构都需要分别增加一维。相应类的c#定义如下：

public struct Int16Triple
{
public int X;
public int Y;
public int Z;
public Int16Triple(int x, int y,int z)
{
X = x;
Y = y;
Z = z;
}
}

public class BitMap3d
{
public const byte WHITE = 255;
public const byte BLACK = 0;
public byte[] data;
public int width;
public int height;
public int depth;public BitMap3d(int width, int height,int depth, byte v)
{
this.width = width;
this.height = height;
this.depth = depth;
data = new byte[width * height * depth];for (int i = 0; i < width * height*depth; i++)
data[i] = v;
}
public BitMap3d(byte[] data, int width, int height,int depth)
{
this.data = data;
this.width = width;
this.height = height;
this.depth = depth;
}
public void SetPixel(int x, int y,int z,byte v)
{
data[x + y * width+z*width*height] = v;
}
public byte GetPixel(int x, int y,int z)
{return data[x + y * width+z*width*height];
}public void ReadRaw(string path)
{
FileStream fs = new FileStream(path, FileMode.Open, FileAccess.Read);
fs.Read(data, 0, width*height*depth);
fs.Close();
}
public void SaveRaw(string path)
{
FileStream fs = new FileStream(path, FileMode.OpenOrCreate, FileAccess.Write);
fs.Write(data, 0, data.Length);
fs.Close();
}
}

public class FlagMap3d
{public int width;
public int height;
public int depth;
BitArray flags;
public FlagMap3d(int width, int height,int depth)
{
this.width = width;
this.height = height;
this.depth = depth;
flags = new BitArray(width * height*depth, false);
}
public void SetFlagOn(int x, int y, int z,bool v)
{
flags[x + y * width+z*width*height] = v;
}
public bool GetFlagOn(int x, int y,int z)
{return flags[x + y * width+z*width*height];
}

}


　　相应实现的c#代码如下：

class FloodFill3d
{

protected BitMap3d bmp;
protected FlagMap3d flagsMap;
protected Container<Int16Triple> container;
protected int count = 0;
public byte min;
public byte max;
protected bool IncludePredicate(Int16Triple p)
{
byte v = bmp.GetPixel(p.X, p.Y, p.Z);
return v > min && v < max;
}
protected void Process(Int16Triple p)
{
count++;
return;
}
protected virtual void ExcuteFloodFill(BitMap3d data, Int16Triple seed)
{
this.bmp = data;
data.ResetVisitCount();
flagsMap = new FlagMap3d(data.width, data.height,data.depth);
Int16Triple[] adjPoints6 = new Int16Triple[6];
flagsMap.SetFlagOn(seed.X, seed.Y, seed.Z,true);
container.Push(seed);
Process(seed);
while (!container.Empty())
{
Int16Triple p = container.Pop();
InitAdj6(ref adjPoints6, ref p);
for (int adjIndex = 0; adjIndex < 6; adjIndex++)
{
Int16Triple t = adjPoints6[adjIndex];
if (t.X < data.width && t.X >= 0 && t.Y < data.height && t.Y >= 0&&t.Z<data.depth&&t.Z>=0)
{
if (!flagsMap.GetFlagOn(t.X, t.Y,t.Z) && IncludePredicate(t))
{
flagsMap.SetFlagOn(t.X, t.Y,t.Z, true);
container.Push(t);
Process(t);
}
}
}
}
return;
}
protected void InitAdj6(ref Int16Triple[] adjPoints6, ref Int16Triple p)
{
adjPoints6[0].X = p.X - 1;
adjPoints6[0].Y = p.Y;
adjPoints6[0].Z = p.Z;

adjPoints6[1].X = p.X + 1;
adjPoints6[1].Y = p.Y;
adjPoints6[1].Z = p.Z;

adjPoints6[2].X = p.X;
adjPoints6[2].Y = p.Y - 1;
adjPoints6[2].Z = p.Z;

adjPoints6[3].X = p.X;
adjPoints6[3].Y = p.Y + 1;
adjPoints6[3].Z = p.Z;

adjPoints6[4].X = p.X;
adjPoints6[4].Y = p.Y;
adjPoints6[4].Z = p.Z-1;

adjPoints6[5].X = p.X;
adjPoints6[5].Y = p.Y;
adjPoints6[5].Z = p.Z+1;
}
}


二、扫描线法

　　二维的扫描线法，是从一个点弹出堆栈后，对其左右扫描出两个端点，然后再对端点范围的两侧即Y-1和Y+1方向进行CheckRange。算法拓展到三维，只需增加检查Z-1和Z+1方向即可。也就是“两侧”变成了“四侧”，即:

（xleft，p.Y-1，p.Z）~（xright，p.Y-1，p.Z）
（xleft，p.Y+1，p.Z）~（xright，p.Y+1，p.Z）
（xleft，p.Y，p.Z-1）~（xright，p.Y，p.Z-1）
（xleft，p.Y，p.Z+1）~（xright，p.Y，p.Z+1）

　　扫描这四条线等价于6向泛洪法，假如扫描下面八条线，就相当于26向泛洪法

（xleft-1，p.Y-1，p.Z）~（xright+1，p.Y-1，p.Z）
（xleft-1，p.Y+1，p.Z）~（xright+1，p.Y+1，p.Z）
（xleft-1，p.Y，p.Z-1）~（xright+1，p.Y，p.Z-1）
（xleft-1，p.Y，p.Z+1）~（xright+1，p.Y，p.Z+1）
（xleft-1，p.Y-1，p.Z-1）~（xright+1，p.Y-1，p.Z-1）
（xleft-1，p.Y+1，p.Z+1）~（xright+1，p.Y+1，p.Z+1）
（xleft-1，p.Y+1，p.Z-1）~（xright+1，p.Y+1，p.Z-1）
（xleft-1，p.Y-1，p.Z+1）~（xright+1，p.Y-1，p.Z+1）

　　扫描线法拓展到三维，基本操作FindXleft，FindXright，CheckRange分别修改至相应的三维即可。其作用没有任何变化。

　　下面是三维扫描线线的c#代码：

class ScanlineFill3d
{
protected int count = 0;
protected Container<Int16Triple> container;//这个容器可以是Queue和Stack中任意一种，这里抽象成一个Container
protected BitMap3d bmp;
public FlagMap3d flagsMap;
protected virtual void ExcuteScanlineFill(BitMap3d data, Int16Triple seed)
{
this.bmp = data;
data.ResetVisitCount();
flagsMap = new FlagMap3d(data.width, data.height, data.depth);
flagsMap.SetFlagOn(seed.X, seed.Y, seed.Z, true);
container.Push(seed);
Process(seed);
while (!container.Empty())
{
Int16Triple p = container.Pop();
int xleft = FindXLeft(p);
int xright = FindXRight(p);
if (p.Y - 1 >= 0)
CheckRange(xleft, xright, p.Y - 1, p.Z);
if (p.Y + 1 < data.height)
CheckRange(xleft, xright, p.Y + 1, p.Z);
if (p.Z - 1 >= 0)
CheckRange(xleft, xright, p.Y, p.Z - 1);
if (p.Z + 1 < data.depth)
CheckRange(xleft, xright, p.Y, p.Z + 1);
}
}//该函数为扫描线法主体
protected void CheckRange(int xleft, int xright, int y, int z)
{
for (int i = xleft; i <= xright; )
{
if ((!flagsMap.GetFlagOn(i, y, z)) && IncludePredicate(i, y, z))
{
int rb = i + 1;
while (rb <= xright && (!flagsMap.GetFlagOn(rb, y, z)) && IncludePredicate(rb, y, z))
{
rb++;
}
rb--;
Int16Triple t = new Int16Triple(rb, y, z);
flagsMap.SetFlagOn(rb, y, z, true);
container.Push(t);
Process(t);
i = rb + 1;
}
else
{
i++;
}
}
}//CheckRange操作
protected int FindXLeft(Int16Triple p)
{
int xleft = p.X - 1;
while (true)
{
if (xleft == -1 || flagsMap.GetFlagOn(xleft, p.Y, p.Z))
{
break;
}
else
{
byte value = bmp.GetPixel(xleft, p.Y, p.Z);
if (IncludePredicate(xleft, p.Y, p.Z))
{
Int16Triple t = new Int16Triple(xleft, p.Y, p.Z);
flagsMap.SetFlagOn(xleft, p.Y, p.Z, true);
Process(t);
xleft--;
}
else
{
break;
}
}
}
return xleft + 1;
}//FindXLeft操作
protected int FindXRight(Int16Triple p)
{
int xright = p.X + 1;
while (true)
{
if (xright == bmp.width || flagsMap.GetFlagOn(xright, p.Y, p.Z))
{
break;
}
else
{
byte value = bmp.GetPixel(xright, p.Y, p.Z);
if (IncludePredicate(xright, p.Y, p.Z))
{
Int16Triple t = new Int16Triple(xright, p.Y, p.Z);
flagsMap.SetFlagOn(xright, p.Y, p.Z, true);
Process(t);
xright++;
}
else
{
break;
}
}
}
return xright - 1;
}//FindXRight操作
public byte min;
public byte max;
protected bool IncludePredicate(int x,int y,int z)
{
byte v = bmp.GetPixel(x, y,z);
return v > min && v < max;
}
protected void Process(Int16Triple p)
{
count++;
}
}


三、区段法

　　区段法扩展到三维也遵循和扫描线法一样的原则，相应的基本操作FindXleft，FindXright，CheckRange需要分别修改到三维的情况。同时在主函数逻辑中，监测相应的相邻区段时要注意ParentDirection所指示的父区段的特殊性。

　　下面是三维区段法的C#代码：

enum ParentDirections
{
Y0 = 1, Y2 = 2,Z0=3,Z2=4, Non = 5
}
enum ExtendTypes
{
LeftRequired = 1, RightRequired = 2, AllRez = 3, UnRez = 4
}
struct Span
{
public int XLeft;
public int XRight;
public int Y;
public int Z;
public ExtendTypes Extended;
public ParentDirections ParentDirection;
}
class SpanFill3d
{
protected int count = 0;
protected BitMap3d bmp;
public FlagMap3d flagsMap;
protected Container<Span> container;//以Span为单位的Queue或Stack容器
protected virtual void ExcuteSpanFill(BitMap3d data, Int16Triple seed)
{
this.bmp = data;
data.ResetVisitCount();
flagsMap = new FlagMap3d(data.width, data.height,data.depth);
Process(seed);
flagsMap.SetFlagOn(seed.X, seed.Y,seed.Z, true);
Span seedspan = new Span();
seedspan.XLeft = seed.X;
seedspan.XRight = seed.X;
seedspan.Y = seed.Y;
seedspan.Z = seed.Z;
seedspan.ParentDirection = ParentDirections.Non;
seedspan.Extended = ExtendTypes.UnRez;
container.Push(seedspan);

while (!container.Empty())
{
Span span = container.Pop();
#region AllRez
if (span.Extended == ExtendTypes.AllRez)
{
if (span.ParentDirection == ParentDirections.Y2)
{
if (span.Y - 1 >= 0)//[spx-spy,y-1,z]
CheckRange(span.XLeft, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Z - 1 >= 0)//[spx-spy,y,z-1]
CheckRange(span.XLeft, span.XRight, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < data.depth)//[spx-spy,y,z+1]
CheckRange(span.XLeft, span.XRight, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Y0)
{
if (span.Y + 1 < bmp.height)//[spx-spy,y+1,z]
CheckRange(span.XLeft, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 > 0)//[spx-spy,y,z-1]
CheckRange(span.XLeft, span.XRight, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < data.depth)//[spx-spy,y,z+1]
CheckRange(span.XLeft, span.XRight, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z2)
{
if (span.Y - 1 >= 0)//[spx-spy,y-1,z]
CheckRange(span.XLeft, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)//[spx-spy,y+1,z]
CheckRange(span.XLeft, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)//[spx-spy,y,z-1]
CheckRange(span.XLeft, span.XRight, span.Y, span.Z - 1, ParentDirections.Z2);
continue;
}
if (span.ParentDirection == ParentDirections.Z0)
{
if (span.Y - 1 >= 0)
CheckRange(span.XLeft, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(span.XLeft, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z + 1 < data.depth)
CheckRange(span.XLeft, span.XRight, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
throw new Exception();
}
#endregion
#region UnRez
if (span.Extended == ExtendTypes.UnRez)
{
int xl = FindXLeft(span.XLeft, span.Y,span.Z);
int xr = FindXRight(span.XRight, span.Y,span.Z);
if (span.ParentDirection == ParentDirections.Y2)
{
if (span.Y - 1 >= 0)
CheckRange(xl, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
{
if (xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.XRight != xr)
CheckRange(span.XRight, xr, span.Y + 1, span.Z, ParentDirections.Y0);
}
if (span.Z - 1 >= 0)
CheckRange(xl, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(xl, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Y0)
{
if (span.Y + 1 < bmp.height)
CheckRange(xl, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Y - 1 >= 0)
{
if (xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.XRight != xr)
CheckRange(span.XRight, xr, span.Y - 1, span.Z, ParentDirections.Y2);
}
if (span.Z - 1 >= 0)
CheckRange(xl, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(xl, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z2)
{
if (span.Y - 1 >= 0)
CheckRange(xl, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(xl, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(xl, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
{
if (xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y, span.Z+1, ParentDirections.Z0);
if (span.XRight != xr)
CheckRange(span.XRight, xr, span.Y, span.Z+1, ParentDirections.Z0);
}
continue;
}
if (span.ParentDirection == ParentDirections.Z0)
{
if (span.Y - 1 >= 0)
CheckRange(xl, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(xl, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
{
if (xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y, span.Z -1, ParentDirections.Z2);
if (span.XRight != xr)
CheckRange(span.XRight, xr, span.Y, span.Z - 1, ParentDirections.Z2);
}
if (span.Z + 1 < bmp.depth)
CheckRange(xl, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Non)
{
if (span.Y + 1 < bmp.height)
CheckRange(xl, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Y - 1 >= 0)
CheckRange(xl, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Z - 1 >= 0)
CheckRange(xl, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(xl, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
throw new Exception();
}
#endregion
#region LeftRequired
if (span.Extended == ExtendTypes.LeftRequired)
{
int xl = FindXLeft(span.XLeft, span.Y,span.Z);
if (span.ParentDirection == ParentDirections.Y2)
{
if (span.Y - 1 >= 0)
CheckRange(xl, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height && xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(xl, span.XRight, span.Y , span.Z-1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(xl, span.XRight, span.Y , span.Z+1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Y0)
{
if (span.Y - 1 >= 0 && xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(xl, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(xl, span.XRight, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(xl, span.XRight, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z2)
{
if (span.Y - 1 >= 0 )
CheckRange(xl, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(xl, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(xl, span.XRight, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth && xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z0)
{
if (span.Y - 1 >= 0)
CheckRange(xl, span.XRight, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(xl, span.XRight, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0 && xl != span.XLeft)
CheckRange(xl, span.XLeft, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth )
CheckRange(xl, span.XRight, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
throw new Exception();
}
#endregion
#region RightRequired
if (span.Extended == ExtendTypes.RightRequired)
{
int xr = FindXRight(span.XRight, span.Y,span.Z);

if (span.ParentDirection == ParentDirections.Y2)
{
if (span.Y - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y - 1,span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height && span.XRight != xr)
CheckRange(span.XRight, xr, span.Y + 1,span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y, span.Z-1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(span.XLeft, xr, span.Y, span.Z+1, ParentDirections.Z0);
continue;
}

if (span.ParentDirection == ParentDirections.Y0)
{
if (span.Y + 1 < bmp.height)
CheckRange(span.XLeft, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Y - 1 >= 0 && span.XRight != xr)
CheckRange(span.XRight, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Z - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(span.XLeft, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z2)
{
if (span.Y - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(span.XLeft, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth && span.XRight != xr)
CheckRange(span.XRight, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
if (span.ParentDirection == ParentDirections.Z0)
{
if (span.Y - 1 >= 0)
CheckRange(span.XLeft, xr, span.Y - 1, span.Z, ParentDirections.Y2);
if (span.Y + 1 < bmp.height)
CheckRange(span.XLeft, xr, span.Y + 1, span.Z, ParentDirections.Y0);
if (span.Z - 1 >= 0 && span.XRight != xr)
CheckRange(span.XRight, xr, span.Y, span.Z - 1, ParentDirections.Z2);
if (span.Z + 1 < bmp.depth)
CheckRange(span.XLeft, xr, span.Y, span.Z + 1, ParentDirections.Z0);
continue;
}
throw new Exception();
}
#endregion
}
}
protected void CheckRange(int xleft, int xright, int y, int z,ParentDirections ptype)
{
for (int i = xleft; i <= xright; )
{
if ((!flagsMap.GetFlagOn(i, y,z)) && IncludePredicate(i, y,z))
{
int lb = i;
int rb = i + 1;
while (rb <= xright && (!flagsMap.GetFlagOn(rb, y,z)) && IncludePredicate(rb, y,z))
{
rb++;
}
rb--;

Span span = new Span();
span.XLeft = lb;
span.XRight = rb;
span.Y = y;
span.Z = z;
if (lb == xleft && rb == xright)
{
span.Extended = ExtendTypes.UnRez;
}
else if (rb == xright)
{
span.Extended = ExtendTypes.RightRequired;
}
else if (lb == xleft)
{
span.Extended = ExtendTypes.LeftRequired;
}
else
{
span.Extended = ExtendTypes.AllRez;
}
span.ParentDirection = ptype;
for (int j = lb; j <= rb; j++)
{
flagsMap.SetFlagOn(j, y,z, true);
Process(new Int16Triple(j, y,z));
}
container.Push(span);

i = rb + 1;
}
else
{
i++;
}
}
}//区段法的CheckRange 注意与扫描线的CheckRange的不同
protected int FindXRight(int x, int y,int z)
{
int xright = x + 1;
while (true)
{
if (xright == bmp.width || flagsMap.GetFlagOn(xright, y,z))
{
break;
}
else
{
if (IncludePredicate(xright, y,z))
{
Int16Triple t = new Int16Triple(xright, y,z);
flagsMap.SetFlagOn(xright, y,z ,true);
Process(t);
xright++;
}
else
{
break;
}
}
}
return xright - 1;
}
protected int FindXLeft(int x, int y,int z)
{
int xleft = x - 1;
while (true)
{
if (xleft == -1 || flagsMap.GetFlagOn(xleft, y, z))
{
break;
}
else
{
if (IncludePredicate(xleft, y,z))
{
Int16Triple t = new Int16Triple(xleft, y,z);
flagsMap.SetFlagOn(xleft, y,z, true);
Process(t);
xleft--;
}
else
{
break;
}
}
}
return xleft + 1;
}
public byte min;
public byte max;
protected bool IncludePredicate(int x, int y, int z)
{
byte v = bmp.GetPixel(x, y, z);
return v > min && v < max;
}
protected void Process(Int16Triple p)
{
count++;
}
}


四、算法测试与实验

　　算法测试采用来自http://www.volvis.org的5组体数据，注意这里的阈值范围是值一个体素值的范围（min,max），在此范围内的体素才纳入区域。也就是说本文的includePredicate不再是跟上一篇一样是判断像素灰度为255的纳入区域，而是判断是否在这个（min,max）的范围内。
图像预览数据名称参数说明

Lobster
大小：301×324×56
文件：Lobster.raw
种子点：(124, 168, 27)
阈值范围：37~255
描述：CT scan of a lobster contained in a block of resin.


engine大小：256×256×128
文件：Engine.raw
种子点：(149, 44, 43)
阈值范围：64~255
描述：CT scan of two cylinders of an engine block.


backpack大小：301×324×56
文件：Backpack.raw
种子点：(53, 390, 160)
阈值范围：46~255
描述：CT scan of a backpack filled with items.


phantom大小：512×512×442
文件：Phantom.raw
种子点：（256，256，227）
阈值范围：42~255
描述：CT scan of a Colon phantom with several different objects and five pedunculated large polyps in the central object.


cube大小：200×200×200
文件：Cube.raw
种子点：（50，50，50）
阈值范围：0~255
描述：全白色，用于测试全部充满的极端情况

　　测试结果如下：
测试数据泛洪法（栈式）泛洪法(队列式)扫描线法（栈式）扫描线法（队列式）区段法（栈式）区段法（队列式）区域点数
lobster86ms76ms56ms64ms40ms47ms277367
engine329ms345ms216ms254ms137ms183ms1154807
backpack596ms638ms426ms470ms319ms351ms19115450
phantom13468ms14118ms6520ms8279ms3701ms4066ms39759257
cube2120ms2259ms1292ms1404ms691ms692ms800000
　　对于其中的backpack数据的单元访问比例统计如下：
算法GetPixel/总点数GetFlagOn/总点数SetFlagOn/总点数Push和Pop/总点数总点数时间花费
泛洪法（栈式）1.9455815.9971031.01.01915450596ms
泛洪法（队列式）1.9455815.9971031.01.01915450638ms
扫描线法（栈式）4.1650885.3374061.00.26681915450426ms
扫描线法（队列式）2.7154245.7745231.00.75891915450470ms
区段法（栈式）1.9114853.5267281.00.21471915450319ms
区段法（队列式）1.9319633.7543331.00.36281915450351ms
　　通过测试可以看出，相比于二维种子点生长，由于多了一维，所以利用X轴数据连续特性所能达到增加效率的效果明显减弱。既可以从时间结果上看，也可以从单元访问比例上看，都能得到这一结论。不过从整体上看，效率上仍然能够得出如下结论：栈式算法略块于队列式算法；区段法快于扫描线法，扫描线法快于泛洪法。
内容来自用户分享和网络整理，不保证内容的准确性，如有侵权内容，可联系管理员处理 点击这里给我发消息
标签： 
相关文章推荐
新的分享
章节导航