FAST PLANNER代码阅读笔记_fast planner github

一、代码运行

首先挂上github链接

https://github.com/HKUST-Aerial-Robotics/Fast-Plannerhttps://github.com/HKUST-Aerial-Robotics/Fast-Planner前端采用Kinodynamic_Astar进行搜索路径，后端采用--优化轨迹。

运行方法：

首先下载功能包，放到src下面，然后catkin_make可以直接进行编译，我用的是Ubuntu18.04。

然后分别在两个终端运行下面两句代码，就可以启动rviz，然后使用2D Nav Goal给它终点坐标即可规划。

source devel/setup.bash && roslaunch plan_manage rviz.launch
source devel/setup.bash && roslaunch plan_manage kino_replan.launch

遇到的报错：pcl-ros没有安装，安装即可

二、代码阅读

这个项目代码量非常大，我是按照功能包的实现顺序来阅读

1、规划管理器plan_manage

这个功能包主要实现前端的kinoastar规划：

fast_planner_node.cpp这个文件主要用来启动ROS节点，然后进行初始化

int main(int argc, char** argv) {
  ros::init(argc, argv, "fast_planner_node");
  ros::NodeHandle nh("~");  //创建一个节点句柄，命名为～
  int planner;
  nh.param("planner_node/planner", planner, -1);
  TopoReplanFSM topo_replan;
  KinoReplanFSM kino_replan;
  if (planner == 1) {//planner默认是1
    kino_replan.init(nh);
  } else if (planner == 2) {
    topo_replan.init(nh);
  }
  ros::Duration(1.0).sleep();
  ros::spin();
  return 0;
}

接下来看KinoReplanFSM类以及其init函数

void KinoReplanFSM::init(ros::NodeHandle& nh) {
  current_wp_  = 0;
  exec_state_  = FSM_EXEC_STATE::INIT;
  have_target_ = false;
  have_odom_   = false;
 
  /*  fsm param  */
  nh.param("fsm/flight_type", target_type_, -1);
  nh.param("fsm/thresh_replan", replan_thresh_, -1.0);
  nh.param("fsm/thresh_no_replan", no_replan_thresh_, -1.0);
 
  nh.param("fsm/waypoint_num", waypoint_num_, -1);
  for (int i = 0; i < waypoint_num_; i++) {
    nh.param("fsm/waypoint" + to_string(i) + "_x", waypoints_[i][0], -1.0);
    nh.param("fsm/waypoint" + to_string(i) + "_y", waypoints_[i][1], -1.0);
    nh.param("fsm/waypoint" + to_string(i) + "_z", waypoints_[i][2], -1.0);
  }
 
  /* initialize main modules */
  planner_manager_.reset(new FastPlannerManager);
  planner_manager_->initPlanModules(nh);
  visualization_.reset(new PlanningVisualization(nh));
 
  /* callback */
  // 创建一个定时器，以某一Duration 秒 的时间间隔执行某一回调函数
  exec_timer_   = nh.createTimer(ros::Duration(0.01), &KinoReplanFSM::execFSMCallback, this);
  safety_timer_ = nh.createTimer(ros::Duration(0.05), &KinoReplanFSM::checkCollisionCallback, this);
 
  waypoint_sub_ =
      nh.subscribe("/waypoint_generator/waypoints", 1, &KinoReplanFSM::waypointCallback, this);
  odom_sub_ = nh.subscribe("/odom_world", 1, &KinoReplanFSM::odometryCallback, this);
 
  replan_pub_  = nh.advertise<std_msgs::Empty>("/planning/replan", 10);
  new_pub_     = nh.advertise<std_msgs::Empty>("/planning/new", 10);
  bspline_pub_ = nh.advertise<plan_manage::Bspline>("/planning/bspline", 10);
}

init函数，以一定的时间间隔调用函数execFSMCallback和checkCollisionCallback。

execFSMCallback函数用来实现对于各种状态的管理

execFSMCallback根据不同状态使用switch语句，对每个状态进行管理，其中重点在GEN_NEW_TRAJ，EXEC_TRAJ和REPLAN_TRAJ

case GEN_NEW_TRAJ: {
      start_pt_  = odom_pos_;
      start_vel_ = odom_vel_;
      start_acc_.setZero();
 
      Eigen::Vector3d rot_x = odom_orient_.toRotationMatrix().block(0, 0, 3, 1);
      start_yaw_(0)         = atan2(rot_x(1), rot_x(0));
      start_yaw_(1) = start_yaw_(2) = 0.0;
 
      bool success = callKinodynamicReplan();
      if (success) {
        changeFSMExecState(EXEC_TRAJ, "FSM"); // 如果成功获取一条轨迹，就改状态为EXEC_TRAJ，去执行
      } else {
        // have_target_ = false;
        // changeFSMExecState(WAIT_TARGET, "FSM");
        changeFSMExecState(GEN_NEW_TRAJ, "FSM");  // 如果没成功，就继续保持生成状态
      }
      break;

2、KinodynamicAstar规划

（1）callKinodynamicReplan函数

GEN_NEW_TRAJ中主要语句是

bool success = callKinodynamicReplan();

这句话调用了callKinodynamicReplan函数

bool KinoReplanFSM::callKinodynamicReplan() {
  bool plan_success =
      planner_manager_->kinodynamicReplan(start_pt_, start_vel_, start_acc_, end_pt_, end_vel_);

然后紧接着调用kinodynamicReplan函数，这个是重点的规划函数

函数输入起点，起始速度，起始加速度，终点位置，终点速度

kinodynamicReplan的主要语句是

  int status = kino_path_finder_->search(start_pt, start_vel, start_acc, end_pt, end_vel, true);

调用了这个search函数，用来进行kinodynamic_Astar规划

（2）搜索轨迹函数KinodynamicAstar::search

重点函数：KinodynamicAstar::search

int KinodynamicAstar::search(Eigen::Vector3d start_pt, Eigen::Vector3d start_v, Eigen::Vector3d start_a,
                             Eigen::Vector3d end_pt, Eigen::Vector3d end_v, bool init, bool dynamic, double time_start)
{
  start_vel_ = start_v;
  start_acc_ = start_a;
 
  PathNodePtr cur_node = path_node_pool_[0];
  cur_node->parent = NULL;
  // 一个state是一个6×1的向量，前面三个存位置，后面三个存速度信息
  cur_node->state.head(3) = start_pt;
  cur_node->state.tail(3) = start_v;
  cur_node->index = posToIndex(start_pt);
  cur_node->g_score = 0.0;
 
  Eigen::VectorXd end_state(6);
  Eigen::Vector3i end_index;
  double time_to_goal;
 
  end_state.head(3) = end_pt;
  end_state.tail(3) = end_v;
  end_index = posToIndex(end_pt);
  cur_node->f_score = lambda_heu_ * estimateHeuristic(cur_node->state, end_state, time_to_goal);
  cur_node->node_state = IN_OPEN_SET;  // open_set_为优先级队列，比较的为 PathNodePtr 的fScore
  open_set_.push(cur_node);
  use_node_num_ += 1;
 
  // 考虑时间实时规划
  if (dynamic)
  {
    time_origin_ = time_start;
    cur_node->time = time_start;
    cur_node->time_idx = timeToIndex(time_start);   // 考虑时间因素的话，就选择带有时间信息的insert
    expanded_nodes_.insert(cur_node->index, cur_node->time_idx, cur_node);
    // cout << "time start: " << time_start << endl;
  }
  else //否则采用普通insert
    expanded_nodes_.insert(cur_node->index, cur_node);
 
  PathNodePtr neighbor = NULL;
  PathNodePtr terminate_node = NULL;
  bool init_search = init;
  const int tolerance = ceil(1 / resolution_); // ceil()函数：求不小于x的最小整数（向上取整）

前面是一些初始化，定义了起点和终点状态，将起点加入open_set

下面进入主循环

while (!open_set_.empty())
  {
    cur_node = open_set_.top();  // 取出cost最小的节点
 
    // Terminate? 判断扩展的节点是否距离终点过近
    bool reach_horizon = (cur_node->state.head(3) - start_pt).norm() >= horizon_;
    bool near_end = abs(cur_node->index(0) - end_index(0)) <= tolerance &&
                    abs(cur_node->index(1) - end_index(1)) <= tolerance &&
                    abs(cur_node->index(2) - end_index(2)) <= tolerance;
 
    if (reach_horizon || near_end)
    {
      terminate_node = cur_node;
      retrievePath(terminate_node);
      if (near_end)
      {
        // Check whether shot traj exist
        estimateHeuristic(cur_node->state, end_state, time_to_goal);
        computeShotTraj(cur_node->state, end_state, time_to_goal);
        if (init_search)
          ROS_ERROR("Shot in first search loop!");
      }
    }
    if (reach_horizon)
    {
      if (is_shot_succ_)
      {
        std::cout << "reach end" << std::endl;
        return REACH_END;
      }
      else
      {
        std::cout << "reach horizon" << std::endl;
        return REACH_HORIZON;
      }
    }
 
    if (near_end)
    {
      if (is_shot_succ_)
      {
        std::cout << "reach end" << std::endl;
        return REACH_END;
      }
      else if (cur_node->parent != NULL)
      {
        std::cout << "near end" << std::endl;
        return NEAR_END;
      }
      else
      {
        std::cout << "no path" << std::endl;
        return NO_PATH;
      }
    }

定义了一些边界条件，如果距离终点或者地图边界过近就返回

// 弹出节点
    open_set_.pop();
    cur_node->node_state = IN_CLOSE_SET;
    iter_num_ += 1;
 
    double res = 1 / 2.0, time_res = 1 / 1.0, time_res_init = 1 / 20.0;
    Eigen::Matrix<double, 6, 1> cur_state = cur_node->state;
    Eigen::Matrix<double, 6, 1> pro_state;
    vector<PathNodePtr> tmp_expand_nodes;
    Eigen::Vector3d um;
    double pro_t;
    vector<Eigen::Vector3d> inputs;
    vector<double> durations;
    if (init_search)
    {
      inputs.push_back(start_acc_);
      for (double tau = time_res_init * init_max_tau_; tau <= init_max_tau_ + 1e-3;
           tau += time_res_init * init_max_tau_)
        durations.push_back(tau);
      init_search = false;
    }
    else
    {
      // inputs里面放入所有可能的三维加速度
      for (double ax = -max_acc_; ax <= max_acc_ + 1e-3; ax += max_acc_ * res)
        for (double ay = -max_acc_; ay <= max_acc_ + 1e-3; ay += max_acc_ * res)
          for (double az = -max_acc_; az <= max_acc_ + 1e-3; az += max_acc_ * res)
          {
            um << ax, ay, az;
            inputs.push_back(um);
          }
      // duration里面放所有可能的tau
      for (double tau = time_res * max_tau_; tau <= max_tau_; tau += time_res * max_tau_)
        durations.push_back(tau);
    }

接下来弹出cost最小的节点，然后对周围所有可能的输入加速度以及时间间隔进行遍历，获取数组inputs以及durations

for (int i = 0; i < inputs.size(); ++i)
      for (int j = 0; j < durations.size(); ++j)
      {
        um = inputs[i];
        double tau = durations[j];
        stateTransit(cur_state, pro_state, um, tau);// 记录状态转移：从当前状态转移到pro_state，这个函数得到pro_state
        pro_t = cur_node->time + tau;  // 时间间隔
 
        Eigen::Vector3d pro_pos = pro_state.head(3);  // 取 pro_state 前3维，作为位置
 
        // Check if in close set
        Eigen::Vector3i pro_id = posToIndex(pro_pos);  // 坐标转换成idx
        int pro_t_id = timeToIndex(pro_t);  // 时间也转换一下
        // 哈希表搜索这个节点，如果是动态的，搜索加上时间
        PathNodePtr pro_node = dynamic ? expanded_nodes_.find(pro_id, pro_t_id) : expanded_nodes_.find(pro_id);
        // 如果没搜索到（==NULL），或者说当前节点在CLOSE LIST内了
        if (pro_node != NULL && pro_node->node_state == IN_CLOSE_SET)
        {
          if (init_search)  // 如果是
            std::cout << "close" << std::endl;
          continue;
        }

然后对所有可能的加速度和时间间隔，使用stateTransit函数计算可能转移到的状态

// 状态转移，从state0转移到state1
void KinodynamicAstar::stateTransit(Eigen::Matrix<double, 6, 1>& state0, Eigen::Matrix<double, 6, 1>& state1,
                                    Eigen::Vector3d um, double tau)
{
  for (int i = 0; i < 3; ++i)
    phi_(i, i + 3) = tau;
 
  Eigen::Matrix<double, 6, 1> integral; // 创建一个向量存增量
  integral.head(3) = 0.5 * pow(tau, 2) * um; //位置为δx = 0.5 * a * t^2 
  integral.tail(3) = tau * um;  // 速度为 δv = a * t
 
  state1 = phi_ * state0 + integral;
}

（3）搜索完毕

搜索完毕后，如果结果为NO_PATH，再次搜索一遍，如果再次为NO_PATH就结束

  if (status == KinodynamicAstar::NO_PATH) {
    cout << "[kino replan]: kinodynamic search fail!" << endl;
 
    // retry searching with discontinuous initial state
    kino_path_finder_->reset();
    status = kino_path_finder_->search(start_pt, start_vel, start_acc, end_pt, end_vel, false);
 
    if (status == KinodynamicAstar::NO_PATH) {
      cout << "[kino replan]: Can't find path." << endl;
      return false;
    } else {
      cout << "[kino replan]: retry search success." << endl;
    }
 
  } else {
    cout << "[kino replan]: kinodynamic search success." << endl;
  }

获取KinoAstar规划的路径，并以0.01的时间间隔采样，存储在kino_path_内

  plan_data_.kino_path_ = kino_path_finder_->getKinoTraj(0.01);
 
  t_search = (ros::Time::now() - t1).toSec();

3、B样条

下面就到了后端利用B样条来优化轨迹

（1）采样getSamples

  // parameterize the path to bspline
 
  double                  ts = pp_.ctrl_pt_dist / pp_.max_vel_;
  vector<Eigen::Vector3d> point_set, start_end_derivatives;
  kino_path_finder_->getSamples(ts, point_set, start_end_derivatives);

首先以最大的控制点距离（参数设置为0.5）除以最大速度，获取间隔时间，然后进入getSamples函数，首先计算所有节点之间的持续时间，加和得到T_sum

void KinodynamicAstar::getSamples(double& ts, vector<Eigen::Vector3d>& point_set,
                                  vector<Eigen::Vector3d>& start_end_derivatives)
{
  /* ---------- path duration ---------- */
  double T_sum = 0.0;
  if (is_shot_succ_)
    T_sum += t_shot_;
  PathNodePtr node = path_nodes_.back();
  while (node->parent != NULL)
  {
    T_sum += node->duration;
    node = node->parent;
  }
  // cout << "duration:" << T_sum << endl;

如果成功到达目标点，那么还要把到达目标点的轨迹时间也加上

下面计算末速度和末加速度

  // Calculate boundary vel and acc
  Eigen::Vector3d end_vel, end_acc;
  double t;
  if (is_shot_succ_)
  {
    t = t_shot_;
    end_vel = end_vel_;
    for (int dim = 0; dim < 3; ++dim)
    {
      Vector4d coe = coef_shot_.row(dim);
      end_acc(dim) = 2 * coe(2) + 6 * coe(3) * t_shot_;
    }
  }
  else
  {
    t = path_nodes_.back()->duration;
    end_vel = node->state.tail(3);
    end_acc = path_nodes_.back()->input;
  }

接下来进行采样

 // Get point samples
  int seg_num = floor(T_sum / ts);
  seg_num = max(8, seg_num);
  ts = T_sum / double(seg_num);
  bool sample_shot_traj = is_shot_succ_;
  node = path_nodes_.back();
 
  for (double ti = T_sum; ti > -1e-5; ti -= ts)
  {
    if (sample_shot_traj)
    {
      // samples on shot traj
      Vector3d coord;
      Vector4d poly1d, time;
 
      for (int j = 0; j < 4; j++)
        time(j) = pow(t, j);
 
      for (int dim = 0; dim < 3; dim++)
      {
        poly1d = coef_shot_.row(dim);
        coord(dim) = poly1d.dot(time);
      }
 
      point_set.push_back(coord);
      t -= ts;
 
      /* end of segment */
      if (t < -1e-5)
      {
        sample_shot_traj = false;
        if (node->parent != NULL)
          t += node->duration;
      }
    }
    else
    {
      // samples on searched traj
      Eigen::Matrix<double, 6, 1> x0 = node->parent->state;
      Eigen::Matrix<double, 6, 1> xt;
      Vector3d ut = node->input;
 
      stateTransit(x0, xt, ut, t);
 
      point_set.push_back(xt.head(3));
      t -= ts;
 
      // cout << "t: " << t << ", t acc: " << T_accumulate << endl;
      if (t < -1e-5 && node->parent->parent != NULL)
      {
        node = node->parent;
        t += node->duration;
      }
    }
  }
  reverse(point_set.begin(), point_set.end());

分别对到达目标点的轨迹和中间的轨迹进行采样，间隔时间为ts，最后存在point_set内

ps：这里的ts重新计算，先根据T_sum和之前的ts算出中间插值样条线的数量seg_num，然后再算出每个seg的平均时间作为ts

最后计算起始加速度，然后获取start_end_derivatives

  // calculate start acc
  Eigen::Vector3d start_acc;
  if (path_nodes_.back()->parent == NULL)
  {
    // no searched traj, calculate by shot traj
    start_acc = 2 * coef_shot_.col(2);
  }
  else
  {
    // input of searched traj
    start_acc = node->input;
  }
 
  start_end_derivatives.push_back(start_vel_);
  start_end_derivatives.push_back(end_vel);
  start_end_derivatives.push_back(start_acc);
  start_end_derivatives.push_back(end_acc);

起始加速度，如果中间没有搜索的轨迹，使用shot轨迹来计算，否则就采用搜索轨迹的第一个点的input作为起始加速度。

start_end_derivatives是一个size为4的vector，存着起始和末速度加速度

（2）B样条控制点生成parameterizeToBspline

创建一个存储B样条控制点的矩阵ctrl_pts，然后运行parameterizeToBspline函数来获取控制点

  Eigen::MatrixXd ctrl_pts;
  NonUniformBspline::parameterizeToBspline(ts, point_set, start_end_derivatives, ctrl_pts);
  NonUniformBspline init(ctrl_pts, 3, ts);

下面是parameterizeToBspline的代码

  int K = point_set.size();
 
  // write A
  Eigen::Vector3d prow(3), vrow(3), arow(3);
  prow << 1, 4, 1;
  vrow << -1, 0, 1;
  arow << 1, -2, 1;
 
  Eigen::MatrixXd A = Eigen::MatrixXd::Zero(K + 4, K + 2);
 
  for (int i = 0; i < K; ++i) A.block(i, i, 1, 3) = (1 / 6.0) * prow.transpose();
 
  A.block(K, 0, 1, 3)         = (1 / 2.0 / ts) * vrow.transpose();
  A.block(K + 1, K - 1, 1, 3) = (1 / 2.0 / ts) * vrow.transpose();
 
  A.block(K + 2, 0, 1, 3)     = (1 / ts / ts) * arow.transpose();
  A.block(K + 3, K - 1, 1, 3) = (1 / ts / ts) * arow.transpose();

计算A矩阵

  // write b
  Eigen::VectorXd bx(K + 4), by(K + 4), bz(K + 4);
  for (int i = 0; i < K; ++i) {
    bx(i) = point_set[i](0);
    by(i) = point_set[i](1);
    bz(i) = point_set[i](2);
  }
 
  for (int i = 0; i < 4; ++i) {
    bx(K + i) = start_end_derivative[i](0);
    by(K + i) = start_end_derivative[i](1);
    bz(K + i) = start_end_derivative[i](2);
  }

计算b矩阵

求解方程Ax=b

  // solve Ax = b
  Eigen::VectorXd px = A.colPivHouseholderQr().solve(bx);
  Eigen::VectorXd py = A.colPivHouseholderQr().solve(by);
  Eigen::VectorXd pz = A.colPivHouseholderQr().solve(bz);

求出来解就是控制点的坐标

  // convert to control pts
  ctrl_pts.resize(K + 2, 3);
  ctrl_pts.col(0) = px;
  ctrl_pts.col(1) = py;
  ctrl_pts.col(2) = pz;

（3）B样条轨迹优化BsplineOptimizeTraj

执行函数BsplineOptimizeTraj来对控制点进行轨迹优化

  // bspline trajectory optimization
 
  t1 = ros::Time::now();
 
  int cost_function = BsplineOptimizer::NORMAL_PHASE;
 
  if (status != KinodynamicAstar::REACH_END) {
    cost_function |= BsplineOptimizer::ENDPOINT;
  }
 
  ctrl_pts = bspline_optimizers_[0]->BsplineOptimizeTraj(ctrl_pts, ts, cost_function, 1, 1);
 
  t_opt = (ros::Time::now() - t1).toSec();

设置 cost_function为BsplineOptimizer::NORMAL_PHASE;

const int BsplineOptimizer::NORMAL_PHASE =
    BsplineOptimizer::SMOOTHNESS | BsplineOptimizer::DISTANCE | BsplineOptimizer::FEASIBILITY;

其中包含了三种优化，一种是平滑性，一种是距离，一种是可行性。

如果没有到达终点，再增加ENDPOINT优化。

接下来执行优化

Eigen::MatrixXd BsplineOptimizer::BsplineOptimizeTraj(const Eigen::MatrixXd& points, const double& ts,
                                                      const int& cost_function, int max_num_id,
                                                      int max_time_id) {
  setControlPoints(points);
  setBsplineInterval(ts);
  setCostFunction(cost_function);
  setTerminateCond(max_num_id, max_time_id);
 
  optimize();
  return this->control_points_;
}

首先设置控制点和时间间隔，维度为3维，间隔时间ts。然后设置CostFunction为目标函数。max_num_id和max_time_id都为1。

void BsplineOptimizer::setControlPoints(const Eigen::MatrixXd& points) {
  control_points_ = points;
  dim_            = control_points_.cols();
}
void BsplineOptimizer::setBsplineInterval(const double& ts) { bspline_interval_ = ts; }
 
void BsplineOptimizer::setCostFunction(const int& cost_code) {
  cost_function_ = cost_code;
 
  // print optimized cost function
  string cost_str;
  if (cost_function_ & SMOOTHNESS) cost_str += "smooth |";
  if (cost_function_ & DISTANCE) cost_str += " dist  |";
  if (cost_function_ & FEASIBILITY) cost_str += " feasi |";
  if (cost_function_ & ENDPOINT) cost_str += " endpt |";
  if (cost_function_ & GUIDE) cost_str += " guide |";
  if (cost_function_ & WAYPOINTS) cost_str += " waypt |";
 
  ROS_INFO_STREAM("cost func: " << cost_str);
}
 
void BsplineOptimizer::setTerminateCond(const int& max_num_id, const int& max_time_id) {
  max_num_id_  = max_num_id;
  max_time_id_ = max_time_id;
}

然后执行优化，返回优化后的控制点坐标。

（4）控制点时间间隔调整reallocateTime

优化完控制点坐标，需要进行时间间隔调整

  // iterative time adjustment
 
  t1                    = ros::Time::now();
  NonUniformBspline pos = NonUniformBspline(ctrl_pts, 3, ts);
 
  double to = pos.getTimeSum();
  pos.setPhysicalLimits(pp_.max_vel_, pp_.max_acc_);
  bool feasible = pos.checkFeasibility(false);
 
  int iter_num = 0;
  while (!feasible && ros::ok()) {
 
    feasible = pos.reallocateTime();
 
    if (++iter_num >= 3) break;
  }
 
  // pos.checkFeasibility(true);
  // cout << "[Main]: iter num: " << iter_num << endl;
 
  double tn = pos.getTimeSum();
 
  cout << "[kino replan]: Reallocate ratio: " << tn / to << endl;
  if (tn / to > 3.0) ROS_ERROR("reallocate error.");
 
  t_adjust = (ros::Time::now() - t1).toSec();

创建一个B样条pos来存储优化后的B样条，进行时间重分配

我们先来看一下这个B样条怎么创建的

NonUniformBspline::NonUniformBspline(const Eigen::MatrixXd& points, const int& order,
                                     const double& interval) {
  setUniformBspline(points, order, interval);
}
void NonUniformBspline::setUniformBspline(const Eigen::MatrixXd& points, const int& order,
                                          const double& interval) {
  control_points_ = points;
  p_              = order;
  interval_       = interval;
 
  n_ = points.rows() - 1;
  m_ = n_ + p_ + 1;
 
  u_ = Eigen::VectorXd::Zero(m_ + 1);
  for (int i = 0; i <= m_; ++i) {
 
    if (i <= p_) {
      u_(i) = double(-p_ + i) * interval_;
    } else if (i > p_ && i <= m_ - p_) {
      u_(i) = u_(i - 1) + interval_;
    } else if (i > m_ - p_) {
      u_(i) = u_(i - 1) + interval_;
    }
  }
}

根据构造函数，可以知道它把维度赋给了p_，然后时间间隔为interval_，

然后使用n_存储控制点的数量，m_存储控制点数量+维度+1，u_来存储每个控制点之间的时间间隔

to来存储控制点之间的时间总和

double NonUniformBspline::getTimeSum() {
  double tm, tmp;
  getTimeSpan(tm, tmp);
  return tmp - tm;
}
void NonUniformBspline::getTimeSpan(double& um, double& um_p) {
  um   = u_(p_);
  um_p = u_(m_ - p_);
}

为pos设置最大的限制   pos.setPhysicalLimits(pp_.max_vel_, pp_.max_acc_);

然后检测pos的可行性

最后进行三次重分配时间的计算

bool NonUniformBspline::reallocateTime(bool show) {
  // SETY << "[Bspline]: total points size: " << control_points_.rows() << endl;
  // cout << "origin knots:\n" << u_.transpose() << endl;
  bool fea = true;
 
  Eigen::MatrixXd P         = control_points_;
  int             dimension = control_points_.cols();
 
  double max_vel, max_acc;

使用P存储控制点，dimension为3

检测速度可行性，然后插入中间点

  /* check vel feasibility and insert points */
  for (int i = 0; i < P.rows() - 1; ++i) {
    Eigen::VectorXd vel = p_ * (P.row(i + 1) - P.row(i)) / (u_(i + p_ + 1) - u_(i + 1));
 
    if (fabs(vel(0)) > limit_vel_ + 1e-4 || fabs(vel(1)) > limit_vel_ + 1e-4 ||
        fabs(vel(2)) > limit_vel_ + 1e-4) {
 
      fea = false;
      if (show) cout << "[Realloc]: Infeasible vel " << i << " :" << vel.transpose() << endl;
 
      max_vel = -1.0;
      for (int j = 0; j < dimension; ++j) {
        max_vel = max(max_vel, fabs(vel(j)));
      }
 
      double ratio = max_vel / limit_vel_ + 1e-4;
      if (ratio > limit_ratio_) ratio = limit_ratio_;
 
      double time_ori = u_(i + p_ + 1) - u_(i + 1);
      double time_new = ratio * time_ori;
      double delta_t  = time_new - time_ori;
      double t_inc    = delta_t / double(p_);
 
      for (int j = i + 2; j <= i + p_ + 1; ++j) {
        u_(j) += double(j - i - 1) * t_inc;
        if (j <= 5 && j >= 1) {
          // cout << "vel j: " << j << endl;
        }
      }
 
      for (int j = i + p_ + 2; j < u_.rows(); ++j) {
        u_(j) += delta_t;
      }
    }
  }

遍历每个控制点，计算两个控制点之间的速度 v = x/t。如果速度没有越limit_vel_界，则以这个速度作为重分配的时间的速度。

然后遍历每个维度，找到三个维度上最大的速度max_vel，计算max_vel和limit_vel_的比值ratio，如果这个重分配的比值ratio大于limit_ratio_，就让它等于limit_ratio_。

接下来计算原来的间隔时间time_ori，给它乘以比率得到time_new，计算差值delta_t，

后面加速度同理

  /* acc feasibility */
  for (int i = 0; i < P.rows() - 2; ++i) {
 
    Eigen::VectorXd acc = p_ * (p_ - 1) *
        ((P.row(i + 2) - P.row(i + 1)) / (u_(i + p_ + 2) - u_(i + 2)) -
         (P.row(i + 1) - P.row(i)) / (u_(i + p_ + 1) - u_(i + 1))) /
        (u_(i + p_ + 1) - u_(i + 2));
 
    if (fabs(acc(0)) > limit_acc_ + 1e-4 || fabs(acc(1)) > limit_acc_ + 1e-4 ||
        fabs(acc(2)) > limit_acc_ + 1e-4) {
 
      fea = false;
      if (show) cout << "[Realloc]: Infeasible acc " << i << " :" << acc.transpose() << endl;
 
      max_acc = -1.0;
      for (int j = 0; j < dimension; ++j) {
        max_acc = max(max_acc, fabs(acc(j)));
      }
 
      double ratio = sqrt(max_acc / limit_acc_) + 1e-4;
      if (ratio > limit_ratio_) ratio = limit_ratio_;
      // cout << "ratio: " << ratio << endl;
 
      double time_ori = u_(i + p_ + 1) - u_(i + 2);
      double time_new = ratio * time_ori;
      double delta_t  = time_new - time_ori;
      double t_inc    = delta_t / double(p_ - 1);
 
      if (i == 1 || i == 2) {
        // cout << "acc i: " << i << endl;
        for (int j = 2; j <= 5; ++j) {
          u_(j) += double(j - 1) * t_inc;
        }
 
        for (int j = 6; j < u_.rows(); ++j) {
          u_(j) += 4.0 * t_inc;
        }
 
      } else {
 
        for (int j = i + 3; j <= i + p_ + 1; ++j) {
          u_(j) += double(j - i - 2) * t_inc;
          if (j <= 5 && j >= 1) {
            // cout << "acc j: " << j << endl;
          }
        }
 
        for (int j = i + p_ + 2; j < u_.rows(); ++j) {
          u_(j) += delta_t;
        }
      }
    }
  }
 
  return fea;

（5）更新轨迹updateTrajInfo

调整完时间间隔，再把新的轨迹更新上去

  double tn = pos.getTimeSum();
 
  cout << "[kino replan]: Reallocate ratio: " << tn / to << endl;
  if (tn / to > 3.0) ROS_ERROR("reallocate error.");
 
  t_adjust = (ros::Time::now() - t1).toSec();
 
  // save planned results
 
  local_data_.position_traj_ = pos;
 
  double t_total = t_search + t_opt + t_adjust;
  cout << "[kino replan]: time: " << t_total << ", search: " << t_search << ", optimize: " << t_opt
       << ", adjust time:" << t_adjust << endl;
 
  pp_.time_search_   = t_search;
  pp_.time_optimize_ = t_opt;
  pp_.time_adjust_   = t_adjust;
 
  updateTrajInfo();
 
  return true;

首先获取新的时间总和tn，然后将pos作为新的轨迹存储在position_traj_

接下来更新轨迹的信息

void FastPlannerManager::updateTrajInfo() {
  local_data_.velocity_traj_     = local_data_.position_traj_.getDerivative();
  local_data_.acceleration_traj_ = local_data_.velocity_traj_.getDerivative();
  local_data_.start_pos_         = local_data_.position_traj_.evaluateDeBoorT(0.0);
  local_data_.duration_          = local_data_.position_traj_.getTimeSum();
  local_data_.traj_id_ += 1;
}

获取速度轨迹，加速度轨迹，起点坐标，持续时间，然后给轨迹id+1

已知位置的轨迹，获取速度的轨迹

NonUniformBspline NonUniformBspline::getDerivative() {
  Eigen::MatrixXd   ctp = getDerivativeControlPoints();
  NonUniformBspline derivative(ctp, p_ - 1, interval_);
 
  /* cut the first and last knot */
  Eigen::VectorXd knot(u_.rows() - 2);
  knot = u_.segment(1, u_.rows() - 2);
  derivative.setKnot(knot);
 
  return derivative;
}