⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 tincore1.cpp

📁 并行TIN生成算法, 基于DeWall算法理论实现
💻 CPP
📖 第 1 页 / 共 4 页
字号:
			c[1]--;	if (c[1] < c[5])	c[1] = c[3] = c[5];	// 说明前一列已经到达搜索边界
			c[2]++;	if (c[2] > c[6])	c[2] = c[4]	= c[6];
			r[1]--;	if (r[1] < r[5])	r[1] = r[3] = r[5];	// 说明前一行已经到达搜索边界
			r[2]++;	if (r[2] > r[6])	r[2] = r[4] = r[6];
			SearchCellNPD(cidx, ppts, ddim, pnt, pid[1], cid[1], dnp, col, c, r);

			c[3] = c[1], c[4] = c[2], r[3] = r[1], r[4] = r[2];
			if (c[3] == c[5] && c[4] == c[6] &&
				r[3] == r[5] && r[4] == r[6])
				break;
		}
	}
	// 检测停止标志
	if (!(*stop))
		return 0;
	// 3.采用空圆法则搜索点pid[2]
	{
		TSS_FLT64 *pt1 = ppts + ddim * pid[0];
		TSS_FLT64 *pt2 = ppts + ddim * pid[1];
		TSS_FLT64  cen[2] = {(pt1[0] + pt2[0]) / 2, (pt1[1] + pt2[1]) / 2};

		// 计算初始搜索范围
		TSS_FLT64  rad = sqrt(EuLen2D(pt1, cen));
		CalcSearchBound(bnd, spn, col, row, cen, rad, c[1], c[2], r[1], r[2]);

		for (i = r[1]; i <= r[2]; i++)
		{
			cin = i * col;
			for (k = c[1]; k <= c[2]; k++)
			{
				cell = cidx + (cin + k);
				if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin + k, pid[2], cid[2], cen))
					return 1;
			}
		}

		// 继续搜索
		c[3] = c[1];	c[4] = c[2];		r[3] = r[1];	r[4] = r[2];
		c[5] = 0;		c[6] = col - 1;		r[5] = 0;		r[6] = row - 1;
		while (1)
		{
			c[1]--;	if (c[1] < c[5])	c[1] = c[5];
			c[2]++;	if (c[2] > c[6])	c[2] = c[6];
			r[1]--;	if (r[1] < r[5])	r[1] = r[5];
			r[2]++;	if (r[2] > r[6])	r[2] = r[6];
			if (SearchCellEMR(uidx, ppts, ddim, pt1, pt2, pid[2], cid[2], cen, col, c, r))
				return 1;	// 已经找到
			
			c[3] = c[1], c[4] = c[2], r[3] = r[1], r[4] = r[2];
			if (c[3] == c[5] && c[4] == c[6] &&
				r[3] == r[5] && r[4] == r[6])
				return 0;
		}
	}

	return 0;
}

// 单元搜索第3点(空圆法则)
inline long SearchCellEMR(IDX_GRID *uidx, IDX_CELL *cell, TSS_FLT64 *ppts, TSS_SNT32 ddim,
						  TSS_FLT64 *pt1, TSS_FLT64 *pt2, TSS_SNT32 cin, TSS_SNT32 &pid,
						  TSS_SNT32 &cid, TSS_FLT64 *cen)
{
	TSS_SNT32 n, nid;
	TSS_FLT64 *pt3, rad;

	for (n = 0; n < cell->cnt; n++)
	{
		nid = cell->pid[n];
		pt3 = ppts + ddim * nid;
//		if (pt3 == pt1 || pt3 == pt2)
//			continue;

		if (!IsLine(pt1, pt2, pt3))	// 不共线
		{
			CalcCircum(pt1, pt2, pt3, cen);		rad = EuLen2D(pt1, cen);	pid = nid;	cid = cin;
			SearchCellEMR(uidx, ppts, ddim, pt1, pt2, pt3, rad, pid, cid, cen);
			return 1;
		}
	}

	return 0;
}

// 检查圆内是否存在点
inline void SearchCellEMR(IDX_GRID *uidx, TSS_FLT64 *ppts, TSS_SNT32 ddim, TSS_FLT64 *pt1,
						  TSS_FLT64 *pt2, TSS_FLT64 *pt3, TSS_FLT64 rad, TSS_SNT32 &pid,
						  TSS_SNT32 &cid, TSS_FLT64 *cen)
{
	TSS_FLT64 *pt4;
	IDX_CELL *cell, *cidx = uidx->idx;
	TSS_SNT32 i, k, n, c[4], r[4], tmp, nid;
	TSS_SNT32 col = uidx->num[0], row = uidx->num[1];

	// 计算搜索边界
	CalcSearchBound(uidx->bnd, uidx->spn, col, row, cen, sqrt(rad), c[1], c[2], r[1], r[2]);

	// 递归搜索
	for (i = r[1]; i <= r[2]; i++)
	{
		tmp = i * col; 
		for (k = c[1]; k <= c[2]; k++)
		{
			cell = cidx + (tmp + k);
			for (n = 0; n < cell->cnt; n++)
			{
				nid = cell->pid[n];
				pt4 = ppts + ddim * nid;
				if (pt4 == pt3 || pt4 == pt2)// || pt4 == pt1)
					continue;

				if (EuLen2D(pt4, cen) < rad) // 有点
				{
					CalcCircum(pt1, pt2, pt4, cen);		rad = EuLen2D(pt1, cen);	pid = nid;	cid = tmp + k;
					SearchCellEMR(uidx, ppts, ddim, pt1, pt2, pt4, rad, pid, cid, cen);
					return;
				}
			}
		}
	}
}

// 批量单元搜索
inline long SearchCellEMR(IDX_GRID *uidx, TSS_FLT64 *ppts, TSS_SNT32 ddim, TSS_FLT64 *pt1,
						  TSS_FLT64 *pt2, TSS_SNT32 &pid, TSS_SNT32 &cid, TSS_FLT64 *cen,
						  TSS_SNT32 col, TSS_SNT32 *c, TSS_SNT32 *r)
{
	TSS_SNT32 i, k, cin;
	IDX_CELL *cell, *cidx = uidx->idx;

	if (r[3] != r[1])	// 即r[3] != r[5]
	{
		for (k = c[1] + 1; k < c[2]; k++)
		{
			cin = r[1] * col + k;		cell = cidx + cin;
			if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
				return 1;
		}
	}
	if (r[4] != r[2])	// 即r[4] != r[6]
	{
		for (k = c[1] + 1; k < c[2]; k++)
		{
			cin = r[2] * col + k;		cell = cidx + cin;
			if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
				return 1;
		}
	}
	if (c[3] != c[1])	// 即c[3] != c[5]
	{
		for (i = r[1] + 1; i < r[2]; i++)
		{
			cin = i * col + c[1];		cell = cidx + cin;
			if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
				return 1;
		}
	}
	if (c[4] != c[2])	// 即c[4] != c[6]
	{
		for (i = r[1] + 1; i < r[2]; i++)
		{
			cin = i * col + c[2];		cell = cidx + cin;
			if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
				return 1;
		}
	}
	// 处理4个角单元
	if (c[3] != c[5] || r[3] != r[5])
	{
		cin = r[1] * col + c[1];	cell = cidx + cin;
		if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
			return 1;
	}
	if (c[3] != c[5] || r[4] != r[6])
	{
		cin = r[2] * col + c[1];	cell = cidx + cin;
		if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
			return 1;
	}
	if (c[4] != c[6] || r[3] != r[5])
	{
		cin = r[1] * col + c[2];	cell = cidx + cin;
		if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
			return 1;
	}
	if (c[4] != c[6] || r[4] != r[6])
	{
		cin = r[2] * col + c[2];	cell = cidx + cin;
		if (SearchCellEMR(uidx, cell, ppts, ddim, pt1, pt2, cin, pid, cid, cen))
			return 1;
	}

	return 0;
}

// 计算外接圆
inline void CalcCircum(TSS_FLT64 *p1, TSS_FLT64 *p2, TSS_FLT64 *p3, TSS_FLT64 *pt)
{
	// 11(+|-), 12(*), 2(/)
	TSS_FLT64 dx21 = p2[0] - p1[0];
	TSS_FLT64 dy21 = p2[1] - p1[1];
	TSS_FLT64 xy21 = dx21 * dx21 + dy21 * dy21;

	TSS_FLT64 dx31 = p3[0] - p1[0];
	TSS_FLT64 dy31 = p3[1] - p1[1];
	TSS_FLT64 xy31 = dx31 * dx31 + dy31 * dy31;

	TSS_FLT64 top1 =   dy21 * xy31 - dy31 * xy21;
	TSS_FLT64 top2 = - dx21 * xy31 + dx31 * xy21;
	TSS_FLT64 deno =   dy21 * dx31 - dy31 * dx21;

	pt[0] = p1[0] + 0.5 * top1 / deno;
	pt[1] = p1[1] + 0.5 * top2 / deno;
}

// 判断是否共线
inline long IsLine(TSS_FLT64 *p1, TSS_FLT64 *p2, TSS_FLT64 *p3)
{
	// 利用三角形面积是否等于0判断
	// (p1[0] * (p2[1] - p3[1]) + p2[0] * (p3[1] - p1[1]) + p3[0] * (p1[1] - p2[1]) : 3(*), 5(+|-)
	// (p1[0] - p2[0]) * (p1[1] - p3[1]) - (p1[0] - p3[0]) * (p1[1] - p2[1])		: 2(*), 5(-)
	if ((p1[0] - p2[0]) * (p1[1] - p3[1]) - (p1[0] - p3[0]) * (p1[1] - p2[1]) == 0.0)
		return 1;
	else
		return 0;
}

// 计算同侧符号
inline long IsSide(TSS_FLT64 *p1, TSS_FLT64 *p2, TSS_FLT64 *p3)
{
	TSS_FLT64 dx21 = p2[0] - p1[0];
//	TSS_FLT64 side = (p3[1] - p1[1]) * (p2[0] - p1[0]) - (p2[1] - p1[1]) * (p3[0] - p1[0]);
	if (dx21 == 0)	// p1, p2的x坐标相等
	{
		TSS_FLT64 side = p3[0] - p1[0];
		if		(side > 0)
			return  1;
		else if (side < 0)
			return -1;
	}
	else	// side / dx21
	{
		TSS_FLT64 side = (p3[1] - p1[1]) * (p2[0] - p1[0]) - (p2[1] - p1[1]) * (p3[0] - p1[0]);// / dx21;
		if		(side > 0)
			return  1;
		else if (side < 0)
			return -1;
	}

	return 0;
}

// 构建三角形网络
long MakeSimplex(IDX_GRID *uidx, TSS_FLT64 *ppts, TSS_SNT32 ddim, TSS_SNT32 *ppid, TSS_SNT32 *pnum,
				 TSS_SNT32 *tnum, TSS_SNT32 *stop, TSS_SNT32 *pid, TSS_SNT32 *cid, CTssList *ATL)
{
	CTssHash *AFH = new CTssHash(128, FACE_EQLH, FACE_HASH);	// 工作链表
	CTssHash *BFH = new CTssHash(128, FACE_EQLH, FACE_HASH);	// 跟踪链表

	TSS_SNT32 *num = uidx->num;
	IDX_CELL *cidx = uidx->idx;

	// 构建初始三角形墙
	TIN_FACE *face = new TIN_FACE;
	SetFace(face, pid[0], pid[1], pid[2], 0, cid[0], cid[1]);
	AFH->PutData(face, face);	ppid[pid[0]]++;	ppid[pid[1]]++;

	face = new TIN_FACE;
	SetFace(face, pid[2], pid[0], pid[1], 0, cid[2], cid[0]);
	AFH->PutData(face, face);	ppid[pid[0]]++;	ppid[pid[2]]++;

	face = new TIN_FACE;
	SetFace(face, pid[1], pid[2], pid[0], 0, cid[1], cid[2]);
	AFH->PutData(face, face);	ppid[pid[1]]++;	ppid[pid[2]]++;

	TIN_TGON *tgon = new TIN_TGON;
	tgon->pid[0] = pid[0];	tgon->pid[1] = pid[1];	tgon->pid[2] = pid[2];
	tgon->tid[0] = tgon->tid[1] = tgon->tid[2] = -1;	ATL->AddTail(tgon);

	// 循环构建
	TIN_FACE *temp;
	TSS_SNT32 sid, eid, mid;
	TSS_FLT64 *pt1, *pt2, *opt;

	*tnum = 1;	face = (TIN_FACE*)(AFH->RemoveHead());
	while (face)
	{
		temp = 0;
		// 搜索第三点
		pid[0] = sid = face->sid[0];	pid[1] = eid = face->eid[0];	cid[0] = face->sid[1];	cid[1] = face->eid[1];
		pt1 = ppts + ddim * sid;	pt2 = ppts + ddim * eid;	opt = ppts + ddim * face->pid;
		if (SearchCell(BFH, uidx, ppts, ddim, stop, pt1, pt2, opt, pid, cid[2]))
		{
			// 加入三角形
			tgon = new TIN_TGON;
			tgon->pid[0] = eid;	tgon->pid[1] = sid;	tgon->pid[2] = pid[2];
			tgon->tid[0] = face->tid;	tgon->tid[1] = tgon->tid[2] = -1;	ATL->AddTail(tgon);
			// 减少计数器
			ppid[sid]--;	ppid[eid]--;
			// 收集边对象
			BFH->PutData(face, face);

			// 更新链表
			face = new TIN_FACE;
			SetFace(face, sid, pid[2], eid, *tnum, cid[0], cid[2]);
			if (PutHash(AFH, BFH, cidx, ppid, face))
			{
				temp = face;
				face = 0;
			}

			if (!face)
				face = new TIN_FACE;
			SetFace(face, pid[2], eid, sid, *tnum, cid[2], cid[1]);
			if (PutHash(AFH, BFH, cidx, ppid, face))
			{
				temp = face;
				face = 0;
			}

			if (face)
				delete face;
			// 删除归零点
			if (ppid[sid] == 0)	{	DelPidEx(cidx + cid[0], sid);	*pnum += 1;	}
			if (ppid[eid] == 0)	{	DelPidEx(cidx + cid[1], eid);	*pnum += 1;	}
			mid = pid[2];
			if (ppid[mid] == 0)	{	DelPidEx(cidx + cid[2], mid);	*pnum += 1;	}

			*tnum += 1;
		}
		else
		{
			// 减少计数器
			ppid[sid]--;	ppid[eid]--;
			// 收集边对象
			BFH->PutData(face, face);
			// 删除归零点
			if (ppid[sid] == 0)	{	DelPidEx(cidx + cid[0], sid);	*pnum += 1;	}
			if (ppid[eid] == 0)	{	DelPidEx(cidx + cid[1], eid);	*pnum += 1;	}
		}

		// 检测停止标志
		if (!(*stop))
		{
			face = (TIN_FACE*)(AFH->RemoveHead());
			while (face)
			{
				delete face;
				face = (TIN_FACE*)(AFH->RemoveNext());
			}
			face = (TIN_FACE*)(BFH->RemoveHead());
			while (face)
			{
				delete face;
				face = (TIN_FACE*)(BFH->RemoveNext());
			}
			break;
		}

		if (temp)
		{
			face = temp;	AFH->DelData(face);
		}
		else
			face = (TIN_FACE*)(AFH->RemoveHead());
	}

	face = (TIN_FACE*)(BFH->RemoveHead());
	while (face)
	{
		delete face;
		face = (TIN_FACE*)(BFH->RemoveNext());
	}
	delete AFH;	delete BFH;

	return 1;
}

// 插入边
inline long PutHash(CTssHash *AFH, CTssHash *BFH, IDX_CELL *cidx, TSS_SNT32 *ppid, TIN_FACE *face)
{
	TSS_SNT32 sid = face->sid[0], eid = face->eid[0];
	TIN_FACE *temp = (TIN_FACE*)(AFH->DelData(face));
	if (temp)	// 已存在
	{
		// 减少计数器
		ppid[sid]--;	ppid[eid]--;
		// 收集边对象
		BFH->PutData(temp, temp);
		return 0;
	}
	else
	{
		// 增加计数器
		ppid[sid]++;	ppid[eid]++;
		// 增加处理边
		AFH->PutData(face, face);
		return 1;
	}

	return 0;
}

// 搜索第3点:空圆法则
inline long SearchCell(CTssHash *BFH, IDX_GRID *uidx, TSS_FLT64 *ppts, TSS_SNT32 ddim, TSS_SNT32 *stop,
					   TSS_FLT64 *pt1, TSS_FLT64 *pt2, TSS_FLT64 *opt, TSS_SNT32 *pid,
					   TSS_SNT32 &cid)
{
	TSS_SNT32 cin;
	TSS_SNT32 c[8], r[8], i, k;
	IDX_CELL *cell, *cidx = uidx->idx;
	TSS_SNT32 col = uidx->num[0], row = uidx->num[1];

	TSS_FLT64 cen[2] = {(pt1[0] + pt2[0]) / 2, (pt1[1] + pt2[1]) / 2};
	// 计算初始搜索范围
	CalcSearchBound(uidx->bnd, uidx->spn, col, row, cen, sqrt(EuLen2D(pt1, cen)), c[1], c[2], r[1], r[2]);
	// 计算搜索方向
	TSS_SNT32 orn = -IsSide(pt1, pt2, opt); // 反向

	pid[2] = -1;
	// 分析是否有必要采用外扩算法优化
	for (i = r[1]; i <= r[2]; i++)
	{
		cin = i * col;
		for (k = c[1]; k <= c[2]; k++)
		{
			cell = cidx + (cin + k);
			if (SearchCellEMR(BFH, uidx, cell, ppts, ddim, pt1, pt2, opt, orn, cin + k, pid, cid, cen))
				return 1;
		}
	}
	////////////////////////////////
	c[3] = c[1];	c[4] = c[2];		r[3] = r[1];	r[4] = r[2];
	c[5] = 0;		c[6] = col - 1;		r[5] = 0;		r[6] = row - 1;
	while (1)
	{
		c[1]--;	if (c[1] < c[5])	c[1] = c[5];
		c[2]++;	if (c[2] > c[6])	c[2] = c[6];
		r[1]--;	if (r[1] < r[5])	r[1] = r[5];
		r[2]++;	if (r[2] > r[6])	r[2] = r[6];
		if (SearchCellEMR(BFH, uidx, ppts, ddim, pt1, pt2, opt, orn, pid, cid, cen, col, c, r))
			return 1;
		
		c[3] = c[1], c[4] = c[2], r[3] = r[1], r[4] = r[2];
		if (c[3] == c[5] && c[4] == c[6] &&
			r[3] == r[5] && r[4] == r[6])
			return 0;

		// 检测停止标志
		if (!(*stop))
			return 0;
	}

	return 0;
}

// 要点:
// 1不形成用过2次的边
// 2不包含归零(闭合)的点
// 方向约束搜索:要求orn = -1 或 1
inline long SearchCellEMR(CTssHash *BFH, IDX_GRID *uidx, IDX_CELL *cell, TSS_FLT64 *ppts, TSS_SNT32 ddim,
						  TSS_FLT64 *pt1, TSS_FLT64 *pt2, TSS_FLT64 *opt, TSS_SNT32 orn,
						  TSS_SNT32 cin, TSS_SNT32 *pid, TSS_SNT32 &cid, TSS_FLT64 *cen)
{
	TIN_FACE fa, fc;
	TSS_FLT64 *pt3, rad;
	TSS_SNT32 n, cnt, *nid;

	cnt = cell->cnt;	nid = cell->pid;
	for (n = 0; n < cnt; n++)
	{
		pt3 = ppts + ddim * nid[n];
		// orn == IsSide(pt1, pt2, pt3)保证pt3不等于pt1或pt2
		if (orn == IsSide(pt1, pt2, pt3))	// 异侧
		{ 
			// 检查是否形成2次边
			SetFace(&fa, pid[0], nid[n]);
			if (BFH->GetData(&fa))
				continue;
			SetFace(&fc, pid[1], nid[n]);

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -