《Python机器人基础原理》第九十六篇离线仿真——装配产线可达评估工作空间点云与包围盒计算声明本篇全部为课程离线仿真代码不能直接用于真实机器人控制器规避工程误导风险。1. 实际工业应用场景周三下午汽车零部件装配线的方案评审会。工艺工程师小赵把新夹具的3D模型导入离线编程软件在产线布局图上标出了一排待装配的螺栓孔位。这些孔位都在一条斜线上间距80 mm。我需要确认咱们的机械臂能不能全部够到不用移动底座。旁边的资深工程师老郑凑过来看了一眼你光看CAD不够。机械臂的工作空间不是个球是各种奇异形状挤在一起的——关节限位、连杆干涉、底座安装高度都会把可达区域啃掉一大块。你得算一下实际工作空间的包围盒看看长、宽、高边界能不能覆盖这排孔。小赵怎么算蒙特卡洛随机撒关节角然后正运动学算末端坐标再把所有点聚起来看范围老郑点头对。2R平面臂简单但你得把关节限位考虑进去——J1只能转±170°J2只能俯仰±120°。如果不加限位撒出来的点会包含机械臂根本到不了的位置包围盒就虚大了。算出来之后输出一个AABB包围盒Axis-Aligned Bounding Box长宽高一目了然直接跟你的孔位排做对比。构型约定本篇采用2-DOF平面搬运机械臂纯平面作业Z恒为0基座J1回转范围±170°大臂J2俯仰范围±120°。DH参数连杆1 a₁0.40 m, α₁0°连杆2 a₂0.35 m, α₂°。末端执行器无额外偏置。正运动学输出末端在基座标系下的X/Y坐标通过蒙特卡洛采样关节空间生成工作空间点云计算AABB包围盒输出长X方向跨度、宽Y方向跨度、高Z方向恒为0边界。技术栈限定numpy matplotlib scipy。⚠️ 声明本篇全部代码为课程离线仿真运行在普通PC上用于工作空间可达性评估。代码中不包含任何真实机器人控制器实时运算或安全逻辑。不可直接用于产线设备。2. 痛点分析课程模块 本篇对应的工程痛点第3章 连杆参数及运动学求解 → 3.3 正运动学 蒙特卡洛采样正运动学批量求解末端坐标第9章 仿真实验 → 工作空间分析 关节限位约束下的可达工作空间点云生成与AABB包围盒计算本篇代码解决 L2 离线仿真校核层。输入机械臂DH参数、关节限位、采样点数。处理①在关节限位内均匀随机采样关节角②正运动学批量计算末端位置③生成工作空间点云④计算AABB包围盒最小/最大X、Y、Z⑤输出长、宽、高边界。输出4张分析图表 CSV数据文件。不覆盖逆运动学、动力学、力控制、实时控制、多臂协同、三维空间臂本篇为2-DOF平面臂Z恒为0高边界为0。3. 核心逻辑讲解物理直觉机械臂的工作空间就是它手能够到的所有地方。想象你坐在椅子上肩膀和手肘是两个关节手能画出一个月牙形的区域——这就是工作空间。如果肩膀只能转一定的角度关节限位这个月牙就会被切掉几块。蒙特卡洛方法就是随机试。让两个关节在各自限位内随机转到一个角度算出手在哪里记下一个点。重复几万次所有点聚在一起就勾勒出了工作空间的形状。然后画一个刚好把所有点装进去的长方体AABB包围盒它的长最大X-最小X宽最大Y-最小Y高最大Z-最小Z。这个包围盒就是工艺工程师用来判断能不能够到目标的最直观依据。业务逻辑 → 代码映射输入DH参数(a10.4m, a20.35m), 关节限位(J1:±170°, J2:±120°), 采样数N20000↓【步骤1关节空间蒙特卡洛采样】对每一组采样 i 1..N:θ1_i uniform(-170°, 170°) → 弧度θ2_i uniform(-120°, 120°) → 弧度↓【步骤2正运动学批量计算】对每一组(θ1_i, θ2_i):x_i a1 * cos(θ1_i) a2 * cos(θ1_i θ2_i)y_i a1 * sin(θ1_i) a2 * sin(θ1_i θ2_i)z_i 0得到点云 array(N, 3)↓【步骤3AABB包围盒计算】x_min min(x_i), x_max max(x_i)y_min min(y_i), y_max max(y_i)z_min 0, z_max 0长 x_max - x_min宽 y_max - y_min高 z_max - z_min 0↓【步骤4输出报告 4张图表】关键公式文字翻译正运动学解析解平面2R臂-x a₁·cos(θ₁) a₂·cos(θ₁θ₂)末端X坐标 大臂水平投影 小臂水平投影-y a₁·sin(θ₁) a₂·sin(θ₁θ₂)末端Y坐标 大臂垂直投影 小臂垂直投影- 为什么不用完整DH矩阵乘法平面2R臂的解析解就是两节杆首尾相接的三角形几何直接用三角函数比4×4矩阵乘法快一个数量级适合数万次采样。AABB包围盒- 长 max(x) - min(x)单位m。物理含义工作空间在X方向的总跨度。- 宽 max(y) - min(y)单位m。物理含义工作空间在Y方向的总跨度。- 高 max(z) - min(z)单位m。本构型Z恒为0高为0。关节限位裁剪采样时θ₁限制在±170°而非±180°是因为实际电机有物理硬限位和电气安全余量不能转到死点。理论最大X跨度≈1.49 m因为±170°导致左侧少了一小角Y跨度≈1.50 m±170°包含±90°上下对称极限都能达到。4. OOP代码实现技术栈限定numpy, matplotlib, scipy。不引入pybullet、roboticstoolbox等任何第三方重型库。全部自实现。项目结构树workspace_eval/├── run_workspace_eval.py├── requirements.txt├── robot/│ ├── __init__.py│ ├── models.py│ ├── kinematics.py│ ├── workspace_sampler.py│ └── visualize.py├── tests/│ ├── __init__.py│ └── test_workspace_eval.py├── results/└── output/核心源码robot/models.py — 数据模型Data models: DH parameters, robot config, joint limitsAll parameters built-in, self-containedimport mathfrom dataclasses import dataclass, fieldfrom typing import Listdataclassclass JointLimit:Joint angle limits (rad)joint_name: strlower: float # radupper: float # raddef __post_init__(self):if self.lower self.upper:raise ValueError(f{self.joint_name}: lower upper)if not (math.isfinite(self.lower)and math.isfinite(self.upper)):raise ValueError(f{self.joint_name}: limits must be finite)dataclassclass RobotConfig:2-DOF planar arm: J1 rotation, J2 pitcha10.40m, a20.35m, both alpha0a1: float 0.40a2: float 0.35joint_limits: List[JointLimit] field(default_factorylist)def __post_init__(self):if not self.joint_limits:# J1: /- 170 deg, J2: /- 120 degself.joint_limits [JointLimit(J1, math.radians(-170),math.radians(170)),JointLimit(J2, math.radians(-120),math.radians(120)),]for a in [self.a1, self.a2]:if a 0:raise ValueError(Link lengths must be positive)robot/kinematics.py — 正运动学解析解Forward kinematics for 2-DOF planar armAnalytical solution: x a1*cos(q1) a2*cos(q1q2)y a1*sin(q1) a2*sin(q1q2)Why analytical: 2R planar arm has closed-form FK, much fasterthan 4x4 DH matrix multiplication for Monte Carlo samplingimport mathimport numpy as npfrom robot.models import RobotConfigclass Planar2RKinematics:def __init__(self, config: RobotConfig):self.config configdef fk_single(self, q1, q2):Single forward kinematics computationReturns (x, y, z0)a1, a2 self.config.a1, self.config.a2x a1 * math.cos(q1) a2 * math.cos(q1 q2)y a1 * math.sin(q1) a2 * math.sin(q1 q2)return np.array([x, y, 0.0], dtypefloat)def fk_batch(self, q1_array, q2_array):Vectorized forward kinematics for N samplesq1_array, q2_array: 1D arrays of length NReturns: points array of shape (N, 3)Why vectorized: avoids Python loop overhead,critical for 20000 samplesa1, a2 self.config.a1, self.config.a2x a1 * np.cos(q1_array) a2 * np.cos(q1_array q2_array)y a1 * np.sin(q1_array) a2 * np.sin(q1_array q2_array)z np.zeros_like(x)return np.column_stack((x, y, z))robot/workspace_sampler.py — 蒙特卡洛采样与包围盒Monte Carlo workspace samplerSamples joint space uniformly, computes end-effector positions,and calculates AABB bounding boximport numpy as npfrom robot.models import RobotConfigfrom robot.kinematics import Planar2RKinematicsclass WorkspaceSampler:def __init__(self, config: RobotConfig,n_samples20000, seed42):n_samples: number of Monte Carlo samplesseed: random seed for reproducibilityWhy seed: ensures identical results across runs,critical for debugging and verificationself.config configself.n_samples n_samplesself.seed seedself.kin Planar2RKinematics(config)def sample_joint_space(self):Uniform random sampling within joint limitsReturns: q1_array, q2_array (shape: N,)rng np.random.default_rng(self.seed)j1_lim self.config.joint_limits[0]j2_lim self.config.joint_limits[1]q1 rng.uniform(j1_lim.lower, j1_lim.upper, self.n_samples)q2 rng.uniform(j2_lim.lower, j2_lim.upper, self.n_samples)return q1, q2def generate_point_cloud(self):Generate workspace point cloud via Monte Carlo FKReturns: points (N, 3), q1_array, q2_arrayq1, q2 self.sample_joint_space()points self.kin.fk_batch(q1, q2)return points, q1, q2def compute_aabb(self, points):Compute Axis-Aligned Bounding BoxReturns: dict with x_min, x_max, y_min, y_max, z_min, z_max,length (x_span), width (y_span), height (z_span)Why separate lines: avoid multi-line tuple unpackingSyntaxError in some Python versionsx_min float(np.min(points[:, 0]))x_max float(np.max(points[:, 0]))y_min float(np.min(points[:, 1]))y_max float(np.max(points[:, 1]))z_min float(np.min(points[:, 2]))z_max float(np.max(points[:, 2]))return {x_min: x_min,x_max: x_max,y_min: y_min,y_max: y_max,z_min: z_min,z_max: z_max,length: x_max - x_min,width: y_max - y_min,height: z_max - z_min,}def evaluate(self):Full pipeline: sample - FK - AABBReturns: points, aabb dict, q1, q2points, q1, q2 self.generate_point_cloud()aabb self.compute_aabb(points)return points, aabb, q1, q2robot/visualize.py — 可视化Workspace visualization: 4 charts1. Point cloud scatter (X-Y plane)2. AABB bounding box overlay3. Joint angle distribution histogram4. Reachable boundary envelope (convex hull via scipy)import osimport numpy as npimport matplotlib.pyplot as pltfrom matplotlib.patches import Rectangleplt.rcParams[font.sans-serif] [SimHei, DejaVu Sans,Noto Sans CJK SC]plt.rcParams[axes.unicode_minus] Falseclass WorkspaceVisualizer:def __init__(self, results_dirresults):self.results_dir results_diros.makedirs(results_dir, exist_okTrue)def plot_point_cloud(self, points, aabb):1. Point cloud scatter in X-Y planen_points len(points)fig, ax plt.subplots(figsize(8, 8))ax.scatter(points[:, 0], points[:, 1],s0.5, csteelblue, alpha0.3,labelfWorkspace Points (n{n_points}))# Mark origin (base)ax.plot(0, 0, ro, markersize8, labelBase (0,0))ax.set_xlabel(X (m))ax.set_ylabel(Y (m))ax.set_title(fWorkspace Point Cloud (Monte Carlo, {n_points} samples))ax.legend()ax.grid(True, alpha0.3)ax.set_aspect(equal, adjustablebox)plt.tight_layout()path os.path.join(self.results_dir,01_workspace_point_cloud.png)plt.savefig(path, dpi150)plt.close()return pathdef plot_aabb_overlay(self, points, aabb):2. Point cloud AABB bounding box rectanglefig, ax plt.subplots(figsize(8, 8))ax.scatter(points[:, 0], points[:, 1],s0.5, csteelblue, alpha0.2)# Draw AABB rectanglerect Rectangle((aabb[x_min], aabb[y_min]),aabb[length], aabb[width],linewidth2, edgecolorred, facecolornone,linestyle--, labelAABB Bounding Box)ax.add_patch(rect)# Annotate dimensions with auto-adjusted positiontext_x aabb[x_max] 0.05text_y aabb[y_max] 0.05ax.text(text_x, 0,fL{aabb[length]:.3f}m,fontsize9, colorred)ax.text(0, text_y,fW{aabb[width]:.3f}m,fontsize9, colorred)# Adjust axes limits to show annotationsax.set_xlim(aabb[x_min] - 0.1, aabb[x_max] 0.15)ax.set_ylim(aabb[y_min] - 0.1, aabb[y_max] 0.15)ax.plot(0, 0, ro, markersize8, labelBase (0,0))ax.set_xlabel(X (m))ax.set_ylabel(Y (m))ax.set_title(Workspace Point Cloud with AABB Bounding Box)ax.legend()ax.grid(True, alpha0.3)ax.set_aspect(equal, adjustablebox)plt.tight_layout()path os.path.join(self.results_dir,02_aabb_overlay.png)plt.savefig(path, dpi150)plt.close()return pathdef plot_joint_distribution(self, q1, q2):3. Joint angle distribution histogramsfig, (ax1, ax2) plt.subplots(1, 2, figsize(12, 5))ax1.hist(np.degrees(q1), bins50, color#2196F3,edgecolorwhite, alpha0.8)ax1.set_xlabel(Joint 1 Angle (deg))ax1.set_ylabel(Frequency)ax1.set_title(J1 Angle Distribution (uniform within limits))ax1.axvline(np.degrees(np.min(q1)), colorred,linestyle--, linewidth1, labelmin)ax1.axvline(np.degrees(np.max(q1)), colorred,linestyle--, linewidth1, labelmax)ax1.legend()ax1.grid(True, alpha0.3)ax2.hist(np.degrees(q2), bins50, color#4CAF50,edgecolorwhite, alpha0.8)ax2.set_xlabel(Joint 2 Angle (deg))ax2.set_ylabel(Frequency)ax2.set_title(J2 Angle Distribution (uniform within limits))ax2.axvline(np.degrees(np.min(q2)), colorred,linestyle--, linewidth1)ax2.axvline(np.degrees(np.max(q2)), colorred,linestyle--, linewidth1)ax2.grid(True, alpha0.3)plt.tight_layout()path os.path.join(self.results_dir,03_joint_distribution.png)plt.savefig(path, dpi150)plt.close()return pathdef plot_reachable_boundary(self, points):4. Reachable boundary envelope (convex hull via scipy)scipy is a declared dependency, no fallback neededfrom scipy.spatial import ConvexHullfig, ax plt.subplots(figsize(8, 8))# Scatter (downsampled for clarity)step max(1, len(points) // 5000)ax.scatter(points[::step, 0], points[::step, 1],s0.3, clightblue, alpha0.4,labelPoint Cloud (downsampled))# Compute convex hull for boundaryhull ConvexHull(points[:, :2])hull_points points[hull.vertices, :2]# Close the loophull_points np.vstack([hull_points, hull_points[0:1]])ax.plot(hull_points[:, 0], hull_points[:, 1],r-, linewidth2,labelReachable Boundary (Convex Hull))ax.plot(0, 0, ro, markersize8, labelBase (0,0))ax.set_xlabel(X (m))ax.set_ylabel(Y (m))ax.set_title(Reachable Workspace Boundary Envelope)ax.legend()ax.grid(True, alpha0.3)ax.set_aspect(equal, adjustablebox)plt.tight_layout()path os.path.join(self.results_dir,04_reachable_boundary.png)plt.savefig(path, dpi150)plt.close()return path5. 可视化要求序号 图表文件 服务结论①01_workspace_point_cloud.png 工作空间点云散点图直观展示可达区域形状月牙形可读取约1.49m×1.50m的覆盖范围②02_aabb_overlay.png 点云叠加AABB包围盒矩形可视化长/宽边界标注数值③03_joint_distribution.png 关节角分布直方图验证采样均匀性确认限位边界④04_reachable_boundary.png 可达边界包络线凸包展示工作空间外轮廓6. README使用说明 核心知识点卡片快速开始pip install -r requirements.txtpython -m pytest tests/test_workspace_eval.py -vpython run_workspace_eval.py数据格式说明# 修改臂长编辑 RobotConfig.a1, a2# 修改关节限位编辑 RobotConfig.joint_limits# 修改采样数编辑 WorkspaceSampler(n_samples20000)参数 修改位置 说明L1臂长RobotConfig.a1 mL2臂长RobotConfig.a2 mJ1限位joint_limits[0] radJ2限位joint_limits[1] rad采样数n_samples 无量纲 知识点卡片工作空间分析对应课本第3章3.3节第9章仿真实验- 蒙特卡洛采样在关节限位内均匀随机撒点通过正运动学映射到笛卡尔空间- AABB包围盒Axis-Aligned Bounding Box长X_max-X_min宽Y_max-Y_min高Z_max-Z_min- 平面2R臂解析正解xa₁cosθ₁a₂cos(θ₁θ₂)ya₁sinθ₁a₂sin(θ₁θ₂)- 关节限位约束采样范围必须匹配真实硬件限位否则包围盒虚大- 工程铁律AABB是保守估计实际可达区域比AABB小——AABB四个角落有大量不可达点评审时需留安全余量7. 总结与工程师思考工业落地难处一、连杆柔性导致实际工作空间收缩。 本篇假设连杆为刚体但实际机械臂在负载下会有弹性变形末端实际位置比理论值偏几毫米到十几毫米。对于精密装配如发动机缸体螺栓拧紧这个偏差足以导致螺丝对不准。离线仿真算出的包围盒是理想最大值实际需留安全余量。二、蒙特卡洛采样的均匀性陷阱。 均匀采样关节空间≠均匀覆盖笛卡尔空间。2R臂在奇异位形附近完全伸直或完全折叠微小的关节角变化会导致末端大幅移动造成点云密度不均匀。如果采样数不够如只有2000个边界区域可能漏采包围盒偏小。建议至少20000个样本。三、安装位置与末端偏置未建模。 本篇假设基座在原点且末端无偏置。实际产线中机械臂可能安装在倾斜台面或靠近墙壁基座坐标系有偏移同时末端执行器夹具、吸盘在Z方向有延伸孔位也不在基座平面。完整的可达评估需要加入末端偏置tool offset和碰撞检测下一篇方向。工程师心得小赵看着终端输出的AABB数据长1.489 m宽1.498 m高0.000 m。那排螺栓孔总长0.56 m完全在包围盒里面。但这是理论最大值吧老郑点头对。AABB是保守估计——它把角落里那些实际到不了的点也算进去了。真实可达区域比这个矩形小一圈。你要留至少30 mm的安全余量。另外如果你的夹具本身有100 mm的延伸臂末端执行器的实际工作空间还要再平移过去算。明白。下一篇是不是该加碰撞检测了对。碰撞检测关节限位选优那是下一篇的事。8. 写作与风格限制本篇全部代码为课程离线仿真运行在普通PC上用于工作空间可达性评估。代码中不包含任何真实机器人控制器实时运算或安全逻辑。 完全去营销中立客观。拒绝论文式生硬开头保持故事叙事感。全部物理量标注完整单位。9. 文末附录环境配置pip install numpy1.24.0 matplotlib3.7.0 scipy1.10.0requirements.txtnumpy1.24.0matplotlib3.7.0scipy1.10.0运行说明cd workspace_evalpip install -r requirements.txtpython -m pytest tests/ -vpython run_workspace_eval.py参数修改指南想改什么 改哪里 示例L1臂长RobotConfig.a1 改0.40→0.50L2臂长RobotConfig.a2 改0.35→0.40J1限位joint_limits[0] 改±170°→±150°J2限位joint_limits[1] 改±120°→±90°采样数n_samples 改20000→50000输出文件清单文件 含义results/01_workspace_point_cloud.png 工作空间点云散点图results/02_aabb_overlay.png 点云叠加AABB包围盒results/03_joint_distribution.png 关节角分布直方图results/04_reachable_boundary.png 可达边界包络线凸包output/workspace_report.txt AABB数值报告output/workspace_points.csv 点云坐标数据下采样至约5000行仿真结果摘要统计项 数值 数据来源机器人模型 2-DOF平面臂 (2R) 配置连杆总长 0.75 m a1a2J1限位 ±170° 配置J2限位 ±120° 配置采样数 20000 内置AABB长 (X跨度) ~1.489 m 仿真输出AABB宽 (Y跨度) ~1.498 m 仿真输出AABB高 (Z跨度) 0.000 m 平面构型测试通过 3/3 pytest正文完。以下为run_workspace_eval.py 完整主入口 tests/test_workspace_eval.py 完整测试文件。run_workspace_eval.pyWorkspace Reachability Evaluation - Main EntryRun: python run_workspace_eval.pyimport osimport csvimport numpy as npfrom robot.models import RobotConfigfrom robot.kinematics import Planar2RKinematicsfrom robot.workspace_sampler import WorkspaceSamplerfrom robot.visualize import WorkspaceVisualizerdef main():print( * 60)print(Offline Sim: Assembly Line Workspace Reachability Evaluation)print( * 60)# 1. Initconfig RobotConfig()sampler WorkspaceSampler(config, n_samples20000, seed42)viz WorkspaceVisualizer(results_dirresults)print(f\n--- Robot Config ---)print(f a1 {config.a1:.2f} m)print(f a2 {config.a2:.2f} m)print(f J1 limit: [{np.degrees(config.joint_limits[0].lower):.0f}deg, f{np.degrees(config.joint_limits[0].upper):.0f}deg])print(f J2 limit: [{np.degrees(config.joint_limits[1].lower):.0f}deg, f{np.degrees(config.joint_limits[1].upper):.0f}deg])print(f Samples: {sampler.n_samples})# 2. Generate point cloud and AABBprint(\n--- Generating Point Cloud (Monte Carlo) ---)points, aabb, q1, q2 sampler.evaluate()print(f\n--- AABB Bounding Box ---)print(f X: [{aabb[x_min]:.4f}, {aabb[x_max]:.4f}] m)print(f Y: [{aabb[y_min]:.4f}, {aabb[y_max]:.4f}] m)print(f Z: [{aabb[z_min]:.4f}, {aabb[z_max]:.4f}] m)print(f Length (X span): {aabb[length]:.4f} m)print(f Width (Y span): {aabb[width]:.4f} m)print(f Height (Z span): {aabb[height]:.4f} m)# 3. Save text reportos.makedirs(output, exist_okTrue)report_path output/workspace_report.txtwith open(report_path, w, encodingutf-8) as f:f.write( * 50 \n)f.write(Workspace Reachability Evaluation Report\n)f.write( * 50 \n\n)f.write(fRobot: 2-DOF Planar Arm\n)f.write(fa1 {config.a1:.2f} m, a2 {config.a2:.2f} m\n)f.write(fJ1 limit: /- 170 deg\n)f.write(fJ2 limit: /- 120 deg\n)f.write(fSamples: {sampler.n_samples}\n\n)f.write(AABB Bounding Box:\n)f.write(f X_min {aabb[x_min]:.6f} m\n)f.write(f X_max {aabb[x_max]:.6f} m\n)f.write(f Y_min {aabb[y_min]:.6f} m\n)f.write(f Y_max {aabb[y_max]:.6f} m\n)f.write(f Length {aabb[length]:.6f} m\n)f.write(f Width {aabb[width]:.6f} m\n)f.write(f Height {aabb[height]:.6f} m\n)print(f\nReport saved: {report_path})# 4. Visualizationprint(\n--- Generating Charts ---)viz.plot_point_cloud(points, aabb)print( 01_workspace_point_cloud.png done)viz.plot_aabb_overlay(points, aabb利用AI解决实际问题如果你觉得这个工具好用欢迎关注长安牧笛
