diff --git a/matlab/addSimulatedConstraint.m b/matlab/addSimulatedConstraint.m
deleted file mode 100644
index a53872b6d..000000000
--- a/matlab/addSimulatedConstraint.m
+++ /dev/null
@@ -1,17 +0,0 @@
-function addSimulatedConstraint(points,angles,sd,id1,id2,graph)
-% addSimulatedConstraint: create a simulated measurement with noise
-% standard deviations sd and add it to graph
-
-key1 = sprintf('x%d', id1);
-key2 = sprintf('x%d', id2);
-
-% ground truth
-delta_x = points(id1,:) - points(id2,:);
-delta_angle = angles(id1) - angles(id2);
-noisy = Pose2(delta_x(1) + sd(1)*randn, delta_x(2) + sd(2)*randn, delta_angle + sd(3)*randn);
-
-% create factor
-factor=Pose2Factor(key1,key2,noisy,diag(sd.*sd));
-
-% add it to the graph
-graph.push_back(factor);
diff --git a/matlab/beijing.mat b/matlab/beijing.mat
deleted file mode 100644
index 6c2ae4b76..000000000
Binary files a/matlab/beijing.mat and /dev/null differ
diff --git a/matlab/beijing_angles.mat b/matlab/beijing_angles.mat
deleted file mode 100644
index 8128a433d..000000000
Binary files a/matlab/beijing_angles.mat and /dev/null differ
diff --git a/matlab/beijing_graph.mat b/matlab/beijing_graph.mat
deleted file mode 100644
index c65162e01..000000000
Binary files a/matlab/beijing_graph.mat and /dev/null differ
diff --git a/matlab/bottom_up_ordering.m b/matlab/bottom_up_ordering.m
deleted file mode 100644
index 0eb58ea1b..000000000
--- a/matlab/bottom_up_ordering.m
+++ /dev/null
@@ -1,19 +0,0 @@
-% find a bottom to up ordering given the tree structure {pred} returned by matlab's graphminspantree
-function [ordering] = bottom_up_ordering(pred)
-
-%% compute the levels of the nodes
-parents = [0];
-node_levels = zeros(length(pred), 1);
-current_level = 1;
-while ~isempty(parents)
-    parents = find(ismember(pred, parents));
-    node_levels(parents) = current_level;
-    current_level = current_level + 1;
-end
-[~, node_order] = sort(node_levels, 'descend');  % the order of the nodes in leaves-to-root order
-
-ordering = Ordering();
-for i = 1:length(node_order)
-    ordering.push_back(sprintf('x%d', node_order(i)));
-end
-end
\ No newline at end of file
diff --git a/matlab/create_config.m b/matlab/example/create_config.m
similarity index 100%
rename from matlab/create_config.m
rename to matlab/example/create_config.m
diff --git a/matlab/create_gaussian_factor_graph.m b/matlab/example/create_gaussian_factor_graph.m
similarity index 100%
rename from matlab/create_gaussian_factor_graph.m
rename to matlab/example/create_gaussian_factor_graph.m
diff --git a/matlab/create_good_ordering.m b/matlab/example/create_good_ordering.m
similarity index 100%
rename from matlab/create_good_ordering.m
rename to matlab/example/create_good_ordering.m
diff --git a/matlab/create_landmarks.m b/matlab/example/create_landmarks.m
similarity index 100%
rename from matlab/create_landmarks.m
rename to matlab/example/create_landmarks.m
diff --git a/matlab/create_linear_factor_graph.m b/matlab/example/create_linear_factor_graph.m
similarity index 100%
rename from matlab/create_linear_factor_graph.m
rename to matlab/example/create_linear_factor_graph.m
diff --git a/matlab/create_ordering.m b/matlab/example/create_ordering.m
similarity index 100%
rename from matlab/create_ordering.m
rename to matlab/example/create_ordering.m
diff --git a/matlab/create_random_landmarks.m b/matlab/example/create_random_landmarks.m
similarity index 100%
rename from matlab/create_random_landmarks.m
rename to matlab/example/create_random_landmarks.m
diff --git a/matlab/create_visibility.m b/matlab/example/create_visibility.m
similarity index 100%
rename from matlab/create_visibility.m
rename to matlab/example/create_visibility.m
diff --git a/matlab/random_walk.m b/matlab/example/random_walk.m
similarity index 100%
rename from matlab/random_walk.m
rename to matlab/example/random_walk.m
diff --git a/matlab/simulate_measurements.m b/matlab/example/simulate_measurements.m
similarity index 100%
rename from matlab/simulate_measurements.m
rename to matlab/example/simulate_measurements.m
diff --git a/matlab/walk.m b/matlab/example/walk.m
similarity index 100%
rename from matlab/walk.m
rename to matlab/example/walk.m
diff --git a/matlab/frank01_CGD.m b/matlab/frank01_CGD.m
deleted file mode 100644
index 0bf9ab5fe..000000000
--- a/matlab/frank01_CGD.m
+++ /dev/null
@@ -1,39 +0,0 @@
-%-----------------------------------------------------------------------
-% frank01.m: try conjugate gradient on our example graph
-%-----------------------------------------------------------------------
-
-% get matrix form H and z
-fg = createGaussianFactorGraph();
-ord = Ordering;
-ord.push_back('x1');
-ord.push_back('x2');
-ord.push_back('l1');
-
-[H,z] = fg.matrix(ord);
-
-% form system of normal equations
-A=H'*H
-b=H'*z
-
-% k=0
-x = zeros(6,1)
-g = A*x-b
-d = -g
-
-for k=1:5
-    alpha = - (d'*g)/(d'*A*d)
-    x = x + alpha*d
-    g = A*x-b
-    beta = (d'*A*g)/(d'*A*d)
-    d = -g + beta*d
-end
-
-% Do gradient descent
-% fg2 = createGaussianFactorGraph();
-% zero = createZeroDelta();
-% actual = fg2.gradientDescent(zero);
-% CHECK(assert_equal(expected,actual,1e-2));
-
-% Do conjugate gradient descent
-% actual2 = fg2.conjugateGradientDescent(zero);
-% CHECK(assert_equal(expected,actual2,1e-2));
diff --git a/matlab/frank02_ToroPose2Graph.m b/matlab/frank02_ToroPose2Graph.m
deleted file mode 100644
index 5704fad2c..000000000
--- a/matlab/frank02_ToroPose2Graph.m
+++ /dev/null
@@ -1 +0,0 @@
-% load the Toro 2D dataset and build pose graph
diff --git a/matlab/frank03_BejingPose2Graph.m b/matlab/frank03_BejingPose2Graph.m
deleted file mode 100644
index 4739ecf71..000000000
--- a/matlab/frank03_BejingPose2Graph.m
+++ /dev/null
@@ -1,44 +0,0 @@
-% frank03: create Beijing matrices in a cleaner way
-
-load beijing.mat;
-load beijing_angles.mat;
-load beijing_graph.mat;
-n=size(points,1);
-
-% create config or load it from file
-if 0
-    load beijing_config.mat;
-else
-    config=Pose2Config();
-    for j=1:n
-        if mod(j,1000) == 0, fprintf(1, 'adding node %d to config\n', j); end   
-        pose=Pose2(points(j,1),points(j,2),angles(j));
-        key = sprintf('x%d', j);
-        config.insert(key,pose);
-    end
-    save('beijing_config.mat','config');
-end
-
-sd = [0.25;0.25;0.01];
-
-% Build factor graph for entire graph
-graph = Pose2Graph;
-
-% First add tree constraints
-[I J] = find(tree);
-for k=length(edge_order):-1:1
-    edge = edge_order(k);
-    if mod(k,1000) == 0, fprintf(1, 'simulating constraint %d\n', k); end   
-    addSimulatedConstraint(points,angles,sd,I(edge),J(edge),graph);
-end
-
-% Then add remaining constraints C
-C=G-tree;
-[I J] = find(C);
-for k=1:length(I)
-    if mod(k,100) == 0, fprintf(1, 'simulating constraint %d\n', k); end   
-    addSimulatedConstraint(points,angles,sd,I(k),J(k),graph);
-end
-
-% generate ordering
-ordering = bottom_up_ordering(pred);
\ No newline at end of file
diff --git a/matlab/frank04_showMatrices.m b/matlab/frank04_showMatrices.m
deleted file mode 100644
index 542838cad..000000000
--- a/matlab/frank04_showMatrices.m
+++ /dev/null
@@ -1,34 +0,0 @@
-% frank04: show matrices associated with Beijing simulation
-
-load beijing.mat;
-load beijing_angles.mat;
-load beijing_graph.mat;
-
-% put spanning tree 'tree' in correct order
-T = tree(node_order,node_order);
-
-% get loop closing constraints
-C=G(node_order,node_order)-T;
-
-close all
-
-% plot on map
-figure
-gplot(C,points(node_order,:),'r')
-hold on
-gplot(T,points(node_order,:),'k')
-axis equal
-
-% show spanning tree, original graph, and loop closures
-figure
-spy(T,'k');
-hold on
-spy(C,'r');
-
-figure
-spy(T);
-
-figure
-spy(C);
-
-
diff --git a/matlab/frank05_rightPreconditioning.m b/matlab/frank05_rightPreconditioning.m
deleted file mode 100644
index 9b547afa8..000000000
--- a/matlab/frank05_rightPreconditioning.m
+++ /dev/null
@@ -1,41 +0,0 @@
-% frank05: show right pre-conditioning
-% first run frank03 or other script to generate graph, config, and ordering
-
-% linearize the non-linear factor graph
-tic
-LFG = graph.linearize_(config);
-toc
-tic
-ijs = LFG.sparse(ordering);
-A = sparse(ijs(1,:),ijs(2,:),ijs(3,:)); 
-toc
-figure(1)
-spy(A);
-
-% isolate the spanning tree part
-A1=A(1:3*nnz(tree),:);
-% add prior
-figure(2)
-spy(A1)
-
-% calculate R1
-tic
-R1 = qr(A1,0);
-toc
-figure(3)
-spy(R1)
-
-% calculate R1
-tic
-R1 = chol(A1'*A1);
-toc
-figure(3)
-spy(R1)
-
-% calculate the entire R factor (expensive)
-tic
-R = qr(A,0);
-toc
-figure(4)
-spy(R)
-
diff --git a/matlab/kai01_BejingPose2Graph.m b/matlab/kai01_BejingPose2Graph.m
deleted file mode 100644
index ab110b414..000000000
--- a/matlab/kai01_BejingPose2Graph.m
+++ /dev/null
@@ -1,79 +0,0 @@
-load beijing.mat;
-load beijing_angles.mat;
-load beijing_graph.mat;
-load beijing_config.mat;
-
-cov = [ 0.25, 0, 0; 0, 0.25, 0; 0, 0, 0.01];
-
-factors = Pose2Graph;
-factors2 = Pose2Graph;
-ord2 = Ordering();
-
-[rows cols] = find(tree);
-for i=length(edge_order):-1:1
-    if mod(i,500) == 0
-        fprintf(1, 'processing edge %d\n', i);
-    end
-    
-    id1 = rows(edge_order(i));
-    id2 = cols(edge_order(i));
-    key1 = sprintf('x%d', id1);
-    key2 = sprintf('x%d', id2);
-    
-    delta_x = points(id1,:) - points(id2,:);
-    delta_angle = angles(id1) - angles(id2);
-    measured = Pose2(delta_x(1), delta_x(2), delta_angle);
-    
-    if pred(id1) == id2 || pred(id2) == id1 %% in the spanning tree
-        factor=Pose2Factor(key1,key2,measured, cov);
-        factors.push_back(factor);
-    else %% not in the spanning tree
-        factors2.push_back(Pose2Factor(key1,key2,measured, cov));
-        ord2.push_back(key1);
-        ord2.push_back(key2);
-    end
-end
-ord2.unique();
-
-if 1
-    config=Pose2Config();
-    n=size(points,1);
-    for j=1:n
-        pose=Pose2(points(j,1),points(j,2),angles(j));
-        key = sprintf('x%d', j);
-        config.insert(key,pose);
-        if mod(j,1000) == 0
-            key
-        end
-    end
-    save('beijing_config.mat','config');
-end
-
-% Spanning tree with bottom-up ordering
-ord = bottom_up_ordering(pred);
-LFG=factors.linearize_(config);
-ijs=LFG.sparse(ord);
-A=sparse(ijs(1,:),ijs(2,:),ijs(3,:)); 
-figure(1)
-spy(A);
-%save('beijing_factors.mat', 'factors');
-
-% LFG2=factors2.linearize_(config);
-% ijs2=LFG2.sparse(ord2);
-% A2=sparse(ijs2(1,:),ijs2(2,:),ijs2(3,:)); 
-% figure(2)
-% spy(A2);
-
-% show R factor
-R = qr(A,0);
-figure(3)
-spy(R)
-
-% show re-ordered R factor
-% P = colamd(A);
-% figure(4)
-% spy(A(:,P))
-% R = qr(A(:,P),0);
-% figure(5)
-% spy(R)
-
diff --git a/matlab/kai02_bottomUp.m b/matlab/kai02_bottomUp.m
deleted file mode 100644
index e3249a105..000000000
--- a/matlab/kai02_bottomUp.m
+++ /dev/null
@@ -1,4 +0,0 @@
-load beijing_graph.mat;
-
-ordering = bottom_up_ordering(pred);
-ordering.print('ordering')
diff --git a/matlab/CHECK.m b/matlab/tests/CHECK.m
similarity index 100%
rename from matlab/CHECK.m
rename to matlab/tests/CHECK.m
diff --git a/matlab/DOUBLES_EQUAL.m b/matlab/tests/DOUBLES_EQUAL.m
similarity index 100%
rename from matlab/DOUBLES_EQUAL.m
rename to matlab/tests/DOUBLES_EQUAL.m
diff --git a/matlab/createGaussianFactorGraph.m b/matlab/tests/createGaussianFactorGraph.m
similarity index 100%
rename from matlab/createGaussianFactorGraph.m
rename to matlab/tests/createGaussianFactorGraph.m
diff --git a/matlab/createNoisyConfig.m b/matlab/tests/createNoisyConfig.m
similarity index 100%
rename from matlab/createNoisyConfig.m
rename to matlab/tests/createNoisyConfig.m
diff --git a/matlab/createZeroDelta.m b/matlab/tests/createZeroDelta.m
similarity index 100%
rename from matlab/createZeroDelta.m
rename to matlab/tests/createZeroDelta.m
diff --git a/matlab/run_tests.m b/matlab/tests/run_tests.m
similarity index 100%
rename from matlab/run_tests.m
rename to matlab/tests/run_tests.m
diff --git a/matlab/testGaussianFactor.m b/matlab/tests/testGaussianFactor.m
similarity index 100%
rename from matlab/testGaussianFactor.m
rename to matlab/tests/testGaussianFactor.m
diff --git a/matlab/testGaussianFactorGraph.m b/matlab/tests/testGaussianFactorGraph.m
similarity index 100%
rename from matlab/testGaussianFactorGraph.m
rename to matlab/tests/testGaussianFactorGraph.m
diff --git a/matlab/testPose2Factor.m b/matlab/tests/testPose2Factor.m
similarity index 100%
rename from matlab/testPose2Factor.m
rename to matlab/tests/testPose2Factor.m
diff --git a/matlab/test_time_average.m b/matlab/timing/test_time_average.m
similarity index 100%
rename from matlab/test_time_average.m
rename to matlab/timing/test_time_average.m
diff --git a/matlab/test_time_new.m b/matlab/timing/test_time_new.m
similarity index 100%
rename from matlab/test_time_new.m
rename to matlab/timing/test_time_new.m
diff --git a/matlab/vila01_BejingPose2Graph.m b/matlab/vila01_BejingPose2Graph.m
deleted file mode 100644
index 4736659dc..000000000
--- a/matlab/vila01_BejingPose2Graph.m
+++ /dev/null
@@ -1,84 +0,0 @@
-load beijing.mat;
-load beijing_angles.mat;
-load beijing_graph.mat;
-%load beijing_config.mat;
-
-cov = [ 0.25, 0, 0; 0, 0.25, 0; 0, 0, 0.01];
-
-
-factors = Pose2Graph;
-factors2 = Pose2Graph;
-ord = Ordering();
-ord2 = Ordering();
-
-for i=1:length(ways)
-    if mod(i,500) == 0
-        fprintf(1, 'processing way %d\n', i);
-    end
-    for j=1:length(ways{i})-1
-        id1 = ways{i}(j);
-        id2 = ways{i}(j+1);
-        key1 = sprintf('x%d', id1);
-        key2 = sprintf('x%d', id2);
-
-        delta_x = points(id1,:) - points(id2,:);
-        delta_angle = angles(id1) - angles(id2);
-        measured = Pose2(delta_x(1), delta_x(2), delta_angle);
-
-        if pred(id1) == id2 || pred(id2) == id1 %% in the spanning tree
-            factor=Pose2Factor(key1,key2,measured, cov);
-            factors.push_back(factor);
-            ord.push_back(key1);
-            ord.push_back(key2);
-        else %% not in the spanning tree
-            factors2.push_back(Pose2Factor(key1,key2,measured, cov));
-            ord2.push_back(key1);
-            ord2.push_back(key2);
-        end
-    end
-end
-ord.unique();
-ord2.unique();
-% ord.reverse();
-% ord2.reverse();
-
-if 0
-    config=Pose2Config();
-    n=size(points,1);
-    for j=1:n
-        pose=Pose2(points(j,1),points(j,2),angles(j));
-        key = sprintf('x%d', j);
-        config.insert(key,pose);
-        if mod(j,1000) == 0
-            key
-        end
-    end
-    save('beijing_config.mat','config');
-end
-
-amd_ord=factors.getOrdering_(); % does not work
-LFG=factors.linearize_(config);
-ijs=LFG.sparse(ord);
-A=sparse(ijs(1,:),ijs(2,:),ijs(3,:)); 
-figure(1)
-spy(A);
-%save('beijing_factors.mat', 'factors');
-
-LFG2=factors2.linearize_(config);
-ijs2=LFG2.sparse(ord2);
-A2=sparse(ijs2(1,:),ijs2(2,:),ijs2(3,:)); 
-figure(2)
-spy(A2);
-
-% show R factor
-R = qr(A,0);
-figure(3)
-spy(R)
-
-% show re-ordered R factor
-P = colamd(A);
-figure(4)
-spy(A(:,P))
-R = qr(A(:,P),0);
-figure(5)
-spy(R)