% FUNCTION: [Ez,r] = MakeWaves( L, f, thetaR, N )
%
% Makes sample fading plots as a function of space.
%
% Author: Greg Durgin,
%
% Inputs:
% L - the total
distance (in meters) to plot
% f - carrier
frequency (in Hz)
% thetaR
- azimuth orientation of movement in space (in degrees)
% N - total number of
plane waves (default is 100)
%
% Outputs:
% Ez
- vector of total z-component of E-field as function of space
% r - vector of
corresponding positions (same size as Ez)
%
% Sample Usage:
% To make a fading profile over 3 meters in space at a
carrier frequency of
% 2.0 GHz as observed with a 45-degree orientation in space:
% >> [Ez,r] = MakeWaves( 3, 2.4e9, 45 );
%
% Notes:
% The code provides a useful example of how to simulate the
small-scale
% fading that occurs when many plane waves superimpose in
space. Feel free
% to adapt this code for your ECE3065
project. The output of this
function
% is the z-component of electric field (in complex phasor form) as a
% function of space. The function automatically graphs Ez as a function of
% r on a dB-scale graph at the end of the routine, so that
you can see the
% types of outputs required for the project.
%
% The function "p" at the end of this file is your
azimuth spectrum of
% multipath
power. The default spectrum is omnidirectional (so the graphs
% will not depend on thetaR at
all!) but this can be changed fairly
% quickly.
function [E1,r,lambda, gamma, thetaMAX]
= MakeWaves( L, f, thetaR,
N )
% Set a default value of 100 for N
if ~exist('N')
N = 100;
end;
% Initialize Variables
wavelength = 3e8/f; % wavelength
M = 16*L/wavelength; % good number of space samples
k = 2*pi/wavelength; % freespace
wavenumber
r = linspace(0,L*(1-1/M),M); % make an array of sequential positions
y2 = -2*ones(size(r)); % make an array of sequential
positions in y component
E1 = zeros(size(r)); % initialize E1 component to
zero
E2 = zeros(size(r)); % initialize E2 component to
zero
theta = linspace(0,2*pi*(1-1/N),N);
% vector azimuth angles
thetaR = thetaR
* pi/180; % converting azimuth
orientation of movement in radian
urban = 120*pi/180; % Azimuth spectrum in urban
area
rural = 3*pi/180; % Azimuth spectrum in rural
area
v = 72.222222222; % maximum speed of the train
F0 = (2*10*urban*exp(j*0*thetaR))/(0^2*urban^2+1)*(1-(-1)^0*urban*exp(-pi/urban));
F1 = (2*10*urban*exp(j*1*thetaR))/(1^2*urban^2+1)*(1-(-1)^1*urban*exp(-pi/urban));
F2 = (2*10*urban*exp(j*2*thetaR))/(2^2*urban^2+1)*(1-(-1)^2*urban*exp(-pi/urban));
lambda = (1-abs(F1)^2/F0^2)^.5;
gamma = abs(F0*F2 - F1^2)/(F0^2 - abs(F1)^2);
thetaMAX =
1/2*angle(F0*F2-F1^2)*180/pi;
rhoRMS = 0:0.005:3;
NR = (sqrt(2*pi)*v*lambda*rhoRMS./wavelength)*sqrt(1+gamma*cos(2*(thetaR-thetaMAX))).*exp(-rhoRMS.^2);
% tau1 = wavelength*(exp(rhoRMS.^2)-1);
% tau2 = sqrt(2*pi)*v*lambda*rhoRMS;
% tau3 = sqrt(1+gamma*cos(2*(thetaR-thetaMAX)));
% tau = tau1./(tau2.*tau3);
tau =
wavelength.*(exp(rhoRMS.^2)-1)./(sqrt(2*pi)*v*lambda*rhoRMS.*sqrt(1+gamma*cos(2*(thetaR-thetaMAX))));
% Step through N azimuth angles and add
z-components of electric
% field from each direction. The expression below that is added to Ez each
% loop is a single plane wave. Its amplitude is related to the multipath
% angle spectrum and its phase is random, but tapers
"k" radians/meter in
% the direction of wave travel. By the end of this loop, Ez
will be the
% sum of N different plane waves, creating the strange
% constructive-destructive interference
pattern as a function of space.
%
for n=theta,
phase =
rand(1)*2*pi;
E1 = E1 + (p(n, thetaR))^.5 * exp( j*( phase + k*r*cos(n-thetaR) ) );
E2 = E2 + (p(n, thetaR))^.5 * exp( j*( phase + k*r*cos(n-thetaR) + k*y2*sin(n-thetaR) ) );
end;
% Normalize E1 and E2 so that average power is 1
Erms1 = (sum(abs(E1).^2)/length(E1))^.5; % root-mean-square
E1 = E1/Erms1; % normalize to average power 1
Erms2 = (sum(abs(E2).^2)/length(E2))^.5;
E2 = E2/Erms2;
figure(1);
subplot(2,1,1), plot(20*log(rhoRMS),
NR);
title(sprintf('Level crossing rate
at rural area %1.0f MHz, \\theta_R=%i',f/1e6,thetaR*180/pi));
xlabel('rho
(dB)');
ylabel('Average level crossing
rate');
subplot(2,1,2), plot(20*log(rhoRMS),
tau);
title(sprintf('Average fading
duration at rural area %1.0f MHz, \\theta_R=%i',f/1e6,thetaR*180/pi));
xlabel('rho
(dB)');
ylabel('Average fading duration');
% Plot a nice graph of power in dB as a function of space
figure(2);
subplot(3,1,1), plot(r, 20*log10(abs(E1)),'r'); %plot dB-scale of E1 as function of r
title(sprintf('Small-Scale Fading
of antenna 1 in Received Power at %1.0f MHz, \\theta_R=%i',f/1e6,thetaR*180/pi));
xlabel( 'Position (m)' );
ylabel( 'Normalized Received Power
(dB)' );
%title(sprintf('Small-Scale Fading
in Received Power at %1.0f MHz, \\theta_R=%i', f/le6, thetaR*180/pi));
subplot(3,1,2), plot(r, 20*log10(abs(E2)), 'b'); %plot
dB-scale of E2 as function of r
title(sprintf('Small-Scale Fading
of antenna 2 in Received Power at %1.0f MHz,
\\theta_R=%i',f/1e6,thetaR*180/pi));
xlabel( 'Position (m)' );
ylabel( 'Normalized Received Power
(dB)' );
% subplot(4,1,3), plot(r, max(20*log10(abs(E1)),
20*log10(abs(E2))), 'm'); % plot
dB-scale Emax as a function of r
% title(sprintf('Small-Scale
Fading of maximum signal in Received Power at %1.0f MHz,
\\theta_R=%i',f/1e6,thetaR*180/pi));
% xlabel( 'Position (m)' );
% ylabel( 'Normalized Received
Power (dB)' );
%plot dB-scale of E1, E2, and Emax
as a function of r
subplot(3,1,3), plot(r, 20*log10(abs(E1)),'r:', r,
20*log10(abs(E2)), 'b--', r, max(20*log10(abs(E1)), 20*log10(abs(E2))),
'k-');
title(sprintf('Comparison of
Small-Scale Fading in Received Power at %1.0f MHz,
\\theta_R=%i',f/1e6,thetaR*180/pi));
xlabel( 'Position (m) Antenna1(....) Antennal(----) Max signal received (solid line)' );
ylabel( 'Normalized Received Power
(dB)' );
return; % finish
% Here is the function defining the azimuth spectrum, the
distribution of
% multipath
power as a function of angle-of-arrival.
%
% The example below is for an omnidirectional
spectrum (a constant).
% Change it to see other types of angle spectra.
%
function P = p(theta, thetaR)
urban =
120*pi/180;
rural = 3*pi/180;
P =
10*exp(-abs((theta - thetaR)/urban));
%P =
10*exp(-abs((theta - thetaR)/rural));
% Here is another example, the cosine-squared spectrum
% P = cos(theta)^2;
return;