From d7432adfe165ab9919ee976f94704f4a21756fe0 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 11:13:17 +0200
Subject: [PATCH 1/9] Add 2D Undecimated Wavelet class and test

---
 mri/numerics/linear.py           |   1 +
 mri/reconstruct/linear.py        | 160 +++++++++++++++++++++++++++++--
 mri/test/test_wavelet_adjoint.py |  20 +++-
 3 files changed, 174 insertions(+), 7 deletions(-)

diff --git a/mri/numerics/linear.py b/mri/numerics/linear.py
index 076b31af..02275f65 100644
--- a/mri/numerics/linear.py
+++ b/mri/numerics/linear.py
@@ -13,3 +13,4 @@
 
 # Package import
 from mri.reconstruct.linear import WaveletN
+from mri.reconstruct.linear import WaveletUD2
\ No newline at end of file
diff --git a/mri/reconstruct/linear.py b/mri/reconstruct/linear.py
index c07d38a8..e2f9689a 100644
--- a/mri/reconstruct/linear.py
+++ b/mri/reconstruct/linear.py
@@ -18,8 +18,8 @@
 from pysap.base.utils import unflatten
 
 # Third party import
-import numpy
-
+import numpy as np
+from modopt.signal.wavelet import get_mr_filters, filter_convolve
 
 class WaveletN(object):
     """ The 2D and 3D wavelet transform class.
@@ -69,7 +69,7 @@ def op(self, data):
         coeffs: ndarray
             the wavelet coefficients.
         """
-        if isinstance(data, numpy.ndarray):
+        if isinstance(data, np.ndarray):
             data = pysap.Image(data=data)
         self.transform.data = data
         self.transform.analysis()
@@ -114,13 +114,161 @@ def l2norm(self, shape):
             the L2 norm.
         """
         # Create fake data
-        shape = numpy.asarray(shape)
+        shape = np.asarray(shape)
+        shape += shape % 2
+        fake_data = np.zeros(shape)
+        fake_data[tuple(zip(shape // 2))] = 1
+
+        # Call mr_transform
+        data = self.op(fake_data)
+
+        # Compute the L2 norm
+        return np.linalg.norm(data)
+
+
+class WaveletUD2(object):
+    """The wavelet undecimated operator using pysap wrapper.
+    """
+    def __init__(self, wavelet_id, nb_scale=4, verbose=0, multichannel=False):
+        """Init function for Undecimated wavelet transform
+
+        Parameters
+        -----------
+        wavelet_id: int
+            ID of wavelet being used
+        nb_scale: int, default 4
+            the number of scales in the decomposition.
+
+        Private Variables:
+            _has_run: Checks if the get_mr_filters was called already
+        """
+        self.wavelet_id = wavelet_id
+        self.multichannel = multichannel
+        self.nb_scale = nb_scale
+        self._opt = [
+            '-t{}'.format(self.wavelet_id),
+            '-n{}'.format(self.nb_scale),
+        ]
+        self._has_run = False
+        self.coeffs_shape = None
+        self.flatten = flatten
+        self.unflatten = unflatten
+
+    def _get_filters(self, shape):
+        """Function to get the Wavelet coefficients of Delta[0][0].
+        This function is called only once and later the
+        wavelet coefficients are obtained by convolving these coefficients
+        with input Data
+        """
+        self.transform = get_mr_filters(
+            tuple(shape),
+            opt=self._opt,
+            coarse=True,
+        )
+        self._has_run = True
+
+    def op(self, data):
+        """ Define the wavelet operator.
+
+        This method returns the input data convolved with the wavelet filter.
+
+        Parameters
+        ----------
+        data: ndarray or Image
+            input 2D data array.
+
+        Returns
+        -------
+        coeffs: ndarray
+            the wavelet coefficients.
+        """
+        if not self._has_run:
+            if self.multichannel:
+                self._get_filters(list(data.shape)[1:])
+            else:
+                self._get_filters(data.shape)
+        if self.multichannel:
+            coeffs = []
+            self.coeffs_shape = []
+            for channel in range(data.shape[0]):
+                coefs_real, coeffs_shape = self.flatten(
+                    filter_convolve(data[channel].real, self.transform))
+                coefs_imag, coeffs_shape = self.flatten(
+                    filter_convolve(data[channel].imag, self.transform))
+                coeffs.append(coefs_real + 1j * coefs_imag)
+                self.coeffs_shape.append(coeffs_shape)
+            return np.asarray(coeffs)
+        else:
+            coefs_real = filter_convolve(data.real, self.transform)
+            coefs_imag = filter_convolve(data.imag, self.transform)
+            coeffs, self.coeffs_shape = self.flatten(
+                coefs_real + 1j * coefs_imag)
+        return coeffs
+
+    def adj_op(self, coefs):
+        """ Define the wavelet adjoint operator.
+
+        This method returns the reconsructed image.
+
+        Parameters
+        ----------
+        coeffs: ndarray
+            the wavelet coefficients.
+        dtype: str, default 'array'
+            if 'array' return the data as a ndarray, otherwise return a
+            pysap.Image.
+
+        Returns
+        -------
+        data: ndarray
+            the reconstructed data.
+        """
+        if not self._has_run:
+            raise RuntimeError(
+                "`op` must be run before `adj_op` to get the data shape",
+            )
+        if self.multichannel:
+            images = []
+            for channel, coeffs_shape in zip(range(coefs.shape[0]),
+                                             self.coeffs_shape):
+                data_real = filter_convolve(
+                    np.squeeze(self.unflatten(coefs.real[channel],
+                                              coeffs_shape)),
+                    self.transform, filter_rot=True)
+                data_imag = filter_convolve(
+                    np.squeeze(self.unflatten(coefs.imag[channel],
+                                              coeffs_shape)),
+                    self.transform, filter_rot=True)
+                images.append(data_real + 1j * data_imag)
+            return np.asarray(images)
+        else:
+            data_real = filter_convolve(
+                np.squeeze(self.unflatten(coefs.real, self.coeffs_shape)),
+                self.transform, filter_rot=True)
+            data_imag = filter_convolve(
+                np.squeeze(self.unflatten(coefs.imag, self.coeffs_shape)),
+                self.transform, filter_rot=True)
+        return data_real + 1j * data_imag
+
+    def l2norm(self, shape):
+        """ Compute the L2 norm.
+        Parameters
+        ----------
+        shape: uplet
+            the data shape.
+        Returns
+        -------
+        norm: float
+            the L2 norm.
+        """
+        # Create fake data
+        shape = np.asarray(shape)
         shape += shape % 2
-        fake_data = numpy.zeros(shape)
+        fake_data = np.zeros(shape)
         fake_data[tuple(zip(shape // 2))] = 1
 
         # Call mr_transform
         data = self.op(fake_data)
 
         # Compute the L2 norm
-        return numpy.linalg.norm(data)
+        return np.linalg.norm(data)
diff --git a/mri/test/test_wavelet_adjoint.py b/mri/test/test_wavelet_adjoint.py
index 39814f5d..c1cc8fce 100644
--- a/mri/test/test_wavelet_adjoint.py
+++ b/mri/test/test_wavelet_adjoint.py
@@ -12,7 +12,7 @@
 import numpy
 
 # Package import
-from mri.reconstruct.linear import WaveletN
+from mri.reconstruct.linear import WaveletN, WaveletUD2
 
 
 class TestAdjointOperatorWaveletTransform(unittest.TestCase):
@@ -80,6 +80,24 @@ def test_Wavelet3D_PyWt(self):
             numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
         print(" Wavelet3 adjoint test passes")
 
+    def test_Wavelet_UD_2D(self):
+        """Test the adjoint operation for Undecimated wavelet
+        """
+        for i in range(self.max_iter):
+            print("Process Wavelet Undecimated test '{0}'...", i)
+            wavelet_op_adj = WaveletUD2(wavelet_id=24,
+                                        nb_scale=4)
+            Img = (numpy.random.randn(self.N, self.N) +
+                   1j * numpy.random.randn(self.N, self.N))
+            f_p = wavelet_op_adj.op(Img)
+            f = (numpy.random.randn(*f_p.shape) +
+                 1j * numpy.random.randn(*f_p.shape))
+            I_p = wavelet_op_adj.adj_op(f)
+            x_d = numpy.vdot(Img, I_p)
+            x_ad = numpy.vdot(f_p, f)
+            numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
+        print("Undecimated Wavelet 2D adjoint test passes")
+
 
 if __name__ == "__main__":
     unittest.main()

From 4a6f87789ab6828c9341a77c457c80ebfcda7680 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 11:17:25 +0200
Subject: [PATCH 2/9] Make default as wavelet 24. Fix PEP8 errors

---
 mri/numerics/linear.py           | 2 +-
 mri/reconstruct/linear.py        | 6 ++++--
 mri/test/test_wavelet_adjoint.py | 3 +--
 3 files changed, 6 insertions(+), 5 deletions(-)

diff --git a/mri/numerics/linear.py b/mri/numerics/linear.py
index 02275f65..4391a9d5 100644
--- a/mri/numerics/linear.py
+++ b/mri/numerics/linear.py
@@ -13,4 +13,4 @@
 
 # Package import
 from mri.reconstruct.linear import WaveletN
-from mri.reconstruct.linear import WaveletUD2
\ No newline at end of file
+from mri.reconstruct.linear import WaveletUD2
diff --git a/mri/reconstruct/linear.py b/mri/reconstruct/linear.py
index e2f9689a..6dea5924 100644
--- a/mri/reconstruct/linear.py
+++ b/mri/reconstruct/linear.py
@@ -21,6 +21,7 @@
 import numpy as np
 from modopt.signal.wavelet import get_mr_filters, filter_convolve
 
+
 class WaveletN(object):
     """ The 2D and 3D wavelet transform class.
     """
@@ -129,12 +130,13 @@ def l2norm(self, shape):
 class WaveletUD2(object):
     """The wavelet undecimated operator using pysap wrapper.
     """
-    def __init__(self, wavelet_id, nb_scale=4, verbose=0, multichannel=False):
+    def __init__(self, wavelet_id=24, nb_scale=4,
+                 verbose=0, multichannel=False):
         """Init function for Undecimated wavelet transform
 
         Parameters
         -----------
-        wavelet_id: int
+        wavelet_id: int, default 24 = undecimated (bi-) orthogonal transform
             ID of wavelet being used
         nb_scale: int, default 4
             the number of scales in the decomposition.
diff --git a/mri/test/test_wavelet_adjoint.py b/mri/test/test_wavelet_adjoint.py
index c1cc8fce..904f2a34 100644
--- a/mri/test/test_wavelet_adjoint.py
+++ b/mri/test/test_wavelet_adjoint.py
@@ -85,8 +85,7 @@ def test_Wavelet_UD_2D(self):
         """
         for i in range(self.max_iter):
             print("Process Wavelet Undecimated test '{0}'...", i)
-            wavelet_op_adj = WaveletUD2(wavelet_id=24,
-                                        nb_scale=4)
+            wavelet_op_adj = WaveletUD2(nb_scale=4)
             Img = (numpy.random.randn(self.N, self.N) +
                    1j * numpy.random.randn(self.N, self.N))
             f_p = wavelet_op_adj.op(Img)

From 055624f3b98d470395535247abea40880829bf25 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 13:22:49 +0200
Subject: [PATCH 3/9] 1) Update to Travis, by adding mr_transform within path
 2) Updates based on Zac's comments 3) Add another test to test multichannel
 with multi-cpu -> This could fail

---
 .travis.yml                      |  21 +++---
 mri/reconstruct/linear.py        | 111 +++++++++++++++++++------------
 mri/test/test_wavelet_adjoint.py |  29 ++++++--
 3 files changed, 105 insertions(+), 56 deletions(-)

diff --git a/.travis.yml b/.travis.yml
index bd52a8bd..ab99b4de 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -44,20 +44,23 @@ install:
     - pip install coverage nose pytest pytest-cov
     - pip install coveralls
     - pip install pycodestyle
-    - pip install git+https://github.com/CEA-COSMIC/pysap@master
-    - pip install git+https://github.com/ghisvail/pyNFFT.git
-    - pip install -b $TRAVIS_BUILD_DIR/build -t $TRAVIS_BUILD_DIR/install --no-clean --upgrade .
-    - ls $TRAVIS_BUILD_DIR/install
-    - export PYTHONPATH=$TRAVIS_BUILD_DIR/install:$PYTHONPATH
-    - if [ $TRAVIS_PYTHON_VERSION == "3.5" ]; then
-        export PATH=$PATH:$TRAVIS_BUILD_DIR/build/temp.linux-x86_64-3.5/extern/bin;
+    - git clone https://github.com/CEA-COSMIC/pysap
+    - ls
+    - pushd pysap
+    - python setup.py install
+    - if [ $TRAVIS_PYTHON_VERSION == "3.5"]; then
+        export PATH=$PATH:$PWD/build/temp.linux-x86_64-3.5/extern/bin;
       fi
     - if [ $TRAVIS_PYTHON_VERSION == "3.6" ]; then
-        export PATH=$PATH:$TRAVIS_BUILD_DIR/build/temp.linux-x86_64-3.6/extern/bin;
+        export PATH=$PATH:$PWD/build/temp.linux-x86_64-3.6/extern/bin;
       fi
     - if [ $TRAVIS_PYTHON_VERSION == "3.7" ]; then
-        export PATH=$PATH:$TRAVIS_BUILD_DIR/build/temp.linux-x86_64-3.7/extern/bin;
+        export PATH=$PATH:$PWD/build/temp.linux-x86_64-3.7/extern/bin;
       fi
+    - popd
+    - export PYTHONPATH=$TRAVIS_BUILD_DIR/install:$PYTHONPATH
+    - pip install git+https://github.com/ghisvail/pyNFFT.git
+    - pip install -b $TRAVIS_BUILD_DIR/build -t $TRAVIS_BUILD_DIR/install --no-clean --upgrade .
 
 script:
     - python setup.py test
diff --git a/mri/reconstruct/linear.py b/mri/reconstruct/linear.py
index 6dea5924..1fbfb268 100644
--- a/mri/reconstruct/linear.py
+++ b/mri/reconstruct/linear.py
@@ -16,11 +16,11 @@
 import pysap
 from pysap.base.utils import flatten
 from pysap.base.utils import unflatten
+from modopt.signal.wavelet import get_mr_filters, filter_convolve
 
 # Third party import
 import numpy as np
-from modopt.signal.wavelet import get_mr_filters, filter_convolve
-
+from joblib import Parallel, delayed
 
 class WaveletN(object):
     """ The 2D and 3D wavelet transform class.
@@ -131,7 +131,7 @@ class WaveletUD2(object):
     """The wavelet undecimated operator using pysap wrapper.
     """
     def __init__(self, wavelet_id=24, nb_scale=4,
-                 verbose=0, multichannel=False):
+                 multichannel=False, n_cpu=1, verbose=0):
         """Init function for Undecimated wavelet transform
 
         Parameters
@@ -140,21 +140,24 @@ def __init__(self, wavelet_id=24, nb_scale=4,
             ID of wavelet being used
         nb_scale: int, default 4
             the number of scales in the decomposition.
-
+        multichannel: bool, default False
+            Boolean value to indicate if the incoming data is from
+            multiple-channels
+        n_cpu: int, default 0
+            Number of CPUs to run on. Only applicable if multichannel=True.
         Private Variables:
             _has_run: Checks if the get_mr_filters was called already
         """
         self.wavelet_id = wavelet_id
         self.multichannel = multichannel
         self.nb_scale = nb_scale
+        self.n_cpu = n_cpu
         self._opt = [
             '-t{}'.format(self.wavelet_id),
             '-n{}'.format(self.nb_scale),
         ]
         self._has_run = False
         self.coeffs_shape = None
-        self.flatten = flatten
-        self.unflatten = unflatten
 
     def _get_filters(self, shape):
         """Function to get the Wavelet coefficients of Delta[0][0].
@@ -169,6 +172,26 @@ def _get_filters(self, shape):
         )
         self._has_run = True
 
+    def _op(self, data):
+        """ Define the wavelet operator for single channel.
+        This is internal function that returns wavelet coefficients for a
+        single channel
+        Parameters
+        ----------
+        data: ndarray or Image
+            input 2D data array.
+
+        Returns
+        -------
+        coeffs: ndarray
+            the wavelet coefficients.
+        """
+        coefs_real = filter_convolve(data.real, self.transform)
+        coefs_imag = filter_convolve(data.imag, self.transform)
+        coeffs, coeffs_shape = flatten(
+            coefs_real + 1j * coefs_imag)
+        return coeffs, coeffs_shape
+
     def op(self, data):
         """ Define the wavelet operator.
 
@@ -190,23 +213,41 @@ def op(self, data):
             else:
                 self._get_filters(data.shape)
         if self.multichannel:
-            coeffs = []
-            self.coeffs_shape = []
-            for channel in range(data.shape[0]):
-                coefs_real, coeffs_shape = self.flatten(
-                    filter_convolve(data[channel].real, self.transform))
-                coefs_imag, coeffs_shape = self.flatten(
-                    filter_convolve(data[channel].imag, self.transform))
-                coeffs.append(coefs_real + 1j * coefs_imag)
-                self.coeffs_shape.append(coeffs_shape)
-            return np.asarray(coeffs)
+            coeffs, self.coeffs_shape = \
+                zip(*Parallel(n_jobs=self.n_cpu)
+                (delayed(self._op)
+                 (data[i])
+                 for i in np.arange(data.shape[0])))
+            coeffs = np.asarray(coeffs)
         else:
-            coefs_real = filter_convolve(data.real, self.transform)
-            coefs_imag = filter_convolve(data.imag, self.transform)
-            coeffs, self.coeffs_shape = self.flatten(
-                coefs_real + 1j * coefs_imag)
+            coeffs, self.coeffs_shape = self._op(data)
         return coeffs
 
+
+    def _adj_op(self, coefs, coeffs_shape):
+        """" Define the wavelet adjoint operator.
+
+        This method returns the reconstructed image for single channel.
+
+        Parameters
+        ----------
+        coeffs: ndarray
+            the wavelet coefficients.
+        coeffs_shape: ndarray
+            The shape of coefficients to unflatten before adjoint operation
+        Returns
+        -------
+        data: ndarray
+            the reconstructed data.
+        """
+        data_real = filter_convolve(
+                np.squeeze(unflatten(coefs.real, coeffs_shape)),
+                self.transform, filter_rot=True)
+        data_imag = filter_convolve(
+                np.squeeze(unflatten(coefs.imag, coeffs_shape)),
+                self.transform, filter_rot=True)
+        return data_real + 1j * data_imag
+
     def adj_op(self, coefs):
         """ Define the wavelet adjoint operator.
 
@@ -216,9 +257,6 @@ def adj_op(self, coefs):
         ----------
         coeffs: ndarray
             the wavelet coefficients.
-        dtype: str, default 'array'
-            if 'array' return the data as a ndarray, otherwise return a
-            pysap.Image.
 
         Returns
         -------
@@ -230,27 +268,14 @@ def adj_op(self, coefs):
                 "`op` must be run before `adj_op` to get the data shape",
             )
         if self.multichannel:
-            images = []
-            for channel, coeffs_shape in zip(range(coefs.shape[0]),
-                                             self.coeffs_shape):
-                data_real = filter_convolve(
-                    np.squeeze(self.unflatten(coefs.real[channel],
-                                              coeffs_shape)),
-                    self.transform, filter_rot=True)
-                data_imag = filter_convolve(
-                    np.squeeze(self.unflatten(coefs.imag[channel],
-                                              coeffs_shape)),
-                    self.transform, filter_rot=True)
-                images.append(data_real + 1j * data_imag)
-            return np.asarray(images)
+            images = Parallel(n_jobs=self.n_cpu)(
+                delayed(self._adj_op)
+                (coefs[i], self.coeffs_shape[i])
+                for i in np.arange(coefs.shape[0]))
+            images = np.asarray(images)
         else:
-            data_real = filter_convolve(
-                np.squeeze(self.unflatten(coefs.real, self.coeffs_shape)),
-                self.transform, filter_rot=True)
-            data_imag = filter_convolve(
-                np.squeeze(self.unflatten(coefs.imag, self.coeffs_shape)),
-                self.transform, filter_rot=True)
-        return data_real + 1j * data_imag
+            images = self._adj_op(coefs, self.coeffs_shape)
+        return images
 
     def l2norm(self, shape):
         """ Compute the L2 norm.
diff --git a/mri/test/test_wavelet_adjoint.py b/mri/test/test_wavelet_adjoint.py
index 904f2a34..2e7fdc2b 100644
--- a/mri/test/test_wavelet_adjoint.py
+++ b/mri/test/test_wavelet_adjoint.py
@@ -20,10 +20,15 @@ class TestAdjointOperatorWaveletTransform(unittest.TestCase):
     """
 
     def setUp(self):
-        """ Set the number of iterations.
+        """ Setup variables:
+        N = Image size
+        max_iter = Number of iterations to test
+        num_channels = Number of channels to be tested with for
+                        multichannel tests
         """
         self.N = 64
         self.max_iter = 10
+        self.num_channels = 10
 
     def test_Wavelet2D_ISAP(self):
         """Test the adjoint operator for the 2D Wavelet transform
@@ -85,18 +90,34 @@ def test_Wavelet_UD_2D(self):
         """
         for i in range(self.max_iter):
             print("Process Wavelet Undecimated test '{0}'...", i)
-            wavelet_op_adj = WaveletUD2(nb_scale=4)
+            wavelet_op = WaveletUD2(nb_scale=4)
             Img = (numpy.random.randn(self.N, self.N) +
                    1j * numpy.random.randn(self.N, self.N))
-            f_p = wavelet_op_adj.op(Img)
+            f_p = wavelet_op.op(Img)
             f = (numpy.random.randn(*f_p.shape) +
                  1j * numpy.random.randn(*f_p.shape))
-            I_p = wavelet_op_adj.adj_op(f)
+            I_p = wavelet_op.adj_op(f)
             x_d = numpy.vdot(Img, I_p)
             x_ad = numpy.vdot(f_p, f)
             numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
         print("Undecimated Wavelet 2D adjoint test passes")
 
+    def test_Wavelet_UD_2D_Multichannel(self):
+        """Test the adjoint operation for Undecmated wavelet Transform in
+        multichannel case"""
+        for i in range(self.max_iter):
+            print("Process Wavelet Undecimated test '{0}'...", i)
+            wavelet_op = WaveletUD2(nb_scale=4, multichannel=True, n_cpu=2)
+            Img = (numpy.random.randn(self.num_channels, self.N, self.N) +
+                   1j * numpy.random.randn(self.num_channels, self.N, self.N))
+            f_p = wavelet_op.op(Img)
+            f = (numpy.random.randn(*f_p.shape) +
+                 1j * numpy.random.randn(*f_p.shape))
+            I_p = wavelet_op.adj_op(f)
+            x_d = numpy.vdot(Img, I_p)
+            x_ad = numpy.vdot(f_p, f)
+            numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
+        print("Undecimated Wavelet 2D adjoint test passes for multichannel")
 
 if __name__ == "__main__":
     unittest.main()

From 0bbd37dc8b014f0b8ca4b4bd9c8ef1d634023a69 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 14:22:37 +0200
Subject: [PATCH 4/9] Push on Travis to check installation with updated Modopt

---
 .travis.yml | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/.travis.yml b/.travis.yml
index ab99b4de..e8b02a75 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -44,9 +44,10 @@ install:
     - pip install coverage nose pytest pytest-cov
     - pip install coveralls
     - pip install pycodestyle
+    - pushd ../
+    - pip install git+https://github.com/chaithyagr/ModOpt.git
     - git clone https://github.com/CEA-COSMIC/pysap
-    - ls
-    - pushd pysap
+    - cd pysap
     - python setup.py install
     - if [ $TRAVIS_PYTHON_VERSION == "3.5"]; then
         export PATH=$PATH:$PWD/build/temp.linux-x86_64-3.5/extern/bin;

From 4fb76d0dbb81eaf49205d1410a9c920bc47e0bd9 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 14:35:27 +0200
Subject: [PATCH 5/9] Upgrade pip for clean python 3.5 installation

---
 .travis.yml | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/.travis.yml b/.travis.yml
index e8b02a75..0bf9fe8a 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -39,6 +39,7 @@ install:
     - git clone https://github.com/CEA-COSMIC/pysap-data.git $HOME/.local/share/pysap/pysap-data
     - ln -s $HOME/.local/share/pysap/pysap-data/pysap-data/* $HOME/.local/share/pysap
     - ls -l $HOME/.local/share/pysap
+    - pip install --upgrade pip
     - pip install numpy
     - pip install cython
     - pip install coverage nose pytest pytest-cov
@@ -49,7 +50,7 @@ install:
     - git clone https://github.com/CEA-COSMIC/pysap
     - cd pysap
     - python setup.py install
-    - if [ $TRAVIS_PYTHON_VERSION == "3.5"]; then
+    - if [ $TRAVIS_PYTHON_VERSION == "3.5" ]; then
         export PATH=$PATH:$PWD/build/temp.linux-x86_64-3.5/extern/bin;
       fi
     - if [ $TRAVIS_PYTHON_VERSION == "3.6" ]; then

From 867c1d30f9e4a8e904436a5c38f426297efe0330 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 14:41:20 +0200
Subject: [PATCH 6/9] Fix PEP8 Errors

---
 .travis.yml                      |  3 +--
 mri/reconstruct/linear.py        | 11 +++++------
 mri/test/test_wavelet_adjoint.py |  1 +
 3 files changed, 7 insertions(+), 8 deletions(-)

diff --git a/.travis.yml b/.travis.yml
index 0bf9fe8a..01281001 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -40,13 +40,12 @@ install:
     - ln -s $HOME/.local/share/pysap/pysap-data/pysap-data/* $HOME/.local/share/pysap
     - ls -l $HOME/.local/share/pysap
     - pip install --upgrade pip
-    - pip install numpy
     - pip install cython
     - pip install coverage nose pytest pytest-cov
     - pip install coveralls
     - pip install pycodestyle
-    - pushd ../
     - pip install git+https://github.com/chaithyagr/ModOpt.git
+    - pushd ../
     - git clone https://github.com/CEA-COSMIC/pysap
     - cd pysap
     - python setup.py install
diff --git a/mri/reconstruct/linear.py b/mri/reconstruct/linear.py
index 1fbfb268..5a3e6f17 100644
--- a/mri/reconstruct/linear.py
+++ b/mri/reconstruct/linear.py
@@ -22,6 +22,7 @@
 import numpy as np
 from joblib import Parallel, delayed
 
+
 class WaveletN(object):
     """ The 2D and 3D wavelet transform class.
     """
@@ -213,17 +214,15 @@ def op(self, data):
             else:
                 self._get_filters(data.shape)
         if self.multichannel:
-            coeffs, self.coeffs_shape = \
-                zip(*Parallel(n_jobs=self.n_cpu)
-                (delayed(self._op)
-                 (data[i])
-                 for i in np.arange(data.shape[0])))
+            coeffs, self.coeffs_shape = zip(*Parallel(n_jobs=self.n_cpu)(
+                delayed(self._op)
+                (data[i])
+                for i in np.arange(data.shape[0])))
             coeffs = np.asarray(coeffs)
         else:
             coeffs, self.coeffs_shape = self._op(data)
         return coeffs
 
-
     def _adj_op(self, coefs, coeffs_shape):
         """" Define the wavelet adjoint operator.
 
diff --git a/mri/test/test_wavelet_adjoint.py b/mri/test/test_wavelet_adjoint.py
index 2e7fdc2b..f760bdbc 100644
--- a/mri/test/test_wavelet_adjoint.py
+++ b/mri/test/test_wavelet_adjoint.py
@@ -119,5 +119,6 @@ def test_Wavelet_UD_2D_Multichannel(self):
             numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
         print("Undecimated Wavelet 2D adjoint test passes for multichannel")
 
+
 if __name__ == "__main__":
     unittest.main()

From bd9757b1771077bcaa5560f9150013f9b46b461f Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Mon, 30 Sep 2019 15:06:47 +0200
Subject: [PATCH 7/9] Add README for Undecimated wavelet filters

---
 .travis.yml | 1 +
 README.rst  | 9 +++++++++
 2 files changed, 10 insertions(+)

diff --git a/.travis.yml b/.travis.yml
index 01281001..3efe88cd 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -40,6 +40,7 @@ install:
     - ln -s $HOME/.local/share/pysap/pysap-data/pysap-data/* $HOME/.local/share/pysap
     - ls -l $HOME/.local/share/pysap
     - pip install --upgrade pip
+    - pip install matplotlib
     - pip install cython
     - pip install coverage nose pytest pytest-cov
     - pip install coveralls
diff --git a/README.rst b/README.rst
index 422da8c3..6d811948 100644
--- a/README.rst
+++ b/README.rst
@@ -15,6 +15,15 @@ This work is made available by a community of people, amoung which the
 CEA Neurospin UNATI and CEA CosmoStat laboratories, in particular A. Grigis,
 J.-L. Starck, P. Ciuciu, and S. Farrens.
 
+Installation instructions
+===============
+
+Install python-pySAP using `pip install python-pySAP`. Later install pysap-mri by calling setup.py
+Note: If you want to use undecimated wavelet transform, please point the `$PATH` environment variable to
+pysap external binary directory:
+
+`export PATH=$PATH:/path-to-pysap/build/temp.linux-x86_64-<PYTHON_VERSION>/extern/bin/`
+
 Important links
 ===============
 

From c6b7a162a1c7512c477e05c27b92452986f37f37 Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Tue, 1 Oct 2019 14:27:07 +0200
Subject: [PATCH 8/9] Fix coding styles and variables cases Add backend support
 for multichannel parallel run Update docs for verbosity level

---
 mri/reconstruct/linear.py        | 23 +++++++++++++++++------
 mri/test/test_wavelet_adjoint.py | 16 ++++++++--------
 2 files changed, 25 insertions(+), 14 deletions(-)

diff --git a/mri/reconstruct/linear.py b/mri/reconstruct/linear.py
index 5a3e6f17..f1a194c5 100644
--- a/mri/reconstruct/linear.py
+++ b/mri/reconstruct/linear.py
@@ -13,14 +13,14 @@
 
 
 # Package import
+from modopt.signal.wavelet import get_mr_filters, filter_convolve
 import pysap
 from pysap.base.utils import flatten
 from pysap.base.utils import unflatten
-from modopt.signal.wavelet import get_mr_filters, filter_convolve
 
 # Third party import
-import numpy as np
 from joblib import Parallel, delayed
+import numpy as np
 
 
 class WaveletN(object):
@@ -131,8 +131,8 @@ def l2norm(self, shape):
 class WaveletUD2(object):
     """The wavelet undecimated operator using pysap wrapper.
     """
-    def __init__(self, wavelet_id=24, nb_scale=4,
-                 multichannel=False, n_cpu=1, verbose=0):
+    def __init__(self, wavelet_id=24, nb_scale=4, multichannel=False,
+                 n_cpu=1, backend='threading', verbose=0):
         """Init function for Undecimated wavelet transform
 
         Parameters
@@ -146,6 +146,11 @@ def __init__(self, wavelet_id=24, nb_scale=4,
             multiple-channels
         n_cpu: int, default 0
             Number of CPUs to run on. Only applicable if multichannel=True.
+        backend: 'threading' | 'multiprocessing', default 'threading'
+            Denotes the backend to use for parallel execution across
+            multiple channels.
+        verbose: int, default 0
+            The verbosity level for Parallel operation from joblib
         Private Variables:
             _has_run: Checks if the get_mr_filters was called already
         """
@@ -153,6 +158,8 @@ def __init__(self, wavelet_id=24, nb_scale=4,
         self.multichannel = multichannel
         self.nb_scale = nb_scale
         self.n_cpu = n_cpu
+        self.backend = backend
+        self.verbose = verbose
         self._opt = [
             '-t{}'.format(self.wavelet_id),
             '-n{}'.format(self.nb_scale),
@@ -214,7 +221,9 @@ def op(self, data):
             else:
                 self._get_filters(data.shape)
         if self.multichannel:
-            coeffs, self.coeffs_shape = zip(*Parallel(n_jobs=self.n_cpu)(
+            coeffs, self.coeffs_shape = zip(*Parallel(n_jobs=self.n_cpu,
+                                                      backend=self.backend,
+                                                      verbose=self.verbose)(
                 delayed(self._op)
                 (data[i])
                 for i in np.arange(data.shape[0])))
@@ -267,7 +276,9 @@ def adj_op(self, coefs):
                 "`op` must be run before `adj_op` to get the data shape",
             )
         if self.multichannel:
-            images = Parallel(n_jobs=self.n_cpu)(
+            images = Parallel(n_jobs=self.n_cpu,
+                              backend=self.backend,
+                              verbose=self.verbose)(
                 delayed(self._adj_op)
                 (coefs[i], self.coeffs_shape[i])
                 for i in np.arange(coefs.shape[0]))
diff --git a/mri/test/test_wavelet_adjoint.py b/mri/test/test_wavelet_adjoint.py
index f760bdbc..f4dcd6fc 100644
--- a/mri/test/test_wavelet_adjoint.py
+++ b/mri/test/test_wavelet_adjoint.py
@@ -91,13 +91,13 @@ def test_Wavelet_UD_2D(self):
         for i in range(self.max_iter):
             print("Process Wavelet Undecimated test '{0}'...", i)
             wavelet_op = WaveletUD2(nb_scale=4)
-            Img = (numpy.random.randn(self.N, self.N) +
+            img = (numpy.random.randn(self.N, self.N) +
                    1j * numpy.random.randn(self.N, self.N))
-            f_p = wavelet_op.op(Img)
+            f_p = wavelet_op.op(img)
             f = (numpy.random.randn(*f_p.shape) +
                  1j * numpy.random.randn(*f_p.shape))
-            I_p = wavelet_op.adj_op(f)
-            x_d = numpy.vdot(Img, I_p)
+            i_p = wavelet_op.adj_op(f)
+            x_d = numpy.vdot(img, i_p)
             x_ad = numpy.vdot(f_p, f)
             numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
         print("Undecimated Wavelet 2D adjoint test passes")
@@ -108,13 +108,13 @@ def test_Wavelet_UD_2D_Multichannel(self):
         for i in range(self.max_iter):
             print("Process Wavelet Undecimated test '{0}'...", i)
             wavelet_op = WaveletUD2(nb_scale=4, multichannel=True, n_cpu=2)
-            Img = (numpy.random.randn(self.num_channels, self.N, self.N) +
+            img = (numpy.random.randn(self.num_channels, self.N, self.N) +
                    1j * numpy.random.randn(self.num_channels, self.N, self.N))
-            f_p = wavelet_op.op(Img)
+            f_p = wavelet_op.op(img)
             f = (numpy.random.randn(*f_p.shape) +
                  1j * numpy.random.randn(*f_p.shape))
-            I_p = wavelet_op.adj_op(f)
-            x_d = numpy.vdot(Img, I_p)
+            i_p = wavelet_op.adj_op(f)
+            x_d = numpy.vdot(img, i_p)
             x_ad = numpy.vdot(f_p, f)
             numpy.testing.assert_allclose(x_d, x_ad, rtol=1e-6)
         print("Undecimated Wavelet 2D adjoint test passes for multichannel")

From 5b4b95ae3f53358c8d43693d260ccb0445717fbe Mon Sep 17 00:00:00 2001
From: chaithyagr <chaithyagr@gmail.com>
Date: Tue, 1 Oct 2019 16:25:59 +0200
Subject: [PATCH 9/9] Travis point to CEA-COSMIC Modopt

---
 .travis.yml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/.travis.yml b/.travis.yml
index 3efe88cd..c5084b5f 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -45,7 +45,7 @@ install:
     - pip install coverage nose pytest pytest-cov
     - pip install coveralls
     - pip install pycodestyle
-    - pip install git+https://github.com/chaithyagr/ModOpt.git
+    - pip install git+https://github.com/CEA-COSMIC/ModOpt.git
     - pushd ../
     - git clone https://github.com/CEA-COSMIC/pysap
     - cd pysap