VelocityConstraint.h

Go to the documentation of this file.

/* ----------------------------------------------------------------------------
 
 * 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
 
 * -------------------------------------------------------------------------- */

 
#pragma once
 
#include <gtsam/config.h>
 
#ifdef GTSAM_ALLOW_DEPRECATED_SINCE_V43
 
#include <gtsam/base/numericalDerivative.h>
#include <gtsam/nonlinear/NonlinearFactor.h>
#include <gtsam/nonlinear/NoiseModelFactorN.h>
#include <gtsam_unstable/dynamics/PoseRTV.h>
 
#include <cassert>
 
namespace gtsam {
 
namespace dynamics {

typedef enum {
  TRAPEZOIDAL, // Constant acceleration
  EULER_START, // Constant velocity, using starting velocity
  EULER_END    // Constant velocity, using ending velocity
} IntegrationMode;
 
}

class VelocityConstraint
    : public gtsam::NoiseModelFactorT<gtsam::Vector3, PoseRTV, PoseRTV> {
public:
  typedef gtsam::NoiseModelFactorT<gtsam::Vector3, PoseRTV, PoseRTV> Base;
 
  // Provide access to the Matrix& version of evaluateError:
  using Base::evaluateError;
 
protected:
 
  double dt_;   
  dynamics::IntegrationMode integration_mode_;  
 
public:

  VelocityConstraint(Key key1, Key key2, const dynamics::IntegrationMode& mode,
      double dt, double mu = 1000)
  : Base(noiseModel::Constrained::All(3, mu), key1, key2), dt_(dt), integration_mode_(mode) {}

    VelocityConstraint(Key key1, Key key2, double dt, double mu = 1000)
    : Base(noiseModel::Constrained::All(3, mu), key1, key2),
      dt_(dt), integration_mode_(dynamics::TRAPEZOIDAL) {}

  VelocityConstraint(Key key1, Key key2, const dynamics::IntegrationMode& mode,
      double dt, const gtsam::SharedNoiseModel& model)
  : Base(model, key1, key2), dt_(dt), integration_mode_(mode) {}

  VelocityConstraint(Key key1, Key key2, double dt, const gtsam::SharedNoiseModel& model)
  : Base(model, key1, key2), dt_(dt), integration_mode_(dynamics::TRAPEZOIDAL) {}
 
  ~VelocityConstraint() override {}

  gtsam::NonlinearFactor::shared_ptr clone() const override {
    return std::static_pointer_cast<gtsam::NonlinearFactor>(
        gtsam::NonlinearFactor::shared_ptr(new VelocityConstraint(*this))); }

  gtsam::Vector3 evaluateError(const PoseRTV& x1, const PoseRTV& x2,
                               OptionalMatrixType H1,
                               OptionalMatrixType H2) const override {
    if (H1) *H1 = gtsam::numericalDerivative21<gtsam::Vector,PoseRTV,PoseRTV>(
        std::bind(VelocityConstraint::evaluateError_, std::placeholders::_1,
            std::placeholders::_2, dt_, integration_mode_), x1, x2, 1e-5);
    if (H2) *H2 = gtsam::numericalDerivative22<gtsam::Vector,PoseRTV,PoseRTV>(
        std::bind(VelocityConstraint::evaluateError_, std::placeholders::_1,
            std::placeholders::_2, dt_, integration_mode_), x1, x2, 1e-5);
    return evaluateError_(x1, x2, dt_, integration_mode_);
  }
 
  void print(const std::string& s = "", const gtsam::KeyFormatter& formatter = gtsam::DefaultKeyFormatter) const override {
    std::string a = "VelocityConstraint: " + s;
    Base::print(a, formatter);
    switch(integration_mode_) {
    case dynamics::TRAPEZOIDAL: std::cout << "Integration: Trapezoidal\n"; break;
    case dynamics::EULER_START: std::cout << "Integration: Euler (start)\n"; break;
    case dynamics::EULER_END: std::cout << "Integration: Euler (end)\n"; break;
    default: std::cout << "Integration: Unknown\n" << std::endl; break;
    }
    std::cout << "dt: " << dt_ << std::endl;
  }
 
private:
  static gtsam::Vector evaluateError_(const PoseRTV& x1, const PoseRTV& x2,
      double dt, const dynamics::IntegrationMode& mode) {
 
    const Velocity3& v1 = x1.v(), v2 = x2.v();
    const Point3& p1 = x1.t(), p2 = x2.t();
    Point3 hx(0,0,0);
    switch(mode) {
    case dynamics::TRAPEZOIDAL: hx = p1 + Point3((v1 + v2) * dt *0.5); break;
    case dynamics::EULER_START: hx = p1 + Point3(v1 * dt); break;
    case dynamics::EULER_END  : hx = p1 + Point3(v2 * dt); break;
    default: assert(false); break;
    }
    return p2 - hx;
  }
};
 
} // \namespace gtsam
 
#endif  // GTSAM_ALLOW_DEPRECATED_SINCE_V43