.. _diagonalize:

==================================
Create linear operator from volume
==================================

It may be interesting to create linear operators that describe a single
contribution for each voxel as in a volume mask. This can be encoded as follows:

.. math::

    \begin{bmatrix}
        w_1 \cdot \mathbf{g} & & & \\
        & w_2 \cdot \mathbf{g} & & \\
        & & \ddots & \\
        & & & w_n \cdot \mathbf{g} \\
    \end{bmatrix}

where :math:`\mathbf{g}` is the generator used for every voxel and
:math:`w_j` is the value of the mask at voxel :math:`j`. Only the voxels
exhibiting non-zero value are considered.

To build such a linear operator, one just needs to provide a three-dimensional
ndarray to the `talon.diagonalize` function.

Example
-------
Let us build a toy volume of dimension 2-by-2-by-2 with values from 0 to 7.

.. code:: python

    import numpy as np
    values = np.arange(2 ** 3).astype(np.float64)

    mask = values.reshape((2, ) * 3)
    print(mask)

Output:

.. code::

    [[[0. 1.]
      [2. 3.]]

     [[4. 5.]
      [6. 7.]]]

To diagonalize the volume, call the corresponding `talon` function.

.. code:: python

    import talon
    indices, weights = talon.diagonalize(mask)

The considered generator is vector :math:`g = [1, 10]^T`.

.. code:: python

    generators = np.array([[1.0, 10.0]])
    linear_operator = talon.operator(generators, indices, weights)

Check the output:

.. code:: python

    print(linear_operator.todense())

    [[ 0.  0.  0.  0.  0.  0.  0.]
     [ 0.  0.  0.  0.  0.  0.  0.]
     [ 1.  0.  0.  0.  0.  0.  0.]
     [10.  0.  0.  0.  0.  0.  0.]
     [ 0.  2.  0.  0.  0.  0.  0.]
     [ 0. 20.  0.  0.  0.  0.  0.]
     [ 0.  0.  3.  0.  0.  0.  0.]
     [ 0.  0. 30.  0.  0.  0.  0.]
     [ 0.  0.  0.  4.  0.  0.  0.]
     [ 0.  0.  0. 40.  0.  0.  0.]
     [ 0.  0.  0.  0.  5.  0.  0.]
     [ 0.  0.  0.  0. 50.  0.  0.]
     [ 0.  0.  0.  0.  0.  6.  0.]
     [ 0.  0.  0.  0.  0. 60.  0.]
     [ 0.  0.  0.  0.  0.  0.  7.]
     [ 0.  0.  0.  0.  0.  0. 70.]]