Merge branch 'soft_demod_volk'

This commit is contained in:
ismagom 2014-07-22 21:09:44 +02:00
commit 85fc5adc01
8 changed files with 366 additions and 64 deletions

View File

@ -27,6 +27,8 @@
#include <complex.h>
#include <stdint.h>
#include "liblte/config.h"
#ifndef CH_AWGN_
@ -37,13 +39,15 @@ typedef _Complex float cf_t;
LIBLTE_API void ch_awgn_c(const cf_t* input,
cf_t* output,
float variance,
int buff_sz);
uint32_t len);
LIBLTE_API void ch_awgn_f(const float* x,
float* y,
float variance,
int buff_sz);
uint32_t len);
LIBLTE_API float ch_awgn_get_variance(float ebno_db,
float rate);
/* High-level API */

View File

@ -35,6 +35,9 @@
typedef _Complex float cf_t;
#define EXPAVERAGE(data, average, nframes) ((data + average * nframes) / (nframes + 1))
/** Return the sum of all the elements */
LIBLTE_API int vec_acc_ii(int *x, uint32_t len);
LIBLTE_API float vec_acc_ff(float *x, uint32_t len);

View File

@ -29,25 +29,31 @@
#include <complex.h>
#include <stdlib.h>
#include <strings.h>
#include <math.h>
#include "gauss.h"
#include "liblte/phy/channel/ch_awgn.h"
void ch_awgn_c(const cf_t* x, cf_t* y, float variance, int buff_sz) {
_Complex float tmp;
int i;
float ch_awgn_get_variance(float ebno_db, float rate) {
float esno_db = ebno_db + 10 * log10f(rate);
return sqrtf(1 / (powf(10, esno_db / 10)));
}
for (i=0;i<buff_sz;i++) {
void ch_awgn_c(const cf_t* x, cf_t* y, float variance, uint32_t len) {
cf_t tmp;
uint32_t i;
for (i=0;i<len;i++) {
__real__ tmp = rand_gauss();
__imag__ tmp = rand_gauss();
tmp *= variance;
y[i] = tmp + x[i];
}
}
void ch_awgn_f(const float* x, float* y, float variance, int buff_sz) {
int i;
void ch_awgn_f(const float* x, float* y, float variance, uint32_t len) {
uint32_t i;
for (i=0;i<buff_sz;i++) {
for (i=0;i<len;i++) {
y[i] = x[i] + variance * rand_gauss();
}
}

View File

@ -36,6 +36,7 @@
void demod_soft_init(demod_soft_t *q) {
bzero((void*)q,sizeof(demod_soft_t));
q->sigma = 1.0;
}
void demod_soft_table_set(demod_soft_t *q, modem_table_t *table) {

View File

@ -34,6 +34,54 @@
#include "soft_algs.h"
#define LLR_APPROX_USE_VOLK
#ifdef LLR_APPROX_USE_VOLK
void
llr_approx(const _Complex float *in, float *out, int N, int M, int B,
_Complex float *symbols, uint32_t(*S)[6][32], float sigma2)
{
int i, s, b;
float num, den;
int change_sign = -1;
float x, y, d[64];
for (s = 0; s < N; s++) { /* recevied symbols */
/* Compute the distances squared d[i] between the received symbol and all constellation points */
for (i = 0; i < M; i++) {
x = __real__ in[s] - __real__ symbols[i];
y = __imag__ in[s] - __imag__ symbols[i];
d[i] = x * x + y * y;
}
for (b = 0; b < B; b++) { /* bits per symbol */
/* initiate num[b] and den[b] */
num = d[S[0][b][0]];
den = d[S[1][b][0]];
/* Minimum distance squared search between recevied symbol and a constellation point with a
'1' and a '0' for each bit position */
for (i = 1; i < M / 2; i++) { /* half the constellation points have '1'|'0' at any given bit position */
if (d[S[0][b][i]] < num) {
num = d[S[0][b][i]];
}
if (d[S[1][b][i]] < den) {
den = d[S[1][b][i]];
}
}
/* Theoretical LLR and approximate LLR values are positive if
* symbol(s) with '0' is/are closer and negative if symbol(s)
* with '1' are closer.
* Change sign if mapping negative to '0' and positive to '1' */
out[s * B + b] = change_sign * (den - num) / sigma2;
}
}
}
#else
/**
* @ingroup Soft Modulation Demapping based on the approximate
* log-likelihood algorithm
@ -50,46 +98,50 @@
* \param S Soft demapping auxiliary matrix
* \param sigma2 Noise vatiance
*/
void llr_approx(const _Complex float *in, float *out, int N, int M, int B,
_Complex float *symbols, uint32_t (*S)[6][32], float sigma2) {
void
llr_approx(const _Complex float *in, float *out, int N, int M, int B,
_Complex float *symbols, uint32_t(*S)[6][32], float sigma2)
{
int i, s, b;
float num, den;
int change_sign = -1;
float x, y, d[64];
for (s=0; s<N; s++) { /* recevied symbols */
/* Compute the distances squared d[i] between the received symbol and all constellation points */
for (i=0; i<M; i++) {
x = __real__ in[s] - __real__ symbols[i];
y = __imag__ in[s] - __imag__ symbols[i];
d[i] = x*x + y*y;
}
for (s = 0; s < N; s++) { /* recevied symbols */
/* Compute the distances squared d[i] between the received symbol and all constellation points */
for (i = 0; i < M; i++) {
x = __real__ in[s] - __real__ symbols[i];
y = __imag__ in[s] - __imag__ symbols[i];
d[i] = x * x + y * y;
}
for (b=0; b<B; b++) {/* bits per symbol*/
/* initiate num[b] and den[b] */
num = d[S[0][b][0]];
den = d[S[1][b][0]];
for (b = 0; b < B; b++) { /* bits per symbol */
/* initiate num[b] and den[b] */
num = d[S[0][b][0]];
den = d[S[1][b][0]];
/* Minimum distance squared search between recevied symbol and a constellation point with a
'1' and a '0' for each bit position */
for (i=1; i<M/2; i++) { /* half the constellation points have '1'|'0' at any given bit position */
if (d[S[0][b][i]] < num) {
num = d[S[0][b][i]];
}
if (d[S[1][b][i]] < den) {
den = d[S[1][b][i]];
}
}
/* Theoretical LLR and approximate LLR values are positive if
* symbol(s) with '0' is/are closer and negative if symbol(s)
* with '1' are closer.
* Change sign if mapping negative to '0' and positive to '1' */
out[s*B+b] = change_sign*(den-num)/sigma2;
}
/* Minimum distance squared search between recevied symbol and a constellation point with a
'1' and a '0' for each bit position */
for (i = 1; i < M / 2; i++) { /* half the constellation points have '1'|'0' at any given bit position */
if (d[S[0][b][i]] < num) {
num = d[S[0][b][i]];
}
if (d[S[1][b][i]] < den) {
den = d[S[1][b][i]];
}
}
/* Theoretical LLR and approximate LLR values are positive if
* symbol(s) with '0' is/are closer and negative if symbol(s)
* with '1' are closer.
* Change sign if mapping negative to '0' and positive to '1' */
out[s * B + b] = change_sign * (den - num) / sigma2;
}
}
}
#endif
/**
* @ingroup Soft Modulation Demapping based on the approximate
* log-likelihood ratio algorithm
@ -107,35 +159,36 @@ void llr_approx(const _Complex float *in, float *out, int N, int M, int B,
* \param sigma2 Noise vatiance
*/
void llr_exact(const _Complex float *in, float *out, int N, int M, int B,
_Complex float *symbols, uint32_t (*S)[6][32], float sigma2) {
_Complex float *symbols, uint32_t(*S)[6][32], float sigma2)
{
int i, s, b;
float num, den;
int change_sign = -1;
float x, y, d[64];
for (s=0; s<N; s++) { /* recevied symbols */
/* Compute exp{·} of the distances squared d[i] between the received symbol and all constellation points */
for (i=0; i<M; i++) {
x = __real__ in[s] - __real__ symbols[i];
y = __imag__ in[s] - __imag__ symbols[i];
d[i] = exp(-1*(x*x + y*y)/sigma2);
}
for (s = 0; s < N; s++) { /* recevied symbols */
/* Compute exp{·} of the distances squared d[i] between the received symbol and all constellation points */
for (i = 0; i < M; i++) {
x = __real__ in[s] - __real__ symbols[i];
y = __imag__ in[s] - __imag__ symbols[i];
d[i] = exp(-1 * (x * x + y * y) / sigma2);
}
/* Sum up the corresponding d[i]'s for each bit position */
for (b=0; b<B; b++) {/* bits per symbol*/
/* initiate num[b] and den[b] */
num = 0;
den = 0;
for (i=0; i<M/2; i++) { /* half the constellation points have '1'|'0' at any given bit position */
num += d[S[0][b][i]];
den += d[S[1][b][i]];
}
/* Theoretical LLR and approximate LLR values are positive if
* symbol(s) with '0' is/are closer and negative if symbol(s)
* with '1' are closer.
* Change sign if mapping negative to '0' and positive to '1' */
out[s*B+b] = change_sign*log(num/den);
}
/* Sum up the corresponding d[i]'s for each bit position */
for (b = 0; b < B; b++) { /* bits per symbol */
/* initiate num[b] and den[b] */
num = 0;
den = 0;
for (i = 0; i < M / 2; i++) { /* half the constellation points have '1'|'0' at any given bit position */
num += d[S[0][b][i]];
den += d[S[1][b][i]];
}
/* Theoretical LLR and approximate LLR values are positive if
* symbol(s) with '0' is/are closer and negative if symbol(s)
* with '1' are closer.
* Change sign if mapping negative to '0' and positive to '1' */
out[s * B + b] = change_sign * log(num / den);
}
}
}

View File

@ -36,5 +36,13 @@ ADD_TEST(modem_qpsk_soft modem_test -n 1020 -m 2 -s)
ADD_TEST(modem_qam16_soft modem_test -n 1020 -m 4 -s)
ADD_TEST(modem_qam64_soft modem_test -n 1020 -m 6 -s)
ADD_EXECUTABLE(soft_demod_test soft_demod_test.c)
TARGET_LINK_LIBRARIES(soft_demod_test lte_phy)
ADD_TEST(modem_bpsk_soft_approx soft_demod_test -n 1020 -m 1)
ADD_TEST(modem_qpsk_soft_approx soft_demod_test -n 1020 -m 2)
ADD_TEST(modem_qam16_soft_approx soft_demod_test -n 1020 -m 4)
ADD_TEST(modem_qam64_soft_approx soft_demod_test -n 1020 -m 6)

View File

@ -0,0 +1,229 @@
/**
*
* \section COPYRIGHT
*
* Copyright 2013-2014 The libLTE Developers. See the
* COPYRIGHT file at the top-level directory of this distribution.
*
* \section LICENSE
*
* This file is part of the libLTE library.
*
* libLTE is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as
* published by the Free Software Foundation, either version 3 of
* the License, or (at your option) any later version.
*
* libLTE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* A copy of the GNU Lesser General Public License can be found in
* the LICENSE file in the top-level directory of this distribution
* and at http://www.gnu.org/licenses/.
*
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <strings.h>
#include <unistd.h>
#include <math.h>
#include <time.h>
#include <stdbool.h>
#include "liblte/phy/phy.h"
int nof_frames = 10;
int num_bits = 1000;
lte_mod_t modulation = 0;
void usage(char *prog) {
printf("Usage: %s [nfv] -m modulation (1: BPSK, 2: QPSK, 4: QAM16, 6: QAM64)\n", prog);
printf("\t-n num_bits [Default %d]\n", num_bits);
printf("\t-f nof_frames [Default %d]\n", nof_frames);
printf("\t-v verbose [Default None]\n");
}
void parse_args(int argc, char **argv) {
int opt;
while ((opt = getopt(argc, argv, "nmvf")) != -1) {
switch (opt) {
case 'n':
num_bits = atoi(argv[optind]);
break;
case 'f':
nof_frames = atoi(argv[optind]);
break;
case 'v':
verbose++;
break;
case 'm':
switch(atoi(argv[optind])) {
case 1:
modulation = LTE_BPSK;
break;
case 2:
modulation = LTE_QPSK;
break;
case 4:
modulation = LTE_QAM16;
break;
case 6:
modulation = LTE_QAM64;
break;
default:
fprintf(stderr, "Invalid modulation %d. Possible values: "
"(1: BPSK, 2: QPSK, 4: QAM16, 6: QAM64)\n", atoi(argv[optind]));
break;
}
break;
default:
usage(argv[0]);
exit(-1);
}
}
if (modulation == 0) {
usage(argv[0]);
exit(-1);
}
}
float mse_threshold() {
switch(modulation) {
case LTE_BPSK:
return 1.0e-6;
case LTE_QPSK:
return 1.0e-6;
case LTE_QAM16:
return 0.11;
case LTE_QAM64:
return 0.18;
default:
return -1.0;
}
}
int main(int argc, char **argv) {
int i;
modem_table_t mod;
demod_soft_t demod_soft;
char *input, *output;
cf_t *symbols;
float *llr_exact, *llr_approx;
parse_args(argc, argv);
/* initialize objects */
if (modem_table_lte(&mod, modulation, true)) {
fprintf(stderr, "Error initializing modem table\n");
exit(-1);
}
/* check that num_bits is multiple of num_bits x symbol */
num_bits = mod.nbits_x_symbol * (num_bits / mod.nbits_x_symbol);
demod_soft_init(&demod_soft);
demod_soft_table_set(&demod_soft, &mod);
demod_soft_sigma_set(&demod_soft, 2.0 / mod.nbits_x_symbol);
/* allocate buffers */
input = malloc(sizeof(char) * num_bits);
if (!input) {
perror("malloc");
exit(-1);
}
output = malloc(sizeof(char) * num_bits);
if (!output) {
perror("malloc");
exit(-1);
}
symbols = malloc(sizeof(cf_t) * num_bits / mod.nbits_x_symbol);
if (!symbols) {
perror("malloc");
exit(-1);
}
llr_exact = malloc(sizeof(float) * num_bits);
if (!llr_exact) {
perror("malloc");
exit(-1);
}
llr_approx = malloc(sizeof(float) * num_bits);
if (!llr_approx) {
perror("malloc");
exit(-1);
}
/* generate random data */
srand(time(NULL));
int ret = -1;
double mse;
struct timeval t[3];
float mean_texec = 0.0;
for (int n=0;n<nof_frames;n++) {
for (i=0;i<num_bits;i++) {
input[i] = rand()%2;
}
/* modulate */
mod_modulate(&mod, input, symbols, num_bits);
/* add noise */
ch_awgn_c(symbols, symbols, ch_awgn_get_variance(5.0, mod.nbits_x_symbol), num_bits / mod.nbits_x_symbol);
/* Compare exact with approximation algorithms */
demod_soft_alg_set(&demod_soft, EXACT);
demod_soft_demodulate(&demod_soft, symbols, llr_exact, num_bits / mod.nbits_x_symbol);
demod_soft_alg_set(&demod_soft, APPROX);
gettimeofday(&t[1], NULL);
demod_soft_demodulate(&demod_soft, symbols, llr_approx, num_bits / mod.nbits_x_symbol);
gettimeofday(&t[2], NULL);
get_time_interval(t);
/* compute exponentially averaged execution time */
mean_texec = EXPAVERAGE((float) t[0].tv_usec, mean_texec, n);
/* check MSE */
mse = 0.0;
for (i=0;i<num_bits;i++) {
float e = llr_exact[i] - llr_approx[i];
mse += e*e;
}
mse/=num_bits;
if (VERBOSE_ISDEBUG()) {
printf("exact=");
vec_fprint_f(stdout, llr_exact, num_bits);
printf("approx=");
vec_fprint_f(stdout, llr_approx, num_bits);
}
if (mse > mse_threshold()) {
goto clean_exit;
}
}
ret = 0;
clean_exit:
free(llr_exact);
free(llr_approx);
free(symbols);
free(output);
free(input);
modem_table_free(&mod);
if (ret == 0) {
printf("Ok Mean Throughput: %.2f. Mbps ExTime: %.2f us\n", num_bits/mean_texec, mean_texec);
} else {
printf("Error: MSE too large (%f > %f)\n", mse, mse_threshold());
}
exit(ret);
}

View File

@ -40,8 +40,6 @@
#define MAX_TIME_OFFSET 128
cf_t dummy[MAX_TIME_OFFSET];
#define EXPAVERAGE(data, average, nframes) ((data + average * nframes) / (nframes + 1))
#define CURRENT_FFTSIZE lte_symbol_sz(q->cell.nof_prb)
#define CURRENT_SFLEN SF_LEN(CURRENT_FFTSIZE, q->cell.cp)