轨迹优化 | 图解欧氏距离场与梯度场算法(附ROS C++/Python实现)

2024-06-17 12:20

本文主要是介绍轨迹优化 | 图解欧氏距离场与梯度场算法(附ROS C++/Python实现),希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

目录

  • 0 专栏介绍
  • 1 什么是距离场?
  • 2 欧氏距离场计算原理
  • 3 双线性插值与欧式梯度场
  • 4 仿真实现
    • 4.1 ROS C++实现
    • 4.2 Python实现

0 专栏介绍

🔥课程设计、毕业设计、创新竞赛、学术研究必备!本专栏涉及更高阶的运动规划算法实战:曲线生成与轨迹优化、碰撞模型与检测、多智能体群控、深度强化学习运动规划、社会性导航、全覆盖路径规划等内容,每个模型都包含代码实现加深理解。

🚀详情:运动规划实战进阶


1 什么是距离场?

距离场(Distance Field)也称为距离变换(Distance Transform),在图像处理和模式识别中是一种重要工具,其核心思想是将图像中每个像素点的值表示为到最近的目标像素的距离度量,更高维的情况依次类推。在计算机视觉领域,距离场是图像分割和配准的基础,在运动规划中则常用于地图建模(广义Voronoi图计算)和轨迹梯度优化。

在这里插入图片描述

给定点集 G \mathcal{G} G和其上的目标点集 Q ⊆ G Q\subseteq \mathcal{G} QG,可以定义函数

D I ( p ) = min ⁡ q ∈ G ( d ( p , q ) + I ( q ) ) \mathcal{D} _I\left( p \right) =\min _{q\in \mathcal{G}}\left( d\left( p,q \right) +I\left( q \right) \right) DI(p)=qGmin(d(p,q)+I(q))

计算任意一点 p p p Q Q Q的最近距离,其中距离度量 d ( ⋅ , ⋅ ) d\left( \cdot ,\cdot \right) d(,)定义了距离场的属性,指示函数

I Q ( q ) = { 0 , i f q ∈ Q ∞ , o t h e r w i s e I_Q\left( q \right) =\begin{cases} 0, \mathrm{if} q\in Q\\ \infty , \mathrm{otherwise}\\\end{cases} IQ(q)={0,ifqQ,otherwise

常用的 d ( ⋅ , ⋅ ) d\left( \cdot ,\cdot \right) d(,)

  • 欧氏距离,此时称为欧氏距离变换(Euclidean Distance Transform, EDT)
  • 曼哈顿距离,此时称为曼哈顿距离变换(Manhattan Distance Transform, MDT)
  • 切比雪夫距离,此时称为切比雪夫距离变换(Chebyshev Distance Transform, CDT)

2 欧氏距离场计算原理

n n n维距离场可以通过一维距离场迭代计算得到,因此只需要讨论一维EDT的计算即可。如下图左侧所示为初始计算轴的EDT计算,右侧所示为更一般的情况,此时障碍物处的采样函数叠加了前轴计算信息。

在这里插入图片描述

注意到 O \mathcal{O} O定义了一系列以障碍物 q ∈ O q\in \mathcal{O} qO为顶点的抛物线,而 ∀ p ∈ G \forall p\in \mathcal{G} pG在系列抛物线形成的下包络的投影组成了距离场 D f ( p ) \mathcal{D} _f\left( p \right) Df(p),而下包络的计算与抛物线交点有关。联立两条抛物线 ( s − q ) 2 + f ( q ) = ( s − r ) 2 + f ( r ) \left( s-q \right) ^2+f\left( q \right) =\left( s-r \right) ^2+f\left( r \right) (sq)2+f(q)=(sr)2+f(r)可得

s = ( f ( r ) + r 2 ) − ( f ( q ) + q 2 ) 2 r − 2 q s=\frac{\left( f\left( r \right) +r^2 \right) -\left( f\left( q \right) +q^2 \right)}{2r-2q} s=2r2q(f(r)+r2)(f(q)+q2)

即任意两条抛物线有且仅有一个交点 s s s。设 K \mathcal{K} K为实际组成下包络的抛物线集合, v ( k ) v\left( k \right) v(k)表示其中第 k k k条抛物线的顶点, z ( k ) z\left( k \right) z(k)表示第 k k k条和第 k − 1 k-1 k1条抛物线的交点,区间 [ z ( k ) , z ( k + 1 ) ) \left[ z\left( k \right) ,z\left( k+1 \right) \right) [z(k),z(k+1))表示第 k k k条抛物线的下包络范围。在遍历求解下包络过程中,对于新的抛物线 e e e,其与 K \mathcal{K} K中最新的一条抛物线 k k k的交点 s s s有两种情况:

  • s > z ( k ) s>z\left( k \right) s>z(k),则将 e e e添加到 K \mathcal{K} K中并更新 v ( k ) v\left( k \right) v(k) z ( k ) z\left( k \right) z(k)
  • s ⩽ z ( k ) s\leqslant z\left( k \right) sz(k),则第 k k k条抛物线不参与构成下包络,应从 K \mathcal{K} K中删除并重新计算 e e e与新的第 k k k条抛物线的交点直至 s s s z ( k ) z(k) z(k)右侧;

算法流程如表所示

在这里插入图片描述

3 双线性插值与欧式梯度场

在求解梯度过程中,需要计算离散距离场的线性插值函数。以二维环境的双线性插值为例,设已知离散空间有四点坐标为 A ( x 1 , y 1 ) A\left( x_1,y_1 \right) A(x1,y1) B ( x 2 , y 1 ) B\left( x_2,y_1 \right) B(x2,y1) C ( x 1 , y 2 ) C\left( x_1,y_2 \right) C(x1,y2) D ( x 2 , y 2 ) D\left( x_2,y_2 \right) D(x2,y2),对其中任意一点 P P P,首先在 x x x方向上插值

{ f ( R 1 ) = f ( x , y 1 ) = x 2 − x x 2 − x 1 f ( A ) + x − x 1 x 2 − x 1 f ( B ) f ( R 2 ) = f ( x , y 2 ) = x 2 − x x 2 − x 1 f ( C ) + x − x 1 x 2 − x 1 f ( D ) \begin{cases} f\left( R_1 \right) =f\left( x,y_1 \right) =\frac{x_2-x}{x_2-x_1}f\left( A \right) +\frac{x-x_1}{x_2-x_1}f\left( B \right)\\ f\left( R_2 \right) =f\left( x,y_2 \right) =\frac{x_2-x}{x_2-x_1}f\left( C \right) +\frac{x-x_1}{x_2-x_1}f\left( D \right)\\\end{cases} {f(R1)=f(x,y1)=x2x1x2xf(A)+x2x1xx1f(B)f(R2)=f(x,y2)=x2x1x2xf(C)+x2x1xx1f(D)

再基于插值点 R 1 R_1 R1 R 2 R_2 R2进行 y y y方向的插值

f ( P ) = f ( x , y ) = y 2 − y y 2 − y 1 f ( R 1 ) + y − y 1 y 2 − y 1 f ( R 2 ) f\left( P \right) =f\left( x,y \right) =\frac{y_2-y}{y_2-y_1}f\left( R_1 \right) +\frac{y-y_1}{y_2-y_1}f\left( R_2 \right) f(P)=f(x,y)=y2y1y2yf(R1)+y2y1yy1f(R2)

展开可得矩阵形式

f ( x , y ) = 1 ( x 2 − x 1 ) ( y 2 − y 1 ) [ x 2 − x x − x 1 ] T [ f ( A ) f ( C ) f ( B ) f ( D ) ] [ y 2 − y y − y 1 ] f\left( x,y \right) =\frac{1}{\left( x_2-x_1 \right) \left( y_2-y_1 \right)}\left[ \begin{array}{c} x_2-x\\ x-x_1\\\end{array} \right] ^T\left[ \begin{matrix} f\left( A \right)& f\left( C \right)\\ f\left( B \right)& f\left( D \right)\\\end{matrix} \right] \left[ \begin{array}{c} y_2-y\\ y-y_1\\\end{array} \right] f(x,y)=(x2x1)(y2y1)1[x2xxx1]T[f(A)f(B)f(C)f(D)][y2yyy1]

在这里插入图片描述

对于离散栅格坐标而言, x 2 − x 1 = y 2 − y 1 = 1 x_2-x_1=y_2-y_1=1 x2x1=y2y1=1,设 Δ x = x − x 1 \varDelta x=x-x_1 Δx=xx1 Δ y = y − y 1 \varDelta y=y-y_1 Δy=yy1,则插值函数简化为

f ( Δ x , Δ y ) = [ 1 − Δ x Δ x ] T [ f ( A ) f ( C ) f ( B ) f ( D ) ] [ 1 − Δ y Δ y ] f\left( \varDelta x,\varDelta y \right) =\left[ \begin{array}{c} 1-\varDelta x\\ \varDelta x\\\end{array} \right] ^T\left[ \begin{matrix} f\left( A \right)& f\left( C \right)\\ f\left( B \right)& f\left( D \right)\\\end{matrix} \right] \left[ \begin{array}{c} 1-\varDelta y\\ \varDelta y\\\end{array} \right] f(Δx,Δy)=[1ΔxΔx]T[f(A)f(B)f(C)f(D)][1ΔyΔy]

则函数 f f f在任意一点的梯度为

∇ f ( Δ x , Δ y ) = [ ∂ f ( Δ x , Δ y ) ∂ Δ x ∂ f ( Δ x , Δ y ) ∂ Δ y ] = [ [ − 1 1 ] T [ f ( A ) f ( C ) f ( B ) f ( D ) ] [ 1 − Δ y Δ y ] [ 1 − Δ x Δ x ] T [ f ( A ) f ( C ) f ( B ) f ( D ) ] [ − 1 1 ] ] \nabla f\left( \varDelta x,\varDelta y \right) =\left[ \begin{array}{c} \frac{\partial f\left( \varDelta x,\varDelta y \right)}{\partial \varDelta x}\\ \frac{\partial f\left( \varDelta x,\varDelta y \right)}{\partial \varDelta y}\\\end{array} \right] =\left[ \begin{array}{c} \left[ \begin{array}{c} -1\\ 1\\\end{array} \right] ^T\left[ \begin{matrix} f\left( A \right)& f\left( C \right)\\ f\left( B \right)& f\left( D \right)\\\end{matrix} \right] \left[ \begin{array}{c} 1-\varDelta y\\ \varDelta y\\\end{array} \right]\\ \left[ \begin{array}{c} 1-\varDelta x\\ \varDelta x\\\end{array} \right] ^T\left[ \begin{matrix} f\left( A \right)& f\left( C \right)\\ f\left( B \right)& f\left( D \right)\\\end{matrix} \right] \left[ \begin{array}{c} -1\\ 1\\\end{array} \right]\\\end{array} \right] f(Δx,Δy)=[Δxf(Δx,Δy)Δyf(Δx,Δy)]=[11]T[f(A)f(B)f(C)f(D)][1ΔyΔy][1ΔxΔx]T[f(A)f(B)f(C)f(D)][11]

4 仿真实现

4.1 ROS C++实现

核心代码如下所示

GradientLayer::updateCosts(nav2_costmap_2d::Costmap2D & master_grid, int min_i, int min_j,int max_i,int max_j)
{if (!enabled_) {return;}unsigned char * master_array = master_grid.getCharMap();unsigned int size_x = master_grid.getSizeInCellsX(), size_y = master_grid.getSizeInCellsY();// Fixing window coordinates with map size if necessary.min_i = std::max(0, min_i);min_j = std::max(0, min_j);max_i = std::min(static_cast<int>(size_x), max_i);max_j = std::min(static_cast<int>(size_y), max_j);// Simply computing one-by-one cost per each cellint gradient_index;for (int j = min_j; j < max_j; j++) {// Reset gradient_index each time when reaching the end of re-calculated window// by OY axis.gradient_index = 0;for (int i = min_i; i < max_i; i++) {int index = master_grid.getIndex(i, j);// setting the gradient costunsigned char cost = (LETHAL_OBSTACLE - gradient_index*GRADIENT_FACTOR)%255;if (gradient_index <= GRADIENT_SIZE) {gradient_index++;} else {gradient_index = 0;}master_array[index] = cost;}}
}

4.2 Python实现

欧氏距离场核心代码

def compute(self, f_get, mat: np.ndarray, dim: int) -> np.array:"""Compute distance field along one-dimension baseon sample funciton f_get.Parameters:f_get (function): sample functionmat (np.ndarray): the matrix to transformdim (int): the dimension to transformReturns:df (np.array): the distance field along dim-dimension"""# initialzationk = 0n = mat.shape[dim]v, z = [0 for _ in range(n)], [0 for _ in range(n + 1)]z[0], z[1] = -self.INF, self.INF# envelopefor q in range(1, n):s = ((f_get(v[k]) + v[k] ** 2) - (f_get(q) + q ** 2)) / (2 * (v[k] - q))while s <= z[k]:k -= 1s = ((f_get(v[k]) + v[k] ** 2) - (f_get(q) + q ** 2)) / (2 * (v[k] - q))k += 1v[k] = qz[k], z[k + 1] = s, self.INF# distance calculationk = 0edf = np.zeros((n, ))for q in range(n):while z[k + 1] < q:k += 1edf[q] = (q - v[k]) ** 2 + f_get(v[k])return edf

欧氏梯度场核心代码:

def gradient(self, df: np.ndarray, x: float, y: float) -> np.array:"""To obtain the gradient at (x, y) in the distance field through bilinear interpolation.Parameters:df (np.ndarray): the distance fieldx/y (float): the query coordinateReturns:g(x, y): the gradient at (x, y)"""m, n = df.shapex, y = max(min(n - 1, x), 0), max(min(m - 1, y), 0)xi, yi = int(x), int(y)dx, dy = x - xi, y - yixi, yi = max(min(n - 1, xi), 0), max(min(m - 1, yi), 0)xp, yp = max(min(n - 1, xi + 1), 0), max(min(m - 1, yi - 1), 0)bl, br = df[yi, xi], df[yi, xp]tl, tr = df[yp, xi], df[yp, xp]return np.array([(1 - dy) * (br - bl) + dy * (tr + tl),-((1 - dx) * (tl - bl) + dx * (tr - br))])

效果如下所示
在这里插入图片描述

一个实际地图的案例如下

在这里插入图片描述

完整工程代码请联系下方博主名片获取


🔥 更多精彩专栏

  • 《ROS从入门到精通》
  • 《Pytorch深度学习实战》
  • 《机器学习强基计划》
  • 《运动规划实战精讲》

👇源码获取 · 技术交流 · 抱团学习 · 咨询分享 请联系👇

这篇关于轨迹优化 | 图解欧氏距离场与梯度场算法(附ROS C++/Python实现)的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/1069424

相关文章

SpringBoot3实现Gzip压缩优化的技术指南

《SpringBoot3实现Gzip压缩优化的技术指南》随着Web应用的用户量和数据量增加,网络带宽和页面加载速度逐渐成为瓶颈,为了减少数据传输量,提高用户体验,我们可以使用Gzip压缩HTTP响应,... 目录1、简述2、配置2.1 添加依赖2.2 配置 Gzip 压缩3、服务端应用4、前端应用4.1 N

SpringBoot实现数据库读写分离的3种方法小结

《SpringBoot实现数据库读写分离的3种方法小结》为了提高系统的读写性能和可用性,读写分离是一种经典的数据库架构模式,在SpringBoot应用中,有多种方式可以实现数据库读写分离,本文将介绍三... 目录一、数据库读写分离概述二、方案一:基于AbstractRoutingDataSource实现动态

Python FastAPI+Celery+RabbitMQ实现分布式图片水印处理系统

《PythonFastAPI+Celery+RabbitMQ实现分布式图片水印处理系统》这篇文章主要为大家详细介绍了PythonFastAPI如何结合Celery以及RabbitMQ实现简单的分布式... 实现思路FastAPI 服务器Celery 任务队列RabbitMQ 作为消息代理定时任务处理完整

Python Websockets库的使用指南

《PythonWebsockets库的使用指南》pythonwebsockets库是一个用于创建WebSocket服务器和客户端的Python库,它提供了一种简单的方式来实现实时通信,支持异步和同步... 目录一、WebSocket 简介二、python 的 websockets 库安装三、完整代码示例1.

揭秘Python Socket网络编程的7种硬核用法

《揭秘PythonSocket网络编程的7种硬核用法》Socket不仅能做聊天室,还能干一大堆硬核操作,这篇文章就带大家看看Python网络编程的7种超实用玩法,感兴趣的小伙伴可以跟随小编一起... 目录1.端口扫描器:探测开放端口2.简易 HTTP 服务器:10 秒搭个网页3.局域网游戏:多人联机对战4.

Java枚举类实现Key-Value映射的多种实现方式

《Java枚举类实现Key-Value映射的多种实现方式》在Java开发中,枚举(Enum)是一种特殊的类,本文将详细介绍Java枚举类实现key-value映射的多种方式,有需要的小伙伴可以根据需要... 目录前言一、基础实现方式1.1 为枚举添加属性和构造方法二、http://www.cppcns.co

使用Python实现快速搭建本地HTTP服务器

《使用Python实现快速搭建本地HTTP服务器》:本文主要介绍如何使用Python快速搭建本地HTTP服务器,轻松实现一键HTTP文件共享,同时结合二维码技术,让访问更简单,感兴趣的小伙伴可以了... 目录1. 概述2. 快速搭建 HTTP 文件共享服务2.1 核心思路2.2 代码实现2.3 代码解读3.

MySQL双主搭建+keepalived高可用的实现

《MySQL双主搭建+keepalived高可用的实现》本文主要介绍了MySQL双主搭建+keepalived高可用的实现,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,... 目录一、测试环境准备二、主从搭建1.创建复制用户2.创建复制关系3.开启复制,确认复制是否成功4.同

Python使用自带的base64库进行base64编码和解码

《Python使用自带的base64库进行base64编码和解码》在Python中,处理数据的编码和解码是数据传输和存储中非常普遍的需求,其中,Base64是一种常用的编码方案,本文我将详细介绍如何使... 目录引言使用python的base64库进行编码和解码编码函数解码函数Base64编码的应用场景注意

Java实现文件图片的预览和下载功能

《Java实现文件图片的预览和下载功能》这篇文章主要为大家详细介绍了如何使用Java实现文件图片的预览和下载功能,文中的示例代码讲解详细,感兴趣的小伙伴可以跟随小编一起学习一下... Java实现文件(图片)的预览和下载 @ApiOperation("访问文件") @GetMapping("