matlab的多用户检测

function [SNR] = simulation_BER_SNR % Rayleigh fading channel (实数) % single path %Shows the BER as a function of the SNR in a Rayleigh fading channel SNR = [0, 2, 4, 6, 8, 10, 12,14,16,18, 20, 22, 24, 26, 28, 30]; N = 32; %the length of spread code K = 10; %the number of users T = 5000; alpha = 1; A=ones(K,T); %for i=1:K % [time, rayl]=rayleigh_fading(10, 5, T/5); % A(i,:) = rayl*alpha; % A is regarded as the rayleigh fading matrix %end Aavg = sum(A.^2,2)/T; SNR_dec = 10.^(SNR / 10); Matcodes = no_orth_spreadcode(N); C = Matcodes(:,1:K);% spread spre matrix % C1=C1(:,1:K); ku = 1; Eb = C(:,ku)' * C(:,ku)*Aavg(ku); mf = C(:,ku); % Matched filter % show BER as performance from different SNR for k=1:length(SNR) disp('loop'); sigma = Eb/SNR_dec(k)/2; Sentbit = sign(randn(K,T));%(10,10000);K is the number of user.考虑多用户（单径） Y=C*(A.*Sentbit) + sqrt(sigma)*randn(N,T);%the received signal(32,10000);C is spread matrix; % A is channel fading matrix; sqrt(sigma)*randn(N,T) is noise. % matched filter Receivedbit1 = sign(mf'*Y);%ku=1,接收第一个用户的信息 Perr1(k) = sum(abs(Receivedbit1-Sentbit(ku,:)))/2/T; % decorrelator Receivedbit2 = sign(inv(C'*C)*C'*Y); Perr2(k) = sum(abs(Receivedbit2(ku,:)-Sentbit(ku,:)))/2/T; % MMSE Receivedbit3 = sign(inv(C'*C+sigma*diag(1./Aavg))*C'*Y); Perr3(k) = sum(abs(Receivedbit3(ku,:)-Sentbit(ku,:)))/2/T; % SIC Z=C'*Y; [Z2,order]=sort(abs(Z)); for j=K:-1:1 for time=1:T Z3=Z(order(j,time),time); sent_bit=sign(Z3); if K>5 R=C'*C; Z(:,time) = Z(:,time) - R(:,order(j,time))*(A(order(j,time),time).*sent_bit);%abs(Z2) end Receivedbit4(order(j,time),time)=sent_bit; end end Perr4(k) = sum(abs(Receivedbit4(ku,:)-Sentbit(ku,:)))/2/T; % DFD F=chol(C'*C); Finv=inv(F'); Yfb = Finv*C'*Y; [Y2,order]=sort(abs(Yfb)); for j=K:-1:1 for time=1:T Y3 = Yfb(order(j,time),time); sent_bit=sign(Y3); Yfb(:,time)=Yfb(:,time)-F(:,order(j,time))*(A(order(j,time),time).*sent_bit);%abs(Yfb(j,:)) Receivedbit5(order(j,time),time)=sent_bit; end end Perr5(k) = sum(abs(Receivedbit5(ku,:)-Sentbit(ku,:)))/2/T; % PIC Z3=C'*Y; Receivedbit6 = sign(Z3); B=C'*C-triu(tril(C'*C)); Z3=C'*Y-B*(A.*Receivedbit6);%triu(tril(Z3)) Receivedbit6 = sign(Z3); Perr6(k) = sum(abs(Receivedbit6(ku,:)-Sentbit(ku,:)))/2/T; end if (nargout==0) semilogy(SNR, Perr1, 'bx-');hold on; semilogy(SNR, Perr2, 'ro-');hold on; semilogy(SNR, Perr3, 'm+-');hold on; semilogy(SNR, Perr4, 'k*-');hold on; semilogy(SNR, Perr5, 'cs-');hold on; semilogy(SNR, Perr6, 'bv-');hold on; xlabel('SNR'); ylabel('BER'); legend('Matched filter', 'Decorrelator','MMSE','SIC','DFD','PIC','Theoretical Single User bound') end hold off; fprintf('Perr calcule.\n'); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function [Ts, z, dB] = rayleigh_fading(f_D, t, f_s) % % function [Ts, z_dB] = rayleigh_fading(f_D, t, f_s) % generates a Rayleigh fading signal for given Doppler frequency f_D % during the time period [0,t], with sampling frequency f_s >= 1000Hz. % % Input(s) % -- f_D : [Hz] [1x1 double] Doppler frequency % -- t : simulation time interval length, time interval [0,t] % -- f_s : [Hz] sampling frequency, set to 1000 if smaller. % Output(s) % -- Ts : [Sec][1xN double] time instances for the Rayleigh signal % -- z_dB : [dB] [1xN double] Rayleigh fading signal % Required parameters if f_s < 1000; f_s=1000; % [Hz] Minimum required sampling rate end; N = ceil(t*f_s); % Number of samples % Ts contains the time instances at which z_dB is specified Ts = linspace(0,t,N); if mod(N,2)==1 N=N+1; % Use even number of samples in calculation end f = linspace(-f_s, f_s, N); % [Hz] Frequency samples used in calculation % Generate complex Gaussian samples with line spectra in frequency domain % Inphase : Gfi_p = randn(2,N/2); CGfi_p = Gfi_p(1,:)+i*Gfi_p(2,:); CGfi = [conj(fliplr(CGfi_p)) CGfi_p ]; % Quadrature : Gfq_p = randn(2,N/2); CGfq_p = Gfq_p(1,:)+i*Gfq_p(2,:); CGfq = [conj(fliplr(CGfq_p)) CGfq_p ]; % Generate fading spectrum, this is used to shape the Gaussian line spectra omega_p = 1; % This makes sure that the average received envelop can be 0dB S_r = omega_p/4/pi./(f_D*sqrt(1-(f/f_D).^2)); % Take care of samples outside the Doppler frequency range, let them be 0 idx1 = find(f>f_D); idx2 = find(f<-f_D); S_r(idx1)=0; S_r(idx2)=0; % Generate r_I(t) and r_Q(t) using inverse FFT: r_I = N*ifft(CGfi.*sqrt(S_r)); r_Q = -i*N*ifft(CGfq.*sqrt(S_r)); % Finally, generate the Rayleigh distributed signal envelope z = sqrt(abs(r_I).^2+abs(r_Q).^2); z_dB = 20*log10(z); % Return correct number of points z_dB = z_dB(1:length(Ts)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function [codeMatrix] = no_orth_spreadcode(SF) % This function returns a matrix whose lignes are code vectors. These codes are % not orthonormal C1=hadamard(SF); % Generates a matrix with orthogonal vectors C2=sign(randn(SF,1)); % Generates a cell code to ensure spreading C3=diag(C2)*C1/sqrt(SF); % Matrix normalized C4=(rand(SF,SF)-.5)/3; % Generates a noise matrix codeMatrix = C4+C3; % Sums the code and the noise to have non-orthogonal codes %end function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%