Merge pull request #63 from r9y9/lpclsp-conversions

Fix lpc2lsp conversion and add lsp2lpc function

docs/changelog.rst

@@ -1,7 +1,12 @@
 Change log
 ==========
 
-v0.1.11 <2018-xx-xx>
+v0.1.12 <2018-xx-xx>
+--------------------
+
+- #63`_: Fix lpc2lsp bug, add lsp2lpc function. Add regression tests for those.
+
+v0.1.11 <2018-02-05>
 --------------------
 
 - `#55`_: Add numpy implementation of cdist
@@ -101,3 +106,4 @@ v0.1.0 <2015-09-05>
 
 .. _#54: https://github.com/r9y9/pysptk/pull/54
 .. _#55: https://github.com/r9y9/pysptk/issues/55
+.. _#63: https://github.com/r9y9/pysptk/pull/63

pysptk/_sptk.pxd

@@ -64,6 +64,7 @@ cdef extern from "SPTK.h":
                            const int maxitr, const double eps)
     int _lpc2par "lpc2par"(double *a, double *k, const int m)
     void _par2lpc "par2lpc"(double *k, double *a, const int m)
+    void _lsp2lpc "lsp2lpc"(double *lsp, double *a, const int m)
     void _lsp2sp "lsp2sp"(double *lsp, const int m, double *x, const int l, const int gain)
 
 

pysptk/_sptk.pyx

@@ -339,30 +339,28 @@ def lpc2c(np.ndarray[np.float64_t, ndim=1, mode="c"] lpc not None,
     return ceps
 
 
-@cython.boundscheck(False)
-@cython.wraparound(False)
 def lpc2lsp(np.ndarray[np.float64_t, ndim=1, mode="c"] lpc not None,
-            numsp=512, maxiter=4, eps=1.0e-6, loggain=False, otype=0,
-            fs=None):
+            numsp=128, maxiter=4, eps=1.0e-6, has_gain=True,
+            loggain=False, otype=0, fs=None):
     cdef np.ndarray[np.float64_t, ndim = 1, mode = "c"] lsp
-    cdef int order = len(lpc) - 1
+    cdef int lpc_start_idx = 1 if has_gain else 0
+    cdef int order = len(lpc) - 1 if has_gain else len(lpc)
     lsp = np.zeros_like(lpc)
-    _lpc2lsp(&lpc[0], &lsp[0], order, numsp, maxiter, eps)
-
+    _lpc2lsp(&lpc[0], &lsp[lpc_start_idx], order, numsp, maxiter, eps)
     if otype == 0:
-        lsp[1:] *= 2 * np.pi
+        lsp[lpc_start_idx:] *= 2 * np.pi
     elif otype == 2 or otype == 3:
         if fs is None:
             raise ValueError("fs must be specified when otype == 2 or 3")
-        lsp[1:] *= fs
+        lsp[lpc_start_idx:] *= fs
 
     if otype == 3:
-        lsp[1:] *= 1000.0
+        lsp[lpc_start_idx:] *= 1000.0
 
-    if loggain:
+    if has_gain:
+      lsp[0] = lpc[0]
+      if loggain:
         lsp[0] = np.log(lpc[0])
-    else:
-        lsp[0] = lpc[0]
 
     return lsp
 
@@ -383,14 +381,66 @@ def par2lpc(np.ndarray[np.float64_t, ndim=1, mode="c"] par not None):
     return lpc
 
 
+def lsp2lpc(np.ndarray[np.float64_t, ndim=1, mode="c"] lsp not None,
+    has_gain=True, loggain=False, fs=None, itype=0):
+    lsp = lsp.copy()
+    cdef int lpc_start_idx = 1 if has_gain else 0
+
+    if loggain and not has_gain:
+        raise ValueError("has_gain must be True if you set loggain=True")
+
+    if itype == 0:
+        lsp[lpc_start_idx:] /= 2 * np.pi
+    elif itype == 2 or itype == 3:
+        if fs is None:
+            raise ValueError("fs must be specified when itype == 2 or 3")
+        lsp[lpc_start_idx:] /= fs
+
+    if itype == 3:
+        lsp[lpc_start_idx:] /= 1000.0
+
+    cdef np.ndarray[np.float64_t, ndim = 1, mode = "c"] lpc
+    cdef int order
+    if has_gain:
+      order = len(lsp) - 1
+    else:
+      order = len(lsp)
+    lpc = np.empty_like(lsp)
+
+    _lsp2lpc(&lsp[lpc_start_idx], &lpc[0], order)
+
+    if has_gain:
+      lpc[0] = lsp[0]
+      if loggain:
+        lpc[0] = np.exp(lpc[0])
+
+    return lpc
+
+
 def lsp2sp(np.ndarray[np.float64_t, ndim=1, mode="c"] lsp not None,
-           fftlen=256):
+           fftlen=256, has_gain=True, loggain=False, fs=None, itype=0):
     assert_fftlen(fftlen)
+    lsp = lsp.copy()
+    cdef int lsp_start_idx = 1 if has_gain else 0
+
+    if itype == 1:
+        lsp[lsp_start_idx:] *= 2 * np.pi
+    elif itype == 2 or itype == 3:
+        if fs is None:
+            raise ValueError("fs must be specified when itype == 2 or 3")
+        lsp[lsp_start_idx:] = lsp[lsp_start_idx:] / fs * 2 * np.pi
+
+    if itype == 3:
+        lsp[lsp_start_idx:] /= 1000.0
+
+    if loggain:
+        lsp[0] = np.log(lsp[0])
+
     cdef np.ndarray[np.float64_t, ndim = 1, mode = "c"] sp
     cdef int sp_length = (fftlen >> 1) + 1
     sp = np.empty(sp_length, dtype=np.float64)
     cdef int order = len(lsp) - 1
-    _lsp2sp(&lsp[0], order, &sp[0], sp_length, 1)
+    _lsp2sp(&lsp[0], order, &sp[0], sp_length, 1 if has_gain else 0)
     return sp
 
 

pysptk/sptk.py

@@ -48,6 +48,7 @@
     lpc2lsp
     lpc2par
     par2lpc
+    lsp2lpc
     lsp2sp
 
 Mel-generalized cepstrum conversions
@@ -968,8 +969,8 @@ def lpc2c(lpc, order=None):
 
 @apply_along_last_axis
 @automatic_type_conversion
-def lpc2lsp(lpc, numsp=512, maxiter=4, eps=1.0e-6, loggain=False, otype=0,
-            fs=None):
+def lpc2lsp(lpc, numsp=128, maxiter=4, eps=1.0e-6, has_gain=True,
+            loggain=False, otype=0, fs=None):
     """LPC to LSP
 
     Parameters
@@ -978,14 +979,17 @@ def lpc2lsp(lpc, numsp=512, maxiter=4, eps=1.0e-6, loggain=False, otype=0,
         LPC
 
     numsp : int, optional
-        Number of unit circle. Default is 512.
+        Number of unit circle. Default is 128.
 
     maxiter : int, optional
         Maximum number of iteration. Default is 4.
 
     eps : float, optional
         End condition for iteration. Default is 1.0e-6.
 
+    has_gain : bool, optional
+        Whether input LPC has gain at the index 0 or not. Default is True.
+
     loggain : bool, optional
         whether the converted lsp should have loggain or not.
         Default is False.
@@ -1019,7 +1023,7 @@ def lpc2lsp(lpc, numsp=512, maxiter=4, eps=1.0e-6, loggain=False, otype=0,
 
     """
 
-    return _sptk.lpc2lsp(lpc, numsp, maxiter, eps, loggain, otype, fs)
+    return _sptk.lpc2lsp(lpc, numsp, maxiter, eps, has_gain, loggain, otype, fs)
 
 
 @apply_along_last_axis
@@ -1075,7 +1079,56 @@ def par2lpc(par):
 
 @apply_along_last_axis
 @automatic_type_conversion
-def lsp2sp(lsp, fftlen=256):
+def lsp2lpc(lsp, has_gain=True, loggain=False, fs=None, itype=0):
+    """LSP to LPC
+
+
+    Parameters
+    ----------
+    lpc : array
+        LPC
+
+    has_gain : bool, optional
+        Whether input LPC has gain at the index 0 or not. Default is True.
+
+    loggain : bool, optional
+        If True, it's assumed that input LPC has loggain and convert it to linear
+        gain. Default is False.
+
+    fs : int, optional
+        Sampling frequency. Default is None and unused.
+
+    itype : int, optional
+        Input LPC format
+            (0)  normalized frequency (0 ~ pi)
+            (1)  normalized frequency (0 ~ 0.5)
+            (2)  frequency (kHz)
+            (3)  frequency (Hz)
+
+        Default is 0.
+
+    Returns
+    -------
+    lsp : array, shape (``order + 1``) if has_gain else (``order``)
+        LPC
+
+    raises
+    ------
+    ValueError
+        if ``fs`` is not specified when itype = 2 or 3.
+        if ``loggain`` and not ``has_gain``.
+
+    See Also
+    --------
+    pysptk.sptk.lpc2lsp
+    """
+
+    return _sptk.lsp2lpc(lsp, has_gain, loggain, fs, itype)
+
+
+@apply_along_last_axis
+@automatic_type_conversion
+def lsp2sp(lsp, fftlen=256, has_gain=True, loggain=False, fs=None, itype=0):
     """LSP to spectrum
 
     Parameters
@@ -1086,17 +1139,30 @@ def lsp2sp(lsp, fftlen=256):
     fftlen : int, optional
         FFT length
 
-    TODO: consider ``otype`` optional argument
+    has_gain : bool, optional
+        Whether input LPC has gain at the index 0 or not. Default is True.
+
+    loggain : bool, optional
+        If True, it's assumed that input LPC has loggain and convert it to linear
+        gain. Default is False.
+
+    fs : int, optional
+        Sampling frequency. Default is None and unused.
+
+    itype : int, optional
+        Input LPC format
+            (0)  normalized frequency (0 ~ pi)
+            (1)  normalized frequency (0 ~ 0.5)
+            (2)  frequency (kHz)
+            (3)  frequency (Hz)
+
+        Default is 0.
 
     Returns
     -------
     sp : array, shape
         Spectrum. ln|H(z)|.
 
-    Notes
-    -----
-    It is asuumed that ``lsp`` has loggain at ``lsp[0]``.
-
     See Also
     --------
     pysptk.sptk.lpc2par

tests/regression/test_mcep.py → tests/regression/test_mgcep.py

@@ -7,6 +7,37 @@
 DATA_DIR = join(dirname(__file__), "..", "data")
 
 
+def test_lpc():
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | dmp +f | awk \
+    # '{print $2}' > test16k_windowed.txt
+    frames = np.loadtxt(
+        join(DATA_DIR, "test16k_windowed.txt")).reshape(759, 512).astype(np.float64)
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | lpc -m 25 -l 512 > test16k.lpc
+    lpc = np.fromfile(join(DATA_DIR, "test16k.lpc"), np.float32).reshape(759, 26).astype(np.float64)
+    lpc_hat = pysptk.lpc(frames, order=25)
+    # yeah may have a bug...
+    assert np.allclose(lpc, lpc_hat, atol=1e-1)
+
+
+def test_lpc2lsp():
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | lpc -m 25 -l 512 > test16k.lpc
+    lpc = np.fromfile(join(DATA_DIR, "test16k.lpc"), np.float32).reshape(759, 26).astype(np.float64)
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | lpc -m 25 -l 512 | lpc2lsp -m 25 > test16k.lsp
+    lsp = np.fromfile(join(DATA_DIR, "test16k.lsp"), np.float32).reshape(759, 26).astype(np.float64)
+    lsp_hat = pysptk.lpc2lsp(lpc)
+    assert np.allclose(lsp, lsp_hat, atol=1e-4)
+
+
+def test_lsp2lpc():
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | lpc -m 25 -l 512 \
+    # | lpc2lsp -m 25 | lsp2lpc -m 25 > test16k.lsp2lpc
+    lpc = np.fromfile(join(DATA_DIR, "test16k.lsp2lpc"), np.float32).reshape(759, 26).astype(np.float64)
+    # frame -l 512 -p 80 < test16k.float | window -l 512 | lpc -m 25 -l 512 | lpc2lsp -m 25 > test16k.lsp
+    lsp = np.fromfile(join(DATA_DIR, "test16k.lsp"), np.float32).reshape(759, 26).astype(np.float64)
+    lpc_hat = pysptk.lsp2lpc(lsp)
+    assert np.allclose(lpc, lpc_hat, atol=1e-4)
+
+
 def test_mcep_from_windowed_frames():
     # frame -l 512 -p 80 < test16k.float | window -l 512 | mcep -q 0 -l 512 -a \
     # 0.41 -m 24 | dmp +f | awk '{print $2}' > test16k_mcep.txt