239 lines
		
	
	
		
			8.3 KiB
		
	
	
	
		
			C++
		
	
	
			
		
		
	
	
			239 lines
		
	
	
		
			8.3 KiB
		
	
	
	
		
			C++
		
	
	
/* ----------------------------------------------------------------------------
 | 
						|
 | 
						|
 * GTSAM Copyright 2010, Georgia Tech Research Corporation,
 | 
						|
 * Atlanta, Georgia 30332-0415
 | 
						|
 * All Rights Reserved
 | 
						|
 * Authors: Frank Dellaert, et al. (see THANKS for the full author list)
 | 
						|
 | 
						|
 * See LICENSE for the license information
 | 
						|
 | 
						|
 * -------------------------------------------------------------------------- */
 | 
						|
 | 
						|
/**
 | 
						|
 * @file RangeISAMExample_plaza2.cpp
 | 
						|
 * @brief A 2D Range SLAM example
 | 
						|
 * @date June 20, 2013
 | 
						|
 * @author Frank Dellaert
 | 
						|
 */
 | 
						|
 | 
						|
// Both relative poses and recovered trajectory poses will be stored as Pose2.
 | 
						|
#include <gtsam/geometry/Pose2.h>
 | 
						|
using gtsam::Pose2;
 | 
						|
 | 
						|
// gtsam::Vectors are dynamic Eigen vectors, Vector3 is statically sized.
 | 
						|
#include <gtsam/base/Vector.h>
 | 
						|
using gtsam::Vector;
 | 
						|
using gtsam::Vector3;
 | 
						|
 | 
						|
// Unknown landmarks are of type Point2 (which is just a 2D Eigen vector).
 | 
						|
#include <gtsam/geometry/Point2.h>
 | 
						|
using gtsam::Point2;
 | 
						|
 | 
						|
// Each variable in the system (poses and landmarks) must be identified with a
 | 
						|
// unique key. We can either use simple integer keys (1, 2, 3, ...) or symbols
 | 
						|
// (X1, X2, L1). Here we will use Symbols.
 | 
						|
#include <gtsam/inference/Symbol.h>
 | 
						|
using gtsam::Symbol;
 | 
						|
 | 
						|
// We want to use iSAM2 to solve the range-SLAM problem incrementally.
 | 
						|
#include <gtsam/nonlinear/ISAM2.h>
 | 
						|
 | 
						|
// iSAM2 requires as input a set set of new factors to be added stored in a
 | 
						|
// factor graph, and initial guesses for any new variables in the added factors.
 | 
						|
#include <gtsam/nonlinear/NonlinearFactorGraph.h>
 | 
						|
#include <gtsam/nonlinear/Values.h>
 | 
						|
 | 
						|
// We will use a non-linear solver to batch-initialize from the first 150 frames
 | 
						|
#include <gtsam/nonlinear/LevenbergMarquardtOptimizer.h>
 | 
						|
 | 
						|
// Measurement functions are represented as 'factors'. Several common factors
 | 
						|
// have been provided with the library for solving robotics SLAM problems:
 | 
						|
#include <gtsam/sam/RangeFactor.h>
 | 
						|
#include <gtsam/slam/BetweenFactor.h>
 | 
						|
#include <gtsam/slam/dataset.h>
 | 
						|
 | 
						|
// Timing, with functions below, provides nice facilities to benchmark.
 | 
						|
#include <gtsam/base/timing.h>
 | 
						|
using gtsam::tictoc_print_;
 | 
						|
 | 
						|
// Standard headers, added last, so we know headers above work on their own.
 | 
						|
#include <fstream>
 | 
						|
#include <iostream>
 | 
						|
#include <random>
 | 
						|
#include <set>
 | 
						|
#include <string>
 | 
						|
#include <utility>
 | 
						|
#include <vector>
 | 
						|
 | 
						|
namespace NM = gtsam::noiseModel;
 | 
						|
 | 
						|
// Data is second UWB ranging dataset, B2 or "plaza 2", from
 | 
						|
// "Navigating with Ranging Radios: Five Data Sets with Ground Truth"
 | 
						|
// by Joseph Djugash, Bradley Hamner, and Stephan Roth
 | 
						|
// https://www.ri.cmu.edu/pub_files/2009/9/Final_5datasetsRangingRadios.pdf
 | 
						|
 | 
						|
// load the odometry
 | 
						|
// DR: Odometry Input (delta distance traveled and delta heading change)
 | 
						|
//    Time (sec)  Delta Distance Traveled (m) Delta Heading (rad)
 | 
						|
using TimedOdometry = std::pair<double, Pose2>;
 | 
						|
std::list<TimedOdometry> readOdometry() {
 | 
						|
  std::list<TimedOdometry> odometryList;
 | 
						|
  std::string data_file = gtsam::findExampleDataFile("Plaza2_DR.txt");
 | 
						|
  std::ifstream is(data_file.c_str());
 | 
						|
 | 
						|
  while (is) {
 | 
						|
    double t, distance_traveled, delta_heading;
 | 
						|
    is >> t >> distance_traveled >> delta_heading;
 | 
						|
    odometryList.emplace_back(t, Pose2(distance_traveled, 0, delta_heading));
 | 
						|
  }
 | 
						|
  is.clear(); /* clears the end-of-file and error flags */
 | 
						|
  return odometryList;
 | 
						|
}
 | 
						|
 | 
						|
// load the ranges from TD
 | 
						|
//    Time (sec)  Sender / Antenna ID Receiver Node ID  Range (m)
 | 
						|
using RangeTriple = std::tuple<double, size_t, double>;
 | 
						|
std::vector<RangeTriple> readTriples() {
 | 
						|
  std::vector<RangeTriple> triples;
 | 
						|
  std::string data_file = gtsam::findExampleDataFile("Plaza2_TD.txt");
 | 
						|
  std::ifstream is(data_file.c_str());
 | 
						|
 | 
						|
  while (is) {
 | 
						|
    double t, range, sender, receiver;
 | 
						|
    is >> t >> sender >> receiver >> range;
 | 
						|
    triples.emplace_back(t, receiver, range);
 | 
						|
  }
 | 
						|
  is.clear(); /* clears the end-of-file and error flags */
 | 
						|
  return triples;
 | 
						|
}
 | 
						|
 | 
						|
// main
 | 
						|
int main(int argc, char** argv) {
 | 
						|
  // load Plaza2 data
 | 
						|
  std::list<TimedOdometry> odometry = readOdometry();
 | 
						|
  size_t M = odometry.size();
 | 
						|
  std::cout << "Read " << M << " odometry entries." << std::endl;
 | 
						|
 | 
						|
  std::vector<RangeTriple> triples = readTriples();
 | 
						|
  size_t K = triples.size();
 | 
						|
  std::cout << "Read " << K << " range triples." << std::endl;
 | 
						|
 | 
						|
  // parameters
 | 
						|
  size_t minK =
 | 
						|
      150;  // minimum number of range measurements to process initially
 | 
						|
  size_t incK = 25;  // minimum number of range measurements to process after
 | 
						|
  bool robust = true;
 | 
						|
 | 
						|
  // Set Noise parameters
 | 
						|
  Vector priorSigmas = Vector3(1, 1, M_PI);
 | 
						|
  Vector odoSigmas = Vector3(0.05, 0.01, 0.1);
 | 
						|
  double sigmaR = 100;        // range standard deviation
 | 
						|
  const NM::Base::shared_ptr  // all same type
 | 
						|
      priorNoise = NM::Diagonal::Sigmas(priorSigmas),  // prior
 | 
						|
      looseNoise = NM::Isotropic::Sigma(2, 1000),      // loose LM prior
 | 
						|
      odoNoise = NM::Diagonal::Sigmas(odoSigmas),      // odometry
 | 
						|
      gaussian = NM::Isotropic::Sigma(1, sigmaR),      // non-robust
 | 
						|
      tukey = NM::Robust::Create(NM::mEstimator::Tukey::Create(15),
 | 
						|
                                 gaussian),  // robust
 | 
						|
      rangeNoise = robust ? tukey : gaussian;
 | 
						|
 | 
						|
  // Initialize iSAM
 | 
						|
  gtsam::ISAM2 isam;
 | 
						|
 | 
						|
  // Add prior on first pose
 | 
						|
  Pose2 pose0 = Pose2(-34.2086489999201, 45.3007639991120, M_PI - 2.021089);
 | 
						|
  gtsam::NonlinearFactorGraph newFactors;
 | 
						|
  newFactors.addPrior(0, pose0, priorNoise);
 | 
						|
  gtsam::Values initial;
 | 
						|
  initial.insert(0, pose0);
 | 
						|
 | 
						|
  // We will initialize landmarks randomly, and keep track of which landmarks we
 | 
						|
  // already added with a set.
 | 
						|
  std::mt19937_64 rng;
 | 
						|
  std::normal_distribution<double> normal(0.0, 100.0);
 | 
						|
  std::set<Symbol> initializedLandmarks;
 | 
						|
 | 
						|
  // set some loop variables
 | 
						|
  size_t i = 1;  // step counter
 | 
						|
  size_t k = 0;  // range measurement counter
 | 
						|
  bool initialized = false;
 | 
						|
  Pose2 lastPose = pose0;
 | 
						|
  size_t countK = 0;
 | 
						|
 | 
						|
  // Loop over odometry
 | 
						|
  gttic_(iSAM);
 | 
						|
  for (const TimedOdometry& timedOdometry : odometry) {
 | 
						|
    //--------------------------------- odometry loop --------------------------
 | 
						|
    double t;
 | 
						|
    Pose2 odometry;
 | 
						|
    std::tie(t, odometry) = timedOdometry;
 | 
						|
 | 
						|
    // add odometry factor
 | 
						|
    newFactors.emplace_shared<gtsam::BetweenFactor<Pose2>>(i - 1, i, odometry,
 | 
						|
                                                           odoNoise);
 | 
						|
 | 
						|
    // predict pose and add as initial estimate
 | 
						|
    Pose2 predictedPose = lastPose.compose(odometry);
 | 
						|
    lastPose = predictedPose;
 | 
						|
    initial.insert(i, predictedPose);
 | 
						|
 | 
						|
    // Check if there are range factors to be added
 | 
						|
    while (k < K && t >= std::get<0>(triples[k])) {
 | 
						|
      size_t j = std::get<1>(triples[k]);
 | 
						|
      Symbol landmark_key('L', j);
 | 
						|
      double range = std::get<2>(triples[k]);
 | 
						|
      newFactors.emplace_shared<gtsam::RangeFactor<Pose2, Point2>>(
 | 
						|
          i, landmark_key, range, rangeNoise);
 | 
						|
      if (initializedLandmarks.count(landmark_key) == 0) {
 | 
						|
        std::cout << "adding landmark " << j << std::endl;
 | 
						|
        double x = normal(rng), y = normal(rng);
 | 
						|
        initial.insert(landmark_key, Point2(x, y));
 | 
						|
        initializedLandmarks.insert(landmark_key);
 | 
						|
        // We also add a very loose prior on the landmark in case there is only
 | 
						|
        // one sighting, which cannot fully determine the landmark.
 | 
						|
        newFactors.emplace_shared<gtsam::PriorFactor<Point2>>(
 | 
						|
            landmark_key, Point2(0, 0), looseNoise);
 | 
						|
      }
 | 
						|
      k = k + 1;
 | 
						|
      countK = countK + 1;
 | 
						|
    }
 | 
						|
 | 
						|
    // Check whether to update iSAM 2
 | 
						|
    if ((k > minK) && (countK > incK)) {
 | 
						|
      if (!initialized) {  // Do a full optimize for first minK ranges
 | 
						|
        std::cout << "Initializing at time " << k << std::endl;
 | 
						|
        gttic_(batchInitialization);
 | 
						|
        gtsam::LevenbergMarquardtOptimizer batchOptimizer(newFactors, initial);
 | 
						|
        initial = batchOptimizer.optimize();
 | 
						|
        gttoc_(batchInitialization);
 | 
						|
        initialized = true;
 | 
						|
      }
 | 
						|
      gttic_(update);
 | 
						|
      isam.update(newFactors, initial);
 | 
						|
      gttoc_(update);
 | 
						|
      gttic_(calculateEstimate);
 | 
						|
      gtsam::Values result = isam.calculateEstimate();
 | 
						|
      gttoc_(calculateEstimate);
 | 
						|
      lastPose = result.at<Pose2>(i);
 | 
						|
      newFactors = gtsam::NonlinearFactorGraph();
 | 
						|
      initial = gtsam::Values();
 | 
						|
      countK = 0;
 | 
						|
    }
 | 
						|
    i += 1;
 | 
						|
    //--------------------------------- odometry loop --------------------------
 | 
						|
  }  // end for
 | 
						|
  gttoc_(iSAM);
 | 
						|
 | 
						|
  // Print timings
 | 
						|
  tictoc_print_();
 | 
						|
 | 
						|
  // Print optimized landmarks:
 | 
						|
  gtsam::Values finalResult = isam.calculateEstimate();
 | 
						|
  for (auto&& landmark_key : initializedLandmarks) {
 | 
						|
    Point2 p = finalResult.at<Point2>(landmark_key);
 | 
						|
    std::cout << landmark_key << ":" << p.transpose() << "\n";
 | 
						|
  }
 | 
						|
 | 
						|
  exit(0);
 | 
						|
}
 |