084eab45c3739f31742361868db67e8daae51cad
[4yp.git] / adaptiveCMA.m
1 function adaptFilterOut = adaptiveCMA(rSampled)
2   %% adaptive filter
3   %% CMA
4   taps = 19; % ODD taps
5   hxx = zeros(taps, 1);
6   %% hxx: real indices  -K, ..., 0, ..., K.   K = floor(taps/2)
7   %%    MATLAB indices   1      1+K     taps
8   %% initialize hxx, hxx[0] = 1, hxx[k] = hxx[-k] = 0
9   hxx(ceil(taps/2)) = 1;
10
11   mu = 1e-3;
12   numSymbs = length(rSampled);
13
14   %% Check average energy of symbols
15   rSampledUnitEnergy = normalizeEnergy(rSampled, numSymbs, 1);
16
17   adaptFilterOut = zeros(numSymbs, 1);
18   converged = 0;
19   convergeCount = 0;
20
21   for it = 1:numSymbs
22     if it <= (taps - 1) / 2;
23       xp = [zeros((taps - 1) / 2 - it + 1, 1); rSampledUnitEnergy(1:it + (taps - 1) / 2)];
24     elseif it + (taps - 1) / 2 > numSymbs
25       xp = [rSampledUnitEnergy(it - (taps - 1) / 2 : end); zeros(it + (taps - 1) / 2 - numSymbs, 1)];
26     else
27       xp = rSampledUnitEnergy(it - (taps - 1) / 2 : it + (taps - 1) / 2);
28     end
29
30     xout = sum(hxx .* xp);
31     ex = 1 - abs(xout) ^ 2;
32
33     if abs(ex) < 1e-3
34       convergeCount = convergeCount + 1;
35     else
36       convergeCount = 0;
37     end
38     if ~converged && convergeCount >= 10
39       converged = 1
40       it
41     end
42
43     adaptFilterOut(it) = xout;
44
45     hxx = hxx + mu * ex * xout * conj(xp);
46   end
47
48 %{
49   %% try MATLAB builtin equalizer
50   alg = cma(mu);
51   eqObj = lineareq(taps, alg);
52   eqObj.Weights((taps + 1) / 2) = 1;
53   rPadded = [rSampledUnitEnergy; zeros((taps - 1) / 2, 1)];
54   matlabEq = equalize(eqObj, rPadded);
55   matlabEq = matlabEq((taps + 1) / 2 : end);
56 %}
57 end