骨骼外轮廓的三维建模
继上一章所得各CT图像的骨骼外轮廓点,我们可以根据二维轮廓线外表面拟合算法进行三维建模。该算法可以将建模过程分解成两两相邻CT图像间的侧表面三角面片生成的子过程。
夏旸师兄在论文中提到的三角面片生成算法和实际代码中应用的算法稍有不同。
首先是论文中提到的算法,该算法的输入是两张CT图像的外轮廓点集,输出是三角面片三个顶点序列的集合。我们假设提取到的轮廓都是单联通的闭合轮廓线,排除拼接过程中的拓扑变化,不考虑分叉的问题,并且两平面轮廓线都处在同一平面上,与Z轴垂直。
现在要解决的是顶点匹配问题,对于轮廓线A、B,有顶点列表ListA
、ListB
。理想状态下,ListA
、ListB
有相同数量的点,这样可以直接得到三角面片网格:
但不幸的是,ListA
、ListB
的数量是可以不一样的。于是对于ListA
中的每一点,需要在ListB
中寻找到与之距离最近的点组成三角形的一条边,假设是Ai
、Bi
。然后在Ai+1
、Bi+1
两个点中选取一点作为三角形的第三个顶点。为了保持网格形状的均匀平滑,会选择更加趋近于锐角三角形的三角面片。
分别计算点Ai
和点Bi+1
的距离S1
,以及点Bi
和点Ai+1
的距离S2
,比较S1
和S2
的长度。如果S1 > S2
,选择Ai+1
作为第三个顶点,反之选择Bi+1
作为第三个顶点。顶点选择后,生成三角面片,会在进行下一轮计算之前更新当前顶点。在这里,系统会选择Bi+1
作为第三个顶点。在进行下一轮计算之前,Bi+1
会替代Bi
进行下一轮的计算,即下一轮会为顶点Ai
、Bi+1
在Ai+1
、Bi+2
中寻找第三个顶点。
实际代码中,并没有计算S1
和S2
,而是通过判断ListA
、ListB
中匹配成功点的比例来决定第三个点的选择。也就是说,假设现在ListA
中已经生成三角面片点的数量所占的比例比ListB
的小,那么即使S1 < S2
,也会选择Ai+1
。咨询师兄为什么不采用论文所述算法,回答是论文算法在实际应用中有的效果并不理想,具体怎么不理想师兄也忘了……
二维轮廓表面拟合算法由工具类MeshDataUtil
负责,函数:
public static List<int> CreateMeshByTwoSlice(List<Vector3> listVectorA, List<Vector3> listVectorB,int offset)
{
if (listVectorA.Count < 2 || listVectorB.Count < 2) return null;
Vector3 firstVectorA = listVectorA[0];
// 在ListB中寻找最近的点
int minIndex = 0;
float minDistance = Vector2.Distance(firstVectorA, listVectorB[0]);
for (int i = 0; i < listVectorB.Count; i++)
{
if (Vector2.Distance(firstVectorA, listVectorB[i]) < minDistance)
{
minDistance = Vector2.Distance(firstVectorA, listVectorB[i]);
minIndex = i;
}
}
int offsetA = offset; // offset是该点序号在所有CT图像轮廓点中的偏移值
int offsetB = listVectorA.Count + offsetA; // ListB中的点,序号偏移量还要加上ListA点的总数
int countA = 0;
int countB = 0;
List<int> triangle = new List<int>(); // 三角形顶点数组
while (countA < listVectorA.Count && countB < listVectorB.Count) // 直到所有点都匹配成功
{
int realIndexA = countA + offsetA;
int realIndexB = (countB + minIndex) % listVectorB.Count + offsetB;
int nextIndexA = (countA + 1) % listVectorA.Count + offsetA;
int nextIndexB = (countB + minIndex + 1) % listVectorB.Count + offsetB;
float percentNextA = (countA + 1) * 1.0f / listVectorA.Count; // ListA中匹配成功点的比例
float percentNextB = (countB + 1) * 1.0f / listVectorB.Count; // ListB中匹配成功点的比例
if (percentNextA == 1 && percentNextB == 1)
{ // 结尾特殊处理
triangle.Add(realIndexA);
triangle.Add(nextIndexA);
triangle.Add(realIndexB);
triangle.Add(realIndexB);
triangle.Add(nextIndexA);
triangle.Add(nextIndexB);
break;
}
if (percentNextA < percentNextB)
{ // 选择Ai+1
triangle.Add(realIndexA);
triangle.Add(nextIndexA);
triangle.Add(realIndexB);
countA++;
}
else
{ // 选择Bi+1
triangle.Add(nextIndexB);
triangle.Add(realIndexB);
triangle.Add(realIndexA);
countB++;
}
}
return triangle;
}
整个的三维建模就是迭代的进行两两平面的三角面片生成,在每次计算完一张CT图像的骨骼外轮廓后计算三角面片,就可以动态地观察到模型的生成过程:
listVectorB = listVertex.GetRange(offset + listVectorA.Count, listContour.Count);
listTriangle.AddRange(MeshDataUtil.CreateMeshByTwoSlice(listVectorA, listVectorB, offset));
surfaceMeshFilter.mesh.SetVertices(listVertex);
surfaceMeshFilter.mesh.SetTriangles(listTriangle, 0);