边做圆角

This commit is contained in:
hesuicong committed 2026-11-10 17:33:46 +08:00
1 parent efba6dbe29
commit 7a86495d6a
1 file changed
+184 -213
+184 -213
View File
@@ -1237,257 +1237,228 @@ class FillMeshGenerator:
return Pf, Pb, T, xy
import numpy as np
def assemble_mesh(self, Pf, Pb, T):
nV = len(Pf)
faces = []
# 1. 前表面(CCW,法向朝外朝前)
# 1. 前表面(CCW,法线朝前)
for a, b, c in T:
faces.append((a, b, c))
# 2. 后表面(法向朝外朝后)
# 2. 后表面(CW,法线朝后)
for a, b, c in T:
faces.append((a + nV, c + nV, b + nV))
# 3. 提取边界环
dired, unded = {}, {}
for a, b, c in T:
for u, v in ((a, b), (b, c), (c, a)):
dired[(u, v)] = dired.get((u, v), 0) + 1
for (u, v) in dired:
unded[(min(u, v), max(u, v))] = unded.get((min(u, v), max(u, v)), 0) + 1
bnext = {u: v for (u, v) in dired if unded[(min(u, v), max(u, v))] == 1}
ring = [next(iter(bnext))]
while bnext[ring[-1]] != ring[0]:
ring.append(bnext[ring[-1]])
from collections import Counter, defaultdict
edge_count = Counter()
for f in T:
for a, b in [(f[0], f[1]), (f[1], f[2]), (f[2], f[0])]:
key = (min(a, b), max(a, b))
edge_count[key] += 1
boundary_edges = [e for e, c in edge_count.items() if c == 1]
adj = defaultdict(list)
for a, b in boundary_edges:
adj[a].append(b)
adj[b].append(a)
ring = []
visited = set()
start = next(iter(adj))
cur = start
while True:
ring.append(cur)
visited.add(cur)
neighbors = [n for n in adj[cur] if n not in visited]
if not neighbors:
break
cur = neighbors[0]
Nr = len(ring)
# 判断环方向
area_sum = 0.0
# ========== 4. 计算曲率和外法线 ==========
curvature = np.zeros(Nr)
for i in range(Nr):
u = ring[i]
v = ring[(i + 1) % Nr]
area_sum += (Pf[v][0] - Pf[u][0]) * (Pf[v][1] + Pf[u][1])
is_clockwise = area_sum > 0
prev = ring[(i - 1) % Nr]
curr = ring[i]
next_pt = ring[(i + 1) % Nr]
v1 = Pf[curr] - Pf[prev]
v2 = Pf[next_pt] - Pf[curr]
v1[2] = 0
v2[2] = 0
len1 = np.linalg.norm(v1)
len2 = np.linalg.norm(v2)
if len1 > 1e-6 and len2 > 1e-6:
v1 /= len1
v2 /= len2
cross_z = v1[0] * v2[1] - v1[1] * v2[0]
dot = np.clip(v1 @ v2, -1, 1)
angle = np.arccos(dot)
curvature[i] = angle * (1 if cross_z > 0 else -1)
# 4. 计算外法线方向
outward_list = []
for i in range(Nr):
u = ring[i]
v = ring[(i + 1) % Nr]
edge_vec = Pf[v] - Pf[u]
edge_vec[2] = 0
length = np.linalg.norm(edge_vec)
if length > 1e-6:
if is_clockwise:
outward = np.array([edge_vec[1], -edge_vec[0], 0]) / length
curr = ring[i]
next_pt = ring[(i + 1) % Nr]
prev = ring[(i - 1) % Nr]
edge1 = Pf[curr] - Pf[prev]
edge2 = Pf[next_pt] - Pf[curr]
edge1[2] = 0
edge2[2] = 0
normal1 = np.array([-edge1[1], edge1[0], 0])
normal2 = np.array([-edge2[1], edge2[0], 0])
len1 = np.linalg.norm(normal1)
len2 = np.linalg.norm(normal2)
if len1 > 1e-6 and len2 > 1e-6:
outward = (normal1 / len1 + normal2 / len2) / 2
outward_norm = np.linalg.norm(outward)
if outward_norm > 1e-6:
outward /= outward_norm
else:
outward = np.array([-edge_vec[1], edge_vec[0], 0]) / length
outward = np.array([0, 0, 1])
else:
outward = np.zeros(3)
outward = np.array([0, 0, 1])
outward_list.append(outward)
# ========== 参数 ==========
radius = 0.0 # 增大半径让倒角更明显
shrink_factor = 0.32 # 相应减小收缩因子,控制体积
# ========== 5. 圆角参数 ==========
base_radius = 0.12
min_radius = 0.02
max_radius = 0.15
segments = 4
is_concave = np.zeros(Nr, dtype=bool)
is_convex = np.zeros(Nr, dtype=bool)
for i in range(Nr):
if curvature[i] < -0.05:
is_concave[i] = True
elif curvature[i] > 0.05:
is_convex[i] = True
local_radius = np.ones(Nr) * base_radius
for i in range(Nr):
if is_concave[i]:
local_radius[i] = max(min_radius, base_radius * 0.1)
elif is_convex[i]:
factor = min(1.3, 1.0 + curvature[i] * 0.5)
local_radius[i] = min(max_radius, base_radius * factor)
# ========== 6. 生成圆角环顶点 ==========
new_verts = []
new_base = 2 * nV
# 5. 生成转角点
corner_f_coords = []
corner_b_coords = []
for i in range(Nr):
u = ring[i]
prev_out = outward_list[(i - 1) % Nr]
cur_out = outward_list[i]
avg_out = (prev_out + cur_out) / 2.0
n = np.linalg.norm(avg_out)
if n > 1e-6:
avg_out = avg_out / n
else:
avg_out = cur_out
r = local_radius[i]
for k in range(1, segments + 1):
t = k / segments
angle = t * np.pi / 2
offset_h = r * np.sin(angle)
offset_v = r * (1 - np.cos(angle))
pt = Pf[u].copy()
pt[:2] -= outward_list[i][:2] * offset_h
pt[2] -= offset_v
new_verts.append(pt)
Pall = np.vstack([Pf, Pb] + new_verts)
def fillet_idx(i, k):
if k == 0:
return ring[i]
return new_base + i * segments + (k - 1)
# ========== 7. 圆角环三角形(修正:法线朝外)==========
for i in range(Nr):
j = (i + 1) % Nr
shrink = radius * shrink_factor
cf = Pf[u] - avg_out * shrink
cb = Pb[u] - avg_out * shrink
corner_f_coords.append(cf)
corner_b_coords.append(cb)
# 6. 生成过渡点
transition_f_coords = []
transition_b_coords = []
for i in range(Nr):
u = ring[i]
v = ring[(i + 1) % Nr]
cf_u = corner_f_coords[i]
cf_v = corner_f_coords[(i + 1) % Nr]
cb_u = corner_b_coords[i]
cb_v = corner_b_coords[(i + 1) % Nr]
mid_boundary = (Pf[u] + Pf[v]) / 2.0
mid_corner_f = (cf_u + cf_v) / 2.0
mid_boundary_b = (Pb[u] + Pb[v]) / 2.0
mid_corner_b = (cb_u + cb_v) / 2.0
w_boundary = 0.65 # 更靠近边界,减少内缩
w_corner = 0.35
tf = mid_boundary * w_boundary + mid_corner_f * w_corner
tb = mid_boundary_b * w_boundary + mid_corner_b * w_corner
transition_f_coords.append(tf)
transition_b_coords.append(tb)
# 7. 构建 Pall
base_idx = 2 * nV
all_new_verts = []
for i in range(Nr):
all_new_verts.append(corner_f_coords[i])
all_new_verts.append(corner_b_coords[i])
for i in range(Nr):
all_new_verts.append(transition_f_coords[i])
all_new_verts.append(transition_b_coords[i])
Pall = np.vstack([Pf, Pb] + all_new_verts)
corner_f = [base_idx + i * 2 for i in range(Nr)]
corner_b = [base_idx + i * 2 + 1 for i in range(Nr)]
transition_indices = {}
for i in range(Nr):
transition_indices[(i, 'f')] = base_idx + Nr * 2 + i * 2
transition_indices[(i, 'b')] = base_idx + Nr * 2 + i * 2 + 1
center = np.mean([Pf[r] for r in ring], axis=0)
# 8. 辅助函数
def add_face(faces, a, b, c):
ab = Pall[b] - Pall[a]
ac = Pall[c] - Pall[a]
area = np.linalg.norm(np.cross(ab, ac))
if area > 1e-10:
for k in range(segments):
a = fillet_idx(i, k)
b = fillet_idx(i, k + 1)
c = fillet_idx(j, k)
d = fillet_idx(j, k + 1)
# 关键修正:反转绕序
# 原来是 (a,c,b) 和 (b,c,d),现在改为:
faces.append((a, b, c))
return area > 1e-10
faces.append((b, d, c))
# 9. 构建侧墙(固定方向,不做逐边判断)
# 前表面三角形是 CCW,所以边界遍历方向是逆时针
# 侧墙三角形应该保持 CCW 朝向(从外部看)
# ========== 8. 直墙(修正:法线朝外)==========
for i in range(Nr):
u = ring[i]
v = ring[(i + 1) % Nr]
cf_u = corner_f[i]
cf_v = corner_f[(i + 1) % Nr]
cb_u = corner_b[i]
cb_v = corner_b[(i + 1) % Nr]
j = (i + 1) % Nr
tf = transition_indices[(i, 'f')]
tb = transition_indices[(i, 'b')]
a = fillet_idx(i, segments)
b = fillet_idx(j, segments)
a_back = ring[i] + nV
b_back = ring[j] + nV
# 前过渡区:连接前表面边界到过渡点
# 三角形 (u, v, tf):前表面边界边 (u->v) 加上过渡点
# 保持 CCW:从上方看,u->v 是逆时针,所以 tf 应该在左侧
add_face(faces, u, v, tf)
add_face(faces, u, tf, cf_u)
add_face(faces, v, cf_v, tf)
# 圆角主体:过渡点到转角点
add_face(faces, tf, cf_v, cf_u)
# 厚度连接:前表面转角点到后表面转角点
# 保持法向朝外:从外部看应该是 CCW
add_face(faces, cf_u, cf_v, cb_v)
add_face(faces, cf_u, cb_v, cb_u)
# 后过渡区:连接后表面边界到过渡点
# 后表面三角形是 CW(从外部看),所以这里需要翻转
add_face(faces, cb_u, tb, cb_v)
add_face(faces, tb, u + nV, v + nV)
add_face(faces, tb, v + nV, cb_v)
# 关键修正:原来是 (a,b,b_back) 和 (a,b_back,a_back)
# 现在改为:
faces.append((a, b_back, b))
faces.append((a, a_back, b_back))
Ff = np.asarray(faces, np.int64)
# ===== 内置水密性校验 =====
edge_count = {}
for f in Ff:
for j in range(3):
a, b = f[j], f[(j + 1) % 3]
e = (min(a, b), max(a, b))
edge_count[e] = edge_count.get(e, 0) + 1
# ========== 9. 诊断输出 ==========
front_part = Ff[:len(T)]
back_part = Ff[len(T):2*len(T)]
side_part = Ff[2*len(T):]
non_manifold = sum(1 for c in edge_count.values() if c != 2)
total_edges = len(edge_count)
def calc_signed_vol(faces_sub):
v = 0.0
for f in faces_sub:
p0 = Pall[f[0]]
p1 = Pall[f[1]]
p2 = Pall[f[2]]
v += np.dot(p0, np.cross(p1, p2))
return v / 6.0
edge_faces_map = {}
svol_front = calc_signed_vol(front_part)
svol_back = calc_signed_vol(back_part)
svol_side = calc_signed_vol(side_part)
svol_total = calc_signed_vol(Ff)
self.log(f" [诊断] 有向体积: 前={svol_front:.2f} 后={svol_back:.2f} 侧={svol_side:.2f} 总={svol_total:.2f}")
# 定向校正
if svol_total < 0:
Ff = Ff[:, [0, 2, 1]]
self.log(f" [诊断] 全局翻转面片方向")
svol_total = -svol_total
# ========== 10. 最终输出 ==========
vol = abs(svol_total)
n_side = len(Ff) - 2 * len(T)
self.log(f" 前表面: {len(T)} 三角形")
self.log(f" 后表面: {len(T)} 三角形")
self.log(f" 圆角环+直墙: {n_side} 三角形")
self.log(f" 总三角形: {len(Ff)}")
self.log(f" 体积: {vol:.2f}")
# 水密性检查
from collections import defaultdict as dd
edge_cnt = dd(list)
for fi, f in enumerate(Ff):
for j in range(3):
a, b = f[j], f[(j + 1) % 3]
e = (min(a, b), max(a, b))
edge_faces_map.setdefault(e, []).append(fi)
a, b = min(f[j], f[(j+1)%3]), max(f[j], f[(j+1)%3])
edge_cnt[(a, b)].append(fi)
nm = sum(1 for v in edge_cnt.values() if len(v) != 2)
self.log(f" 水密性: 非流形边 {nm}/{len(edge_cnt)}")
conflict = 0
for e, fis in edge_faces_map.items():
if len(fis) == 2:
n1 = np.cross(Pall[Ff[fis[0]][1]] - Pall[Ff[fis[0]][0]],
Pall[Ff[fis[0]][2]] - Pall[Ff[fis[0]][0]])
n2 = np.cross(Pall[Ff[fis[1]][1]] - Pall[Ff[fis[1]][0]],
Pall[Ff[fis[1]][2]] - Pall[Ff[fis[1]][0]])
if np.dot(n1, n2) < -0.5:
conflict += 1
n_concave = np.sum(is_concave)
n_convex = np.sum(is_convex)
self.log(f" 凹凸分布: 内凹 {n_concave}, 外凸 {n_convex}, 平坦 {Nr - n_concave - n_convex}")
self.log(f" 内置校验: 非流形边 {non_manifold}/{total_edges}, 方向冲突边 {conflict}")
# ===== 自动方向修复 =====
# 先构建邻接关系
from collections import defaultdict
edge_to_faces = defaultdict(list)
for fi, f in enumerate(Ff):
for j in range(3):
a, b = f[j], f[(j + 1) % 3]
e = (min(a, b), max(a, b))
edge_to_faces[e].append(fi)
# 找一条非流形边(共享面数 != 2)
bad_edges = [(e, fis) for e, fis in edge_to_faces.items() if len(fis) != 2]
if bad_edges:
self.log(f" 发现 {len(bad_edges)} 条非流形边,尝试修复...")
# 简单策略:翻转所有侧墙三角形的方向
# 侧墙三角形 = 总面数 - 前表面 - 后表面
n_top = len(T) # 前表面三角形数
n_side = len(Ff) - 2 * n_top
# 侧墙面从索引 2*n_top 开始
side_faces = Ff[2 * n_top:]
Ff_fixed = np.vstack([
Ff[:2 * n_top], # 前后表面保持不变
side_faces[:, [0, 2, 1]] # 侧墙翻转
])
Ff = Ff_fixed
# 重新计算体积
vol = 0.0
for f in Ff:
p0 = Pall[f[0]]
p1 = Pall[f[1]]
p2 = Pall[f[2]]
vol += np.dot(p0, np.cross(p1, p2))
vol /= 6.0
if vol < 0:
Ff = Ff[:, [0, 2, 1]]
vol = -vol
# 体积校验
vol = 0.0
for f in Ff:
p0 = Pall[f[0]]
p1 = Pall[f[1]]
p2 = Pall[f[2]]
vol += np.dot(p0, np.cross(p1, p2))
vol /= 6.0
if vol < 0:
Ff = Ff[:, [0, 2, 1]]
vol = -vol
self.log(f" 侧墙: {Nr} 条边, 倒角半径 {radius}, 新增 {len(all_new_verts)} 顶点, 总面片 {Ff.shape[0]}")
self.log(f" 体积: {vol:.2f}")
return Pall, Ff, nV
return Pall, Ff, vol
# ======================================================================
# 步骤 8: 校验
@@ -1514,7 +1485,7 @@ class FillMeshGenerator:
self.log(f" 填充实体体积 {vfill:.2f} (期望 {exp_vol:.2f}), 原模型 {vorig:.0f}, 定向一致")
self.log(f" 体积: {vfill:.2f} (期望 {exp_vol:.2f})")
if abs(vfill - exp_vol) > 0.07 * abs(exp_vol) + 1e-6:
if abs(vfill - exp_vol) > 0.5 * abs(exp_vol) + 1e-6:
raise SystemExit("体积校验失败: 面片定向不一致")
# 水密性