Area of Triangles and Polygons (2D)

本文介绍了几种计算二维图形面积的方法,包括三角形和多边形的古代和现代算法。对于三角形,给出了直接公式和向量计算方式;对于多边形,则通过分解为多个三角形并求和的方式给出解决方案。

Area of Triangles and Polygons (2D)

Triangles

Ancient ways

S = b*h/2

   = abSin(θ)/2

   = sqrt[s(s-a)(s-b)(s-c)]

s = (a+b+c)/2

Modern ways

S = |(V1-V0)×(V2-V0)|/2

   = [(x1-x0)(y2-y0)-(x2-x0)(y1-y0)]/2

Polygons Area - 2D

let polygon(Ω) is defined by its vertices Vi=(xi,yi), with V0=Vn. Let P be any point, and for each edge ViVi+1of Ω, form the triangle i=△PViVi+1, then:

S(Ω) = sum{i=0:n-1 | S(i)} and △i=△PViVi+1      (1)

Notice that, for a counterclockwise oriented polygon, use the modern way to calculate the triangle △i area, △PV2V3 and △PVn-1Vn have positive value and △PV0V1 and △PV1V2 have negative value, so the above rule (1) is right.

Let P = (0,0), the area of polygon is:

2S(Ω) = sum{i=0:n-1 | (xiyi+1-xi+1yi)}

C++ Algorithms

point_position

    int point_position(line L, point P)

        // if P is in the left of L , return value > 0

        // if P is in the right of L, return value < 0

        // if P is on the L, reutrn value = 0

int point_position(line L, point P)
{
    return (int)((L.x2-L.x1)*(P.y-L.y1) - (P.x-L.x1)*(L.y2-L.y1));
}

area_triangle(point A, point B, point C)

        double area_triangle(point A, point B, point C)
        {
            return abs((B.x - A.x)*(C.y - A.y) - (C.x - A.x)*(B.y - A.y))/2.0;
        }

area_polygon(vector<point> Poly)

        double area_polygon(vector<point> Poly)
        {
            int i;
            double ret = 0;
            for(i=0;i<Poly.size()-1;i++){
                ret += Poly[i].x*Poly[i+1].y-Poly[i+1].x*Poly[i].y;
            }
            int n = Poly.size()-1;
            ret += Poly[n].x*Poly[0].y-Poly[0].x*Poly[n].y;
            return ret/2.0;
        }

代码有未定义的get_vertices:import math from collections import defaultdict from pathlib import Path def process_folder(input_dir, output_dir): """保持原有文件夹处理逻辑不变""" Path(output_dir).mkdir(parents=True, exist_ok=True) valid_extensions = {'.txt', '.obj', '.geo'} for input_path in Path(input_dir).rglob('*'): if input_path.is_dir() or input_path.suffix.lower() not in valid_extensions: continue relative_path = input_path.relative_to(input_dir) output_path = Path(output_dir) / relative_path.with_name(f"{input_path.stem}_mesh.obj") output_path.parent.mkdir(parents=True, exist_ok=True) try: process_single_file(input_path, output_path) print(f"成功处理:{input_path} -> {output_path}") except Exception as e: print(f"处理失败:{input_path}\n错误详情:{str(e)}\n{'='*50}") def process_single_file(input_path, output_path): """保持原有单文件处理逻辑不变""" with open(input_path, 'r') as f: raw_lines = [line.strip() for line in f if line.strip()] vertices, lines = parse_input_data(raw_lines) obj_lines = enhanced_convert_lines_to_mesh(vertices, lines) with open(output_path, 'w') as f: f.write("\n".join(obj_lines)) def parse_input_data(raw_lines): """保持原有数据解析逻辑不变""" vertices = [] lines = [] for line in raw_lines: if line.startswith('#') or not line.split(): continue parts = line.split() try: if parts[0].lower() == 'v': vertices.append(' '.join(parts[1:4])) elif parts[0].lower() == 'l': start = int(parts[1]) - 1 end = int(parts[2]) - 1 lines.append((start, end)) except (IndexError, ValueError) as e: raise ValueError(f"解析错误:{line}") from e return vertices, lines def enhanced_convert_lines_to_mesh(vertex_data, line_data): """增强型网格生成算法""" # 构建邻接表并缓存坐标 adj, coords = build_adjacency_with_coords(line_data, vertex_data) # 第一步:查找严格闭合的平面多边形 planar_faces = find_planar_faces(adj, coords) # 第二步:处理非闭合结构 edge_faces = process_open_structures(adj, coords) # 合并结果并去重 all_faces = deduplicate_faces(planar_faces + edge_faces) return generate_obj_content(vertex_data, all_faces) def build_adjacency_with_coords(line_data, vertex_data): """构建邻接表并缓存坐标""" adj = defaultdict(set) coord_cache = {} for idx in range(len(vertex_data)): parts = vertex_data[idx].split() coord_cache[idx] = tuple(map(float, parts)) for s, e in line_data: adj[s].add(e) adj[e].add(s) return adj, coord_cache def find_planar_faces(adj, coord_cache, max_length=6): """查找平面多边形""" def dfs(current, path, edges): cycles = [] if len(path) > 2 and current == path[0]: cycle = path[:-1] # 去除重复起点 if validate_planar_cycle(cycle, coord_cache): return [cycle] return [] if len(path) >= max_length: return [] for neighbor in adj[current]: edge = tuple(sorted((current, neighbor))) if edge not in edges: new_edges = edges.copy() new_edges.add(edge) cycles += dfs(neighbor, path + [neighbor], new_edges) return cycles seen = set() planar_faces = [] for start in adj: for cycle in dfs(start, [start], set()): normalized = normalize_cycle(cycle) if normalized not in seen: seen.add(normalized) planar_faces.extend(triangulate_planar_face(normalized, coord_cache)) return planar_faces def validate_planar_cycle(cycle, coord_cache, tolerance=1e-3): """验证多边形平面性(放宽容差)""" if len(cycle) < 3: return False # 计算平均法线 normals = [] for i in range(len(cycle)): a = coord_cache[cycle[i-2]] b = coord_cache[cycle[i-1]] c = coord_cache[cycle[i]] normals.append(calculate_normal(a, b, c)) # 检查法线一致性 ref_normal = normals[0] for n in normals[1:]: if angle_between(ref_normal, n) > math.radians(5): # 允许5度偏差 return False return True def triangulate_planar_face(cycle, coord_cache): """平面多边形三角剖分""" try: # 将三维坐标投影到二维平面 points_3d = [coord_cache[v] for v in cycle] normal = calculate_normal(*points_3d[:3]) basis = create_projection_basis(normal) points_2d = [project_point(p, basis) for p in points_3d] # 执行二维耳切法 triangles = ear_clipping_2d(points_2d) return [(cycle[i], cycle[j], cycle[k]) for i, j, k in triangles] except: return [] def process_open_structures(adj, coord_cache): """处理开放式结构(线段连接)""" faces = [] visited = set() # 查找所有三角形结构 for a in adj: for b in adj[a]: if b in visited: continue for c in adj[b]: if c in adj[a] and c != a and c != b: # 检查是否形成有效三角形 if validate_triangle(a, b, c, coord_cache): faces.append((a, b, c)) visited.update({a, b, c}) return faces def validate_triangle(a, b, c, coord_cache): """验证三角形有效性""" try: p1 = coord_cache[a] p2 = coord_cache[b] p3 = coord_cache[c] # 检查面积 area = 0.5 * math.sqrt(sum(x**2 for x in calculate_normal(p1, p2, p3))) return area > 1e-6 except: return False # ---------- 几何计算辅助函数 ---------- def calculate_normal(a, b, c): """计算法线向量""" v1 = [b[i]-a[i] for i in range(3)] v2 = [c[i]-a[i] for i in range(3)] return cross_product(v1, v2) def cross_product(v1, v2): return [ v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0] ] def angle_between(v1, v2): """计算两向量夹角""" dot = sum(x*y for x, y in zip(v1, v2)) mag1 = math.sqrt(sum(x**2 for x in v1)) mag2 = math.sqrt(sum(y**2 for y in v2)) return math.acos(dot/(mag1*mag2 + 1e-9)) def create_projection_basis(normal): """创建投影坐标系""" axis_x = [1, 0, 0] if abs(normal[1]) > 0.5 else [0, 1, 0] tangent = cross_product(normal, axis_x) binormal = cross_product(normal, tangent) return (tangent, binormal) def project_point(p, basis): """三维坐标投影到二维平面""" return ( p[0]*basis[0][0] + p[1]*basis[0][1] + p[2]*basis[0][2], p[0]*basis[1][0] + p[1]*basis[1][1] + p[2]*basis[1][2] ) def ear_clipping_2d(points): """二维耳切法实现""" indices = list(range(len(points))) triangles = [] while len(indices) > 3: ear_found = False for i in range(len(indices)): a, b, c = get_vertices(indices, i) if is_ear(a, b, c, points, indices): triangles.append((a, b, c)) del indices[i] ear_found = True break if not ear_found: break if len(indices) == 3: triangles.append(tuple(indices)) return triangles def is_ear(a, b, c, points, indices): """检查是否为耳朵""" if not is_convex(points[a], points[b], points[c]): return False return not any(point_in_triangle(points[i], points[a], points[b], points[c]) for i in indices if i not in {a, b, c}) def point_in_triangle(p, a, b, c, tol=1e-6): """精确点包含检查""" area = triangle_area(a, b, c) return abs(area - (triangle_area(a, b, p) + triangle_area(b, c, p) + triangle_area(c, a, p))) < tol def triangle_area(a, b, c): """计算三角形面积""" return 0.5 * abs( (b[0]-a[0])*(c[1]-a[1]) - (c[0]-a[0])*(b[1]-a[1]) ) # ---------- 辅助函数 ---------- def normalize_cycle(cycle): """优化的环标准化""" min_val = min(cycle) rotations = [cycle[i:] + cycle[:i] for i, v in enumerate(cycle) if v == min_val] reversed_rot = [list(reversed(r)) for r in rotations] all_candidates = rotations + reversed_rot return min(all_candidates) def deduplicate_faces(faces): """面片去重""" seen = set() unique = [] for f in faces: key = tuple(sorted(f)) if key not in seen and len(set(key)) == 3: seen.add(key) unique.append(f) return unique def generate_obj_content(vertex_data, faces): """生成OBJ文件内容""" output = ["# Generated Mesh"] output += [f"v {v}" for v in vertex_data] output.append("\n# Faces") output += [f"f {' '.join(str(v+1) for v in f)}" for f in faces] return output if __name__ == "__main__": INPUT_DIR = "D:/idmdownload/Compressed/Tallinn/train/wireframe" OUTPUT_DIR = "D:/idmdownload/Compressed/Tallinn/train/mesh" if not Path(INPUT_DIR).exists(): raise FileNotFoundError(f"输入目录不存在:{INPUT_DIR}") if INPUT_DIR == OUTPUT_DIR: raise ValueError("输入输出目录不能相同") try: process_folder(INPUT_DIR, OUTPUT_DIR) print(f"\n处理完成!输出文件数:{len(list(Path(OUTPUT_DIR).rglob('*.obj')))}") except KeyboardInterrupt: print("\n操作已中断")
最新发布
05-27
评论
成就一亿技术人!
拼手气红包6.0元
还能输入1000个字符
 
红包 添加红包
表情包 插入表情
 条评论被折叠 查看
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值