一、A_star算法简介
0 引言
随着现代技术的发展,飞行器种类不断变多,应用也日趋专一化、完善化,如专门用作植保的大疆PS-X625无人机,用作街景拍摄与监控巡察的宝鸡行翼航空科技的X8无人机,以及用作水下救援的白鲨MIX水下无人机等,决定飞行器性能主要是内部的飞控系统和外部的路径规划问题。就路径问题而言,在具体实施任务时仅靠操作员手中的遥控器控制无人飞行器执行相应的工作,可能会对操作员心理以及技术提出极高的要求,为了避免个人操作失误,进而造成飞行器损坏的危险,一种解决问题的方法就是对飞行器进行航迹规划。
飞行器的测量精度,航迹路径的合理规划,飞行器工作时的稳定性、安全性等这些变化对飞行器的综合控制系统要求越来越高。无人机航路规划是为了保证无人机完成特定的飞行任务,并且能够在完成任务的过程中躲避各种障碍、威胁区域而设计出最优航迹路线的问题。常见的航迹规划算法如图1所示。
图1 常见路径规划算法
文中主要对无人机巡航阶段的航迹规划进行研究,假设无人机在飞行中维持高度与速度不变,那么航迹规划成为一个二维平面的规划问题。在航迹规划算法中,A算法计算简单,容易实现。在改进A算法基础上,提出一种新的、易于理解的改进A算法的无人机航迹规划方法。传统A算法将规划区域栅格化,节点扩展只限于栅格线的交叉点,在栅格线的交叉点与交叉点之间往往存在一定角度的两个运动方向。将存在角度的两段路径无限放大、细化,然后分别用两段上的相应路径规划点作为切点,找到相对应的组成内切圆的圆心,然后作弧,并求出相对应的两切点之间的弧所对应的圆心角,根据下式计算出弧线的长度
式中:R———内切圆的半径;
α———切点之间弧线对应的圆心角。
1 A*算法概述
A算法是在Dijstar算法的基础上引入的启发式函数,通过定义的代价函数来评估代价大小,从而确定最优路径。A算法的代价函数
式中:f(x,y)———初始状态X0(x0,y0)到达目标状态X1(x1,y1)的代价估计;
g(x,y)———状态空间中从初始状态X0(x0,y0)到状态N(x1,y1)的实际代价;
h(x,y)———从状态N(x1,y1)到目标状态X1(x1,y1)最佳路径的估计代价。
要找到最短路径的实质是找到f(x,y)的最小值,其中在式(2)中寻找最短路径的关键在于求估计代价h (x,y)值。设系数λ表示状态N(x1,y1)到X1(x1,y1)最优距离,如果λ<h(x,y),搜索范围小,不能保证得到最优解;λ>h(x,y),搜索范围大,费时,但能找到最优解;λ=h(x,y),搜索到最短路径。其中h(x,y)一般用欧几里德距离(式(3))或者绝对值距离(式(4))计算。
A算法是以起始点为中心,周围8个栅格的中心为下一步预选,并不断地计算预选位置的f(x,y)值,其中f(x,y)值最小的作为当前位置,依次逐层比较,直到当前位置的临近点出现目标点为止,其最小单元如图2所示。
图2 最小单元
A算法的流程如下:
1)创建开始节点START,目标节点TARGET、OPEN列表、CLOSE列表、CLOSE列表初始为空;
2)将START加入到OPEN列表;
3)检查OPEN列表中的节点,若列表为空,则无可行路径;若不为空,选择使f(x,y)值最小的节点k;
4)将节点k从OPEN中去除,并将其添加到CLOSE中,判断节点k是否目标节点TARGET,若是,则说明找到路径;若不是,则继续扩展节点k,生成k节点的子节点集,设q为k的子节点集,对所有节点q计算相应的f(x,y)值,并选择f(x,y)值最小的节点,将该节点放入CLOSE列表中;
5)跳到3),直到算法获得可行路径或无解退出。
二、部分源代码
function varargout = astar_jw(varargin)
if nargin == 0 % Test case, use testCase=2 (Maze) by default
selectedGrid = 2;
[grid, init, goal] = defineGrid(selectedGrid);
% Whether or not to display debug info
printDebugInfo = true;
elseif nargin == 1 % Test case, set testCase to the input
selectedGrid = varargin{1};
[grid, init, goal] = defineGrid(selectedGrid);
% Whether or not to display debug info
printDebugInfo = true;
elseif nargin == 3 % Function call with inputs
grid = varargin{1};
init = varargin{2};
goal = varargin{3};
printDebugInfo = false;
end
% Define all possible moves
delta = [[-1 0]
[ 0 -1]
[ 1 0]
[ 0 1]];
% Add g & f terms to init if necessary
if length(init)==2
init(3) = 0; % g
init(4) = inf; % f
end
% Perform search
[path, directions] = search(grid, init, goal, delta, printDebugInfo);
if nargout==1
varargout{1} = path;
elseif nargout==2
varargout{1} = path;
varargout{2} = directions;
end
end
function [path, directions] = search(grid, init, goal, delta, printDebugInfo)
% This function implements the A* algorithm
% Initialize the open, closed and path lists
open = []; closed = []; path = [];
open = addRow(open, init);
% Initialize directions list
directions = [];
% Initialize expansion timing grid
expanded = -1 * ones(size(grid));
expanded(open(1), open(2)) = 0;
% Compute the heuristic measurement, h
h = computeHeuristic(grid, goal, 'euclidean');
% Open window for graphical debug display if desired
if printDebugInfo; fig = figure; end
% Keep searching through the grid until the goal is found
goalFound = false;
while size(open,1)>0 && ~goalFound
[open, closed, expanded] = expand(grid, open, closed, delta, expanded, h);
% Display debug info if desired
if printDebugInfo
displayDebugInfo(grid, init, goal, open, closed, fig);
end
goalFound = checkForGoal(closed, goal);
end
% If the goal was found lets get the optimal path
if goalFound
% We step from the goal to the start location. At each step we
% select the neighbor with the lowest 'expanded' value and add that
% neighbor to the path list.
path = goal;
% Check to see if the start location is on the path list yet
[~, indx] = ismember(init(1:2), path(:,1:2), 'rows');
while ~indx
% Use our findNeighbors function to find the neighbors of the
% last location on the path list
neighbors = findNeighbors(grid, path, size(path,1), delta);
% Look at the expanded value of all the neighbors, add the one
% with the lowest expanded value to the path
expandedVal = expanded(goal(1),goal(2));
for R = 1:size(neighbors,1)
neighborExpandedVal = expanded(neighbors(R,1), neighbors(R,2));
if neighborExpandedVal<expandedVal && neighborExpandedVal>=0
chosenNeighbor = R;
expandedVal = expanded(neighbors(R,1), neighbors(R,2));
end
end
path(end+1,:) = neighbors(chosenNeighbor,:);
% Check to see if the start location has been added to the path
% list yet
[~, indx] = ismember(init(1:2), path(:,1:2), 'rows');
end
% Reorder the list to go from the starting location to the end loc
path = flipud(path);
% Compute the directions from the path
directions = zeros(size(path)-[1 0]);
for R = 1:size(directions,1)
directions(R,:) = path(R+1,:) - path(R,:);
end
end
% Display results
if printDebugInfo
home
if goalFound
disp(['Goal Found! Distance to goal: ' num2str(closed(goalFound,3))])
disp(' ')
disp('Path: ')
disp(path)
fig = figure;
displayPath(grid, path, fig)
else
disp('Goal not found!')
end
disp(' ')
disp('Expanded: ')
disp(expanded)
disp(' ')
disp([' Search time to target: ' num2str(expanded(goal(1),goal(2)))])
end
end
function A = deleteRows(A, rows)
% The following way to delete rows was taken from the mathworks website
% that compared multiple ways to do it. The following appeared to be the
% fastest.
index = true(1, size(A,1));
index(rows) = false;
A = A(index, :);
end
function A = addRow(A, row)
A(end+1,:) = row;
end
function [open, closed, expanded] = expand(grid, open, closed, delta, expanded, h)
% This function expands the open list by taking the coordinate (row) with
% the smallest f value (path cost) and adds its neighbors to the open list.
% Expand the row with the lowest f
row = find(open(:,4)==min(open(:,4)),1);
% Edit the 'expanded' matrix to note the time in which the current grid
% point was expanded
expanded(open(row,1),open(row,2)) = max(expanded(:))+1;
% Find all the neighbors (potential moves) from the chosen spot
neighbors = findNeighbors(grid, open, row, delta);
% Remove any neighbors that are already on the open or closed lists
neighbors = removeListedNeighbors(neighbors, open, closed);
% Add the neighbors still left to the open list
for R = 1:size(neighbors,1)
g = open(row,3)+1;
f = g + h(neighbors(R,1),neighbors(R,2));
open = addRow(open, [neighbors(R,:) g f] );
end
% Remove the row we just expanded from the open list and add it to the
% closed list
closed(end+1,:) = open(row,:);
open = deleteRows(open, row);
end
function h = computeHeuristic(varargin)
% This function is used to compute the distance heuristic, h. By default
% this function computes the Euclidean distance from each grid space to the
% goal. The calling sequence for this function is as follows:
% h = computeHeuristic(grid, goal[, distanceType])
% where distanceType may be one of the following:
% 'euclidean' (default value)
% 'city-block'
% 'empty' (returns all zeros for heuristic function)
grid = varargin{1};
goal = varargin{2};
if nargin==3
distanceType = varargin{3};
else
distanceType = 'euclidean';
end
[m n] = size(grid);
[x y] = meshgrid(1:n,1:m);
if strcmp(distanceType, 'euclidean')
h = sqrt((x-goal(2)).^2 + (y-goal(1)).^2);
elseif strcmp(distanceType, 'city-block')
h = abs(x-goal(2)) + abs(y-goal(1));
elseif strcmp(distanceType, 'empty')
h = zeros(m,n);
else
warning('Unknown distanceType for determining heuristic, h!')
h = [];
end
end
function neighbors = findNeighbors(grid, open, row, delta)
% This function takes the desired row in the open list to expand and finds
% all potential neighbors (neighbors reachable through legal moves, as
% defined in the delta list).
% Find the borders to the grid
borders = size(grid);
borders = [1 borders(2) 1 borders(1)]; % [xMin xMax yMin yMax]
% Initialize the current location and neighbors list
cenSq = open(row,1:2);
neighbors = [];
% Go through all the possible moves (given in the 'delta' matrix) and
% add moves not on the closed list
for R = 1:size(delta,1)
potMove = cenSq + delta(R,:);
if potMove(1)>=borders(3) && potMove(1)<=borders(4) ...
&& potMove(2)>=borders(1) && potMove(2)<=borders(2)
if grid(potMove(1), potMove(2))==0
neighbors(end+1,:) = potMove;
end
end
end
end
function neighbors = removeListedNeighbors(neighbors, open, closed)
% This function removes any neighbors that are on the open or closed lists
% Check to make sure there's anything even on the closed list
if size(closed,1)==0
return
end
% Find any neighbors that are on the open or closed lists
rowsToRemove = [];
for R = 1:size(neighbors)
% Check to see if the neighbor is on the open list
[~, indx] = ismember(neighbors(R,:),open(:,1:2),'rows');
if indx>0
rowsToRemove(end+1) = R;
else
% If the neighbor isn't on the open list, check to make sure it
% also isn't on the closed list
[~, indx] = ismember(neighbors(R,:),closed(:,1:2),'rows');
if indx>0
rowsToRemove(end+1) = R;
end
end
end
% Remove neighbors that were on either the open or closed lists
if numel(rowsToRemove>0)
neighbors = deleteRows(neighbors, rowsToRemove);
end
end
function goalRow = checkForGoal(closed, goal)
% This function looks for the final goal destination on the closed list.
% Note, you could check the open list instead (and find the goal faster);
% however, we want to have a chance to expand the goal location itself, so
% we wait until it is on the closed list.
[~, goalRow] = ismember(goal, closed(:,1:2), 'rows');
end
function displayDebugInfo(grid, init, goal, open, closed, h)
% Display the open and closed lists in the command window, and display an
% image of the current search of the grid.
home
disp('Open: ')
disp(open)
disp(' ')
disp('Closed: ')
disp(closed)
displaySearchStatus(grid, init, goal, open, closed, h)
pause(0.05)
end
function displaySearchStatus(grid, init, goal, open, closed, h)
% This function displays a graphical grid and search status to make
% visualization easier.
grid = double(~grid);
grid(init(1),init(2)) = 0.66;
grid(goal(1),goal(2)) = 0.33;
figure(h)
imagesc(grid); colormap(gray); axis square; axis off; hold on
plot(open(:,2), open(:,1), 'go', 'LineWidth', 2)
plot(closed(:,2), closed(:,1), 'ro', 'LineWidth', 2)
hold off
end
function displayPath(grid, path, h)
grid = double(~grid);
figure(h)
imagesc(grid); colormap(gray); axis off; hold on
title('Optimal Path', 'FontWeight', 'bold');
plot(path(:,2), path(:,1), 'co', 'LineWidth', 2)
plot(path(1,2), path(1,1), 'gs', 'LineWidth', 4)
plot(path(end,2),path(end,1), 'rs', 'LineWidth', 4)
end
三、运行结果
四、matlab版本及参考文献
1 matlab版本
2014a
2 参考文献
[1] 包子阳,余继周,杨杉.智能优化算法及其MATLAB实例(第2版)[M].电子工业出版社,2016.
[2]张岩,吴水根.MATLAB优化算法源代码[M].清华大学出版社,2017.
[3]巫茜,罗金彪,顾晓群,曾青.基于改进PSO的无人机三维航迹规划优化算法[J].兵器装备工程学报. 2021,42(08)
[4]邓叶,姜香菊.基于改进人工势场法的四旋翼无人机航迹规划算法[J].传感器与微系统. 2021,40(07)
[5]马云红,张恒,齐乐融,贺建良.基于改进A*算法的三维无人机路径规划[J].电光与控制. 2019,26(10)
[6]焦阳.基于改进蚁群算法的无人机三维路径规划研究[J].舰船电子工程. 2019,39(03)