Skip to content

Commit

Permalink
Merge branch 'master' of https://git.lcsr.jhu.edu/mkobila1/gcop
Browse files Browse the repository at this point in the history
  • Loading branch information
msheckells committed Jul 6, 2015
2 parents 79d74d4 + ccffc61 commit ac63468
Show file tree
Hide file tree
Showing 20 changed files with 333 additions and 77 deletions.
4 changes: 2 additions & 2 deletions bin/BulletVehicleEstimation/VehicleDemo.h
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ class btCollisionShape;
#include <stdlib.h>

#include "point3dgps.h"
#include "gndoepv2.h"
#include "deterministicgndoep.h"
#include "params.h"
#include "lqsensorcost.h"
#include "utils.h"
Expand All @@ -41,7 +41,7 @@ using namespace Eigen;
class VehicleDemo : public GlutDemoApplication
{
typedef LqSensorCost<Vector4d, 4, 2, Dynamic, Dynamic, Vector3d, 3> RccarCost;
typedef GnDoep1<Vector4d, 4, 2, Dynamic, Dynamic, Vector3d, 3, Point3dState, 6> RccarDoep;
typedef GnDoep<Vector4d, 4, 2, Dynamic, Dynamic, Vector3d, 3, Point3dState, 6> RccarDoep;

protected:
int iters; //Number of iterations
Expand Down
4 changes: 2 additions & 2 deletions bin/BulletZVehicleTest/VehicleDemo.h
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ class btCollisionShape;
#include <stdlib.h>

#include "point3dgps.h"
#include "gndoepv2.h"
#include "deterministicgndoep.h"
#include "params.h"
#include "lqsensorcost.h"
#include "utils.h"
Expand All @@ -41,7 +41,7 @@ using namespace Eigen;
class VehicleDemo : public GlutDemoApplication
{
typedef LqSensorCost<Vector4d, 4, 2, Dynamic, 10, Vector3d, 3> RccarCost;
typedef GnDoep1<Vector4d, 4, 2, Dynamic, 10, Vector3d, 3, Point3dState, 6> RccarDoep;
typedef GnDoep<Vector4d, 4, 2, Dynamic, 10, Vector3d, 3, Point3dState, 6> RccarDoep;

protected:
int iters; //Number of iterations
Expand Down
3 changes: 3 additions & 0 deletions bin/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -215,6 +215,9 @@ target_link_libraries(gpadd gcop_est gcop_algos gcop_views gcop_systems ${ALL_LI
add_executable(imuekftest imuekftest.cc)
target_link_libraries(imuekftest gcop_est gcop_algos gcop_views gcop_systems ${ALL_LIBS})

add_executable(imumagekftest imumagekftest.cc)
target_link_libraries(imumagekftest gcop_est gcop_algos gcop_views gcop_systems ${ALL_LIBS})

add_executable(insekftest insekftest.cc)
target_link_libraries(insekftest gcop_est gcop_algos gcop_views gcop_systems ${ALL_LIBS})

Expand Down
19 changes: 10 additions & 9 deletions bin/body3dcedemstab.cc
Original file line number Diff line number Diff line change
Expand Up @@ -199,28 +199,29 @@ void solver_process(Viewer* viewer)
params.GetBool("bAuto", ce.ce.bAuto);
params.GetDouble("b", ce.ce.b);

Body3dView<6> dslView(sys, &ctrl.xds);
dslView.lineWidth=5;
dslView.rgba[0] = 0; dslView.rgba[1] = 1; dslView.rgba[2] = 0;
if (viewer)
viewer->Add(dslView);

getchar();


Body3dView<> view(sys, &xs);
view.renderSystem = false;
viewer->Add(view);
view.rgba[0] = 0; view.rgba[1] = 0; view.rgba[2] = 1;

Body3dCeView ceview(ce, view);
viewer->Add(ceview);

Body3dView<6> dslView(sys, &ctrl.xds);
dslView.lineWidth=5;
dslView.rgba[0] = 0; dslView.rgba[1] = 1; dslView.rgba[2] = 0;
if (viewer)
viewer->Add(dslView);

struct timeval timer;

fstream fstr;
fstr.open("logs/body3dcedemstab.txt", std::ios::out | std::ios::trunc);
fstr.precision(20);

getchar();

if (viewer)
viewer->Remove(dslView);

Expand Down
9 changes: 8 additions & 1 deletion bin/body3dforceest.cc
Original file line number Diff line number Diff line change
Expand Up @@ -7,10 +7,17 @@
#include "utils.h"
#include "so3.h"
#include "params.h"
#include "gndoep.h"
#include "sensor.h"
#include "insmanifold.h"

//#define USE_STOCHASTIC_DYNAMICS

#ifdef USE_STOCHASTIC_DYNAMICS
#include "gndoep.h"
#else
#include "deterministicgndoep.h"
#endif

using namespace std;
using namespace Eigen;
using namespace gcop;
Expand Down
1 change: 1 addition & 0 deletions bin/body3dstab.cc
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ using namespace gcop;
//typedef ConstraintCost<Body3dState, 4, 2, Dynamic, 1> DiskCost;



void solver_process(Viewer* viewer)
{
viewer->SetCamera(-7.25, 61, -2.8, -2.55, -10);
Expand Down
138 changes: 138 additions & 0 deletions bin/imumagekftest.cc
Original file line number Diff line number Diff line change
@@ -0,0 +1,138 @@
#include <iostream>
#include "utils.h"
#include "so3.h"
#include "kalmanpredictor.h"
#include "kalmancorrector.h"
#include "unscentedpredictor.h"
#include "unscentedcorrector.h"
#include "imu.h"
#include "imusensor.h"

using namespace std;
using namespace gcop;
using namespace Eigen;

typedef KalmanPredictor<ImuState, 9, 3, Dynamic> ImuKalmanPredictor;
typedef KalmanCorrector<ImuState, 9, 3, Dynamic, Matrix<double, 6, 1>, 6> ImuKalmanCorrector;

typedef UnscentedPredictor<ImuState, 9, 3, Dynamic> ImuUnscentedPredictor;
typedef UnscentedCorrector<ImuState, 9, 3, Dynamic, Matrix<double, 6, 1>, 6> ImuUnscentedCorrector;


int main(int argc, char** argv)
{
Imu imu;

bool mag = true;
// int nz = mag ? 6 : 3;

ImuSensor<6> sensor;
// sensor.sa = 0; // no acceleration bias drift
// sensor.sv = 0; // no acceleration bias drift

ImuKalmanPredictor kp(imu);
ImuKalmanCorrector kc(imu.X, sensor);

// ImuUnscentedPredictor kp(imu);
// ImuUnscentedCorrector kc(imu, sensor);


// ImuKalmanCorrector ukf(imu, sensor);

int N = 1500;

vector<ImuState> xts(N); // true states
xts[0].bg << 0.1, 0, 0; // bias

vector<ImuState> xs(N); // estimated trajectory
xs[0].P.topLeftCorner<3,3>().diagonal().setConstant(.1);
xs[0].P.block<3,3>(3,3).diagonal().setConstant(1e-2);
xs[0].P.bottomRightCorner<3,3>().diagonal().setConstant(1e-10);

// true rotation

double dt = .01;

struct timeval timer;

Vector3d qt;
Vector3d q;

for (int i = 0; i < N-1; ++i) {
double t = i*dt;

// true input (drifted angular velocity)
Vector3d wt(.2, .1, 0); // this is the undrifted
wt = wt + xts[i].bg; // drifted angular velocity

// generate true
imu.Step(xts[i+1], t, xts[i], wt, dt);

// drifted noisy gyro reading
Vector3d w = wt + imu.sv*Vector3d::Ones()*randn();

ImuState x;

timer_start(timer);
kp.Predict(x, t, xs[i], w, dt);
long us = timer_us(timer);
cout << "Predict took " << us << " us." << endl;

SO3::Instance().g2q(qt, xts[i+1].R);
SO3::Instance().g2q(q, x.R);
cout << "rpyt: " << qt.transpose() << endl;
cout << "rpye: " << q.transpose() << endl;

// true measurements of acc and mag
Vector3d at = xts[i+1].R.transpose()*sensor.a0;
Vector3d mt = xts[i+1].R.transpose()*sensor.m0;

// noisy measurements

Vector6d z;
z.head<3>() = at + xts[i+1].ba + sensor.sra*Vector3d::Ones()*randn();
z.head<3>().normalize();
z.tail<3>() = mt + sensor.srm*Vector3d::Ones()*randn();
z.tail<3>().normalize();

/*
Vector3d z;
z = at + xts[i+1].ba + sensor.sra*Vector3d::Ones()*randn();
z.normalize();
*/

cout << "z=" << z.transpose() << endl;

cout << "P-=" << x.P << endl;

// if (i==N-2) {
timer_start(timer);
kc.Correct(xs[i+1], t, x, w, z);
//xs[i+1] = x;
us = timer_us(timer);
cout << "Update took " << us << " us." << endl;
// }

cout << "True attitude: " << xts[i+1].R << endl;
cout << "Estim attitude: " << xs[i+1].R << endl;

cout << "True biases: " << xts[i+1].bg.transpose() << endl;
cout << "Estim biases: " << xs[i+1].bg.transpose() << endl;

cout << "P=: " << xs[i+1].P << endl;

if (std::isnan(xs[i+1].R(0,0)))
break;

SO3::Instance().g2q(qt, xts[i+1].R);
SO3::Instance().g2q(q, xs[i+1].R);
cout << "rpyt: " << qt.transpose() << endl;
cout << "rpye: " << q.transpose() << endl;

Vector3d dR;
SO3::Instance().log(dR, xs[i+1].R.transpose()*xts[i+1].R);
cout << "rotation error: " << dR.norm() << endl;
// cout << "L_2 norm (state) =" << (xt - ukf.x).norm() << endl;
// cout << "L_2 norm (quat) =" << (xt.head(4) - ukf.x.head(4)).norm() << endl;
}
}
4 changes: 2 additions & 2 deletions lib/algos/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ set(headers
doep.h
gndocp.h
gndoep.h
gndoepv2.h
deterministicgndoep.h
pddp.h
ba.h
body2dslam.h
Expand Down Expand Up @@ -41,4 +41,4 @@ set(headers
${headers}
body3ddemcontroller.h
)
ENDIF (DSL_FOUND AND PQP_FOUND)
ENDIF (DSL_FOUND AND PQP_FOUND)
27 changes: 18 additions & 9 deletions lib/algos/gndoepv2.h → lib/algos/deterministicgndoep.h
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,9 @@

#include <unsupported/Eigen/NonLinearOptimization>

/** GnDoep with deterministic dynamics
*/

namespace gcop {


Expand All @@ -32,7 +35,7 @@ namespace gcop {
int _nz = Dynamic,
typename T1 = T,
int _nx1 = _nx>
class GnDoep1 : public Doep<T, _nx, _nu, _np, Tz, _nz, T1, _nx1> {
class GnDoep : public Doep<T, _nx, _nu, _np, Tz, _nz, T1, _nx1> {

typedef Matrix<double, _ng, 1> Vectorgd;
typedef Matrix<double, _ng, _nx> Matrixgxd;
Expand Down Expand Up @@ -90,12 +93,12 @@ namespace gcop {
*/


GnDoep1(System<T, _nx, _nu, _np> &sys, Sensor<T1, _nx1, _nu, _np, Tz, _nz> &sensor,
GnDoep(System<T, _nx, _nu, _np> &sys, Sensor<T1, _nx1, _nu, _np, Tz, _nz> &sensor,
LqSensorCost<T, _nx, _nu, _np, _ng, Tz, _nz> &cost,
vector<double> &ts, vector<T> &xs, vector<Vectorcd> &us, Vectormd &p,
vector<double> &ts1, Func_type _project=NULL, bool update = true);

virtual ~GnDoep1();
virtual ~GnDoep();

/**
* Perform one DOCP iteration. Internally calls:
Expand Down Expand Up @@ -155,7 +158,7 @@ struct Functor
struct GnCost : Functor<double> {
GnCost<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>(int inputs, int values) : Functor<double>(inputs, values) {};

GnDoep1<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1> *doep;
GnDoep<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1> *doep;

typedef Matrix<double, _nx, 1> Vectornd;
typedef Matrix<double, _np, 1> Vectormd;
Expand Down Expand Up @@ -212,9 +215,11 @@ struct Functor
{
std::cout<<"Input: "<<s.transpose()<<endl;
std::cout<<"Fvec: "<<fvec.transpose()<<endl;
std::cout<<"Cost: "<<0.5*(fvec.transpose()*fvec)<<endl;
//std::cout<<"Resp: "<<fvec.tail(np).transpose()<<"\t"<<doep->p<<endl;
}
std::cout<<"Cost: "<<0.5*(fvec.transpose()*fvec)<<endl;
doep->J = 0.5*(fvec.squaredNorm());//#TODO Directly from LevinBergMarquadt if possible
//getchar();
return 0;
}
};
Expand All @@ -223,7 +228,7 @@ struct Functor
using namespace Eigen;

template <typename T, int _nx, int _nu, int _np, int _ng, typename Tz, int _nz, typename T1, int _nx1>
GnDoep1<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::GnDoep1(System<T, _nx, _nu, _np> &sys, Sensor<T1, _nx1, _nu, _np, Tz, _nz> &sensor,
GnDoep<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::GnDoep(System<T, _nx, _nu, _np> &sys, Sensor<T1, _nx1, _nu, _np, Tz, _nz> &sensor,
LqSensorCost<T, _nx, _nu, _np, _ng, Tz, _nz> &cost,
vector<double> &ts,
vector<T> &xs,
Expand All @@ -245,15 +250,15 @@ struct Functor
}

template <typename T, int _nx, int _nu, int _np, int _ng, typename Tz, int _nz, typename T1, int _nx1>
GnDoep1<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::~GnDoep1()
GnDoep<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::~GnDoep()
{
delete lm;
delete numDiff;
delete functor;
}

template <typename T, int _nx, int _nu, int _np, int _ng, typename Tz, int _nz, typename T1, int _nx1>
void GnDoep1<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::Iterate() {
void GnDoep<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>::Iterate() {

if (!lm) {
functor = new GnCost<T, _nx, _nu, _np, _ng, Tz, _nz, T1, _nx1>(inputs, values);
Expand All @@ -267,6 +272,10 @@ struct Functor

const int &np = this->sys.P.n;
s = this->p;//Set the system parameters to initial guess

lm->parameters.maxfev = 1e6;//Maximum nof evaluations is very high
info = lm->minimizeInit(s);
cout <<"info="<<info <<endl;
}

/*
Expand All @@ -277,7 +286,7 @@ struct Functor
*/

// lm.parameters.maxfev=10000;
info = lm->minimize(s);
info = lm->minimizeOneStep(s);

cout <<"info="<<info <<endl;
// check return values
Expand Down
Loading

0 comments on commit ac63468

Please sign in to comment.