1;
+% Gray coded QPSK modulation function
+
+function symbol = qpsk_mod(two_bits)
+ two_bits_decimal = sum(two_bits .* [2 1]);
+ switch(two_bits_decimal)
+ case (0) symbol = 1;
+ case (1) symbol = j;
+ case (2) symbol = -j;
+ case (3) symbol = -1;
+ endswitch
+endfunction
+
+
+% Gray coded QPSK demodulation function
+
+function two_bits = qpsk_demod(symbol)
+ bit0 = real(symbol*exp(j*pi/4)) < 0;
+ bit1 = imag(symbol*exp(j*pi/4)) < 0;
+ two_bits = [bit1 bit0];
+endfunction
+
+
function sim_out = run_sim(sim_in)
- Rs = 100;
+ Rs = 62.5;
Fs = 8000;
M = Fs/Rs;
+ Tcp = 0.004; Ncp = Tcp*M;
+ foffset = 0;
+ woffset = 2*pi*foffset/Fs;
+ bps = sim_in.bps;
EbNodB = sim_in.EbNodB;
verbose = sim_in.verbose;
hf_en = sim_in.hf_en;
Ns = sim_in.Ns; % step size for pilots
Nc = sim_in.Nc; % Number of cols, aka number of carriers
- Nbitsperframe = (Ns-1)*Nc;
- printf("Nbitsperframe: %d\n", Nbitsperframe);
- Nrowsperframe = Nbitsperframe/Nc;
- printf("Nrowsperframe: %d\n", Nrowsperframe);
+ Nbitsperframe = (Ns-1)*Nc*bps;
+ Nrowsperframe = Nbitsperframe/(Nc*bps);
+ if verbose
+ printf("Rs:.........: %4.2f\n", Rs);
+ printf("M:..........: %d\n", M);
+ printf("bps:.........: %d\n", bps);
+ printf("Nbitsperframe: %d\n", Nbitsperframe);
+ printf("Nrowsperframe: %d\n", Nrowsperframe);
+ end
% Important to define run time in seconds so HF model will evolve the same way
% for different pilot insertion rates. So lets work backwards from approx
Nframes = floor((Nrows-1)/Ns);
Nbits = Nframes * Nbitsperframe; % number of payload data bits
- Nr = Nbits/Nc; % Number of data rows to get Nbits total
+ Nr = Nbits/(Nc*bps); % Number of data rows to get Nbits total
if verbose
printf("Nc.....: %d\n", Nc);
% double check if Nbits fit neatly into carriers
- assert(Nbits/Nc == floor(Nbits/Nc), "Nbits/Nc must be an integer");
+ assert(Nbits/(Nc*bps) == floor(Nbits/(Nc*bps)), "Nbits/(Nc*bps) must be an integer");
printf("Nframes: %d\n", Nframes);
% extra row of pilots at end
printf("Nrp....: %d (number of rows including pilots)\n", Nrp);
+ % set up HF model ---------------------------------------------------------------
+
+ if hf_en
+
+ % some typical values, or replace with user supplied
+
+ dopplerSpreadHz = 1.0; path_delay_ms = 1;
+
+ if isfield(sim_in, "dopplerSpreadHz")
+ dopplerSpreadHz = sim_in.dopplerSpreadHz;
+ end
+ if isfield(sim_in, "path_delay_ms")
+ path_delay_ms = sim_in.path_delay_ms;
+ end
+ printf("Doppler Spread: %3.2f Hz Path Delay: %3.2f ms\n", dopplerSpreadHz, path_delay_ms);
+
+ path_delay_samples = path_delay_ms*Fs/1000;
+
+ randn('seed',1);
+ spread1 = doppler_spread(dopplerSpreadHz, Fs, sim_in.Nsec*Fs*1.1);
+ spread2 = doppler_spread(dopplerSpreadHz, Fs, sim_in.Nsec*Fs*1.1);
+
+ % sometimes doppler_spread() doesn't return exactly the number of samples we need
+
+ assert(length(spread1) >= Nrp*M, "not enough doppler spreading samples");
+ assert(length(spread2) >= Nrp*M, "not enough doppler spreading samples");
+
+ hf_gain = 1.0/sqrt(var(spread1)+var(spread2));
+ % printf("nsymb: %d lspread1: %d\n", nsymb, length(spread1));
+ end
+
% simulate for each Eb/No point ------------------------------------
for nn=1:length(EbNodB)
rand('seed',1);
randn('seed',1);
- EbNo = 10 .^ (EbNodB(nn)/10);
+ EbNo = bps * (10 .^ (EbNodB(nn)/10));
variance = 1/(M*EbNo/2);
- % generate tx bits and insert pilot rows and border cols
-
- tx_bits = []; tx_bits_np = [];
- for f=1:Nframes
- tx_bits = [tx_bits; ones(1,Nc+2)];
- for r=1:Nrowsperframe
- arowofbits = rand(1,Nc) > 0.5;
- tx_bits = [tx_bits; [1 arowofbits 1]];
- tx_bits_np = [tx_bits_np; arowofbits];
- end
- end
- tx_bits = [tx_bits; [1 ones(1,Nc) 1]]; % final row of pilots
- [nr nc] = size(tx_bits);
- assert(nr == Nrp);
- % map to BPSK symbols
+ % generate tx bits
- tx_sym = 2*tx_bits - 1;
+ tx_bits = rand(1,Nbits) > 0.5;
-#{
- % upsample to rate Fs using Zero Order Hold (ZOH)
+ % map to symbols in linear array
- tx_sym_os = zeros(Nrp*M, Nc+2);
- for r=1:Nrp
- for rr=(r-1)*M+1:r*M
- tx_sym_os(rr,:) = tx_sym(r,:);
+ if bps == 1
+ tx_sym_lin = 2*tx_bits - 1;
+ end
+ if bps == 2
+ for s=1:Nbits/bps
+ tx_sym_lin(s) = qpsk_mod(tx_bits(2*(s-1)+1:2*s));
end
end
-#}
+
+ % place symbols in multi-carrier frame with pilots and boundary carriers
+
+ tx_sym = []; s = 1;
+ for f=1:Nframes
+ aframe = zeros(Nrowsperframe,Nc+2);
+ aframe(1,:) = 1;
+ for r=1:Nrowsperframe
+ arowofsymbols = tx_sym_lin(s:s+Nc-1);
+ s += Nc;
+ aframe(r+1,2:Nc+1) = arowofsymbols;
+ end
+ tx_sym = [tx_sym; aframe];
+ end
+ tx_sym = [tx_sym; ones(1,Nc+2)]; % final row of pilots
+ [nr nc] = size(tx_sym);
+ assert(nr == Nrp);
% OFDM up conversion and upsampling to rate Fs
w = (0:Nc+1)*2*pi*Rs/Fs;
- W = zeros(M,Nc+2);
+ W = zeros(Nc+2,M);
for c=1:Nc+2
- W(:,c) = exp(j*w(c)*(0:M-1));
+ W(c,:) = exp(j*w(c)*(0:M-1));
end
Nsam = Nrp*M;
- tx = zeros(Nrp*M,1);
+ tx = zeros(1,Nrp*M);
for r=1:Nrp
for c=1:Nc+2
- acarrier = tx_sym(r,c) * W(:,c);
+ acarrier = tx_sym(r,c) * W(c,:);
tx((r-1)*M+1:r*M) += acarrier/M;
end
end
rx = tx;
+
+ if hf_en
+ %rx = hf_gain * tx(1:Nsam) .* spread1(1:Nsam);
+ %rx = hf_gain * [tx(path_delay_samples+1:Nsam) zeros(1,path_delay_samples)] .* spread2(1:Nsam);
+ rx = hf_gain * [tx(path_delay_samples+1:Nsam) zeros(1,path_delay_samples)];
+ end
- noise = sqrt(variance)*(0.5*randn(Nrp*M,1) + j*0.5*randn(Nrp*M,1));
- rx += noise;
+ rx = rx .* exp(j*woffset*(1:Nsam));
+ noise = sqrt(variance)*(0.5*randn(1,Nsam) + j*0.5*randn(1,Nsam));
+ rx += noise;
+
% downconvert, downsample and integrate using OFDM
rx_sym = zeros(Nrp, Nc+2);
for r=1:Nrp
for c=1:Nc+2
- acarrier = rx((r-1)*M+1:r*M) .* conj(W(:,c));
+ acarrier = rx((r-1)*M+1:r*M) .* conj(W(c,:));
rx_sym(r,c) = sum(acarrier);
end
end
-
+
% pilot based phase est, we use known tx symbols as pilots ----------
phase_est_pilot_log = 10*ones(Nrp,Nc+2);
end % r=1:Ns:Nrp-Ns
end % c=2:Nc+1
- % remove pilots to give us just data symbols
-
- rx_np = [];
+
+ % remove pilots to give us just data symbols and demodulate
+
+ rx_bits = []; rx_np = [];
for r=1:Nrp
if mod(r-1,Ns) != 0
- rx_np = [rx_np; rx_corr(r,2:Nc+1)];
+ arowofsymbols = rx_corr(r,2:Nc+1);
+ rx_np = [rx_np arowofsymbols];
+ if bps == 1
+ arowofbits = real(arowofsymbols) > 0;
+ end
+ if bps == 2
+ arowofbits = zeros(1,Nc);
+ for c=1:Nc
+ arowofbits((c-1)*2+1:c*2) = qpsk_demod(arowofsymbols(c));
+ end
+ end
+ rx_bits = [rx_bits arowofbits];
end
end
+ assert(length(rx_bits) == Nbits);
+
% calculate BER stats as a block, after pilots extracted
- rx_bits_np = real(rx_np) > 0;
- errors = xor(tx_bits_np, rx_bits_np);
- Nerrs = sum(sum(errors));
+ errors = xor(tx_bits, rx_bits);
+ Nerrs = sum(errors);
printf("EbNodB: %3.2f BER: %4.3f Nbits: %d Nerrs: %d\n", EbNodB(nn), Nerrs/Nbits, Nbits, Nerrs);
figure(1)
plot(real(tx))
figure(2)
- Tx = abs(fft(tx.*hanning(Nsam)));
+ Tx = abs(fft(tx.*hanning(Nsam)'));
Tx_dB = 20*log10(Tx);
dF = Fs/Nsam;
plot((1:Nsam)*dF, Tx_dB);
plot(rx_np,'+');
axis([-2 2 -2 2]);
+
if hf_en
figure(4); clf;
- plot(abs(hf_model(:,2:Nc+1)));
+ subplot(211)
+ plot(abs(spread1(1:Nsam)));
+ %hold on; plot(abs(spread2(1:Nsam)),'g'); hold off;
+ subplot(212)
+ plot(angle(spread1(1:Nsam)));
end
+
figure(5); clf;
plot(phase_est_log(:,2:Nc+1),'+', 'markersize', 10);
hold on;
plot(phase_est_pilot_log(:,2:Nc+1),'g+', 'markersize', 5);
+
+#{
if sim_in.hf_en
plot(angle(hf_model(:,2:Nc+1)));
end
+#}
+
axis([1 Nrp -pi pi]);
end
function run_single
+ sim_in.bps = 2;
sim_in.Nc = 8;
sim_in.Ns = 8;
sim_in.Nsec = 10;
projects / freetel-svn-tracking.git / commitdiff