-
Notifications
You must be signed in to change notification settings - Fork 6
/
Copy pathDiscreteFourierTransform.m
104 lines (93 loc) · 1.82 KB
/
DiscreteFourierTransform.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
##clc;
##clear all;
##close all;
##
###... Discrete Fourier Transform
##f = 10;
##fs = 100;
##ts = 0:1/fs:1;
##xInput = sin(2*pi*f*ts);
##%xInput = ones(1, 7);
##nPoint = 160;
##xLength = length(xInput);
##
##if (nPoint < xLength)
## error("Can©t perform dft");
##endif
##
##xInput = [xInput zeros(1, nPoint-xLength)];
##
##for k = 0:nPoint-1
## xD(k+1) = 0;
## for n = 0:nPoint-1
## xD(k+1) = xD(k+1)+xInput(n+1)*exp(-i*2*pi*k*n/nPoint);
## endfor
##endfor
##
##n = 0:nPoint-1;
##subplot(3,1,1);
##plot(n,abs(xD));
##axis tight;
##title("Magnitude spectrum");
##
##subplot(3,1,2);
##plot(n,angle(xD));
##axis tight;
##title("Angle spectrum");
##
##for k = 0:nPoint-1
## xInv(k+1) = 0;
## for n = 0:nPoint-1
## xInv(k+1) = xInv(k+1)+xD(n+1)*exp(i*2*pi*n*k/nPoint);
## endfor
## xInv(k+1) = xInv(k+1)/nPoint;
##endfor
##
##n = 0:nPoint-1;
##disp(xInv);
##subplot(3,1,3);
##plot(n, xInv);
##axis tight;
Fs = 8000;
nx = 0:1/Fs:1;
input_signal = [2 3 4 2 1 4 5];
%Perform 8 point DFT
len = 8;
n = 0:len-1;
subplot(2, 2,1);
plot(n, input_signal(1:len));
title('Input Signal');
y = calculate_dft(input_signal, -1, len, 1);
#subplot(2, 2, 1)
#stem(n, abs(y));
#title('Magnitude spectrum');
subplot(2, 2, 2);
real_part = real(y);
stem(n, real_part);
title('Real Part');
subplot(2, 2, 3);
stem(n, imag(y));
title('Imaginary Part');
subplot(2,2,4);
plot(n,angle(y));
axis tight;
title("Angle spectrum");
% dft_output = zeros(1, len);
% for j=1:len
% for k = 1:len
% dft_output(j) = dft_output(j) + input_signal(k)*exp(-1i*2*pi*(k-1)*(j-1)/len);
% end
% end
%
% magnitude = abs(dft_output);
% realP = real(dft_output);
% imgP = imag(dft_output);
%
% subplot(2, 2, 1);
% stem(n, magnitude);
% subplot(2, 2, 2);
% stem(n, realP);
inverse_dft = calculate_dft(y, 1, len, len);
subplot(2, 2, 4)
plot(n, inverse_dft);
title('IDFT');