diff --git a/lte/phy/include/liblte/phy/modem/modem_table.h b/lte/phy/include/liblte/phy/modem/modem_table.h
index cbc4fbee9..794cd0bc4 100644
--- a/lte/phy/include/liblte/phy/modem/modem_table.h
+++ b/lte/phy/include/liblte/phy/modem/modem_table.h
@@ -40,12 +40,15 @@
 typedef _Complex float cf_t;
 typedef struct LIBLTE_API {
   uint32_t idx[2][6][32];
+  uint32_t min_idx[2][64][6];	/* NEW: for each constellation point zone (2, 4, 16, 64 for BPSK, QPSK, 16QAM, 64QAM) the 2x(1, 2, 4, and 6 closest constellation points) for each bit, respectively. */
+  uint32_t d_idx[64][7];	/* NEW: for each constellation point zone (2, 4, 16, 64 for BPSK, QPSK, 16QAM, 64QAM) the 2, 3, 5 and 7 indices to constellation points that need to be computed for any recevied symbol modulated as BPSK, QPSK, 16QAM, and 64QAM, respectively. */
+
 }soft_table_t;
 
 typedef struct LIBLTE_API {
-  cf_t* symbol_table;     // bit-to-symbol mapping
-  soft_table_t soft_table;   // symbol-to-bit mapping (used in soft demodulating)
-  uint32_t nsymbols;        // number of modulation symbols
+  cf_t* symbol_table;     	// bit-to-symbol mapping
+  soft_table_t soft_table;   	// symbol-to-bit mapping (used in soft demodulating)
+  uint32_t nsymbols;        	// number of modulation symbols
   uint32_t nbits_x_symbol;      // number of bits per symbol
 }modem_table_t;
 
diff --git a/lte/phy/lib/modem/src/demod_soft.c b/lte/phy/lib/modem/src/demod_soft.c
index d030a3733..aa9884e71 100644
--- a/lte/phy/lib/modem/src/demod_soft.c
+++ b/lte/phy/lib/modem/src/demod_soft.c
@@ -58,8 +58,10 @@ int demod_soft_demodulate(demod_soft_t *q, const cf_t* symbols, float* llr, int
         q->table->symbol_table, q->table->soft_table.idx, q->sigma);
     break;
   case APPROX:
-    llr_approx(symbols, llr, nsymbols, q->table->nsymbols, q->table->nbits_x_symbol,
+/*    llr_approx(symbols, llr, nsymbols, q->table->nsymbols, q->table->nbits_x_symbol,
         q->table->symbol_table, q->table->soft_table.idx, q->sigma);
+*/    llr_approx(symbols, llr, nsymbols, q->table->nsymbols, q->table->nbits_x_symbol,
+        q->table->symbol_table, q->table->soft_table.idx, q->table->soft_table.d_idx, q->table->soft_table.min_idx, q->sigma);
     break;
   }
   return nsymbols*q->table->nbits_x_symbol;
diff --git a/lte/phy/lib/modem/src/lte_tables.c b/lte/phy/lib/modem/src/lte_tables.c
index e571174c7..b16a6172e 100644
--- a/lte/phy/lib/modem/src/lte_tables.c
+++ b/lte/phy/lib/modem/src/lte_tables.c
@@ -34,6 +34,8 @@
 #include "liblte/phy/modem/modem_table.h"
 #include "lte_tables.h"
 
+void LLR_approx_params(const cf_t* table, soft_table_t *soft_table, int B);
+
 /**
  * Set the BPSK modulation table */
 void set_BPSKtable(cf_t* table, soft_table_t *soft_table, bool compute_soft_demod)
@@ -53,6 +55,18 @@ void set_BPSKtable(cf_t* table, soft_table_t *soft_table, bool compute_soft_demo
   /* BSPK symbols containing a '0' and a '1' (only two symbols, 1 bit) */
   soft_table->idx[0][0][0] = 0;
   soft_table->idx[1][0][0] = 1;
+
+  /* set two matrices for LLR approx. calculation */
+  soft_table->min_idx[0][0][0] = 0; 
+  soft_table->min_idx[0][1][0] = 0;
+  soft_table->min_idx[1][0][0] = 1;
+  soft_table->min_idx[1][1][0] = 1;
+
+  soft_table->d_idx[0][0] = 0; 
+  soft_table->d_idx[0][1] = 1;
+  soft_table->d_idx[1][0] = 0; 
+  soft_table->d_idx[1][1] = 1;
+
 }
 
 /**
@@ -91,6 +105,8 @@ void set_QPSKtable(cf_t* table, soft_table_t *soft_table, bool compute_soft_demo
   soft_table->idx[1][0][1] = 3;
   soft_table->idx[1][1][0] = 1;
   soft_table->idx[1][1][1] = 3;
+
+  LLR_approx_params(table, soft_table, 2);	/* last param indicating B (bits per symbol) */
 }
 
 /**
@@ -156,6 +172,8 @@ void set_16QAMtable(cf_t* table, soft_table_t *soft_table, bool compute_soft_dem
     soft_table->idx[0][3][i] = 2*i;
     soft_table->idx[1][3][i] = 2*i+1;
   }
+
+  LLR_approx_params(table, soft_table, 4);	/* last param indication B (bits per symbol) */
 }
 
 /**
@@ -277,4 +295,96 @@ void set_64QAMtable(cf_t* table, soft_table_t *soft_table, bool compute_soft_dem
     soft_table->idx[0][5][i] = 2*i;
     soft_table->idx[1][5][i] = 2*i+1;
   }
+
+  LLR_approx_params(table, soft_table, 6);	/* last param indication modulation */
+}
+
+/* Precompute two tables for calculating the distances based on the received symbol location relative to the constellation points */
+void LLR_approx_params(const cf_t* table, soft_table_t *soft_table, int B) {
+
+	int i, j, b, k;
+	float x, y, d0, d1, min_d0, min_d1;
+	int M, D;
+	uint32_t min_idx0[64][6], min_idx1[64][6];
+	uint32_t count;
+	int flag;
+
+
+	D = B+1;	/* number of different distances to be computed */
+	//M = pow(2,B);	/* number of constellation points */
+	switch (B) {
+		case 1: 	{M = 2; break;}		/* BPSK */
+		case 2: 	{M = 4; break;}		/* QPSK */
+		case 4: 	{M = 16; break;}	/* 16QAM */
+		case 6: 	{M = 64; break;}	/* 64QAM */
+		default:	{M = 4; break;}		/* QPSK */
+	}
+
+	for (i=0;i<M;i++) {	/* constellation points */
+	        for (b=0;b<B;b++) { /* bits per symbol */
+	        	min_d0 = 100;
+	        	min_d1 = 100;
+
+	        	for (j=0;j<M/2;j++) {	/* half the symbols have a '0', the other half a '1' at any bit position of modulation symbol */
+	        	    x = __real__ table[i] - __real__ table[soft_table->idx[0][b][j]];
+	        	    y = __imag__ table[i] - __imag__ table[soft_table->idx[0][b][j]];
+	        	    d0 = x*x + y*y;
+	        	    if (d0 < min_d0) {
+	        		    min_d0 = d0;
+	        		    min_idx0[i][b] = soft_table->idx[0][b][j];
+	        	    }
+
+	        	    x = __real__ table[i] - __real__ table[soft_table->idx[1][b][j]];
+	        	    y = __imag__ table[i] - __imag__ table[soft_table->idx[1][b][j]];
+	        	    d1 = x*x + y*y;
+	        	    if (d1 < min_d1) {
+	        		    min_d1 = d1;
+	        		    min_idx1[i][b] = soft_table->idx[1][b][j];
+	        	    }
+	        	}
+	        }
+	}
+
+	for (i=0;i<M;i++) {
+		for (j=0;j<D;j++) {
+			soft_table->d_idx[i][j] = -1;   /* intialization */
+		}
+	}
+
+	for (i=0;i<M;i++) {
+		count = 0;
+		for (b=0;b<B;b++) {	/* bit(b) = 0 */
+			flag = 0;
+			for (k=0;k<count;k++) {
+				if (min_idx0[i][b] == soft_table->d_idx[i][k]) {
+					soft_table->min_idx[0][i][b] = k;
+					flag = 1;   /* no new entry to idxdx */
+					break;
+				}
+			}
+
+			if (flag == 0) { /* new entry to min and d_idx */
+				soft_table->d_idx[i][count] = min_idx0[i][b];
+				soft_table->min_idx[0][i][b] = count;
+				count++;
+			}
+		}
+		for (b=0;b<B;b++) {	/* bit(b) = 1 */
+			flag = 0;
+			for (k=0;k<count;k++) {
+				if (min_idx1[i][b] == soft_table->d_idx[i][k]) {
+					soft_table->min_idx[1][i][b] = k;
+					flag = 1;   /* no new entry to d_idx */
+					break;
+				}
+			}
+
+			if (flag == 0) { /* new entry to min and d_idx */
+				soft_table->d_idx[i][count] = min_idx1[i][b];
+				soft_table->min_idx[1][i][b] = count;
+				count++;
+			}
+		}
+	}
 }
+
diff --git a/lte/phy/lib/modem/src/lte_tables.h b/lte/phy/lib/modem/src/lte_tables.h
index 166594fcc..3698945fc 100644
--- a/lte/phy/lib/modem/src/lte_tables.h
+++ b/lte/phy/lib/modem/src/lte_tables.h
@@ -38,6 +38,14 @@
 #define QAM64_LEVEL_3  5/sqrt(42)
 #define QAM64_LEVEL_4  7/sqrt(42)
 
+//////////////// NUEVO //////////////////////
+/* HARD DEMODULATION Thresholds, necessary for obtaining the zone of received symbol for optimized LLR approx implementation */
+#define QAM16_THRESHOLD		2/sqrt(10)
+#define QAM64_THRESHOLD_1	2/sqrt(42)
+#define QAM64_THRESHOLD_2	4/sqrt(42)
+#define QAM64_THRESHOLD_3	6/sqrt(42)
+//=========================================//
+
 #define QAM64_LEVEL_x  2/sqrt(42)
 /* this is not an QAM64 level, but, rather, an auxiliary value that can be used for computing the
  * symbol from the bit sequence */
diff --git a/lte/phy/lib/modem/src/soft_algs.c b/lte/phy/lib/modem/src/soft_algs.c
index a93d2b64a..7a37b19af 100644
--- a/lte/phy/lib/modem/src/soft_algs.c
+++ b/lte/phy/lib/modem/src/soft_algs.c
@@ -31,16 +31,468 @@
 #include <math.h>
 #include <complex.h>
 #include <stdint.h>
+#include <string.h>
 
 #include "soft_algs.h"
 #include "liblte/phy/utils/vector.h"
 
-#define LLR_APPROX_USE_VOLK
+#define QAM16_THRESHOLD		2/sqrt(10)
+#define QAM64_THRESHOLD_1	2/sqrt(42)
+#define QAM64_THRESHOLD_2	4/sqrt(42)
+#define QAM64_THRESHOLD_3	6/sqrt(42)
 
-#ifdef LLR_APPROX_USE_VOLK
 
 typedef _Complex float cf_t; 
 
+//#define LLR_APPROX_USE_VOLK
+#define LLR_APPROX_OPT
+
+#ifdef LLR_APPROX_OPT
+
+float dd[10000][7];	////// NUEVO: 7 distances that are needed to compute LLR approx for 64QAM 
+uint32_t zone[10000];	////// NUEVO: zone of received symbol with respect to grid of QAM constellation diagram
+
+
+
+/**
+ * @ingroup Received modulation symbol zone
+ * Determine location of received modulation symbol
+ *
+ * \param in input symbol (_Complex float)
+ * \param z associated zone in constellation diagram (int)
+ * \param N number of symbols
+ */
+static void zone_QPSK(const cf_t *in, uint32_t *z, int N) {
+
+	int s;
+	float re, im;
+
+	for (s=0;s<N;s++) {
+
+	    re = __real__ in[s];
+	    im = __imag__ in[s];
+
+	    if (re > 0) {
+		if (im > 0) { 	/* 1st Quadrand (upper-right) */
+			z[s] = 0;
+		} else {	/* 4th Quadrand (lower-right) */
+			z[s] = 1;
+		}
+	    } else {
+		if (im > 0) {	/* 2nd Quadrand (upper-left) */
+			z[s] = 2;
+		} else {	/* 3rd Quadrand (lower-left) */
+			z[s] = 3;
+		}
+	    }
+	}
+}
+
+/**
+ * @ingroup Received modulation symbol zone
+ * Determine location of received modulation symbol
+ *
+ * \param in input symbol (_Complex float)
+ * \param z associated zone in constellation diagram (int)
+ * \param N number of symbols
+ */
+static void zone_QAM16(const cf_t *in, uint32_t *z, int N) {
+
+	int s;
+	float re, im;
+
+	for (s=0;s<N;s++) {
+
+	    re = __real__ in[s];
+	    im = __imag__ in[s];
+
+	    if (re > 0) {
+		if (im > 0) { 	/* 1st Quadrand (upper-right) */
+			if (re > QAM16_THRESHOLD) {
+				if (im > QAM16_THRESHOLD) {
+					z[s] = 3;
+				} else {
+					z[s] = 2;
+				}
+			} else {
+				if (im > QAM16_THRESHOLD) {
+					z[s] = 1;
+				} else {
+					z[s] = 0;
+				}
+			}
+		} else {	/* 4th Quadrand (lower-right) */
+			if (re > QAM16_THRESHOLD) {
+				if (im < -QAM16_THRESHOLD) {
+					z[s] = 7;
+				} else {
+					z[s] = 6;
+				}
+			} else {
+				if (im < -QAM16_THRESHOLD) {
+					z[s] = 5;
+				} else {
+					z[s] = 4;
+				}
+			}
+		}
+	    } else {
+		if (im > 0) {	/* 2nd Quadrand (upper-left) */
+			if (re < -QAM16_THRESHOLD) {
+				if (im > QAM16_THRESHOLD) {
+					z[s] = 11;
+				} else {
+					z[s] = 10;
+				}
+			} else {
+				if (im > QAM16_THRESHOLD) {
+					z[s] = 9;
+				} else {
+					z[s] = 8;
+				}
+			}
+		} else {	/* 3rd Quadrand (lower-left) */
+			if (re < -QAM16_THRESHOLD) {
+				if (im < -QAM16_THRESHOLD) {
+					z[s] = 15;
+				} else {
+					z[s] = 14;
+				}
+			} else {
+				if (im < -QAM16_THRESHOLD) {
+					z[s] = 13;
+				} else {
+					z[s] = 12;
+				}
+			}
+		}
+	    }
+	}
+}
+
+/**
+ * @ingroup Received modulation symbol zone
+ * Determine location of received modulation symbol
+ *
+ * \param in input symbol (_Complex float)
+ * \param z associated zone in constellation diagram (int)
+ * \param N number of symbols
+ */
+
+static void zone_QAM64(const cf_t *in, uint32_t *z, int N) {
+
+	int s;
+	float re, im;
+
+	for (s=0;s<N;s++) {
+
+  	    re = __real__ in[s];
+	    im = __imag__ in[s];
+
+	    if (re > 0) {
+
+		if (im > 0) {
+
+			if (re > QAM64_THRESHOLD_2) {
+
+				if (re > QAM64_THRESHOLD_3) {
+
+					if (im > QAM64_THRESHOLD_2) {
+						if (im > QAM64_THRESHOLD_3) {
+							z[s] = 15;
+						} else {
+							z[s] = 14;
+						}
+					} else if (im > QAM64_THRESHOLD_1) {
+						z[s] = 10;
+					} else {
+						z[s] = 11;
+					}
+
+				} else {
+
+					if (im > QAM64_THRESHOLD_2) {
+						if (im > QAM64_THRESHOLD_3) {
+							z[s] = 13;
+						} else {
+							z[s] = 12;
+						}
+					} else if (im > QAM64_THRESHOLD_1) {
+						z[s] = 8;
+					} else {
+						z[s] = 9;
+					}
+				}
+
+			} else if (re > QAM64_THRESHOLD_1) {
+
+				if (im > QAM64_THRESHOLD_2) {
+					if (im > QAM64_THRESHOLD_3) {
+						z[s] = 5;
+					} else {
+						z[s] = 4;
+					}
+				} else if (im > QAM64_THRESHOLD_1) {
+					z[s] = 0;
+				} else {
+					z[s] = 1;
+				}
+
+			} else {
+				if (im > QAM64_THRESHOLD_2) {
+					if (im > QAM64_THRESHOLD_3) {
+						z[s] = 7;
+					} else {
+						z[s] = 6;
+					}
+				} else if (im > QAM64_THRESHOLD_1) {
+					z[s] = 2;
+				} else {
+					z[s] = 3;
+				}
+			}
+
+		} else { /* forth quadrant (lower-right) */
+
+			if (re > QAM64_THRESHOLD_2) {
+
+				if (re > QAM64_THRESHOLD_3) {
+
+					if (im < -QAM64_THRESHOLD_2) {
+						if (im < -QAM64_THRESHOLD_3) {
+							z[s] = 31;
+						} else {
+							z[s] = 30;
+						}
+					} else if (im < -QAM64_THRESHOLD_1) {
+						z[s] = 26;
+					} else {
+						z[s] = 27;
+					}
+
+				} else {
+
+					if (im < -QAM64_THRESHOLD_2) {
+						if (im < -QAM64_THRESHOLD_3) {
+							z[s] = 29;
+						} else {
+							z[s] = 28;
+						}
+					} else if (im < -QAM64_THRESHOLD_1) {
+						z[s] = 24;
+					} else {
+						z[s] = 25;
+					}
+				}
+
+			} else if (re > QAM64_THRESHOLD_1) {
+
+				if (im < -QAM64_THRESHOLD_2) {
+					if (im < -QAM64_THRESHOLD_3) {
+						z[s] = 21;
+					} else {
+						z[s] = 20;
+					}
+				} else if (im < -QAM64_THRESHOLD_1) {
+					z[s] = 16;
+				} else {
+					z[s] = 17;
+				}
+
+			} else {
+				if (im < -QAM64_THRESHOLD_2) {
+					if (im < -QAM64_THRESHOLD_3) {
+						z[s] = 23;
+					} else {
+						z[s] = 22;
+					}
+				} else if (im < -QAM64_THRESHOLD_1) {
+					z[s] = 18;
+				} else {
+					z[s] = 19;
+				}
+			}
+		}
+
+	    } else { /* re < 0 */
+
+		if (im > 0) { /* second quadrant (upper-left) */
+
+			if (re < -QAM64_THRESHOLD_2) {
+
+				if (re < -QAM64_THRESHOLD_3) {
+
+					if (im > QAM64_THRESHOLD_2) {
+						if (im > QAM64_THRESHOLD_3) {
+							z[s] = 47;
+						} else {
+							z[s] = 46;
+						}
+					} else if (im > QAM64_THRESHOLD_1) {
+						z[s] = 42;
+					} else {
+						z[s] = 43;
+					}
+
+				} else {
+
+					if (im > QAM64_THRESHOLD_2) {
+						if (im > QAM64_THRESHOLD_3) {
+							z[s] = 45;
+						} else {
+							z[s] = 44;
+						}
+					} else if (im > QAM64_THRESHOLD_1) {
+						z[s] = 40;
+					} else {
+						z[s] = 41;
+					}
+				}
+
+			} else if (re < -QAM64_THRESHOLD_1) {
+
+				if (im > QAM64_THRESHOLD_2) {
+					if (im > QAM64_THRESHOLD_3) {
+						z[s] = 37;
+					} else {
+						z[s] = 36;
+					}
+				} else if (im > QAM64_THRESHOLD_1) {
+					z[s] = 32;
+				} else {
+					z[s] = 33;
+				}
+
+			} else {
+				if (im > QAM64_THRESHOLD_2) {
+					if (im > QAM64_THRESHOLD_3) {
+						z[s] = 39;
+					} else {
+						z[s] = 38;
+					}
+				} else if (im > QAM64_THRESHOLD_1) {
+					z[s] = 34;
+				} else {
+					z[s] = 35;
+				}
+			}
+		} else { /* third quadrant (lower-left) */
+			if (re < -QAM64_THRESHOLD_2) {
+
+				if (re < -QAM64_THRESHOLD_3) {
+
+					if (im < -QAM64_THRESHOLD_2) {
+						if (im < -QAM64_THRESHOLD_3) {
+							z[s] = 63;
+						} else {
+							z[s] = 62;
+						}
+					} else if (im < -QAM64_THRESHOLD_1) {
+						z[s] = 58;
+					} else {
+						z[s] = 59;
+					}
+
+				} else {
+
+					if (im < -QAM64_THRESHOLD_2) {
+						if (im < -QAM64_THRESHOLD_3) {
+							z[s] = 61;
+						} else {
+							z[s] = 60;
+						}
+					} else if (im < -QAM64_THRESHOLD_1) {
+						z[s] = 56;
+					} else {
+						z[s] = 57;
+					}
+				}
+
+			} else if (re < -QAM64_THRESHOLD_1) {
+
+				if (im < -QAM64_THRESHOLD_2) {
+					if (im < -QAM64_THRESHOLD_3) {
+						z[s] = 53;
+					} else {
+						z[s] = 52;
+					}
+				} else if (im < -QAM64_THRESHOLD_1) {
+					z[s] = 48;
+				} else {
+					z[s] = 49;
+				}
+
+			} else {
+				if (im < -QAM64_THRESHOLD_2) {
+					if (im < -QAM64_THRESHOLD_3) {
+						z[s] = 55;
+					} else {
+						z[s] = 54;
+					}
+				} else if (im < -QAM64_THRESHOLD_1) {
+					z[s] = 50;
+				} else {
+					z[s] = 51;
+				}
+			}
+
+		}
+	    }
+	}
+}
+
+static void compute_zone(const cf_t *in, uint32_t *z, int N, int B)
+{
+    switch(B) {
+	case 1:	{memset(zone, 0, N*sizeof(int)); break;}/* BPSK */
+	case 2: {zone_QPSK(in, z, N); break;}		/* QPSK */
+	case 4: {zone_QAM16(in, z, N); break;}		/* 16QAM */
+	case 6: {zone_QAM64(in, z, N); break;}		/* 64QAM */
+    }
+}
+
+static void compute_square_dist(const cf_t *in, cf_t *symbols, uint32_t (*idx)[7], int N, int B) {
+  int s, b;
+  float *d_ptr;
+	float x, y;
+  cf_t symbols_extract[7];
+
+  for (s=0;s<N;s++) {	/* N: number of received symbols */
+    d_ptr = dd[s];
+    for (b=0;b<B+1;b++) {
+	symbols_extract[b] = symbols[idx[zone[s]][b]];	/* only subset of distances to constellation points needed for LLR approx */
+	//x = __real__ in[s] - __real__ symbols[idx[zone[s]][b]];
+	//y = __imag__ in[s] - __imag__ symbols[idx[zone[s]][b]];
+	//dd[s][b] = x*x + y*y;
+	//printf("\n%f + j %f", __real__ symbols_extract[b], __imag__ symbols_extract[b]);
+    }
+    vec_square_dist(in[s], symbols_extract, d_ptr, B+1);/* B+1 distances to be computed */
+  }
+}
+
+static void compute_llr(int N, int B, uint32_t (*min)[64][6], float sigma2, float *out) {
+  int s, b;
+  for (s = 0; s < N; s++) {     
+    for (b = 0; b < B; b++) {   /* bits per symbol */
+      out[s * B + b] = (dd[s][min[0][zone[s]][b]] - dd[s][min[1][zone[s]][b]]) / sigma2;
+    }
+  }
+}
+
+void llr_approx(const _Complex float *in, float *out, int N, int M, int B,
+           _Complex float *symbols, uint32_t(*S)[6][32], uint32_t (*idx)[7], uint32_t (*min)[64][6], float sigma2)
+{
+  if ((M == 2) || (M == 4) || (M == 16) || (M == 64)) {
+	compute_zone(in, zone, N, B);
+    	compute_square_dist(in, symbols, idx, N, B);
+    	compute_llr(N, B, min, sigma2, out);
+  }
+}
+
+#endif 
+
+#ifdef LLR_APPROX_USE_VOLK
+
 float d[10000][64];
 float num[10000], den[10000];
 
@@ -76,14 +528,14 @@ static void compute_min_dist(uint32_t (*S)[6][32], int N, int B, int M) {
 static void compute_llr(int N, int B, float sigma2, float *out) {
   int s, b;
   for (s = 0; s < N; s++) {     
-    for (b = 0; b < B; b++) {   /* bits per symbol */
+    for (b = 0; b < B; b++) {  /* bits per symbol */
       out[s * B + b] = (num[s * B + b]-den[s * B + b]) / 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)
+           _Complex float *symbols, uint32_t(*S)[6][32], uint32_t (*idx)[7], uint32_t (*min)[64][6], float sigma2)
 {
   
   if (M <= 64) {
@@ -95,7 +547,9 @@ void llr_approx(const _Complex float *in, float *out, int N, int M, int B,
   }
 }
 
-#else 
+#endif
+
+#ifdef LLR_APPROX_OLD
 
 /**
  * @ingroup Soft Modulation Demapping based on the approximate
@@ -113,9 +567,8 @@ void llr_approx(const _Complex float *in, float *out, int N, int M, int B,
  * \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], uint32_t (*idx)[7], uint32_t (*min)[64][6], float sigma2)
 {
   int i, s, b;
   float num, den;
@@ -150,9 +603,9 @@ llr_approx(const _Complex float *in, float *out, int N, int M, int B,
        * with '1' are closer.
        * Change sign if mapping negative to '0' and positive to '1' */
       out[s * B + b] = change_sign * (den - num) / sigma2;
-      if (s<10)
+/*      if (s<10)
         printf("out[%d]=%f=%f/%f\n",s*B+b,out[s*B+b], num,den);
-    }
+  */  }
   }
 
 }
diff --git a/lte/phy/lib/modem/src/soft_algs.h b/lte/phy/lib/modem/src/soft_algs.h
index c683700a6..91a5e93aa 100644
--- a/lte/phy/lib/modem/src/soft_algs.h
+++ b/lte/phy/lib/modem/src/soft_algs.h
@@ -26,7 +26,7 @@
  */
 
 
-void llr_approx(const _Complex float *in, 
+/*void llr_approx(const _Complex float *in, 
                 float *out, 
                 int N, 
                 int M, 
@@ -34,6 +34,17 @@ void llr_approx(const _Complex float *in,
                 _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], 
+                uint32_t (*idx)[7],	/*64x7 table of integers [0..63], indices to 7 distances to be computed */
+		uint32_t (*min)[64][6],	/*2x64x6 table of integers [0..6], indices to 2x6 nearest symbols */
+                float sigma2);
 
 void llr_exact(const _Complex float *in, 
                float *out, 
@@ -43,3 +54,4 @@ void llr_exact(const _Complex float *in,
                _Complex float *symbols, 
                uint32_t (*S)[6][32], 
                float sigma2);
+
diff --git a/lte/phy/lib/modem/test/modem_test.c b/lte/phy/lib/modem/test/modem_test.c
index 0fc09a642..c173edfed 100644
--- a/lte/phy/lib/modem/test/modem_test.c
+++ b/lte/phy/lib/modem/test/modem_test.c
@@ -36,6 +36,9 @@
 
 #include "liblte/phy/phy.h"
 
+time_t start, finish;
+struct timeval x, y;
+
 int num_bits = 1000;
 lte_mod_t modulation;
 bool soft_output = false, soft_exact = false;
@@ -98,6 +101,10 @@ int main(int argc, char **argv) {
   cf_t *symbols;
   float *llr;
 
+//  unsigned long strt, fin;
+//  strt = x->tv_usec;
+//  fin = y->tv_usec;
+  
   parse_args(argc, argv);
 
   /* initialize objects */
@@ -156,7 +163,12 @@ int main(int argc, char **argv) {
 
   /* demodulate */
   if (soft_output) {
+
+    gettimeofday(&x, NULL);
     demod_soft_demodulate(&demod_soft, symbols, llr, num_bits / mod.nbits_x_symbol);
+    gettimeofday(&y, NULL);
+    printf("\nElapsed time [ns]: %u\n", y.tv_usec - x.tv_usec);
+    
     for (i=0;i<num_bits;i++) {
       output[i] = llr[i]>=0 ? 1 : 0;
     }
diff --git a/lte/phy/lib/modem/test/soft_demod_test.c b/lte/phy/lib/modem/test/soft_demod_test.c
index 12389bfd9..96b76c5b0 100644
--- a/lte/phy/lib/modem/test/soft_demod_test.c
+++ b/lte/phy/lib/modem/test/soft_demod_test.c
@@ -129,6 +129,7 @@ int main(int argc, char **argv) {
   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) {