tcston2.cpp

并行TIN生成算法, 基于DeWall算法理论实现

					continue;
				// 增加共圆点
				TIN_X *X = new TIN_X;
				X->cid = ci;
				X->pid = pi;
				L->AddTail(X);
			}
			else if (dn <  d2 && p4 != *p3)	// 有点
			{
				// 紧缩搜索
				TIN_X *X = (TIN_X*)L->GetHead();
				X->cid = ci;
				X->pid = pi;
				// 清理其候选
				X = (TIN_X*)L->RemoveAt(1);
				while (X)
				{
					delete X;
					X = (TIN_X*)L->RemoveAt(1);
				}
				// 计算外接圆圆心
				ComCen2(p1, p2, p4, ce);
				// 计算外接圆半径平方
				d2 = EuLen2D(p4, ce);
				// 空圆检测和紧缩搜索
				*p3 = p4;

				return 1;
			}
		}
	}

	return 0;
}
// 空圆检测:有异侧检测
// 优化开关
#define OPTI_OU1
inline TCS_S32 RevEmp2(TIN_G *G, TCS_F64 *P, TCS_F64 *p1, TCS_F64 *p2, TCS_S32 on, TCS_F64 *p3, TCS_F64 *ce, TCS_F64 &d2, CTcsBL *L)
{
	TCS_S32 i, k, c[6], r[6], o[2];
	TCS_S32 in, ci, co = G->num[0], ro = G->num[1];

REV_EMP2:
	// 计算搜索边界
	ComBnd2(G, ce, sqrt(d2), c[4], c[5], r[4], r[5]);

#ifdef	OPTI_OU1
	// 采用外扩算法优化
	// 从圆心单元开始
	ComIdx2(G, ce, o);
	// 修正开始单元
	if		(o[0] < 0)		o[0] = 0;
	else if (o[0] > co - 1)	o[0] = co - 1;
	if		(o[1] < 0)		o[1] = 0;
	else if (o[1] > ro - 1)	o[1] = ro - 1;

	// 中心单元空圆搜索
	k = o[0], i = o[1];	ci = i * co + k;
	in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
	if		(in == 1)	// 紧缩搜索
		goto REV_EMP2;
	else if (in == 2)	// 包含点
		return 1;

	// 继续空圆检测
	if (c[5] - c[4] > 0 || r[5] - r[4] > 0)
	{
		c[0] = c[1] = c[2] = c[3] = k;
		r[0] = r[1] = r[2] = r[3] = i;
		while (1)
		{
			c[0]--;	if (c[0] < c[4])	c[0] = c[4];
			c[1]++;	if (c[1] > c[5])	c[1] = c[5];
			r[0]--;	if (r[0] < r[4])	r[0] = r[4];
			r[1]++;	if (r[1] > r[5])	r[1] = r[5];
			/////////////////////////////////////////////////////////////
			if (r[2] != r[0])	// 即r[2] != r[4]
			{
				ci = r[0] * co;
				for (k = c[0] + 1; k < c[1]; k++)
				{
					in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci + k, L);
					if		(in == 1)
						goto REV_EMP2;
					else if (in == 2)
						return 1;
				}
			}
			if (r[3] != r[1])	// 即r[3] != r[5]
			{
				ci = r[1] * co;
				for (k = c[0] + 1; k < c[1]; k++)
				{
					in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci + k, L);
					if		(in == 1)
						goto REV_EMP2;
					else if (in == 2)
						return 1;
				}
			}
			if (c[2] != c[0])	// 即c[2] != c[4]
			{
				for (i = r[0] + 1; i < r[1]; i++)
				{
					ci = i * co + c[0];
					in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
					if		(in == 1)
						goto REV_EMP2;
					else if (in == 2)
						return 1;
				}
			}
			if (c[3] != c[1])	// 即c[3] != c[5]
			{
				for (i = r[0] + 1; i < r[1]; i++)
				{
					ci = i * co + c[1];
					in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
					if		(in == 1)
						goto REV_EMP2;
					else if (in == 2)
						return 1;
				}
			}
			// 处理4个角单元
			if (c[2] != c[4] || r[2] != r[4])
			{
				ci = r[0] * co + c[0];
				in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
				if		(in == 1)
					goto REV_EMP2;
				else if (in == 2)
					return 1;
			}
			if (c[2] != c[4] || r[3] != r[5])
			{
				ci = r[1] * co + c[0];
				in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
				if		(in == 1)
					goto REV_EMP2;
				else if (in == 2)
					return 1;
			}
			if (c[3] != c[5] || r[2] != r[4])
			{
				ci = r[0] * co + c[1];
				in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
				if		(in == 1)
					goto REV_EMP2;
				else if (in == 2)
					return 1;
			}
			if (c[3] != c[5] || r[3] != r[5])
			{
				ci = r[1] * co + c[1];
				in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci, L);
				if		(in == 1)
					goto REV_EMP2;
				else if (in == 2)
					return 1;
			}				
			/////////////////////////////////////////////////////////////
			
			c[2] = c[0], c[3] = c[1], r[2] = r[0], r[3] = r[1];
			if (c[2] == c[4] && c[3] == c[5] &&
				r[2] == r[4] && r[3] == r[5])
				break;
		}
	}

#else
	// 先列后行顺序
	for (i = r[4]; i <= r[5]; i++)
	{
		ci = i * co;
		for (k = c[4]; k <= c[5]; k++)
		{
			in = RevEmp2(G, P, p1, p2, on, &p3, ce, d2, ci + k, L);
			if		(in == 1)
				goto REV_EMP2;
			else if (in == 2)
				return 1;
		}
	}

#endif

	return 0;
}
// 空圆检测:无异侧检测

// 网格单元搜索:有异侧检测
// 1:成功;0失败
inline TCS_S32 FndEmp2(TIN_G *G, TCS_F64 *P, TCS_F64 *p1, TCS_F64 *p2, TCS_F64 *op, TCS_S32 on, TCS_S32 ci, TCS_F64 *ce, TCS_F64 &d2, CTcsBL *L)
{
	TCS_S32 n, pi;
	TCS_F64 *p3, dn;
	TIN_C *C = G->idx + ci;

	for (n = 0; n < C->ocn; n++)
	{
		pi = C->pid[n];
		p3 = P + 2 * pi;
		// 异侧判断
		// 下面可以精确处理p1 == p3、p2 == p3和共线问题
		if (on == IsSide2(p1, p2, p3))
		{
			// 计算外接圆圆心
			ComCen2(p1, p2, p3, ce);
			// 计算外接圆半径平方
			d2 = EuLen2D(p3, ce);
			// 首先检测op是否在外接圆内
			dn = EuLen2D(op, ce);
			// 精确处理浮点数相等问题
			if (d2 > dn)
				continue;

			// 增加候选点
			TIN_X *X = new TIN_X;
			X->cid = ci;
			X->pid = pi;
			L->AddTail(X);
			// 空圆检测和紧缩搜索
			if (RevEmp2(G, P, p1, p2, on, p3, ce, d2, L))
				continue;

			return 1;
		}
	}

	return 0;
}
// 网格单元搜索:无异侧检测

// 批量单元搜索
// 1:成功;0:失败
inline TCS_S32 FndEmp2(TIN_G *G, TCS_F64 *P, TCS_F64 *p1, TCS_F64 *p2, TCS_F64 *op, TCS_S32 on,
					   TCS_F64 *ce, TCS_F64 &d2, CTcsBL *L, TCS_S32 co, TCS_S32 *c, TCS_S32 *r)
{
	TCS_S32 i, k, ci;

	if (r[2] != r[0])	// 即r[2] != r[4]
	{
		ci = r[0] * co;
		for (k = c[0] + 1; k < c[1]; k++)
		{
			if (FndEmp2(G, P, p1, p2, op, on, ci + k, ce, d2, L))
				return 1;
		}
	}
	if (r[3] != r[1])	// 即r[3] != r[5]
	{
		ci = r[1] * co;
		for (k = c[0] + 1; k < c[1]; k++)
		{
			if (FndEmp2(G, P, p1, p2, op, on, ci + k, ce, d2, L))
				return 1;
		}
	}
	if (c[2] != c[0])	// 即c[2] != c[4]
	{
		for (i = r[0] + 1; i < r[1]; i++)
		{
			ci = i * co + c[0];
			if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
				return 1;
		}
	}
	if (c[3] != c[1])	// 即c[3] != c[5]
	{
		for (i = r[0] + 1; i < r[1]; i++)
		{
			ci = i * co + c[1];
			if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
				return 1;
		}
	}
	// 处理4个角单元
	if (c[2] != c[4] || r[2] != r[4])
	{
		ci = r[0] * co + c[0];
		if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
			return 1;
	}
	if (c[2] != c[4] || r[3] != r[5])
	{
		ci = r[1] * co + c[0];
		if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
			return 1;
	}
	if (c[3] != c[5] || r[2] != r[4])
	{
		ci = r[0] * co + c[1];
		if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
			return 1;
	}
	if (c[3] != c[5] || r[3] != r[5])
	{
		ci = r[1] * co + c[1];
		if (FndEmp2(G, P, p1, p2, op, on, ci, ce, d2, L))
			return 1;
	}

	return 0;
}

// 搜索第3点:空圆法则
// 优化开关
#define OPTI_OU2
inline TCS_V00 FndEmp2(TIN_G *G, TIN_S *S, TCS_F64 *P, TCS_F64 *ce, TCS_F64 &d2, CTcsBL *L)
{
	TCS_S32 i, k, ci, c[6], r[6], o[2];

	TCS_F64 *sp = P + 2 * S->sid[0];
	TCS_F64 *ep = P + 2 * S->eid[0];
	TCS_F64 *op = P + 2 * S->pid[0];
	TCS_S32 co = G->num[0], ro = G->num[1];
	// 计算搜索方向
	TCS_S32 on = -IsSide2(sp, ep, op); // 反向
	// 计算初始搜索范围
	ce[0] = (sp[0] + ep[0]) / 2;
	ce[1] = (sp[1] + ep[1]) / 2;d2 = EuLen2D(sp, ce);
	ComBnd2(G, ce, ::sqrt(d2), c[4], c[5], r[4], r[5]);

#ifdef OPTI_OU2
	// 采用外扩算法优化
#ifdef OPTI_OU4
	// 计算垂直平分线与圆的交点
	TCS_F64 p1[2], p2[2];

	if		(sp[0] == ep[0])	// 垂直
	{
		TCS_F64 dd = ::sqrt(d2);
		p1[0] = ce[0] - dd;
		p2[0] = ce[0] + dd;
		p1[1] = p2[1] = ce[1];
	}
	else if (sp[1] == ep[1])
	{
		TCS_F64 dd = ::sqrt(d2);
		p1[0] = p2[0] = ce[0];
		p1[1] = ce[1] - dd;
		p2[1] = ce[1] + dd;
	}
	else
	{
		TCS_FLT64 s = (ep[0] - ce[0]) / (ep[1] - ce[1]);
		TCS_FLT64 d = sqrt(d2 / ( 1.0 + s * s));
		
		p1[0] = ce[0] + s * d;p1[1] = ce[1] + d;
		p2[0] = ce[0] - s * d;p2[1] = ce[1] - d;
	}
	if (IsSide2(sp, ep, p1) == on)
		ComIdx2(G, p1, o);
	else
		ComIdx2(G, p2, o);
	if		(o[0] < 0)		o[0] = 0;
	else if (o[0] > co - 1)	o[0] = co - 1;
	if		(o[1] < 0)		o[1] = 0;
	else if (o[1] > ro - 1)	o[1] = ro - 1;
#else
	// 从圆心单元开始
	ComIdx2(G, ce, o);
	// 修正开始单元
	if		(o[0] < 0)		o[0] = 0;
	else if (o[0] > co - 1)	o[0] = co - 1;
	if		(o[1] < 0)		o[1] = 0;
	else if (o[1] > ro - 1)	o[1] = ro - 1;
#endif

	k = o[0], i = o[1];	ci = i * co + k;
	if (FndEmp2(G, P, sp, ep, op, on, ci, ce, d2, L))
		return;

	if (c[5] - c[4] > 0 || r[5] - r[4] > 0)
	{
		c[0] = c[1] = c[2] = c[3] = k;
		r[0] = r[1] = r[2] = r[3] = i;
		while (1)
		{
			c[0]--;	if (c[0] < c[4])	c[0] = c[4];
			c[1]++;	if (c[1] > c[5])	c[1] = c[5];
			r[0]--;	if (r[0] < r[4])	r[0] = r[4];
			r[1]++;	if (r[1] > r[5])	r[1] = r[5];
			if (FndEmp2(G, P, sp, ep, op, on, ce, d2, L, co, c, r))
				return;
			
			c[2] = c[0], c[3] = c[1], r[2] = r[0], r[3] = r[1];
			if (c[2] == c[4] && c[3] == c[5] &&
				r[2] == r[4] && r[3] == r[5])
				break;
		}
	}

#else
	// 先列后行搜索
	for (i = r[4]; i <= r[5]; i++)
	{
		ci = i * co;
		for (k = c[4]; k <= c[5]; k++)
		{
			if (FndEmp2(G, P, sp, ep, op, on, ci + k, ce, d2, L))
				return;
		}
	}

#endif

	// 继续外扩搜索
	c[0] = c[4];	c[1] = c[5];		r[0] = r[4];	r[1] = r[5];
	c[2] = c[0];	c[3] = c[1];		r[2] = r[0];	r[3] = r[1];
	c[4] = 0;		c[5] = co - 1;		r[4] = 0;		r[5] = ro - 1;
	while (1)
	{
		c[0]--;	if (c[0] < c[4])	c[0] = c[4];
		c[1]++;	if (c[1] > c[5])	c[1] = c[5];
		r[0]--;	if (r[0] < r[4])	r[0] = r[4];
		r[1]++;	if (r[1] > r[5])	r[1] = r[5];
		if (FndEmp2(G, P, sp, ep, op, on, ce, d2, L, co, c, r))
			return;
		
		c[2] = c[0], c[3] = c[1], r[2] = r[0], r[3] = r[1];
		if (c[2] == c[4] && c[3] == c[5] &&
			r[2] == r[4] && r[3] == r[5])
			break;
	}
}
// 判断边与X分隔线相交
inline TCS_S32 LxSide2(TCS_F64 X, TCS_F64 *P, TIN_S *S)
{
	TCS_F64 *sp = P + 2 * S->sid[0];
	TCS_F64 *ep = P + 2 * S->eid[0];

	if ((sp[0] <= X && X <= ep[0]) ||
		(sp[0] >= X && X >= ep[0]))
		return 1;

	return 0;
}
// 判断边与X分隔线重合
inline TCS_S32 SxSide2(TCS_F64 X, TCS_F64 *P, TIN_S *S)
{
	TCS_F64 *sp = P + 2 * S->sid[0];
	TCS_F64 *ep = P + 2 * S->eid[0];

	if (sp[0] == X && ep[0] == X)
		return 1;

	return 0;
}
// 判断边与Y分隔线相交
inline TCS_S32 LySide2(TCS_F64 Y, TCS_F64 *P, TIN_S *S)
{
	TCS_F64 *sp = P + 2 * S->sid[0];
	TCS_F64 *ep = P + 2 * S->eid[0];

	if ((sp[1] <= Y && Y <= ep[1]) ||
		(sp[1] >= Y && Y >= ep[1]))
		return 1;

	return 0;
}
// 判断边与Y分隔线重合
inline TCS_S32 SySide2(TCS_F64 Y, TCS_F64 *P, TIN_S *S)
{
	TCS_F64 *sp = P + 2 * S->sid[0];
	TCS_F64 *ep = P + 2 * S->eid[0];

	if (sp[1] == Y && ep[1] == Y)
		return 1;