% this handles the exact-gradient-descent CMA algorithms for the% bergulator, and is called by ui_traj.m%% if not available for your platform, you can create a .mex version % by typing "mex -O berg_cmegd.c" at the system prompt % Copyright 1997-1998 Phil Schniterfunction [f,y,srcme] = berg_cmegd(N,N_f,D_f,D_e,is_FSE,kappa,mu,r,f_init,alg,sigma2_n,C); f_cur = f_init; f = zeros(N_f,floor((N-1)/D_f)+1);    	% equalizer history y = zeros(1,floor((N-1)/D_e)+1);      	% soft decisions srcme = zeros(1,floor((N-1)/D_e)+1);  	% sqrt(CM) error history if alg==0,                		% PAM source with real noise   for i = 1:N,     reg = r(1, N_f+(1+is_FSE)*i:-1:(1+is_FSE)*i+1);     y_cur = reg*f_cur;      		% for constellation plot...     h = C*f_cur;     hTh = norm(h)^2;     fTf = norm(f_cur)^2;     srcme_cur = sqrt( kappa^2 - 2*kappa*(hTh + sigma2_n*fTf)...        + 3*hTh^2 + (kappa-3)*norm(h.^2).^2 ...        + 6*sigma2_n*hTh*fTf + 3*sigma2_n^2*fTf^2 );     dJdf = (kappa-3)*C'*(h.^3) + (3*hTh + 3*sigma2_n*fTf - kappa)*...        (C'*h + sigma2_n*f_cur);     f_cur = f_cur  - mu*dJdf;     f(:,floor((i-1)/D_f)+1) = f_cur;     y(floor((i-1)/D_e)+1) = y_cur;     srcme(floor((i-1)/D_e)+1) = srcme_cur;   end elseif alg==1,                 	% PAM source with rot.inv. noise   for i = 1:N,     reg = r(1, N_f+(1+is_FSE)*i:-1:(1+is_FSE)*i+1);     y_cur = reg*f_cur;      % for constellation plot...     h = C*f_cur;     hTh = norm(h)^2;     fTf = norm(f_cur)^2;     srcme_cur = sqrt( (kappa-3)*sum(abs(h).^4) + 2*hTh^2 + abs(h.'*h)^2 ...       + 2*sigma2_n^2*fTf^2 + 4*sigma2_n*hTh*fTf...       - 2*kappa*(hTh + sigma2_n*fTf) + kappa^2 );     dJdf = (kappa-3)*C'*(h.*(abs(h).^2)) + (h.'*h)*C'*conj(h)...       + (2*hTh + 2*sigma2_n*fTf - kappa)*(C'*h + sigma2_n*f_cur);     f_cur = f_cur - mu*dJdf;     f(:,floor((i-1)/D_f)+1) = f_cur;     y(floor((i-1)/D_e)+1) = y_cur;     srcme(floor((i-1)/D_e)+1) = srcme_cur;   end elseif alg==2,                          % QAM source with rot.inv. noise   for i = 1:N,     reg = r(1, N_f+(1+is_FSE)*i:-1:(1+is_FSE)*i+1);     y_cur = reg*f_cur;      % for constellation plot...     h = C*f_cur;     hTh = norm(h)^2;     fTf = norm(f_cur)^2;     srcme_cur = sqrt( (kappa-2)*sum(abs(h).^4) + 2*hTh^2 ...       + 2*sigma2_n^2*fTf^2 + 4*sigma2_n*hTh*fTf...       - 2*kappa*(hTh + sigma2_n*fTf) + kappa^2 );     dJdf = (kappa-2)*C'*(h.*(abs(h).^2))...       + (2*hTh + 2*sigma2_n*fTf - kappa)*(C'*h + sigma2_n*f_cur);     f_cur = f_cur - mu*dJdf;     f(:,floor((i-1)/D_f)+1) = f_cur;     y(floor((i-1)/D_e)+1) = y_cur;     srcme(floor((i-1)/D_e)+1) = srcme_cur;   end else   error('That algorithm is not supported by berg_cmegd.m'); end;