# Welcome to Intelligent Additive Manufacturing Lab

We're aiming to create a knowledge hub for 3D printing of the future.

## Decomposition Method

### Singular Value Decomposition

A singular value decomposition (SVD) is a factorization of a matrix into several constitutive components all of which have a specific meaning in applications. Also, the SVD is a transformation that stretches/compresses and rotates a vector. SVD decomposition is given by Equation (1).

(1)

Where  and  are unitary matrices.  is diagonal, which components are non-negative and ordered from the largest to smallest.

In this page, we will use it to combine with the Principal Component Analysis.

### Principal Component Analysis

When we have a matrix , by implement SVD of the matrix , the data is projected into in terms of the new orthonormal basis . So, in different measurements, the measurements have zeros covariance thus, the principal components were picking out. As a result, the principle information was extracted, and the noise was eliminated.

The calculating procedure is below, it can show how to erase noise and get the new diagonal matrix.

First, we project matrix in terms of basis :

As a result, get Equation (2).

(2)

Since is a unitary matrix, it should have .

(3)

Then, follows that to get the covariance of .

(4)

(5)

As a result, get equation (6).

(6)

clear all; clc; close all;

%%
% Tracking for camera 1 test 1
[x11, y11] = trackPG(vidFrames1_1,[321,228],[20,20,20,20]);
% Tracking for camera 2 test 1
[x21, y21] = trackPG(vidFrames2_1,[274,274],[20,20,20,20]);
% Tracking for camera 3 test 1
[x31, y31] = trackPG(vidFrames3_1,[318,271],[15,15,15,15]);

% SVD test 1
A1(1,:) = x11 - mean(x11);
A1(2,:) = y11 - mean(y11);
A1(3,:) = x21(10:235) - mean(x21(10:235));
A1(4,:) = y21(10:235) - mean(y21(10:235));
A1(5,:) = x31(1:226) - mean(x31(1:226));
A1(6,:) = y31(1:226) - mean(y31(1:226));
[U1,s1,V1] = svd(A1.',0);

sigment1 = diag(s1);
en11 = sigment1(1)/sum(sigment1);
en21 = sum(sigment1(1:2))/sum(sigment1);
en31 = sum(sigment1(1:3))/sum(sigment1);

figure(1)
subplot(2,1,1)
plot(1:length(sigment1),sigment1,'ro')
title('Singular Values')

subplot(2,1,2)
plot(1:6,V1(:,1),'.-',1:6,V1(:,2),'.-',1:6,V1(:,3),'.-')
legend('Mode 1','Mode 2','Mode 3','Location','southeast')
xticks([1 2 3 4 5 6])
title('Linear POD Modes')

%%
% Tracking for camera 1 test 2
[x12, y12] = trackP(vidFrames1_2,[325,308],[17,17*0.9,15,15*0.9]);
% Tracking for camera 2 test 2
[x22, y22] = trackP(vidFrames2_2,[314,357],[34,34,34,35]);
% Tracking Camera 2 test 2
% Tracking for camera 3 test 2
[x32, y32] = trackPG(vidFrames3_2,[349,245],[44,45,44,44]);

% SVD test 2
A2 = zeros(6,314);
A2(1,:) = x12 - mean(x12);
A2(2,:) = y12 - mean(y12);
A2(3,:) = x22(15:328) - mean(x22(15:328));
A2(4,:) = y22(15:328) - mean(y22(15:328));
A2(5,:) = x32(1:314) - mean(x32(1:314));
A2(6,:) = y32(1:314) - mean(y32(1:314));
[U2,s2,V2] = svd(A2.',0);

sigment2 = diag(s2);
en12 = sigment2(1)/sum(sigment2);
en22 = sum(sigment2(1:2))/sum(sigment2);
en32 = sum(sigment2(1:3))/sum(sigment2);

figure(2)
subplot(2,1,1)
plot(1:length(sigment2),sigment2,'ro')
title('Singular Values')

subplot(2,1,2)
plot(1:6,V2(:,1),'.-',1:6,V2(:,2),'.-',1:6,V2(:,3),'.-')
legend('Mode 1','Mode 2','Mode 3','Location','southeast')
xticks([1 2 3 4 5 6])
title('Linear POD Modes')

% Tracking for camera 1 test 3, tracking pink

[x13, y13] = trackPP(vidFrames1_3,[330,290],[20,20,20,20]);
% Tracking for camera 2 test 3, tracking pink
[x23, y23] = trackPP(vidFrames2_3,[252,294],[20,20,20,20]);
% Tracking for camera 3 test 3
[x33, y33] = trackPG(vidFrames3_3,[353,229],[17,17,17,17]);

% SVD test 3
A3 = zeros(6,206);
A3(1,:) = x13(20:225) - mean(x13(20:225));
A3(2,:) = y13(20:225) - mean(y13(20:225));
A3(3,:) = x23(5:210) - mean(x23(5:210));
A3(4,:) = y23(5:210) - mean(y23(5:210));
A3(5,:) = x33(10:215) - mean(x33(10:215));
A3(6,:) = y33(10:215) - mean(y33(10:215));
[U3,s3,V3] = svd(A3.',0);
sigment3 = diag(s3);

en13 = sigment3(1)/sum(sigment3);
en23 = sum(sigment3(1:2))/sum(sigment3);
en33 = sum(sigment3(1:3))/sum(sigment3);

figure(3)
subplot(2,1,1)
plot(1:length(sigment3),sigment3,'ro')
title('Singular Values')

subplot(2,1,2)
plot(1:6,V3(:,1),'.-',1:6,V3(:,2),'.-',1:6,V3(:,3),'.-')
legend('Mode 1','Mode 2','Mode 3','Location','southeast')
xticks([1 2 3 4 5 6])
title('Linear POD Modes')

% Tracking for camera 1 test 4
[x14, y14] = trackPP(vidFrames1_4,[402,263],[44,45,44,44]);
% Tracking for camera 2 test 4
[x24, y24] = trackPP(vidFrames2_4,[244,243],[44,45,44,44]);
% Tracking for camera 3 test 4
[x34, y34] = trackP(vidFrames3_4,[363,244],[40,40,30,30]);

% SVD test 4
A4 = zeros(6,392);
A4(1,:) = x14 - mean(x14);
A4(2,:) = y14 - mean(y14);
A4(3,:) = x24(1:392) - mean(x24(1:392));
A4(4,:) = y24(1:392) - mean(y24(1:392));
A4(5,:) = x34(1:392) - mean(x34(1:392));
A4(6,:) = y34(1:392) - mean(y34(1:392));
[U4,s4,V4] = svd(A4.',0);

sigment4 = diag(s4);
en14 = sigment4(1)/sum(sigment4);
en24 = sum(sigment4(1:2))/sum(sigment4);
en34 = sum(sigment4(1:3))/sum(sigment4);

figure(4)
subplot(2,1,1)
plot(1:length(sigment4),sigment4,'ro')
title('Singular Values')

subplot(2,1,2)
plot(1:6,V4(:,1),'.-',1:6,V4(:,2),'.-',1:6,V4(:,3),'.-')
legend('Mode 1','Mode 2','Mode 3','Location','southeast')
title('Linear POD Modes')

% Plot all position vectors

figure(5)
plot(1:226,y11,'r')
hold on
plot(1:226,y21(10:235),'k')
plot(1:226,x31(1:226),'b')
hold off
title('Test 1: Ideal Case')
lgd = legend('Camera 1', 'Camera 2', 'Camera 3','Location','southeast');
lgd.FontSize =7;
ylabel('Displacement')
xlabel('Frame')

figure(6)
plot(1:314,y12,'r')
hold on
plot(1:314  ,y22(15:328),'k')
plot(1:314,x32(1:314),'b')
lgd = legend('Camera 1', 'Camera 2', 'Camera 3','Location','southeast');
lgd.FontSize =7;
hold off
title('Test 2: Noisy Case')
ylabel('Displacement')
xlabel('Frame')

figure(7)
plot(1:206,y13(20:225),'r')
hold on
plot(1:206  ,y23(5:210),'k')
plot(1:206,x33(10:215),'b')
hold off
title('Test 3: Horizontal Displacement')
ylabel('Displacement')
xlabel('Frame')
lgd = legend('Camera 1', 'Camera 2', 'Camera 3','Location','southeast');
lgd.FontSize =7;

figure(8)
plot(1:392,y14,'r')
hold on
plot(1:392  ,y24(1:392),'k')
plot(1:392,x34(1:392),'b')
hold off
title('Test 4: Horizontal Displacement and Rotation')
ylabel('Displacement')
xlabel('Frame')
lgd = legend('Camera 1', 'Camera 2', 'Camera 3','Location','southeast');
lgd.FontSize =7;

%------------------------------------------------------------------------------

function [xp,yp] = trackPG(vidFrames,guess,width)
dim = size(vidFrames);
L = dim(4); %length of video in frames
xp = zeros(L,1);
yp = zeros(L,1);
maxloc = [0,0];
maxval = 0;

init_guess = guess;
xp(1) = init_guess(1);
yp(1) = init_guess(2);

%center the first search window to the initial point
center = init_guess;

for f = 2:L
frame = double(rgb2gray(vidFrames(:,:,:,f)));
for x = center(1)-width(1):center(1)+width(2)
for y = center(2)-width(3):center(2)+width(4)
point = frame(y,x,:);
if point > maxval
maxloc = [x,y];
maxval = point;
end
end
end
xp(f) = maxloc(1);
yp(f) = maxloc(2);
maxloc = [0,0];
maxval = 0;

center(1) = xp(f);
center(2) = yp(f);
end
end
%-------------------------------------------------------------
function [xp,yp] = trackP(vidFrames,guess,width)
dim = size(vidFrames);
L = dim(4);     %length of video in frames
xp = zeros(L,1);
yp = zeros(L,1);
maxloc = [0,0];
maxval = 0;

init_guess = guess;
xp(1) = init_guess(1);
yp(1) = init_guess(2);

%center the first search window to the initial point
center = init_guess;

for f = 2:L
frame = double(vidFrames(:,:,:,f));
for x = center(1)-width(1):center(1)+width(2)
for y = center(2)-width(3):center(2)+width(4)
point = frame(y,x,:);
colorsum = point(1,1,1)+point(1,1,2)+point(1,1,3);
if colorsum > maxval
maxloc = [x,y];
maxval = colorsum;
end
end
end
xp(f) = maxloc(1);
yp(f) = maxloc(2);
maxloc = [0,0];
maxval = 0;

center(1) = xp(f);
center(2) = yp(f);
end
end
%-------------------------------------------------------------
function [xp,yp] = trackPP(vidFrames,guess,width)
dim = size(vidFrames);
L = dim(4); %length of video in frames
xp = zeros(L,1);
yp = zeros(L,1);
maxloc = [0,0];
maxred = 0;

init_guess = guess;
xp(1) = init_guess(1);
yp(1) = init_guess(2);

%center the first search window to the initial point
center = init_guess;

for f = 2:L
frame = double(vidFrames(:,:,:,f));
for x = center(1)-width(1):center(1)+width(2)
for y = center(2)-width(3):center(2)+width(4)
point = frame(y,x,:);
if ((point(1,1,1) > maxred)&(point(1,1,2)<150)&(point(1,1,3)<150))
maxloc = [x,y];
maxred = point;
end
end
end
xp(f) = maxloc(1);
yp(f) = maxloc(2);
maxloc = [0,0];
maxred = 0;

center(1) = xp(f);
center(2) = yp(f);
end
end