Merge branch 'release/v3.1.1'

Conflicts:
	setup.py

prox_tv/_prox_tv.py

@@ -183,6 +183,43 @@ def _call(fn, *args):
     return fn(*args_m)
 
 
+def force_float_scalar(x):
+    r"""Forces an scalar value into float format
+    
+    Parameters
+    ----------
+    x: scalar value to check
+    
+    Returns
+    -------
+    float
+        Float representation of the provided value
+    """    
+    if not isinstance(x, float):
+        return float(x)
+    else:
+        return x
+
+
+def force_float_matrix(x):
+    r"""Forces a numpy matrix into float format
+    
+    Parameters
+    ----------
+    x: numpy array
+        matrix to check
+    
+    Returns
+    -------
+    numpy array
+        Float representation of the provided matrix
+    """    
+    if x.dtype != np.dtype('float64'):
+        return x.astype('float')
+    else:
+        return x
+
+
 def tv1_1d(x, w, sigma=0.05, method='tautstring'):
     r"""1D proximal operator for :math:`\ell_1`.
 
@@ -194,7 +231,7 @@ def tv1_1d(x, w, sigma=0.05, method='tautstring'):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The signal we are approximating.
     w : float
         The non-negative weight in the optimization problem.
@@ -216,6 +253,8 @@ def tv1_1d(x, w, sigma=0.05, method='tautstring'):
     """
     assert method in ('tautstring', 'pn', 'condat', 'dp')
     assert w >= 0
+    w = force_float_scalar(w)
+    x = force_float_matrix(x)
     y = np.zeros(np.size(x))
     if method == 'tautstring':
         _call(_lib.tautString_TV1, x, w, y, np.size(x))
@@ -240,7 +279,7 @@ def tv1w_1d(x, w, method='tautstring', sigma=0.05):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The signal we are approximating.
     w : numpy array
         The non-negative weights in the optimization problem.
@@ -256,6 +295,8 @@ def tv1w_1d(x, w, method='tautstring', sigma=0.05):
     """
     assert np.all(w >= 0)
     assert np.size(x)-1 == np.size(w)
+    w = force_float_matrix(w)
+    x = force_float_matrix(x)
     y = np.zeros(np.size(x))
     if method == 'tautstring':
         _call(_lib.tautString_TV1_Weighted, x, w, y, np.size(x))
@@ -277,7 +318,7 @@ def tv2_1d(x, w, method='mspg'):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The signal we are approximating.
     w : float
         The non-negative weight in the optimization problem.
@@ -295,6 +336,8 @@ def tv2_1d(x, w, method='mspg'):
     """
     assert w >= 0
     assert method in ('ms', 'pg', 'mspg')
+    w = force_float_scalar(w)
+    x = force_float_matrix(x)
     info = np.zeros(_N_INFO)
     y = np.zeros(np.size(x), order='F')
     if method == 'ms':
@@ -318,7 +361,7 @@ def tvp_1d(x, w, p, method='gpfw', max_iters=0):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The signal we are approximating.
     w : float
         The non-negative weight in the optimization problem.
@@ -342,6 +385,9 @@ def tvp_1d(x, w, p, method='gpfw', max_iters=0):
     assert method in methods
     assert w >= 0
     assert p >= 1
+    w = force_float_scalar(w)
+    p = force_float_scalar(p)
+    x = force_float_matrix(x)
     info = np.zeros(_N_INFO)
     y = np.zeros(np.size(x), order='F')
     _call(methods[method], x, w, y, info, np.size(x), p, _ffi.NULL)
@@ -360,7 +406,7 @@ def tv1_2d(x, w, n_threads=1, max_iters=0, method='dr'):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The signal we are approximating.
     w : float
         The non-negative weight in the optimization problem.
@@ -384,7 +430,8 @@ def tv1_2d(x, w, n_threads=1, max_iters=0, method='dr'):
     """
     assert w >= 0
     assert method in ('dr', 'pd', 'yang', 'condat', 'chambolle-pock')
-    x = np.asfortranarray(x)
+    x = np.asfortranarray(x, dtype='float64')
+    w = force_float_scalar(w)
     y = np.asfortranarray(np.zeros(x.shape))
     info = np.zeros(_N_INFO)
     if method == 'yang':
@@ -414,7 +461,7 @@ def tv1w_2d(x, w_col, w_row, max_iters=0, n_threads=1):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The MxN matrix we are approximating.
     w_col : numpy array
         The (M-1) x N matrix of column weights :math:`w^c`.
@@ -431,10 +478,10 @@ def tv1w_2d(x, w_col, w_row, max_iters=0, n_threads=1):
     M, N = x.shape
     assert w_col.shape == (M-1, N)
     assert w_row.shape == (M, N-1)
-    x = np.asfortranarray(x)
+    x = np.asfortranarray(x, dtype='float64')
     y = np.zeros(x.shape, order='F')
-    w_col = np.asfortranarray(w_col)
-    w_row = np.asfortranarray(w_row)
+    w_col = np.asfortranarray(w_col, dtype='float64')
+    w_row = np.asfortranarray(w_row, dtype='float64')
     info = np.zeros(_N_INFO)
     _call(_lib.DR2L1W_TV, M, N, x, w_col, w_row, y, n_threads, max_iters, info)
     return y
@@ -477,7 +524,11 @@ def tvp_2d(x, w_col, w_row, p_col, p_row, n_threads=1, max_iters=0):
     assert p_row >= 1
 
     info = np.zeros(_N_INFO)
-    x = np.asfortranarray(x)
+    x = np.asfortranarray(x, dtype='float64')
+    w_col = force_float_scalar(w_col)
+    w_row = force_float_scalar(w_row)
+    p_col = force_float_scalar(p_col)
+    p_row = force_float_scalar(p_row)
     y = np.zeros(np.shape(x), order='F')
     _call(_lib.DR2_TV,
           x.shape[0], x.shape[1], x, w_col, w_row, p_col, p_row, y,
@@ -500,7 +551,7 @@ def tvgen(x, ws, ds, ps, n_threads=1, max_iters=0):
 
     Parameters
     ----------
-    y : numpy array
+    x : numpy array
         The matrix signal we are approximating.
     ws : list
         Weights to apply in each penalty term.
@@ -526,7 +577,9 @@ def tvgen(x, ws, ds, ps, n_threads=1, max_iters=0):
     assert max_iters >= 0
 
     info = np.zeros(_N_INFO)
-    x = np.asfortranarray(x)
+    x = np.asfortranarray(x, dtype='float64')
+    ws = force_float_matrix(ws)
+    ps = force_float_matrix(ps)
     y = np.zeros(np.shape(x), order='F')
 
     # Run algorithm depending on the structure of the data and the requested penalties

prox_tv/prox_tv_test.py

@@ -13,6 +13,28 @@ def test_tv1w_1d():
         solutions = [tv1w_1d(x, w, method=method) for method in methods]
         for i in range(len(solutions)-1):
             assert np.allclose(solutions[i], solutions[i+1])
+
+
+def test_tv1w_1d_uniform_weights():
+    for _ in range(20):
+        dimension = np.random.randint(1e1, 1e3)
+        x = 100*np.random.randn(dimension)
+        w1 = np.random.rand()
+        w = np.ones(dimension-1) * w1 
+        solw = tv1w_1d(x, w)
+        sol = tv1_1d(x, w1)
+        assert np.allclose(solw, sol)
+
+
+def test_tv1w_1d_uniform_weights_small_input():
+    for _ in range(1000):
+        dimension = np.random.randint(2, 4)
+        x = 100*np.random.randn(dimension)
+        w1 = np.random.rand()
+        w = np.ones(dimension-1) * w1 
+        solw = tv1w_1d(x, w)
+        sol = tv1_1d(x, w1)
+        assert np.allclose(solw, sol)
 
 
 def test_tv1_1d():
@@ -65,3 +87,42 @@ def test_tv1_2d():
         solutions.append(tv1w_2d(x, w_cols, w_rows, max_iters=5000))
         for i in range(1, len(solutions)):
             assert np.allclose(solutions[i], solutions[0], atol=1e-3)
+
+
+def test_tv1w_2d_uniform_weights():
+    for _ in range(20):
+        rows = np.random.randint(1e1, 3e1)
+        cols = np.random.randint(1e1, 3e1)
+        x = 100*np.random.randn(rows, cols)
+        w1 = np.random.rand()
+        w_rows = np.ones([rows-1, cols]) * w1 
+        w_cols = np.ones([rows, cols-1]) * w1 
+        solw = tv1w_2d(x, w_rows, w_cols, max_iters=5000)
+        solw1 = tv1_2d(x, w1, max_iters=5000)
+        assert np.allclose(solw, solw1, atol=1e-3)
+
+
+def test_tv1w_2d_uniform_weights_small_input():
+    for _ in range(1000):
+        rows = np.random.randint(2, 4)
+        cols = np.random.randint(2, 4)
+        x = 100*np.random.randn(rows, cols)
+        w1 = np.random.rand()
+        w_rows = np.ones([rows-1, cols]) * w1 
+        w_cols = np.ones([rows, cols-1]) * w1 
+        solw = tv1w_2d(x, w_rows, w_cols, max_iters=5000)
+        solw1 = tv1_2d(x, w1, max_iters=5000)
+        assert np.allclose(solw, solw1, atol=1e-3)
+
+
+def test_tv1w_2d_emengd():
+    r"""Issue reported by emengd
+    
+    Make the solver fail due to missing checks on integer arguments
+    """        
+    a = -np.array([[1,2,3],[4,5,6],[7,8,9]])/10.
+    sol1 = tv1w_2d(a, np.array([[1,1,1],[1,1,1]]), 
+                   np.array([[1,1],[1,1],[1,1]]), max_iters=100)
+    sol2 = tv1_2d(a, 1)
+    assert np.allclose(sol1, sol2, atol=1e-3)
+

setup.py

@@ -27,7 +27,7 @@ def finalize_options(self):
 
 setup(
     name="prox_tv",
-    version="3.1.0",
+    version="3.1.1",
     description="Toolbox for fast Total Variation proximity operators",
     long_description="proxTV is a toolbox implementing blazing fast implementations of Total Variation proximity operators. While the core algorithms are implemented in C to achieve high efficiency, Matlab and Python interfaces are provided for ease of use. The library provides efficient solvers for a variety of Total Variation proximity problems, with address input signals of any dimensionality (1d, images, video, ...) and different norms to apply in the Total Variation term.",
     packages=['prox_tv'],
@@ -50,7 +50,7 @@ def finalize_options(self):
     package_data={
         'prox_tv': ['src/*.h', 'src/demos/*']
     },
-    author="Alvaro Barbero, Suvrit Sra",
+    author="Alvaro Barbero, Suvrit Sra, Josip Djolonga (python bindings)",
     author_email="alvaro.barbero@uam.es",
     url='https://github.com/albarji/proxTV',
     license='BSD',