# Compressive Sensing Operators¶

Let’s import necessary libraries

```import jax.numpy as jnp
# Random numbers
from jax import random
# For plotting diagrams
import matplotlib.pyplot as plt
## CR-Sparse modules
import cr.nimble as cnb
# Linear operators
from cr.sparse import lop
# Some random number seeds
key = random.PRNGKey(0)
keys = random.split(key, 10)
```

Size of the sensing matrix

```m, n = 4, 8
```

### Unnormalized atoms / columns¶

Construct the operator wrapping the sensing matrix

```Phi = lop.rademacher_dict(keys[0], m, n, normalize_atoms=False)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-1.  1. -1. -1.  1. -1. -1.  1.]
[ 1. -1. -1. -1.  1. -1.  1.  1.]
[-1. -1. -1. -1. -1.  1. -1. -1.]
[ 1. -1. -1.  1.  1. -1.  1. -1.]]
```

Let’s print the contents of the adjoint matrix

```print(lop.to_adjoint_matrix(Phi))
```
```[[-1.  1. -1.  1.]
[ 1. -1. -1. -1.]
[-1. -1. -1. -1.]
[-1. -1. -1.  1.]
[ 1.  1. -1.  1.]
[-1. -1.  1. -1.]
[-1.  1. -1.  1.]
[ 1.  1. -1. -1.]]
```

Let’s apply the forward operator

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[-24  -2 -18  16]
```

```z = Phi.trans(y)
print(z)
```
```[ 56 -20  28  60   8  -8  56 -24]
```

Let’s JIT compile the trans and times functions

```Phi = lop.jit(Phi)
```

Let’s verify the forward operator after JIT compile

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[-24  -2 -18  16]
```

Let’s verify the adjoint operator after JIT compile

```z = Phi.trans(y)
print(z)
```
```[ 56 -20  28  60   8  -8  56 -24]
```

### Normalized atoms / columns¶

Construct the operator wrapping the sensing matrix

```Phi = lop.rademacher_dict(keys[0], m, n, normalize_atoms=True)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-0.5  0.5 -0.5 -0.5  0.5 -0.5 -0.5  0.5]
[ 0.5 -0.5 -0.5 -0.5  0.5 -0.5  0.5  0.5]
[-0.5 -0.5 -0.5 -0.5 -0.5  0.5 -0.5 -0.5]
[ 0.5 -0.5 -0.5  0.5  0.5 -0.5  0.5 -0.5]]
```

Check the column wise norms

```print(cnb.norms_l2_cw(lop.to_matrix(Phi)))
```
```[1. 1. 1. 1. 1. 1. 1. 1.]
```

Let’s print the contents of the adjoint matrix

```print(lop.to_adjoint_matrix(Phi))
```
```[[-0.5  0.5 -0.5  0.5]
[ 0.5 -0.5 -0.5 -0.5]
[-0.5 -0.5 -0.5 -0.5]
[-0.5 -0.5 -0.5  0.5]
[ 0.5  0.5 -0.5  0.5]
[-0.5 -0.5  0.5 -0.5]
[-0.5  0.5 -0.5  0.5]
[ 0.5  0.5 -0.5 -0.5]]
```

Let’s apply the forward operator

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[-12.  -1.  -9.   8.]
```

```z = Phi.trans(y)
print(z)
```
```[14. -5.  7. 15.  2. -2. 14. -6.]
```

Let’s JIT compile the trans and times functions

```Phi = lop.jit(Phi)
```

Let’s verify the forward operator after JIT compile

```y = Phi.times(x)
print(y)
```
```[-12.  -1.  -9.   8.]
```

Let’s verify the adjoint operator after JIT compile

```z = Phi.trans(y)
print(z)
```
```[14. -5.  7. 15.  2. -2. 14. -6.]
```

### Column wise application on input¶

Number of signals

```k = 5
```

Prepare a signal matrix

```X = random.randint(keys[1], (n, k), -5, 5)
print(X)
```
```[[-4 -5 -2 -4  2]
[ 0  2 -2  1  3]
[ 2 -2  4 -4 -2]
[ 2  4 -5  1 -1]
[ 2 -5 -2 -2 -4]
[-4  2 -3 -1 -4]
[-1 -1 -1  4  1]
[-3 -4 -1 -4  1]]
```

Let’s apply the forward operator column wise

```Y = Phi.times(X)
print(Y)
```
```[[  2.   -2.5   1.   -0.5   2. ]
[ -3.  -10.5   0.   -1.5   2. ]
[ -1.    6.5   3.    3.5  -2. ]
[  2.   -2.5  -4.    3.5   0. ]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[ -1.   -8.5  -4.   -0.5   1. ]
[  2.    2.    1.   -3.    1. ]
[  0.    4.5   0.   -2.5  -1. ]
[  2.    2.   -4.    1.   -1. ]
[  1.  -11.   -3.   -1.    3. ]
[ -1.   11.    3.    1.   -3. ]
[ -1.   -8.5  -4.   -0.5   1. ]
[ -1.   -8.5   1.   -4.5   3. ]]
```

### Row wise application on input¶

Construct the operator wrapping the sensing matrix

```Phi = lop.rademacher_dict(keys[0], m, n, normalize_atoms=True, axis=1)
```

Prepare a signal matrix with each signal along a row

```XT = X.T
```

Let’s apply the forward operator column wise

```Y = Phi.times(XT)
print(Y)
```
```[[  2.   -3.   -1.    2. ]
[ -2.5 -10.5   6.5  -2.5]
[  1.    0.    3.   -4. ]
[ -0.5  -1.5   3.5   3.5]
[  2.    2.   -2.    0. ]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[ -1.    2.    0.    2.    1.   -1.   -1.   -1. ]
[ -8.5   2.    4.5   2.  -11.   11.   -8.5  -8.5]
[ -4.    1.    0.   -4.   -3.    3.   -4.    1. ]
[ -0.5  -3.   -2.5   1.   -1.    1.   -0.5  -4.5]
[  1.    1.   -1.   -1.    3.   -3.    1.    3. ]]
```

### Square sensing matrix¶

Construct the operator wrapping the sensing matrix

```Phi = lop.rademacher_dict(keys[0], n, normalize_atoms=False)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[ 1.  1. -1. -1. -1. -1.  1.  1.]
[-1.  1.  1.  1. -1. -1. -1. -1.]
[-1. -1.  1.  1.  1. -1.  1. -1.]
[ 1. -1.  1. -1.  1.  1.  1. -1.]
[-1.  1. -1.  1.  1.  1. -1. -1.]
[-1.  1. -1.  1.  1. -1.  1.  1.]
[ 1. -1.  1. -1. -1. -1. -1.  1.]
[-1. -1.  1. -1.  1. -1.  1. -1.]]
```

## Gaussian sensing operators¶

### Unnormalized atoms / columns¶

Construct the operator wrapping the sensing matrix

```Phi = lop.gaussian_dict(keys[0], m, n, normalize_atoms=False)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-0.76261741 -0.54906996  0.35375743  0.44785326  0.20306852  0.4002937
-0.10082346 -1.38919403]
[ 0.68214811 -0.73785154 -0.31989581 -0.37379108  0.02430423  0.60747778
0.24837508  0.39562076]
[ 0.08833386  0.49151262 -0.45341975 -0.63881669 -0.23097626  0.10892717
-0.75912382  0.09899961]
[-0.97121466  0.38958316 -0.38142231  1.05170598 -0.59028129 -0.18775354
-0.78933777  0.06830431]]
```

Let’s print the contents of the adjoint matrix

```print(lop.to_adjoint_matrix(Phi))
```
```[[-0.76261741  0.68214811  0.08833386 -0.97121466]
[-0.54906996 -0.73785154  0.49151262  0.38958316]
[ 0.35375743 -0.31989581 -0.45341975 -0.38142231]
[ 0.44785326 -0.37379108 -0.63881669  1.05170598]
[ 0.20306852  0.02430423 -0.23097626 -0.59028129]
[ 0.4002937   0.60747778  0.10892717 -0.18775354]
[-0.10082346  0.24837508 -0.75912382 -0.78933777]
[-1.38919403  0.39562076  0.09899961  0.06830431]]
```

Let’s apply the forward operator

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[  5.79912122   1.61042927 -12.54007687  -3.16525231]
```

```z = Phi.trans(y)
print(z)
```
```[ -1.3575335  -11.76911601   8.42952898   6.67708692   5.98160846
2.52798389  11.83322722  -8.87664847]
```

Let’s JIT compile the trans and times functions

```Phi = lop.jit(Phi)
```

Let’s verify the forward operator after JIT compile

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[  5.79912122   1.61042927 -12.54007687  -3.16525231]
```

Let’s verify the adjoint operator after JIT compile

```z = Phi.trans(y)
print(z)
```
```[ -1.3575335  -11.76911601   8.42952898   6.67708692   5.98160846
2.52798389  11.83322722  -8.87664847]
```

### Normalized atoms / columns¶

Construct the operator wrapping the sensing matrix

```Phi = lop.gaussian_dict(keys[0], m, n, normalize_atoms=True)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-0.53952553 -0.4932261   0.46508795  0.32887168  0.30488928  0.52725899
-0.08942487 -0.95844277]
[ 0.48259628 -0.6628074  -0.42056978 -0.27448567  0.03649064  0.80015779
0.22029504  0.27294953]
[ 0.06249316  0.4415227  -0.59611486 -0.46910168 -0.34679027  0.14347673
-0.6733011   0.06830253]
[-0.6871009   0.34996011 -0.50145921  0.77229829 -0.88625475 -0.24730527
-0.70009922  0.047125  ]]
```

Check the column wise norms

```print(cnb.norms_l2_cw(lop.to_matrix(Phi)))
```
```[1. 1. 1. 1. 1. 1. 1. 1.]
```

Let’s print the contents of the adjoint matrix

```print(lop.to_adjoint_matrix(Phi))
```
```[[-0.53952553  0.48259628  0.06249316 -0.6871009 ]
[-0.4932261  -0.6628074   0.4415227   0.34996011]
[ 0.46508795 -0.42056978 -0.59611486 -0.50145921]
[ 0.32887168 -0.27448567 -0.46910168  0.77229829]
[ 0.30488928  0.03649064 -0.34679027 -0.88625475]
[ 0.52725899  0.80015779  0.14347673 -0.24730527]
[-0.08942487  0.22029504 -0.6733011  -0.70009922]
[-0.95844277  0.27294953  0.06830253  0.047125  ]]
```

Let’s apply the forward operator

```x = jnp.array([3, -2, 4, 8, 2, 2, 6, -1])
y = Phi.times(x)
print(y)
```
```[  5.94539092   1.61735679 -11.34757498  -5.10351374]
```

```z = Phi.trans(y)
print(z)
```
```[  0.37032318 -10.80065612  11.40858028   2.89306098  10.32994614
4.06291433  11.0379306   -6.27243136]
```

Let’s JIT compile the trans and times functions

```Phi = lop.jit(Phi)
```

Let’s verify the forward operator after JIT compile

```y = Phi.times(x)
print(y)
```
```[  5.94539092   1.61735679 -11.34757498  -5.10351374]
```

Let’s verify the adjoint operator after JIT compile

```z = Phi.trans(y)
print(z)
```
```[  0.37032318 -10.80065612  11.40858028   2.89306098  10.32994614
4.06291433  11.0379306   -6.27243136]
```

### Column wise application on input¶

Number of signals

```k = 5
```

Prepare a signal matrix

```X = random.randint(keys[1], (n, k), -5, 5)
print(X)
```
```[[-4 -5 -2 -4  2]
[ 0  2 -2  1  3]
[ 2 -2  4 -4 -2]
[ 2  4 -5  1 -1]
[ 2 -5 -2 -2 -4]
[-4  2 -3 -1 -4]
[-1 -1 -1  4  1]
[-3 -4 -1 -4  1]]
```

Let’s apply the forward operator column wise

```Y = Phi.times(X)
print(Y)
```
```[[ 5.21151717  5.54975384  1.13780878  2.47243    -8.19423762]
[-7.48728951 -3.88963009 -2.91612777 -2.26915609 -2.76095358]
[-3.17949966  2.30739835 -0.17883383 -0.30940293  3.31914138]
[ 3.06551754 12.67579867 -2.02564721  4.9074167   3.78756451]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[ -8.73010144 -13.43671088  -0.64054117  -5.82025136   0.69355917]
[  2.06116221   5.29558478   0.58377917   1.86533811   8.66256504]
[  5.93085649  -3.51488594   2.87800132  -0.17219376  -6.52776613]
[  7.62807664  11.59989779  -0.30588428   5.371094    -0.5688957 ]
[ -0.29849024 -10.48404468   2.0977516   -3.57090899  -7.10687709]
[ -4.45749663  -2.98989408  -1.25814771  -1.77009412  -6.99014907]
[ -2.12085781 -11.7810428    0.79440614  -3.94833782  -4.76190827]
[ -7.11129849  -5.62584597  -1.99415376  -2.77891873   7.50530159]]
```

### Row wise application on input¶

Construct the operator wrapping the sensing matrix

```Phi = lop.gaussian_dict(keys[0], m, n, normalize_atoms=True, axis=1)
```

Prepare a signal matrix with each signal along a row

```XT = X.T
```

Let’s apply the forward operator column wise

```Y = Phi.times(XT)
print(Y)
```
```[[ 5.21151717 -7.48728951 -3.17949966  3.06551754]
[ 5.54975384 -3.88963009  2.30739835 12.67579867]
[ 1.13780878 -2.91612777 -0.17883383 -2.02564721]
[ 2.47243    -2.26915609 -0.30940293  4.9074167 ]
[-8.19423762 -2.76095358  3.31914138  3.78756451]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[ -8.73010144   2.06116221   5.93085649   7.62807664  -0.29849024
-4.45749663  -2.12085781  -7.11129849]
[-13.43671088   5.29558478  -3.51488594  11.59989779 -10.48404468
-2.98989408 -11.7810428   -5.62584597]
[ -0.64054117   0.58377917   2.87800132  -0.30588428   2.0977516
-1.25814771   0.79440614  -1.99415376]
[ -5.82025136   1.86533811  -0.17219376   5.371094    -3.57090899
-1.77009412  -3.94833782  -2.77891873]
[  0.69355917   8.66256504  -6.52776613  -0.5688957   -7.10687709
-6.99014907  -4.76190827   7.50530159]]
```

### Square sensing matrix¶

Construct the operator wrapping the sensing matrix

```Phi = lop.gaussian_dict(keys[0], n, normalize_atoms=False)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-5.06479786e-01 -1.76665935e-01  2.54534851e-01 -5.06364763e-01
-1.96507134e-01  1.79333415e-01  6.42379919e-01  2.10040808e-01]
[ 1.07240051e-01 -5.94812765e-01  1.25249038e-01  1.14459574e-01
-2.02012833e-01 -4.09751511e-01 -4.33936884e-01 -1.34516542e-01]
[ 1.11053768e-01  3.48672071e-01  3.68653905e-01  9.45349248e-02
-3.29041791e-01 -1.77899332e-01 -4.04006639e-01 -1.00726911e-01]
[ 1.33743289e-01 -5.08041866e-01 -4.99556060e-01  1.33236537e-01
-1.68340852e-01 -2.78095509e-01 -9.48804677e-02 -3.28437869e-01]
[-1.03963687e-03 -5.37674760e-02 -2.64978470e-01 -5.69247421e-01
-2.17436962e-01 -3.91220044e-01 -2.19387753e-01 -2.30803981e-01]
[ 1.71456228e-01  1.59848991e-01 -2.55230677e-01 -2.51277665e-01
2.37702097e-02 -3.31785957e-01 -1.86295948e-01  9.79684620e-01]
[-4.17172092e-02  2.52724150e-01  2.02487423e-01 -3.80860117e-01
-2.98349811e-01 -1.20580149e-01  8.39526917e-01 -2.28584283e-01]
[ 2.06869121e-01  1.04520297e+00  8.87791963e-01  2.20836838e-01
1.17742766e-01 -1.85396223e-01  4.36565908e-03 -2.27983571e-01]]
```

## Random matrix with orthonormal rows¶

Construct the operator wrapping the sensing matrix

```Phi = lop.random_orthonormal_rows_dict(keys[0], m, n)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-0.48805761  0.12995918  0.4365591   0.01555415  0.05653164 -0.14781949
-0.62155505 -0.37776645]
[-0.46559454  0.31376094 -0.50148534  0.45550458  0.37354653  0.06455508
0.2203263  -0.18318237]
[ 0.34953016 -0.08652247 -0.37115078  0.21179339 -0.33253457 -0.53588747
-0.13369294 -0.52163662]
[-0.34327893 -0.77342544  0.15727696  0.3162613  -0.21117691  0.13210122
0.2701341  -0.15555246]]
```

Check the row wise norms

```print(cnb.norms_l2_rw(lop.to_matrix(Phi)))
```
```[1. 1. 1. 1.]
```

Let’s apply the forward operator column wise

```Y = Phi.times(X)
print(Y)
```
```[[ 5.3156526   3.44362845  3.71437922 -0.58889084 -2.10908012]
[ 2.58851023  4.55426428 -4.95769928  5.43997158 -1.15770268]
[ 1.46024812  2.47991672 -0.14152435  2.96448472  3.78835967]
[ 1.56595683  2.79219775  1.19267863  2.27984255 -3.206765  ]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[-3.82670014 -3.8928266   0.03655821 -1.9918535   3.99333346]
[ 0.16549613 -0.49764393 -1.98301766 -0.38946395  1.51507918]
[ 0.72681332 -1.26182404  4.3478673  -3.72685574 -2.2505701 ]
[ 2.06628118  3.53634503 -1.85325656  3.81765653 -0.77196994]
[ 0.45115468  0.48159786 -1.84675587  0.53154779 -1.13425141]
[-1.19431947 -1.17513811 -0.63570689 -0.84923291 -2.21672446]
[-2.50586037 -0.71426    -3.05989872  1.7841294  -0.31689645]
[-3.48755201 -3.86309685 -0.60670455 -2.67506262 -0.46851652]]
```

The frame operator

```F = lop.frame(Phi)
```

Check that it’s close to identity matrix

```print(jnp.round(lop.to_matrix(F),2))
```
```[[ 1. -0.  0.  0.]
[-0.  1.  0.  0.]
[ 0.  0.  1. -0.]
[ 0.  0. -0.  1.]]
```

## Random orthonormal basis operator¶

Construct the operator wrapping the sensing matrix

```Phi = lop.random_onb_dict(keys[0], n)
```

Let’s print the contents of the matrix

```print(lop.to_matrix(Phi))
```
```[[-0.83021008  0.10099603  0.1924291  -0.18709837 -0.05440769 -0.42034524
0.05004106  0.21531801]
[ 0.17578544 -0.49159467  0.73866229 -0.34880127  0.22399819  0.03967139
-0.03076634  0.08633818]
[ 0.1820368   0.2083679   0.2342978  -0.08887811 -0.75177769 -0.00617792
-0.49250229  0.22996581]
[ 0.21922894 -0.43921983 -0.12485129  0.04649406 -0.44739928 -0.55798428
0.46709786 -0.10762012]
[-0.00170415 -0.03951229 -0.294542   -0.6724019   0.03554208 -0.17193851
-0.34738298 -0.55502752]
[ 0.28104713  0.04021688 -0.36937749 -0.38856286  0.23274741 -0.19465774
0.00554408  0.73476564]
[-0.06838189  0.20710371  0.01707164 -0.47829639 -0.29978828  0.51519234
0.60668761 -0.00266108]
[ 0.33909513  0.6824461   0.35125433 -0.0373899   0.19031866 -0.4212764
0.21687527 -0.18377793]]
```

Check the row wise norms

```print(cnb.norms_l2_rw(lop.to_matrix(Phi)))
```
```[1. 1. 1. 1. 1. 1. 1. 1.]
```

Check the column wise norms

```print(cnb.norms_l2_cw(lop.to_matrix(Phi)))
```
```[1. 1. 1. 1. 1. 1. 1. 1.]
```

Let’s apply the forward operator column wise

```Y = Phi.times(X)
print(Y)
```
```[[ 4.20807229  1.73982564  4.26812837  2.33307439  0.60917887]
[ 0.13764291 -6.08988076  4.70769159 -5.45427267 -3.25084291]
[-2.11354663  2.00161547  2.38539797 -2.92598768  3.37874579]
[ 0.15927087 -0.45449504  1.91737998  2.9814184   3.7050189 ]
[ 0.84423256 -0.12512179  3.61341624  1.4045029   0.78271592]
[-3.60578444 -6.63795754 -0.79905496 -3.18274785  2.3980138 ]
[-3.90797221  0.74207081  0.63229793  2.4558273   0.67111076]
[ 1.67154977 -2.45856137  0.39897891 -0.4730898   3.0173381 ]]
```

Let’s apply the adjoint operator column wise

```Z = Phi.trans(Y)
print(Z)
```
```[[-4. -5. -2. -4.  2.]
[ 0.  2. -2.  1.  3.]
[ 2. -2.  4. -4. -2.]
[ 2.  4. -5.  1. -1.]
[ 2. -5. -2. -2. -4.]
[-4.  2. -3. -1. -4.]
[-1. -1. -1.  4.  1.]
[-3. -4. -1. -4.  1.]]
```

The gram operator

```G = lop.gram(Phi)
```

Check that it’s close to identity matrix

```print(jnp.round(lop.to_matrix(G),2))
```
```[[ 1.  0.  0.  0.  0.  0.  0. -0.]
[ 0.  1.  0. -0.  0.  0. -0. -0.]
[ 0.  0.  1. -0. -0.  0.  0.  0.]
[ 0. -0. -0.  1.  0. -0. -0.  0.]
[ 0.  0. -0.  0.  1. -0. -0. -0.]
[ 0.  0.  0. -0. -0.  1. -0.  0.]
[ 0. -0.  0. -0. -0. -0.  1.  0.]
[-0. -0.  0.  0. -0.  0.  0.  1.]]
```

The frame operator

```F = lop.frame(Phi)
```

Check that it’s close to identity matrix

```print(jnp.round(lop.to_matrix(F),2))
```
```[[ 1. -0. -0.  0.  0.  0.  0.  0.]
[-0.  1.  0.  0.  0. -0.  0. -0.]
[-0.  0.  1. -0.  0.  0. -0. -0.]
[ 0.  0. -0.  1.  0.  0. -0.  0.]
[ 0.  0.  0.  0.  1. -0.  0.  0.]
[ 0. -0.  0.  0. -0.  1.  0.  0.]
[ 0.  0. -0. -0.  0.  0.  1.  0.]
[ 0. -0. -0.  0.  0.  0.  0.  1.]]
```

Total running time of the script: (0 minutes 2.243 seconds)

Gallery generated by Sphinx-Gallery