clear 

file = "C:\Bachelorarbeit\Messung 16.03\Messung1_500Hz.csv";

test = readtable(file, 'VariableNamingRule', 'preserve');


datenneu = table2array(test(:,2:7)); % data without timestamp

l1 = length(datenneu);

timestamp = 0:0.002:(l1*0.002); %create timestamp, 0.002 wegen 500 Hz

data = [timestamp(2:end)', datenneu]; % data with timestamp

Time = data(:,1);

ax1 = data(:,2)*1000;% *1000 auf mm/s^2
ax2 = data(:,5)*1000;

ay1 = data(:,3)*1000;
ay2 = data(:,6)*1000;

az1 = data(:,4)*1000;
az2 = data(:,7)*1000;
% + Zahl für Offset



fs=500; %Sampling rate



figure(1) 

fcax = 2; %Cut-Off Frequency Beschleunigung in X Richtung
[cax,dax]=butter(2,fcax/(fs/2),'high');

%Beschleunigung X Richtung
ax1filter = filter(cax,dax,ax1);
ax2filter = filter(cax,dax,ax2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
ax1filter = ax1filter - mean(ax1filter);
ax2filter = ax2filter - mean(ax2filter);


subplot(3,1,1) %3 Zeilen, 1 Spalte, 1. Plot
plot(Time, ax1, Time, ax1filter,Time,ax2,Time,ax2filter)
ylim([-100 100])
xlim([0 16.5])

legend({'ax1', 'ax1filter','ax2','ax2filter'})

xlabel('time [s]')
ylabel('Beschleunigung [mm/s^2]')
title('X Beschleunigung-Zeit Verlauf')

grid minor


fcay = 6; %Cut-Off Frequency Beschleunigung in Y Richtung
[cay,day]=butter(2,fcay/(fs/2),'high');

%Beschleunigung Y Richtung
ay1filter = filtfilt(cay,day,ay1);
ay2filter = filtfilt(cay,day,ay2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
ay1filter = ay1filter - mean(ay1filter);
ay2filter = ay2filter - mean(ay2filter);

subplot(3,1,2)
plot(Time, ay1, Time, ay1filter,Time,ay2,Time,ay2filter)
ylim([-200 200])
xlim([0 16.5])

legend({'ay1', 'ay1filter','ay2','ay2filter'})

xlabel('time [s]')
ylabel('Beschleunigung [mm/s^2]')
title('Y Beschleunigung-Zeit Verlauf')

grid minor


fcaz1 = 6; %Cut-Off Frequency Beschleunigung in X Richtung
[caz1,daz1]=butter(2,fcaz1/(fs/2),'high');

fcaz2 = 3; %Cut-Off Frequency Beschleunigung in X Richtung
[caz2,daz2]=butter(2,fcaz2/(fs/2),'high');

%Beschleunigung Z Richtung

az1filter = filtfilt(caz1,daz1,az1);
az2filter = filtfilt(caz2,daz2,az2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
az1filter = az1filter - mean(az1filter);
az2filter = az2filter - mean(az2filter);

subplot(3,1,3)
plot(Time, az1, Time, az1filter,Time,az2,Time,az2filter)
ylim([-500 500])
xlim([0 16.5])

legend({'az1', 'az1filter','az2','az2filter'})

xlabel('time [s]')
ylabel('Beschleunigung [mm/s^2]')
title('Z Beschleunigung-Zeit Verlauf')

grid minor


%Integrieren mit der Trapezregel
vx1 = cumtrapz(Time,ax1filter)*1000; % *1000 auf mym/s
vx2 = cumtrapz(Time,ax2filter)*1000; 

vy1 = cumtrapz(Time,ay1filter)*1000;
vy2 = cumtrapz(Time,ay2filter)*1000;

vz1 = cumtrapz(Time,az1filter)*1000;
vz2 = cumtrapz(Time,az2filter)*1000;



figure(2)

%Geschwindigkeit X Richtung

fcvx1 = 10; %Cut-Off Frequency Beschleunigung in X Richtung
[cvx1,dvx1]=butter(1,fcvx1/(fs/2),'high');

fcvx2 = 10; %Cut-Off Frequency Beschleunigung in X Richtung
[cvx2,dvx2]=butter(1,fcvx2/(fs/2),'high');

vx1filter = filtfilt(cvx1,dvx1,vx1);
vx2filter = filtfilt(cvx2,dvx2,vx2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
vx1filter = vx1filter - mean(vx1filter);
vx2filter = vx2filter - mean(vx2filter);

subplot(3,1,1)
plot(Time, vx1, Time, vx1filter ,Time,vx2, Time, vx2filter)
ylim([-100 100])
xlim([0 16.5])

legend({'vx1','vx1filter','vx2','vx2filter'})

xlabel('time [s]')
ylabel('Geschwindigkeit [μm/s]')
title('X Geschwindigkeit-Zeit Verlauf')

grid minor



%Geschwindigkeit Y Richtung
fcvy1 = 5; %Cut-Off Frequency Beschleunigung in X Richtung
[cvy1,dvy1]=butter(1,fcvy1/(fs/2),'high');

fcvy2 = 5; %Cut-Off Frequency Beschleunigung in X Richtung
[cvy2,dvy2]=butter(1,fcvy2/(fs/2),'high');

vy1filter = filtfilt(cvy1,dvy1,vy1);
vy2filter = filtfilt(cvy2,dvy2,vy2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
vy1filter = vy1filter - mean(vy1filter);
vy2filter = vy2filter - mean(vy2filter);

subplot(3,1,2)
plot(Time, vy1, Time, vy1filter ,Time,vy2, Time, vy2filter)
ylim([-100 100])
xlim([0 16.5])

legend({'vy1','vy1filter','vy2','vy2filter'})

xlabel('time [s]')
ylabel('Geschwindigkeit [μm/s]')
title('Y Geschwindigkeit-Zeit Verlauf')

grid minor

%Geschwindigkeit Z Richtung
fcvz1 = 5; %Cut-Off Frequency Beschleunigung in X Richtung
[cvz1,dvz1]=butter(1,fcvz1/(fs/2),'high');

fcvz2 = 5; %Cut-Off Frequency Beschleunigung in X Richtung
[cvz2,dvz2]=butter(1,fcvz2/(fs/2),'high');

vz1filter = filtfilt(cvz1,dvz1,vz1);
vz2filter = filtfilt(cvz2,dvz2,vz2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
vz1filter = vz1filter - mean(vz1filter);
vz2filter = vz2filter - mean(vz2filter);

subplot(3,1,3)
plot(Time, vz1, Time, vz1filter ,Time,vz2, Time, vz2filter)
ylim([-200 200])
xlim([0 16.5])

legend({'vz1','vz1filter','vz2','vz2filter'})

xlabel('time [s]')
ylabel('Geschwindigkeit [μm/s]')
title('Z Geschwindigkeit-Zeit Verlauf')

grid minor


%Integrieren mit der Trapezregel
sx1 = cumtrapz(Time,vx1filter)*1000; % *1000 auf mym/s
sx2 = cumtrapz(Time,vx2filter)*1000; 

sy1 = cumtrapz(Time,vy1filter)*1000;
sy2 = cumtrapz(Time,vy2filter)*1000;

sz1 = cumtrapz(Time,vz1filter)*1000;
sz2 = cumtrapz(Time,vz2filter)*1000;


figure(3)
fc = 150; %Cut-Off Frequency Beschleunigung in X Richtung
[b,a]=butter(1,fc/(fs/2),'high');

%Weg X Richtung

sx1filter = filtfilt(b,a,sx1);
sx2filter = filtfilt(b,a,sx2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
sx1filter = sx1filter - mean(sx1filter);
sx2filter = sx2filter - mean(sx2filter);

subplot(3,1,1)
plot(Time, sx1, Time, sx1filter ,Time,sx2, Time, sx2filter)
ylim([-100 100])
xlim([0 16.5])

legend({'sx1','sx1filter','sx2','sx2filter'})

xlabel('time [s]')
ylabel('Weg [nm]')
title('X Weg-Zeit Verlauf')

grid minor


%Weg Y Richtung

sy1filter = filtfilt(b,a,sy1);
sy2filter = filtfilt(b,a,sy2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
sy1filter = sy1filter - mean(sy1filter);
sy2filter = sy2filter - mean(sy2filter);

subplot(3,1,2)
plot(Time, sy1, Time, sy1filter ,Time,sy2, Time, sy2filter)
ylim([-100 100])
xlim([0 16.5])

legend({'sy1','sy1filter','sy2','sy2filter'})

xlabel('time [s]')
ylabel('Weg [nm]')
title('Y Weg-Zeit Verlauf')

grid minor


%Weg Z Richtung

sz1filter = filtfilt(b,a,sz1);
sz2filter = filtfilt(b,a,sz2);

%Mittelwert abziehen, damit Funktion um 0 schwingt
sz1filter = sz1filter - mean(sz1filter);
sz2filter = sz2filter - mean(sz2filter);

subplot(3,1,3)
plot(Time, sz1, Time, sz1filter ,Time,sz2, Time, sz2filter)
ylim([-200 200])
xlim([0 16.5])

legend({'sz1','sz1filter','sz2','sz2filter'})

xlabel('time [s]')
ylabel('Weg [nm]')
title('Z Weg-Zeit Verlauf')

grid minor


figure(4) %Relativweg X und Y und Z

%Relativweg muss im Genauigkeitsbereich +/- 15 nm sein
Relx = sx1filter-sx2filter;
Rely = sy1filter-sy2filter;
Relz = sz1filter-sz2filter;

subplot(3,1,1)
plot(Time, Relx)
ylim([-25 25])
xlim([0 16.5])

xlabel('time [s]')
ylabel('Relativweg [nm]')
title('Relativweg X')

grid minor


subplot(3,1,2)
plot(Time, Rely)
ylim([-25 25])
xlim([0 16.5])

xlabel('time [s]')
ylabel('Relativweg [nm]')
title('Relativweg Y')

grid minor


subplot(3,1,3)
plot(Time, Relz)
ylim([-40 40])
xlim([0 16.5])


xlabel('time [s]')
ylabel('Relativweg [nm]')
title('Relativweg Z')

grid minor


%%FFT - Analyse wird mit den ungefilterten Daten durchgeführt

AX1=fft(ax1filter);
AX2=fft(ax2filter);

AY1=fft(ay1filter);
AY2=fft(ay2filter);

AZ1=fft(az1filter);
AZ2=fft(az2filter);


Fs=500; %sampling frequency [Hz]
N=8192; %number of samples

f= (0:N-1)*Fs/N; %frequency [Hz] samples/sec/samples


figure(5) %FFT Analyse Sensor 1

AX1_norm=2/N*AX1;

subplot(3,1,1)
plot(f,abs(AX1_norm))

ylim([0 1])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum X1-Beschleunigung')

grid minor

AY1_norm=2/N*AY1;

subplot(3,1,2)
plot(f,abs(AY1_norm))

ylim([0 1])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum Y1-Beschleunigung')

grid minor

AZ1_norm=2/N*AZ1;

subplot(3,1,3)
plot(f,abs(AZ1_norm))

ylim([0 1])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum Z1-Beschleunigung')

grid minor


figure(6) %FFT Analyse Sensor 2

AX2_norm=2/N*AX2;

subplot(3,1,1)
plot(f,abs(AX2_norm))

ylim([0 2])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum X2-Beschleunigung')

grid minor

AY2_norm=2/N*AY2;

subplot(3,1,2)
plot(f,abs(AY2_norm))

ylim([0 2])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum Y2-Beschleunigung')

grid minor

AZ2_norm=2/N*AZ2;

subplot(3,1,3)
plot(f,abs(AZ2_norm))

ylim([0 5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Beschleunigung [mm/s^2]')
title('Frequenzspektrum Z2-Beschleunigung')

grid minor

VX1=fft(vx1filter);
VX2=fft(vx2filter);

VY1=fft(vy1filter);
VY2=fft(vy2filter);

VZ1=fft(vz1filter);
VZ2=fft(vz2filter);


figure(7) %FFT Analyse Sensor 1

VX1_norm=2/N*VX1;

subplot(3,1,1)
plot(f,abs(VX1_norm))

ylim([0 5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum X1-Geschwindigkeit')

grid minor

VY1_norm=2/N*VY1;

subplot(3,1,2)
plot(f,abs(VY1_norm))

ylim([0 5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum Y1-Geschwindigkeit')

grid minor

VZ1_norm=2/N*VZ1;

subplot(3,1,3)
plot(f,abs(VZ1_norm))

ylim([0 4])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum Z1-Geschwindigkeit')

grid minor


figure(8) %FFT Analyse Sensor 2

VX2_norm=2/N*VX2;

subplot(3,1,1)
plot(f,abs(VX2_norm))

ylim([0 5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum X2-Geschwindigkeit')

grid minor

VY2_norm=2/N*VY2;

subplot(3,1,2)
plot(f,abs(VY2_norm))

ylim([0 5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum Y2-Geschwindigkeit')

grid minor

VZ2_norm=2/N*VZ2;

subplot(3,1,3)
plot(f,abs(VZ2_norm))

ylim([0 10])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Geschwindigkeit [μm/s]')
title('Frequenzspektrum Z2-Geschwindigkeit')

grid minor


SX1=fft(sx1filter);
SX2=fft(sx2filter);

SY1=fft(sy1filter);
SY2=fft(sy2filter);

SZ1=fft(sz1filter);
SZ2=fft(sz2filter);


figure(9) %FFT Analyse Sensor 1

SX1_norm=2/N*SX1;

subplot(3,1,1)
plot(f,abs(SX1_norm))

ylim([0 0.5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum X1-Weg')

grid minor

SY1_norm=2/N*SY1;

subplot(3,1,2)
plot(f,abs(SY1_norm))

ylim([0 0.5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum Y1-Weg')

grid minor

SZ1_norm=2/N*SZ1;

subplot(3,1,3)
plot(f,abs(SZ1_norm))

ylim([0 0.5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum Z1-Weg')

grid minor


figure(10) %FFT Analyse Sensor 2

SX2_norm=2/N*SX2;

subplot(3,1,1)
plot(f,abs(SX2_norm))

ylim([0 0.5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum X2-Weg')

grid minor

SY2_norm=2/N*SY2;

subplot(3,1,2)
plot(f,abs(SY2_norm))

ylim([0 0.5])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum Y2-Weg')

grid minor

SZ2_norm=2/N*SZ2;

subplot(3,1,3)
plot(f,abs(SZ2_norm))

ylim([0 1])
xlim([0 250])
xlabel('Frequenz [Hz]')
ylabel('Weg [nm]')
title('Frequenzspektrum Z2-Weg')

grid minor