calc_dsegamma_qam.m

有关信道估计和信道均衡的仿真程序

% this chooses DSE-CMA gamma for the bergulator under uniform dither
% and an M-QAM source
 
% Copyright 1997-1998 Phil Schniter

function gam_opt = calc_dsegamma_QAM(dse_alpha,alphabet,kappa);

 if isreal(alphabet), error('expecting CMPLX-valued alphabet!'); end;
 alf = alphabet(1:length(alphabet(:))/2);
 M_a = length(alf);
 M = 10;

 % check if kappa will work
 if dse_alpha >...
	 max([(abs(alf).^2-kappa).*real(alf),(abs(alf).^2-kappa).*imag(alf)]),
   gam_opt = kappa;
   return;
 end;

 % estimate best gamma
 [gam,inds] = sort([kappa,linspace(0.98,2,2*M)]);
 [dum,kap_ind] = min(inds);
 J = zeros(size(gam));
 for i=1:length(gam),
   J(i) = sum(sum( (min(dse_alpha*ones(1,M_a),max(-dse_alpha*ones(1,M_a),...
     	(abs(alf).^2-gam(i)).*real(alf))).*conj(alf)) +...
     	j*(min(dse_alpha*ones(1,M_a),max(-dse_alpha*ones(1,M_a),...
	(abs(alf).^2-gam(i)).*imag(alf))).*conj(alf)) ));
 end;
 [J_opt,ind] = min(abs(J));
 gam_opt = gam(ind);
 J_kap = abs(J(kap_ind));

%figure(7)
%plot(gam,abs(J))

 % refine estimate
 for k=1:ceil(16/log10(M)),
   gam = linspace(gam_opt-1/M^k,gam_opt+1/M^k,2*M);
   J = zeros(size(gam));
   for i=1:length(gam),
     J(i) = sum(sum( (min(dse_alpha*ones(1,M_a),max(-dse_alpha*ones(1,M_a),...
       	(abs(alf).^2-gam(i)).*real(alf))).*conj(alf))+...
	j*(min(dse_alpha*ones(1,M_a),max(-dse_alpha*ones(1,M_a),...
	(abs(alf).^2-gam(i)).*imag(alf))).*conj(alf)) ));
   end;
   [J_opt,ind] = min(abs(J));
   gam_opt = gam(ind);
 end;

return

 % output results
 format long
 gam_opt
 kappa
 format short

 figure(7)
 hold on; plot(gam_opt,J_opt,'ro'); hold off;
 hold on; plot(kappa, abs( sum(sum( (min(dse_alpha*ones(1,M_a),...
	max(-dse_alpha*ones(1,M_a), (abs(alf).^2-kappa).*real(alf))...
	).*conj(alf)) + j*(min(dse_alpha*ones(1,M_a), ...
	max(-dse_alpha*ones(1,M_a),(abs(alf).^2-kappa).*imag(alf))...
	).*conj(alf)) )) ),'gd'); hold off;
 zoom on; grid on;

 figure(8)
 x = j*linspace(0,2,200).'*ones(1,200)+ones(200,1)*linspace(0,2,200);
 surf(real(x(1,:)),real(x(1,:)),abs( min(dse_alpha*ones(size(x)),...
	max(-dse_alpha*ones(size(x)),(abs(x).^2-gam_opt).*...
	real(x)))+j*min(dse_alpha*ones(size(x)),max(-dse_alpha*ones(...
	size(x)),(abs(x).^2-gam_opt).*imag(x))) ));
 hold on; 
 sqrt_M_s = sqrt(length(alphabet(:)));
 a_pos = alphabet(1:sqrt_M_s/2, sqrt_M_s/2+1:sqrt_M_s); a_pos = a_pos(:);
 plot3(real(a_pos),imag(a_pos), abs( min(dse_alpha*ones(size(a_pos)),...
	max(-dse_alpha*ones(size(a_pos)),(abs(a_pos).^2-gam_opt).*...
	real(a_pos)))+j*min(dse_alpha*ones(size(a_pos)),max(-dse_alpha*...
	ones(size(a_pos)),(abs(a_pos).^2-gam_opt).*imag(a_pos))) ),'r*'); 
 hold off;
 xlabel('real');ylabel('imag');
 title('abs(DSE-CMA error function) for \gamma_{opt}')
 grid on;