-
Notifications
You must be signed in to change notification settings - Fork 0
/
singleSample.m
115 lines (106 loc) · 3.89 KB
/
singleSample.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
%Please copy-paste the "Samples" folder in the same directory as this Octave file as such:
% |
% |---sappho_analyser.m
% |---Samples
% |--Sappho_00XYZ.txt
clear;
clc;
pkg load signal;
%=====BASIC DEFINITIONS=====
fileName = 'Samples/Sappho_00068.txt';
pixelsNum = 1500;
pixels = 1:1:pixelsNum;
distanceFromPDA = 3.5; %centimetres (this is the distance between the PDA and the object that is being measured - not the distance between the laser and the PDA)
sensorLength = 1.1; %centimetres
angleMax = rad2deg(atan(sensorLength/distanceFromPDA)); %find the maximum angle (in degrees) between a pixel and the laser pointer
step = 2*angleMax/(pixelsNum); %Δθ per pixel
angles = -(angleMax-step):step:angleMax; %aq array of 1500 angles;
%=====FILE HANDLER====
file = fopen(fileName,'r');
data = dlmread(file, ' ', 4,0); %Ignore the first lines of the file that contain text. Read ONLY numeric data.
fid = fclose(file);
%=====DATA INVERTING & FILTERING=====
data = abs(data-4096);
n = 10; % filter order
dataFiltered = filter(ones(n, 1)/n, 1, data); %this is a moving average filter
%=====RESHAPE MATRIX=====
framesNum = size(data,1) / pixelsNum; %Calculate the number of the frames, based on the array size
data = reshape(data,[pixelsNum, framesNum]); %reshape the Nx1 array to a matrix with dimensions "PIXELS x FRAMES"
dataFiltered = reshape(dataFiltered,[pixelsNum, framesNum]);
%=====OVERLAY EVERY FRAME ON A PLOT (UNFILTERED)=====%
figure(1)
hold
title([fileName, ': Repeatability Test: Unfiltered Frames'],'Interpreter', 'none');
xlabel('Pixels');
ylabel('Brightness (a.u.)');
set(gca, "linewidth", 2, "fontsize", 12)
for i=1:framesNum
temp = data(:,i);
plot(temp)
legend
i++;
endfor
%=====OVERLAY EVERY FRAME ON A PLOT (FILTERED)=====%
figure(2)
hold
title([fileName, ': Repeatability Test: Frames + Moving Average Filter'],'Interpreter', 'none');
xlabel('Pixels');
ylabel('Brightness (a.u.)');
set(gca, "linewidth", 2, "fontsize", 12)
for i=1:framesNum
temp = dataFiltered(:,i);
plot(temp)
legend
i++;
endfor
%set(gca, "linewidth", 4, "fontsize", 12)
%=====CALCULATE AVERAGE AND MEDIAN SIGNAL=====%
figure(3)
hold
title([fileName, ': Frame Average vs Median: Linear Plot'],'Interpreter', 'none');
xlabel('Pixels');
ylabel('Brightness (a.u.)');
set(gca, "linewidth", 2, "fontsize", 12)
for i=1:pixelsNum
avgData(i) = mean(dataFiltered(i,:)); %If the original signal is too noisy (e.g. because of the light source) replace "data" with "dataFiltered" here
mdnData(i) = median(dataFiltered(i,:));
endfor
plot(avgData)
plot(mdnData)
legend({"Average", "Median"});
%=====NORMALISE DATA=====
figure(4)
hold
title([fileName, ': Frame Average vs Median: Normalised Plot'],'Interpreter', 'none');
xlabel('Pixels');
ylabel('Brightness (normalised)');
set(gca, "linewidth", 2, "fontsize", 12)
normAvg = (avgData - min(avgData)) / ( max(avgData) - min(avgData)); %max value -> 1, min value -> 0
normMdn = (mdnData - min(mdnData)) / ( max(mdnData) - min(mdnData));
plot(normAvg)
plot(normMdn)
legend ({"Norm. Average", "Norm. Median"});
%=====MAX-MIN % ERROR CALCULATION=====
figure(5)
hold
title([fileName, ': Max-Min Relative Error: Unfiltered vs Filtered'],'Interpreter', 'none');
xlabel('Pixels');
ylabel('Brightness (a.u.)');
set(gca, "linewidth", 2, "fontsize", 12)
for i=1:pixelsNum %Ignore the first N pixels because the filtering is quirky there
relError(i) = (max(data(i,:)) - min(data(i,:)) ) / max(data(i,:));
i++; %calculate the relative min-max error per pixel
endfor
for i=1:pixelsNum-1 %Ignore the first N pixels because the filtering is quirky there
if i<=n
relErrorFiltered(i) = 0;
elseif i>n
relErrorFiltered(i) = (max(dataFiltered(i,:)) - min(dataFiltered(i,:)) ) / max(dataFiltered(i,:));
endif
i++; %calculate the relative min-max error per pixel
endfor
relError = relError * 100;
relErrorFiltered = relErrorFiltered * 100;
plot(relError)
plot(relErrorFiltered)
legend ({"Unfilt. Rel. Error", "Filt. Rel. Error"});