%  This program simulates a rectangular pulse shape signal.
clear;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Simulation Parameters                  %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Simulation Parameters                  %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%fmax=input('Enter simulation bandwidth= ');
fmax=128
%N=input('Enter number of samples= ');
N=16384
fs=2*fmax;                % Sampling Frequency
df=2*fmax/N               % Frequency spacing
dt=1./(df.*N)             % Time spacing
t=(1:N)*dt-dt;            % Time samples
Tmax=N*dt                 % Simulation time
f=(1:N)*df-df;            % Frequency samples
f2=f-N/2*df;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                System Parameters                      %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
L=64                       % Number of data bits
Tb=1;                      % Duration of signal
Nspb=Tb./dt                % Number of samples per bit
N0=input('Noise power spectral density, N_0= ');
%N0=0.004
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Generate Pulse Shape                   %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

pf(1:N)=0.0*ones(1,N)*1e-160; % Make signal very close to zero
pf(1:N/2)=pf(1:N/2)+Tb*sinc(f(1:N/2)*Tb).*exp(-j*pi*f(1:N/2)*Tb);
pf(N/2+2:N)=conj(fliplr(pf(2:N/2)));
pt=ifft(pf)./dt;

b=sign(randn(1,L));
data=b(1)*ones(1,N/2);
for l=2:L
data=data+b(l)*exp(-j*2*pi*f(1:N/2)*(l-1)*Tb);
end
xf(1:N/2)=pf(1:N/2).*data;

xf(N/2+2:N)=conj(fliplr(xf(2:N/2)));
xt=real(ifft(xf)./dt);    % Real part throws away the imaginary part which is 
                          % nonzero only because of roundoff

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Plot the Signal                        %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure(1)
subplot(2,1,1), plot(t,xt)
grid on
axis([0 10 -1.5 1.5])
xlabel('time')
ylabel('x(t)')
subplot(2,1,2), plot(f2,10*log10(abs(fftshift(xf))))
grid on
axis([-10 10 -50 10])
xlabel('frequency')
ylabel('|X(f)| (dB)')


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Generate Noise                         %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

sigma=sqrt(N0*fmax)
nt=sigma*randn(1,N);
nf=fft(nt)*dt;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Plot the Noise                         %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

figure(2)
subplot(3,1,1), plot(t,nt)
grid
axis([0 10 -1 1])
xlabel('time')
ylabel('n(t)')
subplot(3,1,2), plot(f2,10*log10(abs(fftshift(nf))))
grid
axis([-30 30 -30 0])
xlabel('frequency')
ylabel('N(f)')
a=-1:.01:1;
subplot(3,1,3), hist(nt,a)
grid


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Add signal and noise                   %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
yt=xt+nt;
yf=fft(yt)*dt;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Plot the Signal and Noise              %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure(3)
subplot(2,1,1), plot(t,yt)
axis([0 10 -2 2])
grid
xlabel('time')
ylabel('y(t)')
subplot(2,1,2), plot(f2,10*log10(abs(fftshift(yf))))
grid on
axis([-30 30 -30 10])
xlabel('frequency')
ylabel('|X(f)| (dB)')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%           Filter the Received Signal                  %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

zf=yf.*pf;
zt=real(ifft(zf)./dt);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                       %
%                Plot the Filter Output                 %
%                                                       %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure(4)
subplot(2,1,1), plot(t,zt)
axis([0 10 -2 2])
grid
xlabel('time')
ylabel('z(t)')
subplot(2,1,2), plot(f2,10*log10(abs(fftshift(zf))))
grid on
axis([-30 30 -30 10])
xlabel('frequency')

eyediagram(zt,Nspb)
grid on
