% Code to demonstrate the implementation of a narrowband fading channel
% with memory by use of frequency based filtering detailed in section 4.1.2
% Tim Brown, January 2012

N = 1000; % Arbritrary number of samples, 1000 is useful to see clear comparison.
Delta_t = 0.1; % Arbritrary change in time in seconds
P = 1; % Arbritrary power level
G = 1; % Arbritrary antenna gain

% Generate bath tub filter
freq = [-(0.5/Delta_t):((1/Delta_t)/(N-1)):(0.5/Delta_t)];
f_d_max = 0.25 / Delta_t; % This allows for a sampling rate twice the maximum Doppler shift
s = P .* G ./ (pi .* sqrt((f_d_max^2) - (freq.^2)));

% Generate N frequency samples
H = (randn(N, 1) + j*randn(N, 1))/ sqrt(2);

% Derive time domain samples
H = H.*s.';
h = ifft(H);

% Take the Fourier transform of both the channel with and without memory
% after normalisation
h = h ./ mean(abs(h));
fft_h = fftshift(fft(h));
h_nomem = (randn(N, 1) + j*randn(N, 1))/ sqrt(2);
h_nomem = h_nomem ./ mean(abs(h));
fft_w = fftshift(fft(h_nomem));

% Plot a comparison of the channels with and without memory in the
% frequency domain
f = [-(0.5/Delta_t):((1/Delta_t)/(N-1)):(0.5/Delta_t)];
figure
plot(f,20.*log10(abs(fft_h)),f,20.*log10(abs(fft_w)))
set(gca,'FontSize',10)
set(gcf,'color','white')
xlabel('Frequency (Hz)')
ylabel('Magnitude (dB)')
legend('Without Memory','With Memory')

% Plot in the time domain, the first 100 for clarity, though this can be
% changed by adjusting 'range' below
range = 1:100;
t = (0:(N-1))*Delta_t;
figure
plot(t(range),20.*log10(abs(h(range))),t(range),20.*log10(abs(h_nomem(range))))
set(gca,'FontSize',10)
set(gcf,'color','white')
xlabel('Time (s)')
ylabel('Magnitude (dB)')
legend('Without Memory','With Memory')