From 792f9c50f873eead98c343eaaadd7eb108ebf587 Mon Sep 17 00:00:00 2001
From: Jens Ahrensfeld <jens@jayfield.org>
Date: Sun, 19 Apr 2015 18:47:46 +0000
Subject: [PATCH] - come and go

git-svn-id: http://moon:8086/svn/matlab/trunk@33 801c6759-fa7c-4059-a304-17956f83a07c
---
 ofdm/cfo_coarse_fract_eval.m | 20 +++++++++
 ofdm/cfo_fine_eval.m         | 39 +++++++++++++++++
 ofdm/cfoml_eval.m            | 33 --------------
 ofdm/fdet_eval.m             | 34 ---------------
 ofdm/fdet_wav.m              | 83 ------------------------------------
 ofdm/timing_coarse.m         | 54 +++++++++++++++++++++++
 ofdm/timing_coarse_eval.m    | 38 +++++++++++++++++
 ofdm/timing_frame.m          | 23 ++++++++++
 8 files changed, 174 insertions(+), 150 deletions(-)
 create mode 100644 ofdm/cfo_coarse_fract_eval.m
 create mode 100644 ofdm/cfo_fine_eval.m
 delete mode 100644 ofdm/cfoml_eval.m
 delete mode 100644 ofdm/fdet_eval.m
 delete mode 100644 ofdm/fdet_wav.m
 create mode 100644 ofdm/timing_coarse.m
 create mode 100644 ofdm/timing_coarse_eval.m
 create mode 100644 ofdm/timing_frame.m

diff --git a/ofdm/cfo_coarse_fract_eval.m b/ofdm/cfo_coarse_fract_eval.m
new file mode 100644
index 0000000..9b4d281
--- /dev/null
+++ b/ofdm/cfo_coarse_fract_eval.m
@@ -0,0 +1,20 @@
+function cfo_coarse_fract_eval(x, cpoff)
+
+N = 1024;
+Ng = 256;
+L = Ng;
+
+for d=1:3,
+    k = 1;
+    Nd = N-2+d;
+    for nt=cpoff+1:N+Ng:length(x)-2*N-Ng,
+        x1 = x(nt+1:nt+L);
+        x2 = x(nt+Nd+1:nt+Nd+L);
+        c(k, d) = sum(x1 .* conj(x2));
+        k = k + 1;
+    end
+end
+domega_f = 1/(2*pi*N) * angle(c)
+
+figure(1);
+plot(0:length(c)-1, domega_f*48000); grid;
diff --git a/ofdm/cfo_fine_eval.m b/ofdm/cfo_fine_eval.m
new file mode 100644
index 0000000..b75b95b
--- /dev/null
+++ b/ofdm/cfo_fine_eval.m
@@ -0,0 +1,39 @@
+function cfo_fine_eval(f, kn)
+
+Kdown = 4;
+fs = 48000/Kdown;
+
+N = 1024/Kdown;
+K = 100;
+n = (0:K*N-1)';
+
+phi = pi/2;
+x = cos(f/fs*2*pi.*n + phi) + kn*randn(K*N, 1);
+
+Xlast = [];
+phi_x_ = [];
+phi_z_ = [];
+Zm = ones(N, 1);
+for k=1:K-1,
+    xp = x((k-1)*N+1:k*N);
+    X = fft(xp);
+    if (~isempty(Xlast))
+        Z = Xlast .* conj(X);
+        [v, bin] = max(abs(Z(1:N/2)));
+        phi_x_ = [phi_x_ 0.5*angle(X(bin))/pi/N*fs];
+        phi_z_ = [phi_z_ 0.5*angle(Z(bin))/pi/N*fs];
+    end
+    Xlast = X;
+    
+end
+df = -mean(phi_z_)
+f_est = (bin-1)/N*fs + df
+
+figure(1)
+plot(abs(Z)); grid;
+
+figure(2)
+plot(phi_x_); grid; ylabel('phi_x');
+
+figure(3)
+plot((bin-1)/N*fs - phi_z_); grid; ylabel('d_f');
diff --git a/ofdm/cfoml_eval.m b/ofdm/cfoml_eval.m
deleted file mode 100644
index 4756b88..0000000
--- a/ofdm/cfoml_eval.m
+++ /dev/null
@@ -1,33 +0,0 @@
-function cfoml_eval(x)
-
-N = 1024;       % number of subcarriers
-M = 103;        % number of subcarriers that carry information
-omega = 2*pi/N; % carrier spacing
-
-fa = 48000;
-fo = 0.0;
-
-x = 0.5*exp(i*2*pi.*(0:length(x)-1)'*fo/fa) .* x;
-
-U = exp(i*omega*(0:N-1)'*(0:N-1));
-
-W = U(:, 1:M);
-V = U(:, (M+1):N);
-
-phi = 0;
-dphi = pi/100;
-
-for j=1:100,
-phik(j) = phi;
-P = calcP(phi, N);
-phi = phi + dphi;
-
-S(j) = sum(abs(x'*P*U(:, M+1:N)).^2);
-end
-
-close all;
-
-plot(phik, real(S)); grid
-
-function P = calcP(phi, N)
-P = diag(exp(i*(0:N-1)*phi));
diff --git a/ofdm/fdet_eval.m b/ofdm/fdet_eval.m
deleted file mode 100644
index 37265e2..0000000
--- a/ofdm/fdet_eval.m
+++ /dev/null
@@ -1,34 +0,0 @@
-function fdet_eval(f_true, phi_true, knoise_dB, p)
-%
-% Example: fdet_eval(1200, pi/4, -10, [1100 1350 1])
-
-fs = 48000;
-N = 1024;
-K = 40;
-
-omega_min = p(1)/fs;
-omega_max = p(2)/fs;
-omega_step = p(3)/fs;
-
-% Generate signal
-x = exp(j*(2*pi*f_true/fs*(0:K*N-1)'+ phi_true)) + sqrt(10^(knoise_dB/10))*randn(K*N,1);
-
-[omega, err, dmod, omega_coarse, omega_est] = fdet(x, N, 0, omega_min, omega_step, omega_max);
-
-f_coarse = omega_coarse*fs
-f_est = omega_est*fs
-error_ppm = 1E6*(1-f_est/f_true)
-
-close all;
-subplot(4,1,1)
-plot(1:lge(omega), real(dmod)/N, 1:lge(omega), imag(dmod)/N, 1:lge(omega), abs(dmod)/N); grid; legend('Re', 'Im', 'abs()');
-
-subplot(4,1,2)
-plot(1:lge(omega), err*fs, '-'); grid; legend('error');
-
-subplot(4,1,3)
-plot(1:lge(omega), omega*fs); grid; legend('f');
-
-subplot(4,1,4)
-plot(1:lge(omega), angle(dmod), '-*', 1:lge(omega), diff(angle([0 dmod])), 'r-*'); grid; legend('phi', 'dphi');
-
diff --git a/ofdm/fdet_wav.m b/ofdm/fdet_wav.m
deleted file mode 100644
index 853c1be..0000000
--- a/ofdm/fdet_wav.m
+++ /dev/null
@@ -1,83 +0,0 @@
-function fdet_wav(filename, foff, p)
-%
-% Example: fdet_wav('drm-15435-if.wav', -174.65, [700 850 1])
-% Example: ofdm_tx(40, 1)
-%          fdet_wav('xv.wav', 0, [700 800 1])
-%       or fdet_wav('yv.wav', 0, [700 800 1])
-
-close all;
-fs = 48000;
-fc = fs/4;
-N = 1024;
-Ng = 256;
-
-omega_min = p(1)/fs;
-omega_max = p(2)/fs;
-omega_step = p(3)/fs;
-
-% Generate signal
-[xwav fs nbits] = wavread(filename);
-
-[len ndim] = size(xwav);
-
-if (ndim > 1)
-    fprintf('Mode: Baseband\n');
-    x = xwav(:, 1) + j*xwav(:, 2);
-else
-    fprintf('Mode: Passband\n');
-    x = xwav.*exp(-j*(2*pi*fc/fs*(0:len-1)'));
-end
-
-fprintf('Perform frequency correction\n');
-x = x.*exp(j*(2*pi*foff/fs*(0:len-1)'));
-
-% Remove cyclic prefix
-fprintf('Perform CP calculation\n');
-x0 = x(1:2*(N+Ng));
-x1 = x(N+Ng+1:3*(N+Ng));
-
-c = xcorr(x0, x1);
-[v, ii] = max(abs(c));
-cp_off = mod(ii, Ng+N)
-
-% Detect CFO
-fprintf('Perform frequency detection\n');
-[omega, err, dmod, omega_coarse, omega_est, Xm] = fdet(x, N, Ng, cp_off, omega_min, omega_step, omega_max);
-
-f_coarse = omega_coarse*fs
-f_est = omega_est*fs
-
-fprintf('Receive Symbols\n');
-K = fix((length(x)-Ng)/(N+Ng));
-for k=10:30
-    xp = x((k-1)*(N+Ng)+cp_off+1:k*(N+Ng)+cp_off);
-    Xp(:,k) = fft(xp(Ng+1:N+Ng), N);
-end;
-
-close all;
-
-plot(abs(c)); grid;
-figure;
-
-bin = fix(N*750/48000)+1
-plot(real(Xp(bin, :)), imag(Xp(bin, :)), '+'); grid;
-
-figure;
-plot(0:N-1, abs(Xm)); grid;
-
-figure;
-subplot(5,1,1)
-plot(1:lge(omega), real(dmod)/N, 1:lge(omega), imag(dmod)/N, 1:lge(omega), abs(dmod)/N); grid; legend('Re', 'Im', 'abs()');
-
-subplot(5,1,2)
-plot(1:lge(omega), err*fs, '-'); grid; legend('error');
-
-subplot(5,1,3)
-plot(1:lge(omega), omega*fs); grid; legend('f');
-
-subplot(5,1,4)
-plot(1:lge(omega), angle(dmod), '-', 1:lge(omega), diff(angle([0 dmod])), 'r-'); grid; legend('phi', 'dphi');
-
-subplot(5,1,5)
-plot(1:lge(omega), log10(abs(dmod) ./ (1E-6+abs(real(dmod)))), '-'); grid; legend('diff');
-
diff --git a/ofdm/timing_coarse.m b/ofdm/timing_coarse.m
new file mode 100644
index 0000000..07599b6
--- /dev/null
+++ b/ofdm/timing_coarse.m
@@ -0,0 +1,54 @@
+function [n, Nd, valid, Zf] = timing_coarse(x, N, Ng, L, Zi)
+
+if (~isfield(Zi, 'valid'))
+    Zi = struct('valid', 1, 'n', 0, 'xp', [], 'p', 1.0);
+end
+Zf = Zi;
+Zf.xp = [Zf.xp' x']';
+
+n = Zf.n;
+Nd = N;
+valid = 0;
+
+if (length(Zf.xp) < 3*N+Ng)
+    return
+end
+
+RCWIN_ROLLOFF = 0.125;
+[rcwin Nt wng] = txwin_eval(N, Ng, RCWIN_ROLLOFF);
+
+% RC Windowing at TX
+% Calculate merged CP'
+% -----------------------------------------------
+% | CP |    :    :    | CP |    :    :    | CP |     
+% ----------------------------------------------------> nt
+%        x0             x1             x2
+%
+% CP' = x0 * (1-w) + x2 * w
+
+k = 1;
+for nt=N+Ng:N+Ng+N-2,
+    for d=1:3,
+        Nd = N-2+d;
+        x0 = Zf.xp(nt-Nd+1:nt-Nd+L);
+        x1 = Zf.xp(nt+1:nt+L);
+        x2 = Zf.xp(nt+Nd+1:nt+Nd+L);
+        xcps = x0 .* (1-wng') + x2 .* (wng');
+        n = 0.5*(abs(x1).^2 + abs(xcps).^2);
+        c(k, d) = sum(abs(x1 .* conj(xcps)) - 1.0*n);
+    end
+    k = k + 1;
+end
+[v1, n1] = max(c);
+[v2, n2] = max(v1);
+Nd = N-2+n2;
+n = mod(n1(n2)-1, N+Ng);
+
+Zf.xp = Zf.xp(N+Ng+1:length(Zf.xp));
+Zf.n = n;
+valid = 1;
+if (0)
+    plot(d); grid;
+    pause(0.01);
+end
+
diff --git a/ofdm/timing_coarse_eval.m b/ofdm/timing_coarse_eval.m
new file mode 100644
index 0000000..314612c
--- /dev/null
+++ b/ofdm/timing_coarse_eval.m
@@ -0,0 +1,38 @@
+function timing_coarse_eval(x)
+
+N = 1024;
+Ng = 256;
+L = Ng;
+RCWIN_ROLLOFF = 0.125;
+
+[rcwin Nt wng] = txwin_eval(N, Ng, RCWIN_ROLLOFF);
+
+% RC Windowing at TX
+% Calculate merged CP'
+% -----------------------------------------------
+% | CP |    :    :    | CP |    :    :    | CP |     
+% ----------------------------------------------------> nt
+%        x0             x1             x2
+%
+% CP' = x0 * (1-w) + x2 * w
+
+k = 1;
+for nt=N+Ng:8*(N+Ng),
+    for d=1:3,
+        Nd = N-2+d;
+        x0 = x(nt-Nd+1:nt-Nd+L);
+        x1 = x(nt+1:nt+L);
+        x2 = x(nt+Nd+1:nt+Nd+L);
+        xcps = x0 .* (1-wng') + x2 .* (wng');
+        n = 0.5*(abs(x1).^2 + abs(xcps).^2);
+        c(k, d) = sum(abs(x1 .* conj(xcps)) - 1.0*n);
+    end
+    k = k + 1;
+end
+
+[v1, n1] = max(c);
+[v2, n2] = max(v1);
+n = mod(n1(n2)-1, N+Ng)
+
+figure(1);
+plot(0:length(c)-1, 1-abs(c)); grid;
diff --git a/ofdm/timing_frame.m b/ofdm/timing_frame.m
new file mode 100644
index 0000000..09fdb4c
--- /dev/null
+++ b/ofdm/timing_frame.m
@@ -0,0 +1,23 @@
+function [kc, valid, Zf] = frame_timing(Z, bins_timing_pilots, sym_per_frame, Zi)
+
+N = length(Z);
+
+if (~isfield(Zi, 'valid'))
+    Zi = struct('valid', 1, 'Zk', []);
+end
+
+if 1
+    valid = 0;
+    kc = 0;
+    if ~isempty(Zi.Zk)
+        Xk = Z .* conj(Zi.Zk);
+        if 0
+            plot(0:length(Xk)-1, abs(Xk), '-*'); grid;
+            pause (0.01); 
+        end
+        valid = 1;
+        kc = mean(abs(Xk));
+    end
+    Zi.Zk = Z;
+end
+Zf = Zi;