CMA Equaliser and FSK I'm looking into equalisation for FSK over HF channels for Digital Voice applications.  As a starting point, I am trying to understand the CMA algorithm.

I found a nice CMA Matlab example for BPSK in a previous post here on this site.

I have modified the simulation to use 2FSK but it doesn't converge to the same result.  The first plot is for BPSK modulation, the 2nd plot for FSK.  The combination of the channel and equaliser impulse response is converging to 1 for BPSK.





Any suggestions of how I can get this simulation working for FSK would be much appreciated!

Here is the source code.  the variable "tx_type" switches the tx modulation source.

N = 20000;          % # symbols
h = [1 0.45 -0.2];  % channel impulse response
h = h/norm(h);
Le = 20;            % equalizer length
mu = 1E-3;          % step size
snr = 30;           % snr in dB
M = 10;             % oversample rate

tx_type = "bpsk";   % select modulation type here "bpsk" or "fsk"

if strcmp(tx_type, "bpsk")
  s0 = round( rand(N,1) )*2 - 1;     % BPSK signal
  s0M = zeros(N*M,1);                % oversampled BPSK signal
  k = 1;
  for i=1:M:N*M
   s0M(i:i+M-1) = s0(k);
    k ++;
  end
end

if strcmp(tx_type, "fsk")
  tx_bits = round(rand(1,N));

  % continuous phase FSK modulator

  w1 = pi/4;
  w2 = pi/2;
  tx_phase = 0;
  tx = zeros(M*N,1);

  for i=1:N
    for k=1:M
      if tx_bits(i)
        tx_phase += w2;
      else
        tx_phase += w1;
      end
      tx((i-1)*M+k) = exp(j*tx_phase);
    end
  end

  s0M = real(tx);
end

s = filter(h,1,s0M);                % filtered signal

% add Gaussian noise at desired snr

n = randn(N*M,1);
vs = var(s);
vn = vs*10^(-snr/10);
n = sqrt(vn)*n;
r = s + n;          % received signal

e = zeros(N*M,1);   % error
w = zeros(Le,1);    % equalizer coefficients
w(Le)=1;            % actual filter taps are flipud(w)!

yd = zeros(N*M,1);

for i = 1:N*M-Le,
    x = r(i:Le+i-1);
    y = w'*x;
    yd(i)=y;
    e(i) = y^2 - 1;
    w = w - mu * e(i) * y * x;
end

np = 100;           % # sybmols to plot (last np will be plotted); np &lt; N!

figure(1); clf;
subplot(311), plot(e.*e), title('error')
subplot(312), stem(conv(flipud(w),h)), title('equalized channel impulse response')
subplot(313);
plot(1:np,s0M(N-np+1:N),1:np,yd(N-np+1-Le+1:N-Le+1))
title('transmitted and equalized signal'), legend('transmitted','equalized'), axis([0,np,-1.5,2])


Thanks,

David
The confusion comes from the fact that what is tagged as "transmitted" isn't the real transmitted waveform but its baseband representation, which are $\left\{+1,-1\right\}$ symbols in 1 dimension. For 2-FSK, that representation would be $\left\{\left(+1,0\right), \left(0,+1\right)\right\}$ symbols in 2 dimensions.

A common convenient way to handle such 2 dimensional symbols is to represent them with complex number, mapping the coordinate $\left(+1,0\right)$ to the purely real value $+1$ and the coordinate $\left(0,+1\right)$ to the purely imaginary value $j=\sqrt{-1}$. Using that mapping, the symbol generation code for 2-FSK becomes:

if strcmp(tx_type, "fsk")
  tx_bits = round(rand(1,N));
  s0 = zeros(N,1);
  s0M = zeros(N*M,1);

  for i=1:N
    if tx_bits(i)
      s0(i) = 1i;
    else
      s0(i) = 1;
    end
  end

  k = 1;
  for i=1:M:N*M
   s0M(i:i+M-1) = s0(k);
    k ++;
  end
end


The problem then is that the CMA weight update rule from the original post is given for a purely real signal. The adjusted weight update rule for complex signals is given by:

$$
\begin{align}
  y &amp;= \vec w_{n}^H x\\
  e_n &amp;= \left|y\right|^2 - 1 \\
  \vec w_{n+1} &amp;= \vec w_{n} - \mu \vec \nabla e_n \\
           &amp;= \vec w_{n} - \mu' e_n \Re\left\{y^\star x\right\}
\end{align}
$$

This can be implemented with:

for i = 1:N*M-Le,
    x = r(i:Le+i-1);
    y = w'*x;
    yd(i)=y;
    e(i) = abs(y)^2 - 1;
    w = w - mu * e(i) * real(conj(y) * x);
end


Plotting the results (showing the real part and imaginary parts in separate graphs) should give you something like the following: