Plane Thermoelastic Waves in Infinite Half-Space Caused


FACTA UNIVERSITATIS   

Series: Electronics and Energetics Vol. 30, No 4, December 2017, pp. 477 - 510 

DOI: 10.2298/FUEE1704477D 

ON SOME COMMON COMPRESSIVE SENSING RECOVERY 

ALGORITHMS AND APPLICATIONS

 

Anđela Draganić, Irena Orović, Srđan Stanković  

University of Montenegro, Faculty of Electrical Engineering, Podgorica, Montenegro 

Abstract. Compressive Sensing, as an emerging technique in signal processing is reviewed 

in this paper together with its common applications. As an alternative to the traditional 

signal sampling, Compressive Sensing allows a new acquisition strategy with significantly 

reduced number of samples needed for accurate signal reconstruction. The basic ideas and 

motivation behind this approach are provided in the theoretical part of the paper. The 

commonly used algorithms for missing data reconstruction are presented. The Compressive 

Sensing applications have gained significant attention leading to an intensive growth of 

signal processing possibilities. Hence, some of the existing practical applications assuming 

different types of signals in real-world scenarios are described and analyzed as well. 

Key words: compressive sensing, optimization algorithms, sampling theorem, under-

sampled data  

1. THE BASIC COMPRESSIVE SENSING CONCEPTS 

Continuous time, bandlimited signals, sampled according to the Shannon-Nyquist 

sampling theorem, may produce a large number of samples to be further processed. 

Having in mind that the signal samples are acquired at the rate being at least twice the 

maximal signal frequency, the conventional sampling may be inefficient, especially in 

applications dealing with the high-frequency signals. Also, a large number of sensors 

required for acquisition may lead to large power consumption. Hence, the compression 

arises as a necessary step in conventional signal processing. Most of the signals we are 

dealing with contain redundant information, and this fact is exploited during the 

compression step. The compression discards certain percent of the samples in the sparse 

transformation domain, assuming that the majority of samples are insignificant for signal 

analysis. The fact that  most signals exhibit sparsity in certain transformation domain is 

used in the Compressive Sensing (CS) theory [1]-[12]. Namely, one of the ideas behind 

the CS was to avoid compression after acquisition and to directly acquire data in the 

compressed form. In other words, the CS offers the possibility to acquire less data then it 

                                                           
Received April 19, 2017 

Corresponding author: Anđela Draganić  

University of Montenegro, Faculty of Electrical Engineering, Džordža Vašingtona bb, 81000 Podgorica Montenegro 

(e-mail: andjelad@ac.me) 



478 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

is commonly done, but still to be able to reconstruct the entire information afterwards. 

The missing signal information can also appear as a consequence of omitting samples that 

are exposed to different kinds of noise or losing some parts of the signal during 

transmission. These missing samples can be recovered using the CS reconstruction 

algorithms [7]-[27].  

Some of the concepts that are nowadays used within the CS approaches date from the 

early seventies. The least square solutions, based on the norm minimization, are used by 

Claerbout and Muir in 1973 [11]. In 1986, Santosa and Symes proposed an application of 

the ℓ1-norm in recovering sparse spike trains. The ℓ1-minimization of the image gradient - 

total variation minimization was proposed in 1990s by Rudin, Osher and Fatemi [17] for 

removing noise from images. In the early 2000s Blu, Marziliano and Vetterli showed that 

the K-sparse signals can be sampled and recovered by using only 2K parameters. The idea 

of the CS starts to grow from the moment when it was shown that a small set of non-

adaptive measurements can provide exact signal reconstruction, which proved the basic 

idea behind the data acquisition in the compressed form [3],[4]. Later, in [19], the CS is 

analyzed in terms of signal recovery when the missing samples are result of signal 

degradation due to the presence of noise. The influence of the number of missing samples 

on the spectral signal representation is examined, and the reconstruction procedure is 

proposed.  

Our focus in this paper is on the practical applications of the CS approaches. 

However, there are some specific requirements that are imposed to the measurements in 

order to be able to apply CS signal reconstruction. Signal sparsity is one of the conditions 

required in CS approach, and can be satisfied in different domains: time, frequency or 

time-frequency domains [28]-[48]. This condition is valid for a variety of real-world 

signals. The other condition is incoherence which will be explained later in the text.  

There is a large number of CS applications, from those assuming one-dimensional 

signals to various image processing and video applications. Some of the CS applications 

are adopted to work in  real-time. The constant growth and development in the field of the 

CS applications aims to reduce the complexity of devices, to speed up the acquisition and 

transmission procedure and to decrease power consumption. Since CS is used to extract 

as much as possible information from minimal available data, it is important to highlight 

its use in biomedicine, especially in Magnetic Resonance Imaging – MRI. By lowering the 

number of coefficients required for MR image reconstruction, the time of patient exposure to 

the MR device is reduced, and consequently, the negative influence of the MR device is lower. 

Another useful application is in radar imaging, where the CS exploits the sparsity in the 

frequency domain.  Furthermore, it is used in communication and network systems, sparse 

channel estimation, wireless sensor networks (WSNs), in cognitive radios for spectrum sensing, 

etc. Some of the applications will be described later in the text. 

2. MATHEMATICAL BACKGROUND OF A COMPRESSIVE SENSING CONCEPT 

Assume that we are dealing with the signal x of length N, that is sparse in the 

transform domain (defined by the direct transformation  ). Then, the vector of acquired 

samples y can be defined as [1]-[12]: 

 1 ,y X
   (1) 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 479 

where matrix Ω is used to randomly under-sample the observed signal and 1 is the 

inverse transform matrix. The signal sparsity is K, meaning that only K out of N 

coefficients from the transform domain are non-zero and we assume that only M out of N 

samples are acquired in the vector y (M<<N and M>2K). The vector X is the vector of the 

transform domain coefficients, i.e.: 

 .X x    (2) 

Different transform domains can be used: discrete Fourier transform domain – DFT 

[3],[6],[12], discrete cosine transform domain – DCT [3],[6],[12],[46], wavelet domain, 

Hermite transform domain [48]-[55], time-frequency domain [56], etc. Apart from the 

sparsity, another important property is incoherence, which enables successful signal 

reconstruction from a small set of acquired samples. Namely, the measurement matrix Ω 

should be incoherent with the transform domain matrix  . The coherence between the 

two matrices represents the highest correlation between any two column/row vectors of 

the matrices. A measure of correlation between the two matrices is defined as follows 

[11],[12]: 

 
1,

( ) max , ,k j
k j N

N  
 

   (3) 

where N is a signal length, Ωk  and j  are row and column vectors of the matrices Ω and 

 , respectively. The coherence takes values from the interval: 

 1 ( ) .N    (4) 

The value of the coherence is greater if the two matrices are more correlated. In the CS 

scenario the value of the coherence should be as low as possible.  

The system of equations (1) can be written as follows: 

 1 1,M M N Ny A X     (5) 

where A denotes CS matrix: 1A
 . The system is under-determined since it has M 

equations and N unknowns. Therefore, the optimization techniques are used in order to 

find an optimal solution for this system. Optimal solution is related to the signal sparsity – 

the sparsest solution is the optimal one.  

There are a number of algorithms used to obtain a sparse solution of the system. 

Some of them are based on the convex optimization [1]-[7],[11],[12]: basis pursuit, 

Dantzig selector, and gradient-based algorithms. They provide high reconstruction accuracy, 

but they are computationally demanding. The commonly used and less computationally 

demanding are greedy algorithms – Matching Pursuit and Orthogonal Matching Pursuit 

[1]-[6],[8],[12]. Also, recently proposed threshold based algorithms provide high reconstruction 

accuracy with low computational complexity: (e.g. Iterative Hard Thresholding - IHT, 

Iterative Soft Thresholding - IST) [1],[6],[13], Automated Threshold Based Iterative 

Solution [12],[57], Adaptive Gradient-Based Algorithm, [3],[12],[26],[27], etc. 



480 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

3. COMPRESSIVE SENSING ALGORITHMS 

The sparsity of the signal can be defined as a number of nonzero elements within a 

vector. It can be described by using the ℓ0-norm [13]: 

 
0 0

1 1, 0

lim 1

i

N N
p

i
p

i i x

x K


  

   x

,  
(6)

 

and represents the cardinality of the support of x: 

 0
card{supp( )} K x x

.  (7) 

Therefore, the solution of the undetermined system of equations (5), in the cases when the 

signal x is sparse in the transform domain, can be reduced to the minimization of the ℓ0-

norm, i.e.: 

 
0

min subject to X y AX . (8) 

The ℓ0-norm is not feasible in practice, since small noise in the signal will be assumed 

as a non-zero sample. Therefore, the ℓ1-norm is commonly used. The optimization 

problem based on the ℓ1-norm is recast as follows [12],[13]:  

 
1

min subject to X y AX . (9) 

In the sequel, some of the commonly used algorithms for sparse reconstruction are described. 

3.1. Convex optimizations 

Basis Pursuit and Basis Pursuit Denoising 

The equation (8) represents a non-convex combinatorial optimization problem. 

Solution of this problem requires exhaustive searches over subsets of columns of the 

matrix A. For a K-sparse signal of length N, the total number of K-position subsets is 
N

K

 
 
  , 

which is not computationally feasible. Other approach solves a convex optimization 

problem through linear programming, which is computationally more efficient. Commonly 

used convex optimization algorithms are: Basis Pursuit, Basis Pursuit De-Noising (BPDN), 

Least Absolute Shrinkage and Selection Operator (LASSO), Least Angle Regression 

(LARS), etc.  

The approach based on the convex ℓ1-minimization that provides near optimal 

solution, can be defined as:  

 
1

min subject to X y AX . (10) 

This approach is known as a Basis Pursuit (BP) [6],[12]. It aims at decomposing a signal 

into a superposition of dictionary elements that have the smallest ℓ1-norm of the 

coefficients. BP can be solved by using a primal-dual interior point method. The problem 

(10) can be recast as follows, in the case of real y, A and X [12]: 

 min , subject to ,t t t tX y AX    ,  (11) 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 481 

where variable t is introduced to avoid absolute value in 1 1

N

ii
X X . The steps of 

the primal-dual interior point method are described within the Algorithm 1. In the cases of 

noisy measurements, y=AX+n, where n denotes noise and 2
n  , the optimization 

problem is known as Basis Pusuit Denoising (BPD) and is defined as [6]: 

 1 2
min subject to - X y AX 

 (12) 

Algorithm 1: Primal-dual interior point method 

 Set 0
T

X X A y  , for the known measurement vector y. 

 Set  max0t   X X0 0 . Parameters γ and λ are user-defined. 

 The next step is forming a Lagrangian function:  

1 1
, , , , ( ) ( ) , where

0 0 0 0 0 0 0 0

t t
t g f t g

t t t t

X X
X AX y

X X X X

   
       

    

1 1

0 0 0 0
g

t t
A

X X

 
   

    .  

 Update each argument of the Lagrangian function by step direction (Δ) and step length (u).  

Step directions for the Algorithm 1 are obtained by finding the first derivatives of Λ in 

terms of its arguments. Step lengths are calculated using the backtracking line search [12]. 

For example, a new value for X is obtained as X=X+u(ΔX). 

Adaptive gradient based algorithm 

Adaptive gradient based algorithm proposed in [26], belongs to the group of convex 

optimization approaches. It starts from the chosen initial values of the available signal 

samples. The initial value is iteratively changed for +Δ and –Δ, and the concentration 

improving is measured in the sparsity domain. The gradient vector, used to update the signal 

values, is obtained as a difference between the ℓ1-norms of the vectors changed for +Δ and 

changed for –Δ. This gradient value is used to update the values of the missing samples.  

The performance of this algorithm can be efficient even for the signals that are not 

strictly sparse. The algorithm for both 1D and 2D cases is summarized in the Algorithm 2. 

3.2. Greedy algorithms 

The greedy algorithms represent the second group of algorithms used to obtain the 

sparsest solution of the system (5). These algorithms are less computationally complex 

and therefore much faster compared to the ℓ1-norm based optimization techniques, but are 

also less precise. The greedy algorithms are based on finding the elements of the 

transform matrix called dictionary that best matches the signal through iterations. 

Commonly used greedy algorithms are Matching Pursuit (MP), Orthogonal Matching 

Pursuit (OMP), Compressive Sampling Matching Pursuit (CoSaMP), etc.  

The procedure for the OMP algorithm is described within the Algorithm 3. 



482 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

Algorithm 2: Adaptive gradient based algorithm  

Input: set of the positions of the available samples Ωa and set of the missing samples position: 

Ωm=N\Ωa; measurement vector y; in the 2D signal case n is: n=(nx,ny)  

 Set (0)
( ), for 

( )
0, for 

a

m

n n
n

n


 



y
y ,  and 0k   

 Set (0)max ( )n y  

 repeat 

   set ( )( ) ( )k
p n ny y  

     repeat 

      1k k   
        for tN  do 

          if mt  then 

                         
 

 

( )
1 2

( )

( ) { ( ) }, in the 2D  case ( , )

( ) { ( ) }, denotes 1D DFT  or  2D DFT

k

k

X f n f f f

X f n





   

   

y

y
;                  

                            
( )

1 1
( )k t    X X  , 

             else  

                   
( )( ) 0k t 

 
           end if 

                
( 1) ( ) ( )( ) ( ) ( )k k kt t t  y y

 
       end for 

              

1

2 2
1

2 2

arccos

k k

k
k k







 


 
  

       until 170k
   

      / 10  

       2 2
( ) ( )

1010log ( ) ( ) / ( )
m m

k k
pn n

R n n n
 

  y y y  until max required precisionR R    

 return ( ) ( )k ny  

Output: reconstructed signal ( ) ( )k ny  

Algorithm 3: Orthogonal matching pursuit 

 Input: Compressive sensing matrix A=Ωℑ , measurement vector y 

 Initialization of the variables:  

 initial residual r0=y; initial solution x0=0; matrix of chosen atoms ϒ0=[]. 

 Do following steps until the stopping criterion is met: 

 arg max ,1
1,...,

n n i
i M

r A  


                - finding maximum correlation column 

 1 nn n  A   
                       - update matrix of chosen atoms 

  
2

arg min 1 1 2n n n nx r x

x

        - solving least square problem 

 1n n n nr r x                          - residual update 

  n=n+1 

 Output: XP and rP, where P denotes number of iterations. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 483 

3.3. Threshold based algorithms 

Iterative hard and soft thresholding  

Thresholding algorithms are based on an adaptive threshold applied within several 

iterations. They are much faster than algorithms based on convex relaxation. An iteration 

can be described in terms of threshold function as [1],[6],[11],[13]: 

 1T ( ( ))i ix f X  .  (13) 

The thresholding function is denoted as Tε, while f is the function that modifies the output 

of the previous iterate and X is a sparse vector. The signal can be recovered from its 

measurement by using hard or soft thresholding. Therefore, there are two types of 

iterative thresholding algorithms: iterative hard thresholding (IHT) and iterative soft 

thresholding (IST).  

 IHT algorithm sets all but the K largest components, in terms of signal X magnitudes, 

to zero. The hard thresholding function HK is defined as [6],[11],[13]: 

 
,

( )
0, otherwise

i i
K

X X
H

 
 


X


 . (14) 

The ε is the K largest component of X [6]. The algorithm is summarized within the 

Algorithm 4 [11]. Soft thresholding function is applied to each element of the vector X 

and is defined as [6]: 

 

,

( ) 0, .

,

i

i i

i

X

S X X

X



 



 

 


 


  

X

X

  (15) 

Algorithm 4: Iterative hard thresholding 

Input: signal sparsity K, transform matrix ℑ , measurement matrix Ω, CS matrix 1 A  , 

measurement vector y 

Output: an approximation of the signal X  

             0 0X   

for i=1,…, until stopping criterion is met do 

              1 1( )T
i K i iH    X X A y AX  

end for 

return iX X  

Automated threshold based solution 

A non-iterative and iterative threshold based solutions for sparse signal reconstruction 

are proposed in [98].  

The proposed solutions are based on the model of noise appearing as a consequence 

of missing samples. By using a predefined probability of error P, a general threshold T 

that separates signal components from spectral noise in the transform domain is defined.  



484 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

Algorithm 5: Automated threshold based iterative solution 

Input: M, N,  y=x(Ωa), Ωa={n1,…,nM},  , Φ,  A , 
N

 .    

o Set k = ; 

for i=1 : i=i+1: until all components are detected 

o Calculate variance: 
2

2

1

( )

1

M

i

N M y i
M

N M








 ;  

o For a given P calculate:  2 1/log 1 ( ) NT P T    ; 

o Calculate the initial DFT vector Xi:  
1( )i
 X y  ; 

o Update set k: arg{ / }i T N  k k X  ; 

o Calculate  
1

H


F = A A Ay ;                                                 (CS matrix A contains rows defined 

by  the  

                                                                                                  set k, and M columns of the DFT 

matrix) 

o Update y:  

2

: = ( ) - ( )
ak

j
N

ak X k e

 

  k y x  ; 

o Update the initial DFT vector X according to the new vector y; 

o Update 
22 /A M y   and 2 2

1

N M
A M

N






  ; 

o If 2 2

N
A    break ; 

end for 

The algorithm uses DFT as a domain of sparsity but the same concept can be applied to 

other transform domains. This approach can provide successful signal x reconstruction 

within a single iteration of the reconstruction algorithm. However, if the number of 

available samples M is very low, the iterative version of the algorithm is derived as well, 

updating the threshold value. If the inputs of the algorithm are vector of the M available 

samples y, signal length N, set of the available samples positions Ωa={n1,…,nM}, 

transform and measurement matrices   and Φ, CS matrix  A  and Gaussian noise 

variance ζN, then the iterative version can be described using the Algorithm 5. 

4. CS APPLICATIONS 

 The CS theory stating that the compressible signals can be efficiently reconstructed 

using a small set of incoherent measurements, motivated the researchers to explore 

possible fields of applications. Having in mind that many real-world signals satisfy the 

sparsity property, the applications ranges from the speech and audio signals [58]-[61], 

radar and communications [62]-[83], underwater, acoustic and linear frequency 

modulated signals [84]-[86], image reconstruction [87]-[98], biomedical applications 

[98]-[101], etc. The review of CS applications for different 1D signals, images and video 

data will be addressed in the sequel. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 485 

4.1. CS devices (analog to information, single pixel camera, random lens imager) 

Let us firstly consider some of the hardware devices that are based on the CS principles. 

(A) Duarte et al. in [102] proposed single pixel camera concept. This CS camera 

architecture is an optical computer, composed of a Digital Micromirror Device - DMD, 

two lenses and a single photon detector. It also contains an analog-to-digital (A/D) 

converter that computes random linear measurements of the scene under view. The image 

is recovered from the acquired measurements by a digital computer. Compared to the 

conventional silicon-based cameras, single pixel camera is a simpler, smaller, and cheaper 

and can operate efficiently across a much broader spectral range.  

(B) Fergus et al. in [103] developed a random lens imaging technique. The technique uses a 

normal Digital single-lens reflex - DSLR camera, whose lens is replaced with a transparent 

material. The mirrors in this material are randomly distributed. The authors modified a Pentax 

stereo adapter in order to make one of the mirrors have a random reflective surface. The CS 

measurements are in the form of images, obtained using this system. The new camera set-up has 

to be calibrated, in order to reconstruct the original image. 

(C) Trakimas et al. in [104] proposed the design and implementation of an analog-to-

information converter (AIC). The presented AIC is designed in a way that can sample at the 

Nyquist rate, but has also the CS operation mode. This design shows minimal complexity 

compared to conventional Nyquist rate sampling architectures. When dealing with signals 

with sparse frequency representation, this design has increased power efficiency of the 

sampling operation. To generate the pseudorandom sequence, a PN clock generator is used 

and it can be configured to provide a synchronous clock signal when Nyquist sampling is 

required. 

4.2. CS in biomedical applications 

CS finds usage in numerous biomedical applications, such as in Magnetic Resonance 

Imaging (MRI) [105]-[108], then electroencephalography (EEG), electrocardiography 

(ECG), electrooculography (EOG) and electromyography (EMG) signals, [109]-[120], 

etc. Some of the specific biomedical applications are given in the sequel. 

(A) Lowering the time of patient exposition to the harmful MR waves was the primary 

motivation of CS usage in MRI. However, the MR acquisition time is proportional to the 

dimensionality of the MR dataset, i.e., the number of spatial frequencies acquired. Scan 

time can be reduced by lowering the amount of data acquired, but still it has to be able to 

recover the whole information.  

(B)  Lustig et al. in [98] implemented CS approach for rapid MRI imaging. Sparsity of 

the MRI in the transform domain is exploited for achieving two goals: reducing the scan 

time and improving the resolution of the observed fast spin-echo brain images and 3D 

contrast enhanced angiographs. The non-linear conjugate gradient solution is used for the 

optimization problem solving. The problem is in the form: 

 
2

12
arg min ,u    

x

x y x   (16) 

where x is the image of interest,   is an operator that transforms signal from pixel 

representation into sparse representation, u  is an undersampled Fourier transform, y 



486 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

denotes measured k-space (e.g. frequency space) data from the scanner and   is a 

regularization parameter. The conjugate gradient procedure is described in detail in [98]. 

(C) Bioucas-Dias et al. introduced TwIST: Two-Step Iterative Shrinkage/ Thresholding 

Algorithms for Image Restoration [101]. This algorithm is introduced as an improved 

version of the Iterative Shrinkage Thresholding Algorithm (IST) – to overcome the 

problem of its slow convergence in the cases when the measurement matrix A is ill-posed 

or ill-conditioned. The optimization problem is well-posed if ? has a solution, if a solution 

is unique and the solution changes continuously on the data [121]. Otherwise, the problem 

is ill-posed. If small perturbation of the y in the problem y=AX leads to large perturbation 

of the solution, the problem is ill-conditioned [121].  

The algorithm is successfully applied on image deconvolution problems, as well as 

reconstruction of the images with missing samples. Considering the system of equations 

y=AX, the t-th iteration of the TwIST algorithm can be defined as follows: 

 
1 0

1 1

( ),

(1 ) ( ) ( ),t t t t

G

G

X X

X X X X



    



    
  (17) 

where μ and δ are nonzero parameters. The starting value for the vector X, X0, can be 

user-defined or X0=A
-1

y. Function Gη is defined by using denoising operator Ψη as: 

 ( ) ( ( )),TG    X X A y AX   (18) 

where 
21( ) arg min ( ) / 2regv

X

X y AX       and Φreg(X) is a regularization function. 

 The application of the TwIST algorithm in MRI reconstruction is shown. Fig. 1 shows 

an example of MRI reconstruction when only 2% of the image samples are available. The 

samples are acquired from the 2D DFT domain using a mask. The mask is formed of 

radial lines and placed around the origin. The TV regularization is done according to 

[101]. The original image, mask and the reconstructed image are shown in Fig. 1. 

  Original

  

Estimate

100 200 300 400 500

100

200

300

400

500  
a) b) c) 

Fig. 1 a) Original image; b) Mask in the 2D DFT domain; c) Image reconstructed from 

available samples (2% of the total number of samples) 

(D) Trzasko and Manduca in [122] proposed a method for under-sampled MR images 

recovering by using homotopic approximation of the ℓ0-norm. It is shown that the 

computed local minima of the homotopic ℓ0-minimization problem allows very highly 

undersampled K-space image reconstruction. The optimization problem can be defined 

starting from the relation: 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 487 

 
0

arg min subject to = ,
u

u u u f      (19) 

where Ψ is wavelet, curvelet, etc. operator, Φ is Fourier sampling operator and f is the 

continuous signal. If the ℓ0 semi norm is replaced with the ℓ1 norm, as proposed by 

Candes and Donoho [122], the optimization problem can be recast as: 

 
2

1 2
arg min subject to ,n

u

u u u f        (20) 

where measured data fn is noisy and ε denotes the statistic of the noise process. Chartrand 

[122] proposed an alternative to the ℓ0 semi norm that provides better sampling bounds 

compared to the ℓ1 and that is computationally feasible. He proposed the usage of the ℓp 

semi norms (0<p<1].  

The zero semi-norm of the signal can be defined as: 

  0
( ) 0 ,u u n


  1   (21) 

with Ω that denotes the image domain and 1 is the indicator function. Any semimetric 

functional ρ that satisfies following relation [122]: 

    
0

lim ( ) , ( ) 0u n u n


 
 

   1  , (22) 

and if the σ is sufficiently small, can be used as a sparsity prior. Based on the previous 

relations, a new reconstruction paradigm can be defined as [122]: 

  
2

20
min lim ( ) , subject to n

u
u n u f


  


     , (23) 

where Ψ is sparsifying basis and Ф is the CS matrix. The class of functionals that satisfy 

(22) are homotopic with ℓ0 semi norm. The proposed method is tested on the MR images, 

but the application to other medical images, such as the x-ray CT will be investigated in 

the future.  

(E) Abdulghani et al. in [111] analyzed performance of different reconstruction 

algorithms, applied for reconstruction of the EEG signals. They show that the best 

reconstruction results, with minimal error, are obtained by using the BP algorithm. This 

approach provides better reconstruction quality compared to the greedy approaches, but 

has greater computational cost. Six different EEG dictionaries have been observed and it 

is proved that the B-Spline dictionaries are the most suitable for CS of the EEG signals. 

(F) The approach proposed in [118]-[120] is used for the reconstruction of the under-

sampled ECG signals, more precisely, the QRS complexes within the ECG signals. It is 

shown that the reconstruction of ECG signals can be done using just a few coefficients 

from the Hermite transform domain. The reconstruction is done by using gradient-based 

algorithm (Algorithm 2), where DFT is replaced by the Hermite transform.  

 A continuous signal x(t) can be represented in terms of Hermite functions as follows: 

 
0

( ) ( )p pp
x t c t




   , (24) 

where cp denotes Hermite expansion coefficients, while ( )p t  denotes the p-th Hermite 

function. Signal x(t) needs to be sampled at the non-uniform points that correspond to the 



488 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

roots of the Hermite polynomial or, uniformly sampled signal should be interpolated to 

obtain the requested signal values. If c and x are Hermite coefficients and signal vector, 

respectively, and H denotes the transform matrix, then the signal expansion in the matrix 

form can be written as c=Hx: 

 

   

   

   

2 2
0 1 1 1 0 1

0 1
2 2

1 21 1 1 1 1 1

2 21
1 1 1 1 1 1

( ) / ( ) ... ( ) / ( ) ( )

( )1 ( ) / ( ) ... ( ) / ( )

... ...... ... ...
( )

( ) / ( ) ... ( ) / ( )

M M M M

M M M M

P M
M M M M M M

t t t tc x t

c x tt t t t

M

c x t
t t t t

 

 


   

        
       
     
    
            

  (25) 

To obtain the values at the non-uniform points, the signal sampled according to the 

sampling theorem is interpolated by using the sinc interpolation formula: 

 
 sin ( ) /

( ) ( ) ,
( ) /

K
m

m

n K m

t nT T
x t x nT

t nT T

 


 





   (26) 

where m = 1,…, M, and T is the sampling period. The time-axis scaling parameter λ is 

introduced instead of ζ, in order to avoid stretching and compressing of the basis 

functions. The parameter ζ is fixed (ζ=1) and an optimal value of the parameter λ 

(producing the best possible concentration - sparsity) is found according to the l1-norm 

optimization: 

  
1

1 0 1
min min min HT ( )

M

opt p mp
c x tc

  
 




    , (27) 

where the operator HT{.} is used to denote the Hermite transform. The results obtained 

by applying the proposed procedure are shown in Fig. 2. The observed QRS complex is 

reconstructed by using only 55% of the available samples, with MSE of order 10
-3

.  
 

10 20 30 40 50

-200

-100

0

100

200

 10 20 30 40 50
-0.6

-0.5

-0.4

-0.3

-0.2

-0.1

0

0.1

 
a) b) 

Fig. 2 Reconstruction results for QRS complex: a) desired (solid line) and reconstructed 

signal (dashed line), b) Hermite coefficients of desired (blue) and reconstructed 

signal (green)   

4.3. Applications in communications and radar signal processing 

The CS applications in radars [62]-[72] and communications [123]-[126] are widely 

studied in the literature. Hence, we describe some of these applications in wireless, ultra 

wideband (UWB) and Inverse synthetic aperture radar (ISAR) systems. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 489 

(A) Zhang et al. in [123] proposed CS application in UWB communications. CS is used 

to reduce the high data-rate of ADC at receiver. The receiver design is simpler with only 

one low-rate A/D. Here, the CS exploits the time sparsity of the signal through a filter-

based CS approach applied on continuous time signals. 

(B) Weiss in [125] used the CS in distributed radar network, composed of Wireless 

Local Area Network - WLAN routers as transmitters. The corresponding receivers are 

widely separated in the context of sparse modeling. The CS is applied in order to reduce 

the number of samples transferred to the central processing stage and it is used for the 

estimation of multiple targets positions and velocities. The paper uses IEEE 802.11b 

signal. The proposed approach is compared with the matched filter technique that is the 

traditional approach for target determination.  

(C) Particular interest of using the CS approach in communications raised for signals 

belonging to the two interfering standards in wireless communications are: Bluetooth and 

IEEE 802.11b that share Industrial, Scientific and Medical (ISM) frequency band 

[75],[82]. The first standard uses frequency hopping spread spectrum modulated signals 

(FHSS), while the second standard applies direct sequence spread spectrum (DSSS) 

modulation. Different sparsification basis are considered for different signals. Original 

signal, consisting of FHSS and DSSS components, does not satisfy sparsity property, 

which is not the case with the separated components. Therefore, the first step is 

components separation, done by using an eigenvalue decomposition method. The EVD of 

the covariance matrix C is defined as: 

 *1
1 ( ) ( )i i i

N
i u n u nT

C = U U 
  ,  (28) 

where U is an eigenvectors matrix,   is a diagonal eigenvalues matrix where eigenvalues 

are sorted in decreasing order, λi are eigenvalues and ui are eigenvectors. The covariance 

matrix is defined based on the TF representation – the inverse form of the S-method is 

used for definition of the covariance matrix: 

 

4/2
1

/2

1
( , ) ( , )

1

N j mk
N

K
k N

C n m n m SM n k e
N







  


  . (29) 

The EVD is applied to the matrix CK (where K denotes the number of signal components) 

according to the relation (28), resulting in the eigenvectors that correspond to the signal 

components. The choice of the sparsity domain is made between the DFT and HT domain, 

by measuring the concentration using the ℓ1-norm and choosing the one producing the 

minimal concentration. It is shown that the HT is chosen as a sparsifying basis for the 

FHSS, while the DFT domain is used for the DSSS components. In the case of HT, the 

eigenvectors should be resampled at non-uniform points being proportional to the roots of 

the Hermite polynomial [127].  

Signal having 3 FHSS and 4 DSSS components is observed. The first 3 eigenvectors 

correspond to the FHSS components and therefore the optimized HT is calculated (Fig. 

3). It can be seen that the optimized HT provides more compact representation compared 

to the standard HT. The rest of the eigenvectors correspond to the DSSS components and 

the DFT is used as a sparsification basis in this case.  

 



490 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

 

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Standard HT, component 1

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Optimized HT, component 1

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Standard HT, component 2

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Optimized HT, component 2

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Standard HT, component 3

20 40 60 80 100 120
-0.2

-0.1

0

0.1

0.2

Optimized HT, component 3

 

Fig. 3 HT and the optimized HT of the first 3 eigenvectors  

that correspond to the FHSS signal components 

50 100 150 200 250
-0.5

0

0.5

 

50 100 150 200 250
-0.5

0

0.5

 

50 100 150 200 250
-0.5

0

0.5

 
 a) 

50 100 150 200 250
-0.2

0

0.2

 

50 100 150 200 250
-0.2

0

0.2

 

50 100 150 200 250
-0.2

0

0.2

 

50 100 150 200 250
-0.2

0

0.2

 
b) 

Fig. 4 Separated and reconstructed components of the a) FHSS and b) IEEE 802.11b 

signals. Blue is the original component, brown is the CS reconstructed component. 

The percent of the randomly selected samples is 50% from each eigenvector 

After sparsification domain is chosen, the eigenvectors are randomly under-sampled 

and reconstructed from only 50% of the acquired samples. Under-sampled vectors are 

reconstructed by using the ℓ1-norm minimization according to (9). 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 491 

The reconstructed components from the 50 % of measurements are shown in Fig. 4 (blue 

denotes the original while brown denotes the reconstructed component). 

(D) Lu et al. in [128] proposed a novel Distributed Streaming Compressive Spectrum 

Sensing (DSCSS) algorithm. This approach is used for wideband spectrum sensing under 

decentralized cognitive radio network (CRN) scenario. The algorithm has low computational 

complexity that makes it suitable for the on-line applications. DSCSS uses an analog-to-

information (AIC) convertor in acquisition process. An AIC showed to be suitable for 

streaming framework. A sliding time-window of length equal to NT is applied to the sensed 

signal, while M sub-Nyquist samples are produced according to one time-windowed signal.  

DSCSS approach does not require sparsity of the signal to be a priori known. In practical 

applications, sparsity of the wide-band spectrum is not available. Therefore, DSCSS 

estimates sparsity and the support set of the spectrum. This support set is then exchanged 

within the cognitive radio network as an a priori information, in order to obtain a cooperative 

sensing gain. This is done with a goal to overcome wireless fading effects.  

(E) The CS application in recovery of the ISAR signal rigid body, by exploiting the 

concept of sparsity, is proposed in [58], [72]. The ISAR signal contains both the rigid 

body and micro-Doppler segment caused by fast-moving parts of a target. The observed 

segments of the signal partially overlap in the time-frequency plane. The method finds 

and removes the overlapping values and recovers a rigid body signal. The separation of 

the rigid body stationary and the micro-Doppler nonstationary signal segments is done by 

using the sorted Short-Time Fourier transform (STFT) values along the time axis. Then a 

certain percent of the strongest STFT values, for each frequency, is removed, which 

results in reduction (or elimination) of the micro-Doppler nonstationary components. The 

remaining part of the STFT will contain the rigid body components. 

A radar return signal after coherent processing and filtering  consists  of two 

components, rigid body xr(t) and micro-Doppler xm(t): 

 

2

1

( ) ( ) ( ), ( )

f toiK j
N

r m r i
i

x t x t x t x t r e





   , (30) 

where ri denotes components’ amplitudes and foi are components’ frequencies. The STFT 

is used for the TF representation of x(t): 

 

21

0

( , ) ( ) , or ( ) ( ),

lfL j
L

L L

l

S t f x t l e t t

 



   S E x   (31) 

where 

 2

( ) [ ( ,0),..., ( , 1)] , ( ) = [ ( ), ( 1),..., ( 1)] ,

( , ) ,

T T
L

fl
j

L

t S t S t L t x t x t x t M

E l f e




    



S x

  (32) 

and EL is the L×L DFT matrix. Non-overlapping STFT is used in TF analysis, calculated 

with step L in time t, and we can write the following relations: 

 ,[ (0) , ( ) ,..., ( ) ] .T T T T
L L L L NM N - M S = S S S E x   (33) 

Note that the S(t,f) denotes scalar STFT value at time t and frequency f, SL(t) denotes 

STFT vector that contains L frequencies at instant t, while S is a vector of all STFT values 



492 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

for all frequencies f and all time instants t. The N×N matrix EL,N  and the signal vector x 

are formed as: 

 
,

...

...
,

... ... ... ...

...

[ (0) , ( ) ,..., ( - ) ] = [ (0), (1),..., ( -1)] .

M M M

M M M
L N

M M M

T T T TM N M x x x N

E 0 0

0 E 0
E

0 0 E

x =

 
 
 
 
 
 

x x x

  (34) 

Vector x can be written using the inverse DFT matrix 1
N


E  and DFT vector X, and 

following relations are obtained: 

 1 1
, ,,N L N L N N

   x E X S = E x E E X = X ,  (35) 

where    is the N×N transformation matrix. After calculating the STFT according to (35), 

the sorting operation along is performed along the frequency (for each f): 

  ( , ) ( , ) ,sorts t f sort S t f   (36) 

where t=0, …., L-1. The sorting operation is performed in ascending order and certain 

percentage P of low value coefficients and percentage Q of high value coefficients are 

removed from the sorted STFT: 

 ( ) { ( , ), , 1,..., }cs SORTs f S n f n P P L Q      (37) 

where scs(f) denotes the vector of available STFT coefficients at frequency f, while SCS 

denotes  vector for all frequencies.  

frequency

ti
m

e

 

20 40 60

10

20

30

40

50

60

frequency

ti
m

e

 

20 40 60

10

20

30

40

50

60

frequency

ti
m

e

 

20 40 60

10

20

30

40

50

60

 
a) b) c) 

200 400 600 800 1000 1200

50

100

150

200 400 600 800 1000 1200

50

100

150

 
d) 

Fig. 5 Real radar signal: a) STFT, b) sorted STFT, c) STFT that remains after discarding 

certain region from the sorted STFT, d) the original DFT transform - blue, and the 

reconstructed DFT 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 493 

The vector of all available STFT points is SCS=AX, where the matrix A
 
is formed by 

omitting the rows which corresponds to the removed STFT values. The CS problem can 

be formulated as follows: 

 
1

min subject to CSX S AX  . (38) 

The results of the procedure applied on real radar signal, are shown in Fig. 5. The signal 

consists of the rigid body and three corner reflectors rotating at ∼ 60 RPM. The 80% of 

the samples from the sorted STFT are omitted and rigid body is reconstructed from the 

rest of the samples. 

(F) The same approach is used in [75] for separation of signals belonging to two 

different interfering wireless standards - Bluetooth and IEEE 802.11b standard. The 

proposed procedure works even in the case of overlapping signal components. Based on 

the a priori knowledge on the signals’ nature, it is possible to select and extract a small set 

of time-frequency points that entirely belong to the IEEE 802.11b signal. The extracted 

points are used to recover the full signal by using the described CS-based approach.  

     frequency

ti
m

e

 

 frequency

ti
m

e

 

 frequency

ti
m

e

 

 
a) b) c) 

Fig. 6 a) STFT of the original signal, b) sorted STFT, c) remaining STFT values  

after discarding certain region from the sorted STFT 

Once the components of IEEE 802.11b signal are extracted, the remaining components in 

the time-frequency plane belong to the Bluetooth signal. The example of signal consisting  

of Bluetooth and IEEE 802.11b components that overlap is shown in Fig. 6. Six components 

belong to the IEEE 802.11b signal, while twelve components (of short duration) belong to 

the Bluetooth signal. The STFT is sorted and a certain percent of the lowest and highest 

energy samples are removed (sorted and remaining STFT are shown in Fig. 6b and 5c). The 

recovered STFT is shown in Fig. 7a, and it corresponds to the IEEE 802.11b signal, while 

the remaining part of the TF coefficients correspond to the Bluetooth signal.  

frequency

ti
m

e

 

 frequency

ti
m

e

 

 
a) b) 

Fig. 7 Reconstructed STFT of: a) IEEE 802.11b, b) Bluetooth signal 

(G) Zhang et al. in [76] presented a high-resolution inverse synthetic aperture radar 

imaging and showed that the proposed CS imaging outperforms the conventional range-



494 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

Doppler approach, regarding the image resolution. More details regarding the proposed 

procedure follow.  

If the signal of interest is defined as: 

   2

1
( ) / ( )

K j f tk
k ak

x t E t T e t
 


  ,  (39) 

where K denotes the number of the strongest scattering centers, γ(t) is the synthetic 

additive noise in the range cell, Ek is the scattering amplitude, ρ(t/Ta) is the unit 

rectangular function, fk is the carrier frequency and Ta is the time width of chirp pulse. If 

we denote the transform matrix as  , then the previous relation can be rewritten as: 

 ( ) ( ),x t t      (40) 

where θ is the vector whose non-zero components correspond to the complex amplitudes 

of the K strongest scattering centers. The optimization problem is recast to: 

 
1 2

min subject to ( )x t    ,  (41) 

where ξ is the noise level, and measurement vector y is: y=Ωx with ΩM×N (M<N) 

measurement matrix. Here, the assumption is that the precise motion compensation is 

done. It is not always the case in large-size maneuvering targets, which is one of the 

drawbacks of the proposed method.  

4.4. Compressive sensing image reconstruction 

(A) Musić et al. in [87] proposed a method that combines CS based image reconstruction and 

object detection algorithm. The method is used in search and rescue application and allows 

image recovering if up to 80% of pixels is unavailable or corrupted by noise. Preserving 

good quality of the reconstructed image is of particular importance for object detection 

algorithm. The solution is based on the 2D DCT domain and gradient descent recovery 

method.  

Through the iterations, a two version of image I are formed: I
+
(n,m) and I

-
(n,m), where 

(n,m) denotes pixels positions: 

 
( , ) ( , ) ( , )

( , ) ( , ) ( , )

i k k

i k k

I n m I n m n n m m

I n m I n m n n m m









    

    
 , (42) 

where the correction factor Δ is Δ=mean{I(n,m)} and δ is 2D discrete delta function. The 

next step is gradient vector calculation: 

 
1 1

1
( ) ( ) , for ( , )

( , )

0, otherwise

k k
i k k

I I n m
G n m N

 
   

 



 , (43) 

where N is the total number of pixels,  denotes 2D DCT and Ɵ is the set of positions 

corresponding to the missing pixels: Ɵ={(nk,mk), k=1,…,NM}, NM-number of missing pixels. 

In each iteration, the image is updated: 1( , ) ( , ) ( , )i i iI n m I n m G n m   . By minimizing the 

ℓ1-norm of the sparsity measure the missing samples are changed toward the exact values.  

Finding a potential targets on the CS-reconstructed image is done by observing parts 

of the image. The image preprocessing assumes conversion from the RGB to YCbCr color 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 495 

model.  The image is firstly divided into the set of C non-overlapping clusters denoted as 

K, using means shift clustering algorithm. The clusters correspond to significant image 

features, i.e. dominant colors. If an image I consists of set of clusters K, then: 

 1
C
c jI K   . (44) 

Clusters are further divided into the sets of segments Sj,k. One segment represents a 

spatially connected component region. If the total number of segments in a cluster is 

denoted as Q, then the cluster can be defined as: 

 ,1
Q

j j kkK S   . (45) 

   
a) b) c) 

Fig. 8 Reconstruction results for search-and-rescue image. 70% of pixels at random positions 

are missing: a) Original image; b) Image with missing samples; c) Reconstructed image  

The segment Sj,k will be considered as a target if the following conditions are satisfied: 

 1 , 2( ) ( )j j k As K T s S T Q N      , (46) 

where s denotes size of the set, T1 and T2 are threshold values and NA is the maximum 

allowed number of candidate segments in a cluster.  

The algorithm’s performance in case of 70% of missing samples is shown in Fig. 8. 

The reconstruction is done on 16 × 16 blocks in 50 iterations. Green squares represent 

correct detections (objects found in both, original and reconstructed images), red squares 

represent object found on the original, but not found on the reconstructed image and 

orange squares represent object found on the reconstructed but not found on the original 

image. 

(B) An application of the adaptive threshold based algorithm for the ISAR images 

reconstruction is done within the paper [96]. The algorithm, initially proposed in [98] for 

1D data, is adapted for the 2D data and exhibits sparsity in 2D DFT domain. It is based 

on the analytically derived threshold that precisely separates signal and non-signal 

components. The approach provides efficient results even if less than 10% of the samples 

are available.  

If we denote the full and partial sets of signal samples as SF and SP: 

 
{ ( , ) : {1,..., }, {1,..., }},

{ ( , ) : {1,..., }, {1,..., }}, ,

F

P

S x y x I y J

S x y x M y L

  

   

F

P P F

S

S S S
 (47) 



496 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

where I×J denotes the total number of samples, M and L number of samples along x and y 

directions respectively (M<I, L<J), the variance of noise, that appears as a consequence of 

missing signal samples, in 2D case can be calculated as follows: 

 2 2

1 1

K

i

i

IJ ML
E ML

IJ








  , (48) 

where Ei denotes amplitude of the i-th component. The above relation is valid for K signal 

components. If the probability that all DFT values at noisy components positions are 

below the signal components is P=0.99, then the threshold that separates signal and non-

signal components in the 2D DFT domain can be calculated as follows: 

 1/( )log(1 )IJT P    . (49) 

Algorithm 7: Iterative SFAR-2D algorithm 

Input: vector of the available samples y, set of the available samples positions Ωa, 2D DFT and 

inverse 2D DFT matrices    and 1 ,  and k =     

o Calculate variance: 2 2 2 2

1 1
,

1

K K

i ii i

IJ ML
E ML E E

IJ


 


 


   ; 

o Calculate threshold T for a given P=0.99; 

o Calculate the initial DFT vector Y that corresponds to y:
a

 Y y  ;  

o for i=1 : i=i+1: until all components are detected 

                 arg{ }T  k k Y   ; 

                  Calculate  
1

1, ( , )H H
a


  F = A A A y A k   ; 

                         Update y, Y, E2 and σ2: y=y-AY(k);
aY = y  ;  22 / ( )E ML y  and        

                        2 2

1 1

K

i

i

IJ ML
E ML

IJ








  ; 

                       end for  

Output: Reconstructed DFT F 

Assume that the set of the available samples positions is denoted as Ωa, 1 1{( , ),..., ( , )}a M Lx y x y  , 
1  denotes inverse 2D DFT, y is vector of measurements and k denotes set of positions. The 

vector of initial DFT is calculated as:  

 
aY = y , (50) 

where   enotes 2D DFT matrix obtained as a Kronecker product of two I×J DFT 

matrices: DFT DFT   and matrix a contains all rows of the full 2D DFT matrix   

and only those columns defined by the set Ωa.  

 The iterative version of the algorithm is summarized through the Algorithm 7 (Simple 

and Fast Algorithm for Reconstruction of 2D signals - SFAR 2D [96]). The algorithm is 

tested on the simulated radar signal. It is assumed that only 25% of the samples are 

available. The 2D DFT of the original signal is shown in Fig. 9a. The components of the 

initial DFT, selected by the threshold, and the reconstructed 2D DFT of the signal, are 

shown in Fig. 9b and c, respectively. 

 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 497 

  

 
a) b) c) 

Fig. 9 Procedure applied on simulated radar signal: a) 2D DFT of the original signal;  

b) components of 2D signal detected in a single step, c) 2D DFT of the 

reconstructed signal  

(C) Wu et al. proposed method for Multivariate Compressive Sensing for Image 

Reconstruction in the Wavelet Domain - multivariate pursuit algorithm (MPA) [130]. 

Statistical structures of the wavelet coefficients are exploited, having in mind significantly 

statistical dependency that exists among the wavelet coefficients of images. The observed 

images are sparse or compressive in wavelet domain.  

Unlike the traditional CS where the measurements are taken by firstly rearranging the 

coefficients of images to a vector and then projecting it with sensing matrix, here the 

coefficients are rearranged according to the partitioned neighborhoods.  

Algorithm 8: MPA 

Input: CS matrix 'K MA , measurement matrix Ω.  

Output: Approximation X    to the true coefficients matrix X 

1. Set n=0, 0X , 0 rY Y   

2. Set 
00,

n n X r rY  . Compute the vector of residual correlations c: 0Tc A r   

3. Form the index set I: 

{ (1) ,..., ( )},

...
(1) ( )1 1

I I M'

q qj j
c c
I I M'j j

 
  

I =

 ; 

4. Estimate successively the components of 
n

X  according to I: 

a) i=1, ri=r0; 

b) Solve the following sub-problem to estimate 
( )

n
I i

x :   ( ) ( )

T
i n i

I i I i
 r a x  ; 

c) Update the residual:  1
( ) ( )

T
i i n

I i I i
  r r a x  ,  ( )I ia  - I(i)-th column of A; 

d) i=i+1;  if 'i M   go to step b); 

5. Update
n

 X X X  , 1 'n M rY r  ; 

6. n=n+1 and go to the step 2. Stop if n=N+1 or 
2

1 / ( ) 0,1 1

K q n Kq
k jk j


   
  
 
  rY  .                        



498 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

The CS measurements consist of two parts: the scaling coefficients measured directly, and 

the compressed samples of the wavelet coefficients, obtained by projecting the wavelet 

coefficients with the sensing matrix. If we denote the CS matrix as 'K MA , measurement 

matrix Ω, ' /M M q    , where M is the total number of wavelet coefficients and .     

denotes ceiling function (case q=1 corresponds to the coefficients vector of traditional 

CS), then the vector of multivariate measurements YK×q (K<<M’) is: 

 Y AX ,  (51) 

where K×q is the number of measurements. Parameter q denotes the size of the wavelet 

neighborhood. The Multivariate Pursuit Algorithm – MPA steps can be summarized 

through the Algorithm 8. By exploiting the statistical dependences among wavelet 

coefficients in multivariate algorithms, the reconstruction performance is much improved 

compared to the state-of-the art methods. The multivariate algorithms also have higher 

computational efficiency. 

(D) Bobin et al. in [131], [132] exploited the CS approach in astronomical images and 

astronomical instrument design. There is a growing interest in astronomical data 

compression, having in mind that the conventional data compression cannot be used in many 

cases. The system that encodes into the analog domain can be designed with an optical 

system that directly measures incoherent projections of the input image. There are three main 

properties that have to be under control: resolution (point spread function), sensitivity 

(ability to detect low level signals) and photometry. 

The aim is to recover the original signal x from the compressible signal y= Ωx. If the 

measurement matrix is denoted Ω, transform matrix is denoted as   denotes the 

transform matrix, and X= x, then the optimal X can be found as a solution of: 

 
1

min subject to =
X

X y X  . (52) 

A more realistic case is when observations are corrupted by noise: y=Ω(x+n), where n is a 

white Gaussian noise with variance 2
n . The projected noise nΩ=Ωn is still white Gaussian 

noise with the same variance. The optimization problem can now be recast as: 

 
1 2

min subject to - e 
X

X y X  , (53) 

where e denotes an upper bound of  2
n . Suppose that N observations of the same sky 

area are available yi, i=1,…,N:  

 = , 1,...,
ii i i N   y X +n ,  (54) 

and Ωλi, i=1,…,N are N independent random submatrices of Ω, card(λi)=M. Recovering x 

from N compressed observations, i.e. the decompression problem (53) can now be 

defined as: 

 
1 2

1

min subject to -
i

N

i

i

Ne


  
X

X y X ,  (55) 

Or can be recast into the Lagrangian form: 

 
2

1 2
1

1
min -

2 i

N

i

i




  
X

X y X  . (56) 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 499 

That can be solved by using projected Landweber iterative algorithm [131]. The choice of 

the parameter μ is important issue, due to its role as a balance between sparsity constraint 

and data fit of the solution. More details can be found in [131].  

4.5. Compressive sensing in video data 

(A) Beside the application in still image recovering, there are various situations when we 

need to apply the CS in video data [133]. However, the CS application in video has higher 

demands in terms of the complexity of imaging architectures, as well as reconstruction 

algorithms, in comparison  with the still-image CS [133]-[135]. Also, it can happen that a 

compressive camera does not capture a sufficient number of measurements to recover the 

frames of the video. Video recovery problems have high memory requirements and 

algorithm implementations require large dense matrix systems. Therefore, it requires 

high-performance hardware and fast iterative algorithms, in order to provide sufficiently 

high throughout. 

Several compressive imaging architectures are designed: 

 Spatial multiplexing cameras - apply CS multiplexing in space in order to improve 

the spatial resolution of videos obtained from sensor arrays whose spatial 

resolution is low; 

 Temporal multiplexing cameras - apply CS multiplexing in time with an aim to 

improve the temporal resolution of videos obtained from sensor arrays with low 

temporal resolution; 

 Spectral and angular multiplexing cameras - apply CS multiplexing to sense 

variations of light in a scene beyond the spatial and temporal dimensions. 

Several methods for CS video recovery are developed, among them, variational (constrained 

and unconstrained) and greedy methods.  

Constrained variational approach deals with the problems in the form: 

 
,

arg min ( , , ) ( ) subject tof g    
x X

x x y X X x  . (57) 

Function f is used to model the video acquisition process – optics, modulation and 

sampling, while the function g is a regularization function.  One example of the f and g 

functions definitions is in the case of frame-by-frame recovery by using 2D wavelet 

transform: 

 
2

2 1
( , , ) , ( ) ,f g    x y y x X X   (58) 

where x denotes vectorized image frame, Ω is the sensing matrix, Ψ is the 2D wavelet 

transform matrix and μ>0 is a regularization parameter.  

If we are dealing with an invertible transform Ψ then the problem (57) can be defined 

as follows: 

 1arg min ( , , ) ( )f g  
X

X X y X  (59) 

where X denotes the video or single video frame in the transform domain.  

 Greedy algorithms are used for the unconstrained problems. They are based on 

iterative constructing a sparse set of non-zero transform coefficients and finding solution 

of the minimization problem 
2

1

2
- X y . The minimization problem solution can be 



500 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

found using the following greedy algorithms: OMP, regularized OMP (ROMP) and  

stagewise OMP (StOMP), CoSaMP [133]. 

(B) A new approach for estimation of the motion parameters in compressive sensed 

video sequences under a reduced number of randomly chosen video frames is proposed in 

[134]. The method focuses on the velocity estimation and combines sparse reconstruction 

algorithms with time-frequency analysis, applied to μ-propagation signal. The μ-

propagation maps the video frames sequence into the frequency modulated signal, or into 

the high nonlinear phase signal. If a video frame at the instant t is defined as: 

 ( , , ) ( , ) ( , )F x y t p x y o x y     , (60) 

where Δx=x-x0-bxt, Δy=y-y0-byt, o(x,y) denotes the moving object, p is background, (x0,y0) 

denote an initial object position and (bx,by) is the velocity. The projection of the frame 

onto the x-axis is defined as: 

 ( , ) ( , ) ( , ) ( , ) ( ) ( )
y y y

R x t F x y p x y o x y P x O x           . (61) 

Finding derivative of the R(x,t) with respect to t and assuming the constant background, 

the following signal is obtained: 

 
( , ) ( )

( ) ( , 1) ( , )x

R x t O x
R x b R x t R x t

t x


     

 

 
.  (62) 

The velocity estimation is done by applying the TF analysis to the signal in the form: 

 ( ) ( ) j x

x

m t R x e    , (63) 

having in mind that the instantaneous frequency corresponds to the moving object 

velocity. As TF representation, the S-method can be used since it provides cross-terms 

free representation and is more suitable in the noisy signal cases. It is defined based on 

the STFT as: 

 *( , ) ( , ) ( , )
L

M

i L

S t f STFT t f j STFT t f j 


   ,  (64) 

where L is the S-method window width, while the STFT(t,f) is defined as the FT of the 

windowed signal m(t), with window function w(η): ( , ) ( ) ( ) jSTFT t f w m t e 


     . 

The CS is employed to reduce the number of frames required for the IF estimation. In 

other words, the CS is used to assure motion parameters estimation from an incomplete 

set of frames. If the subset of frame is denoted as S, S(x,y,ts)⸦F(x,y,t), where only M 

frames are acquired, ts={t1,…,tM}, then the μ propagation vector will contain small 

number of samples, i.e. we will have signal m(ts). For each windowed signal part used for 

the STFT calculation, we have the measurement vector y(tsi): 

 ( ) ( ) ( ),si si si sy t w m t t t     ,  (65) 

instead of desired vector x(t)=w(η)m(t+η). The FT of the vector y(tsi) will produce low 

resolution in the STFT, and therefore, the CS is used in this step to recover missing 

samples in the vector y(tsi) and improve the resolution. If we denote the desired signal as 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 501 

x, measurement vector as y, the measurement and transform matrices as Φ and   

respectively, then the relation follows:  

   y x = X = AX , (66) 

where X corresponds to the STFT coefficients at certain available time instant tsi. To find 

x or its spectral representation X from an incomplete measurement vector y, the following 

optimization problem can be used: 

 
1

min subject to X y AX  , (67) 

performed for each available time instant.  

 
a) 

frames

Initial SM

50 100 150

50

100

150

 frames

CS-based SM

50 100 150

50

100

150

 
 b) c) 

   0 50 100 150 200
60

80

100

120

    0 50 100 150 200
60

80

100

120

 
 d) e) 

Fig. 10  a) Several frames from the observed video sequence; b) Initial S-method of 

variable µ-propagation vector; c) CS based S-method of variable µ-propagation 

vector; Velocity estimation using: d) initial S-method and e) CS based S-method  

The results obtained by using real video sequence are shown in Fig. 10. The percentage of 

the available frames is 40%, due to the compressive acquisition. The moving of metronome’s 

pendulum is observed and some of the frames from the video sequence are shown in Fig. 

10a. The S-method of the μ-propagation vector calculated using the available samples, is 

shown in Fig. 10b, while the CS-based S-method is shown in Fig. 10c. The corresponding 

velocity estimations graphs are shown in Fig. 10e and f. It is shown that the initial S-method 

produces error in velocity estimation, while precise results are obtained by using the CS-

based method.  



502 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

4.6. CS in watermarking 

(A) Data protection in terms of CS has been discussed in [136]-[141]. Fakhr in [137] 

proposed a watermark embedding and recovery technique based on the CS framework, 

tested under MP3 compression. The sparsity of both, the host and the watermark signal is 

assumed. The watermark is embedded into the measurement vector y. If we denote signal 

with x, transform domain matrix as  , a sparse signal b as watermark of length L, then 

the random watermark creation is described as: 

 , w b  (68) 

where ΩM×L is the random Gaussian matrix, and M is the measurement vector length. 

Matrix Ω is used for random expansion of the sparse vector b. The embedding is done as 

follows: 

    y x a b , (69) 

resulting in watermarked measurement vector. Embedding strength a is adapted for each 

frame of the audio signal as: 2
1

0.04 ,
M

i
i


  a X X x . The advantage of the proposed 

method is that, in order to recover the clean signal, the optimization problem has to be 

solved and thus, matrix Ω has to be known.  In this paper, for the optimization problem 

solving three methods are used: Direct Justice Pursuit, Multiplying by the Inverse of Ω 

and Multiplying by the annihilator of Ω. 

(B) An image watermarking procedure in the CS scenario is proposed in [139]. The 

randomly chosen pixels that serve as CS measurements are used to bring the watermark. 

The image is firstly divided into the blocks and measurements are selected from each 

block. Samples are taken from the space domain, while the image sparsity is assumed in 

the DFT domain.  

If we denote the N×N image block as Ij, vector of measurements for j-th block as yj, Tj 

vector of transform domain coefficients (DFT) of the block Ij,   as the Fourier transform 

matrix and Ωj as the measurement matrix for the block j, then the measurement vector is 

defined as: 

 j j j j j    y I T .  (70) 

The watermarked measurement vector jy   is obtained as follows: 

 j jj  y y  , (71) 

where μ denotes watermark strength and ω is M×1 watermark vector (M denotes the 

number of measurements). The vector of watermarked coefficients is used to recover the 

image according to the total variation optimization: 

 min ( ) subject toj j j jTV   
T

T y T  . (72) 

The reconstructed image block IRj is obtained as Rj j I T . Watermark detection is based 

on using the standard correlator that requires measurement matrix Ωj to be known:  

 ( ) i ii
D   y  . (73) 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 503 

The procedure is tested on 256×256 image, divided into 16×16 block. From each block, 

50% of the pixels is randomly chosen and serve as a measurement in the reconstruction 

process, and carry watermark as well. The results are shown in Fig. 11. PSNR between 

original and watermarked/reconstructed image is 31.79 dB. 

Original image Reconstructed image

 

             
500 1000 1500 2000 2500

-2

0

2

4

6

8

10
x 10

-3

 

 

Right keys

Wrong trials

 
Fig. 11 a) Original image; b) Watermarked and reconstructed image; c) Detector responses  

for 25 right keys and 2500 wrong trials (100 wrong trials for each right key) 

5. CONCLUSION 

The paper focuses on the Compressive Sensing, as an approach that records an intensive 

development in signal processing in recent years. An overview of the Compressive Sensing 

applications and commonly used algorithms for reconstruction of the signals with missing 

data is given. Algorithms for the reconstruction of both, 1D and 2D signals, are described in 

the paper. The paper covers the applications starting from the radar signal processing, 

communications, biomedical signals and image reconstruction, through natural image 

reconstruction, velocity estimation in video signal processing, CS-based protection of the 

digital data and hardware devices designed based on the CS principles. Experimental results 

are provided in order to show the performance of the presented algorithms and approaches.  

Acknowledgement: The paper is a part of the research supported by the Montenegrin Ministry of 

Science, project grant: “New ICT Compressive sensing based trends applied to: multimedia, 

biomedicine and communications (CS-ICT)” (Montenegro Ministry of Science, Grant No. 01-

1002). The authors are thankful to Dr Josip Musić for testing the performance of the object 

detection in Search&Rescue images. 



504 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

REFERENCES  

[1] G. Pope, “Compressive Sensing: a Summary of Reconstruction Algorithms”, Eidgenossische Technische 

Hochschule, Zurich, Switzerland, 2008. 

[2] E. Candes, J. Romberg, “l1-magic: Recovery of Sparse Signals via Convex Programming”, October 

2005. 

[3] D. Donoho, “Compressed sensing,” IEEE Transactions on IT, vol. 52, no.4, 2006, pp. 1289 - 1306. 

[4] E. J. Candes, J. Romberg, T. Tao, "Robust uncertainty principles: exact signal reconstruction from 

highly incomplete frequency information," IEEE Transactions on Information Theory, vol. 52, no. 2, pp. 

489-509, Feb. 2006. 

[5] LJ. Stankovic, M. Dakovic, S. Stankovic, I. Orovic, "Sparse Signal Processing," in the Book: Digital 

Signal Processing, L. Stankovic, CreateSpace, Amazon, 2015. 

[6]  I. Orovic, V. Papic, C. Ioana, X. Li, S. Stankovic, "Compressive Sensing in Signal Processing: 

Algorithms and Transform Domain Formulations," Mathematical Problems in Engineering, Review 

paper, 2016. 

[7] E. J. Candes and T. Tao, “Decoding by linear programming,” Information Theory, IEEE Transactions 

on, vol. 51, no. 12, pp. 4203–4215, 2005.  

[8] G. Davis, S. Mallat, and M. Avellaneda, “Adaptive greedy approximations,” Constructive approximation, vol. 

13, no. 1, pp. 57–98, 1997. 

[9] Y. Arjoune, N. Kaabouch, H. El Ghazi, A. Tamtaoui, "Compressive sensing: Performance comparison of 

sparse recovery algorithms," In 2017 IEEE 7th Annual Computing and Communication Workshop and 

Conference (CCWC), Las Vegas, NV, 2017, pp. 1-7. 

[10] S. Stankovic, I. Orovic, M. Amin, "L-statistics based Modification of Reconstruction Algorithms for 

Compressive Sensing in the Presence of Impulse Noise," Signal Processing, vol.93, no.11, November 

2013, pp. 2927-2931, 2013. 

[11] Y. C. Eldar and G. Kutyniok, "Compressed Sensing: Theory and Applications", Cambridge University 

Press, May 2012. 

[12] S. Stankovic, I. Orovic, E. Sejdic, "Multimedia Signals and Systems: Basic and Advance Algorithms for 

Signal Processing," Springer-Verlag, New York, 2015. 

[13] V.M. Patel and R. Chellappa, “Sparse Representations and Compressive Sensing for Imaging and 

Vision,” SpringerBriefs in Electrical and Computer Engineering, 2013. 

[14] T. Blumensath, M. E. Davies, “Iterative Thresholding for Sparse Approximations”, Journal of Fourier Analysis 

and Applications, vol. 14, no. 5-6, pp 629-654, December 2008. 

[15] T. Blumensath, M. E. Davies, "Gradient Pursuits," IEEE Transactions on Signal Processing, vol.56, 

no.6, pp.2370-2382, June 2008. 

[16] R. Mihajlovic, M. Scekic, A. Draganic, S. Stankovic, "An Analysis of CS Algorithms Efficiency for 

Sparse Communication Signals Reconstruction," In Proceedings of the 3rd Mediterranean Conference on 

Embedded Computing, MECO, 2014. 

[17] L. I. Rudin, S. Osher, E. Fatemi, “Nonlinear total variation based noise removal algorithms”, Physica D: 

Nonlinear Phenomena, vol. 60, Issues 1–4, 1 November 1992, pp. 259-268  

[18] S. Stankovic, I. Orovic, "Robust Complex-Time Distributions based on Reconstruction Algorithms," In 

Proceedings of the 2nd Mediterranean Conference on Embedded Computing MECO - 2013, Budva, 

Montenegro, 2013, pp. 105-108. 

[19] LJ. Stankovic, S. Stankovic, M. Amin, "Missing Samples Analysis in Signals for Applications to L-

estimation and Compressive Sensing," Signal Processing, vol. 94, Jan 2014, pp. 401-408, 2014. 

[20] S. Stankovic, LJ. Stankovic, I. Orovic, "A relationship between the Robust Statistics Theory and Sparse 

Compressive Sensed Signals Reconstruction," IET Signal Processing, Special issue on Compressive 

Sensing and Robust Transforms, vol. 8, Issue 3, pp. 223 - 229, May, 2014 

[21] S. Bahmani, “Algorithms for Sparsity-Constrained Optimization”, Springer Theses, Series Volume 261, 

ISBN 978-3-319-01880-5, 2014. M.  

[22] T. Zhang, “Sparse Recovery with Orthogonal Matching Pursuit Under RIP,” IEEE Trans. on 

Information Theory, vol. 57, no. 9, pp. 6215-6221, 2011. 

[23] S. Stankovic, I. Orovic, LJ. Stankovic, A. Draganic, "Single-Iteration Algorithm for Compressive Sensing 

Reconstruction," Telfor Journal, vol. 6, no. 1, pp. 36-41, 2014. 

[24] A. Draganic, I. Orovic, N. Lekic, M. Dakovic, S. Stankovic, "Architecture for Single Iteration 

Reconstruction Algorithm," In Proceedings of the 4th Mediterranean Conference on Embedded 

Computing. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 505 

[25] J. A. Tropp, A. C. Gilbert, “Signal Recovery From Random Measurements Via Orthogonal Matching 

Pursuit,” IEEE Transaction on Information Theory, vol. 53, no.12, 2007. 

[26] LJ. Stanković, M. Daković, and S. Vujović, “Adaptive Variable Step Algorithm for Missing Samples 

Recovery in Sparse Signals,” IET Signal Processing, vol. 8, no. 3, pp. 246 -256, 2014. 

[27] S. Vujovic, M. Dakovic, I. Orovic, S. Stankovic, "An Architecture for Hardware Realization of 

Compressive Sensing Gradient Algorithm," In Proceedings of the 4th Mediterranean Conference on 

Embedded Computing MECO - 2015,  Budva, Montenegro. 

[28] Y. Wang, J. Xiang, Q. Mo, S. He, “Compressed sparse time–frequency feature representation via 

compressive sensing and its applications in fault diagnosis”, Measurement, vol. 68, pp. 70–81, May 

2015. 

[29] S. Stankovic, I. Orovic, "An Ideal OMP based Complex-Time Distribution," 2nd Mediterranean 

Conference on Embedded Computing MECO - 2013, pp. 109-112, June 2013, Budva, Montenegro. 

[30] Y. C. Eldar "Sampling Theory: Beyond Bandlimited Systems", Cambridge University Press, April 2015. 

[31] P. Flandrin, P. Borgnat, "Time-Frequency Energy Distributions Meet Compressed Sensing," IEEE Transactions 

on Signal Processing, vol.58, no.6, pp.2974, 2982, June 2010.  

[32] I. Orovic, S. Stankovic, T. Thayaparan, "Time-Frequency Based Instantaneous Frequency Estimation of 

Sparse Signals from an Incomplete Set of Samples," IET Signal Processing, Special issue on 

Compressive Sensing and Robust Transforms, vol. 8, Issue 3, pp. 239 - 245, May, 2014. 

[33] I. Orovic, S. Stankovic, M. Amin, "Compressive Sensing for Sparse Time-Frequency Representation of 

Nonstationary Signals in the Presence of Impulsive Noise," SPIE Defense, Security and Sensing, 

Baltimore, Maryland, United States, 2013. 

[34] P. Borgnat, and P. Flandrin, "Time-frequency localization from sparsity constraints," In Proceedings of 

the IEEE International Conference on Acoustics, Speech and Signal Proceesing ICASSP-08, Las Vegas 

(NV), 2008, pp. 3785–3788. 

[35] M. Brajović, B. Lutovac, I. Orović, M. Daković, S. Stanković, “Sparse Signal Recovery Based on 

Concentration Measures and Genetic Algorithm,” In Proceedings of the 13th Symposium on Neural 

Networks and Applications NEUREL 2016, Belgrade, Serbia, November 2016. 

[36] X. Li, G. Bi, “Time-frequency representation reconstruction based on the compressive sensing”, In 

Proceedings of the 9th IEEE Conference on Industrial Electronics and Applications, Hangzhou, 2014, 

pp. 1158-1162. 

[37] S. Stankovic, I. Orovic, M. Amin, "Compressed Sensing Based Robust Time-Frequency Representation 

for Signals in Heavy-Tailed Noise," In Proceedings of the Information Sciences, Signal Processing and 

their Applications, ISSPA 2012, Canada, 2012. 

[38] P. K. Mishra, R. Bharath, P. Rajalakshmi, U. B. Desai, "Compressive sensing ultrasound beamformed imaging 

in time and frequency domain," In Proceedings of the 17th International Conference on E-health Networking, 

Application & Services (HealthCom), Boston, MA, 2015, pp. 523-527. 

[39] I. Orovic, S. Stankovic, T. Thayaparan, LJ. Stankovic, "Multiwindow S-method for Instantaneous 

Frequency Estimation and its Application in Radar Signal Analysis," IET Signal Processing, vol. 4, no. 

4, pp. 363-370, 2010 

[40] I. Orovic, S. Stankovic, "A Class of Highly Concentrated Time-Frequency Distributions Based on the 

Ambiguity Domain Representation and Complex-Lag Moment," EURASIP Journal on Advances in 

Signal Processing, vol. 2009, Article ID 935314, 9 pages, 2009. 

[41] S. Stankovic, I. Orovic, LJ. Stankovic, “Polynomial Fourier Domain as a Domain of Signal Sparsity”, 

Signal Processing, vol. 130, Issue C, pp. 243-253, January 2017. 

[42] S. Stankovic, I. Orovic, T. Pejakovic, M. Orovic, "Compressive sensing reconstruction of signals with 

sinusoidal phase modulation: application to radar micro-Doppler," In Proceedings of the 22nd 

Telecommunications Forum , TELFOR, 2014. 

[43] H. Su, Y. Zhang, "Time-frequency analysis based on Compressive Sensing," In Proceedings of the 2nd 

International Conference on Cloud Computing and Internet of Things (CCIOT), Dalian, 2016, pp. 138-

142. 

[44] I. Volaric, V. Sucic, Z. Car, "A compressive sensing based method for cross-terms suppression in the 

time-frequency plane," In Proceedings of the IEEE 15th International Conference on Bioinformatics and 

Bioengineering (BIBE), Belgrade, 2015, pp. 1-4. 

[45] LJ. Stankovic, S. Stankovic, I. Orovic, M. Amin, "Robust Time-Frequency Analysis based on the L-

estimation and Compressive Sensing," IEEE Signal Processing Letters, vol. 20, no. 5, pp. 499-502, 

2013. 



506 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

[46] G. Hua, Y. Hiang, G. Bi, “When Compressive Sensing meets Data Hiding”, IEEE Signal Processing 

Letters, vol. 23, no. 4, April 2016. 

[47] A. Draganic, M. Brajovic, I. Orovic, S. Stankovic, "A Software Tool for Compressive Sensing based 

Time-Frequency Analysis," In Proceedings of the 57th International Symposium, ELMAR-2015, Zadar, 

Croatia, 2015. 

[48] I. Orovic, S. Stankovic, T. Chau, C. M. Steele, E. Sejdic, "Time-frequency analysis and Hermite 

projection method applied to swallowing accelerometry signals," EURASIP Journal on Advances in 

Signal Processing, vol. 2010, Article ID 323125, 7 pages, 2010. 

[49] A. Krylov, D. Korchagin, “Fast Hermite projection method,” In Proceedings of the 3rd International Conference 

on Image Analysis and Recognition (ICIAR ’06), vol. 1, pp. 329–338, Povoa de Varzim, Portugal, September 

2006. 

[50] S. Stankovic, I. Orovic, A. Krylov, "The Two-Dimensional Hermite S-method for High Resolution 

Inverse Synthetic Aperture Radar Imaging Applications," IET Signal Processing, vol. 4, no. 4, pp. 352-

362, 2010. 

[51] M. Brajović, I. Orović, M. Daković, S. Stanković, “Compressive Sensing of Signals Sparse in 2D 

Hermite Transform Domain,” 58th International Symposium ELMAR-2016, Zadar, Croatia, September 

2016. 

[52] A. Sandryhaila, S. Saba, M. Püschel, J. Kovačević, “Efficient Compression of QRS Complexes Using 

Hermite Expansion,” IEEE Transactions on Signal Processing, vol. 60, no. 2, pp. 947-955, February 

2012. 

[53] A. Draganić, I. Orović, S. Stanković, “Robust Hermite transform based on the L-estimate principle,” In 

Proceedings of the 23rd Telecommunications Forum, TELFOR 2015. 

[54] S. Stankovic, LJ. Stankovic, I. Orovic, "Compressive sensing approach in the Hermite transform 

domain," Mathematical Problems in Engineering, vol. 2015 (2015), Article ID 286590, 9 pages. 

[55] M. Brajovic, I. Orovic, M. Dakovic, S. Stankovic, "The Analysis of Missing Samples in Signals Sparse 

in the Hermite Transform Domain," In Proceedings of the 23rd Telecommunications Forum TELFOR, 

2015, Belgrade, Serbia, 2015. 

[56] I. Orovic, S. Stankovic, "Improved Higher Order Robust Distributions based on Compressive Sensing 

Reconstruction," IET Signal Processing, vol. 8, Issue: 7, pp. 738 - 748, May 2014. 

[57] S. Stankovic, I. Orovic, LJ. Stankovic, "An Automated Signal Reconstruction Method based on Analysis 

of Compressive Sensed Signals in Noisy Environment," Signal Processing, vol. 104, Nov 2014, pp. 43 - 

50, 2014. 

[58] M. G. Christensen, J. Østergaard, S. H. Jensen, "On compressed sensing and its application to speech 

and audio signals," Conference Record of the Forty-Third Asilomar Conference on Signals, Systems and 

Computers, Pacific Grove, CA, 2009, pp. 356-360. 

[59] M. Scekic, R. Mihajlovic, I. Orovic, S. Stankovic, "CS Performance Analysis for the Musical Signals 

Reconstruction," In Proceedings of the 3rd Mediterranean Conference on Embedded Computing, 

MECO, 2014. 

[60] D. Wu, W. P. Zhu, M. N. S. Swamy, "A compressive sensing method for noise reduction of speech and 

audio signals," In Proceedings of the IEEE 54th International Midwest Symposium on Circuits and 

Systems (MWSCAS), Seoul, 2011, pp. 1-4. 

[61] L. Sun, X. Shao, Z. Yang, "An Adaptive Multiscale Framework for Compressed Sensing of Speech 

Signal," In Proceedings of the 6th International Conference on Wireless Communications Networking 

and Mobile Computing (WiCOM), Chengdu, 2010, pp. 1-4. 

[62] M. Dakovic, LJ. Stankovic, S. Stankovic, "A Procedure for Optimal Pulse Selection Strategy in Radar 

Imaging Systems," In Proceedings of the International Workshop on Compressed Sensing Theory and its 

Applications to Radar, Sonar and Remote Sensing (CoSeRa), 19-22 September, Aachen, Germany, 2016. 

[63] A. Bacci, E. Giusti, D. Cataldo, S. Tomei, M. Martorella, "ISAR resolution enhancement via 

compressive sensing: A comparison with state of the art SR techniques," In Proceedings of the 4th 

International Workshop on Compressed Sensing Theory and its Applications to Radar, Sonar and 

Remote Sensing (CoSeRa), Aachen, 2016, pp. 227-231. 

[64] S. Costanzo, A. Rocha, M. D. Migliore, “Compressed Sensing: Applications in Radar and Communications”, 

The Scientific World Journal, vol. 2016 (2016), Article ID 5407415, 2 pages, Editorial. 

[65] LJ. Stankovic, S. Stankovic, T. Thayaparan, M. Dakovic, I. Orovic, "Separation and Reconstruction of 

the Rigid Body and Micro-Doppler Signal in ISAR Part I-Theory ," IET Radar, Sonar & Navigation, 

vol. 9, no. 9, pp. 1147-1154, 2015. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 507 

[66] LJ. Stankovic, S. Stankovic, T. Thayaparan, M. Dakovic, I. Orovic, "Separation and Reconstruction of 

the Rigid Body and Micro-Doppler Signal in ISAR Part II-Statistical Analysis," IET Radar, Sonar & 

Navigation, vol. 9, no. 9, pp. 1155-1161, 2015. 

[67] A. Draganic, I. Orovic, S. Stankovic, X. Li, "ISAR Reconstruction from Incomplete Data using Total 

Variation Optimization," In Proceedings of the 5th Mediterranean Conference on Embedded Computing, 

(MECO 2016). 

[68] L. C. Potter, E. Ertin, J. T. Parker, M. Cetin, "Sparsity and Compressed Sensing in Radar Imaging," In 

Proceedings of the IEEE, vol. 98, no.6, pp.1006-1020, June 2010. 

[69] M. Dakovic, LJ. Stankovic, S. Stankovic, "Gradient Algorithm Based ISAR Image Reconstruction From 

the Incomplete Dataset," In Proceedings of the 3rd International Workshop on Compressed Sensing 

Theory and its Applications to Radar, Sonar and Remote Sensing, CoSeRa, 2015. 

[70] J. Ender, “On compressive sensing applied to radar”, Signal Processing, vol. 90, Issue 5, May 2010, pp. 

1402–1414. 

[71] LJ. Stankovic, S. Stankovic, I. Orovic, Y. Zhang, "Time-Frequency Analysis of Micro-Doppler Signals Based on 

Compressive Sensing," Compressive Sensing for Urban Radar, Ed. M. Amin, CRC-Press, 2014. 

[72] LJ. Stankovic, I. Orovic, S. Stankovic, M. Amin, "Compressive Sensing Based Separation of Non-

Stationary and Stationary Signals Overlapping in Time-Frequency," IEEE Transactions on Signal 

Processing, vol. 61, no. 18, pp. 4562-4572, Sept. 2013. 

[73] S. Stankovic, LJ. Stankovic, I. Orovic, "L-statistics combined with compressive sensing," SPIE Defense, 

Security and Sensing,  Baltimore, Maryland, United States, 2013. 

[74] M. A. Hadi, S. Alshebeili, K. Jamil, F. E. Abd El-Samie, “Compressive sensing applied to radar systems: 

an overview”, Signal, Image and Video Processing, December 2015, Volume 9, Supplement 1, pp 25–

39. 

[75] A. Draganic, I. Orovic, S. Stankovic, "Blind Signals Separation in wireless communications based on 

Compressive Sensing," In Proceedings of the 22nd Telecommunications Forum, TELFOR, 2014. 

[76] L. Zhang, M. Xing, C. W. Qiu J. Li, Z. Bao, “Achieving higher resolution ISAR imaging with limited 

pulses via compressed sampling,” IEEE Geoscience and Remote Sensing Letters, vol.6, no.3, pp.567–

571, 2009. 

[77] I. Orovic, A. Draganic, S. Stankovic, "Sparse Time-Frequency Representation for Signals with Fast 

Varying Instantaneous Frequency," IET Radar, Sonar & Navigation, vol. 9, Issue 9, pp. 1260 – 1267. 

[78] S. Li, G. Zhao, W. Zhang, Q. Qiu, H. Sun, "ISAR Imaging by Two-Dimensional Convex Optimization-

Based Compressive Sensing," IEEE Sensors Journal, vol. 16, no. 19, pp. 7088-7093, Oct.1, 2016. 

[79] P. Zhang, Z. Hu, R. C. Qiu, B. M. Sadler, “A Compressed Sensing Based UltraWideband Communication 

System,” In Proceedings of the IEEE International Conference on Communications, 14-18 June 2009. 

[80] A. Draganic, I. Orovic, S. Stankovic, M. Amin, "Rekonstrukcija FHSS signala zasnovana na principu 

kompresivnog odabiranja," In Proceedings of the TELFOR 2012, Belgrade, 2012 

[81] J. Meng, J. Ahmadi-Shokouh, H. Li, E. J. Charlson, Z. Han, S. Noghanian, E. Hossain, “Sampling Rate 

Reduction for 60 GHz UWB Communication using Compressive Sensing, ” In Proceedings of the 

Asilomar Conf. on Signals, Systems, and Computers, Monterey, California, November 2009. 

[82] A. Draganic, I. Orovic, S. Stankovic, X. Li, Z. Wang, "Reconstruction and classification of wireless 

signals based on Compressive Sensing approach," In Proceedings of the 5th Mediterranean Conference 

on Embedded Computing, (MECO 2016). 

[83]  B. Jokanovic, M. Amin, S. Stankovic, "Instantaneous frequency and time-frequency signature estimation using 

compressive sensing," SPIE Defense, Security and Sensing, Baltimore, Maryland, United States, 2013, 

http://dx.doi.org/10.1117/12.2016636 

[84] C. Bernard, C. Ioana, I. Orovic, S. Stankovic, "Analysis of underwater signals with nonlinear time-frequency 

structures using warping based compressive sensing algorithm," In Proceedings of the MTS/IEEE North 

American OCEANS conference, October 2015, Washington, DC, United States, 2015. 

[85] I. Murgan, A. Digulescu, I. Candel, C. Ioana, “Compensation of position offset of acoustic transducers 

using compressive sensing concept”, In Proceedings of the OCEANS 2016 MTS/IEEE Monterey, Sep 

2016, Monterey, United States. pp. 1-4. 

[86] I. Orovic, S. Stankovic, LJ. Stankovic, "Compressive Sensing Based Separation of LFM Signals," In 

Proceedings of the 56th International Symposium ELMAR 2014, Zadar, Croatia, 2014. 

[87] J. Musić, T. Marasović, V. Papić, I. Orović, S. Stanković, "Performance of Compressive Sensing Image 

Reconstruction for Search and Rescue," IEEE Geoscience and Remote Sensing Letters, vol. 13, no. 11, 

pp. 1739-1743, Nov. 2016.  

http://dx.doi.org/10.1117/12.2016636


508 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

[88] J. Music, I. Orovic, T. Marasovic, V. Papic, S. Stankovic, "Gradient Compressive Sensing for Image 

Data Reduction in UAV based Search and Rescue in the Wild," Mathematical Problems in Engineering, 

November,  2016 

[89] A. Akbari, D. Mandache, M. Trocan and B. Granado, "Adaptive saliency-based compressive sensing 

image reconstruction," In Proceedings of the IEEE International Conference on Multimedia & Expo 

Workshops (ICMEW), Seattle, WA, 2016, pp. 1-6. 

[90] N. Eslahi, A. Aghagolzadeh, "Compressive Sensing Image Restoration Using Adaptive Curvelet 

Thresholding and Nonlocal Sparse Regularization," IEEE Transactions on Image Processing, vol. 25, 

no. 7, pp. 3126-3140, July 2016. 

[91] J. Wen, Z. Chen, Y. Han, J. D. Villasenor, S. Yang, "A compressive sensing image compression 

algorithm using quantized DCT and noiselet information," In Proceedings of the IEEE International 

Conference on Acoustics, Speech and Signal Processing, Dallas, TX, 2010, pp. 1294-1297. 

[92] I. Stankovic, I. Orovic, S. Stankovic, M. Dakovic, "Iterative Denoising of Sparse Images," In 

Proceedings of the 39th International Convention on Information and Communication Technology, 

Electronics and Microelectronics, (MIPRO 2016), 2016. 

[93] M. Medenica, S. Zukovic, A. Draganic, I. Orovic, S. Stankovic, "Comparison of the algorithms for CS 

image reconstruction," ETF Journal of Electrical Engineering 2014, 09/2014; vol. 20, no. 1, pp. 29-39. 

[94] C.-S. Lu, H.-W. Chen, “Compressive image sensing for fast recovery from limited samples: A variation 

on compressive sensing”, Information Sciences, vol. 325, 20 December 2015, Pages 33–47. 

[95] M. Maric, I. Orovic, S. Stankovic, "Compressive Sensing based image processing in TrapView pest 

monitoring system," In Proceedings of the 39th International Convention on Information and 

Communication Technology, Electronics and Microelectronics, (MIPRO 2016). 

[96] S. Stankovic, I. Orovic, "An Approach to 2D Signals Recovering in Compressive Sensing Context," 

Circuits Systems and Signal Processing, 2016. 

[97] Z. Zhu, K. Wahid, P. Babyn, D. Cooper, I. Pratt, Y. Carter, “Improved Compressed Sensing-Based 

Algorithm for Sparse-View CT Image Reconstruction”, Computational and Mathematical Methods in 

Medicine, vol. 2013 (2013), Article ID 185750, 15 pages. 

[98] I. Stankovic, I. Orovic, S. Stankovic, "Image Reconstruction from a Reduced Set of Pixels using a 

Simplified Gradient Algorithm," In Proceedings of the 22nd Telecommunications Forum TELFOR 2014, 

Belgrade, Serbia, 2014. 

[99] M. Lustig, D. Donoho, J. Pauly, “Sparse MRI: The application of compressed sensing for rapid MR imaging,” 

Magn. Reson. Med., vol. 58, no. 6, pp. 1182–1195, 2007 

[100] C. G. Graff, E. Y. Sidky, “Compressive sensing in medical imaging”, Applied Optics, 2015 Mar 10; vol. 

54, no. 8, C23–C44. 

[101] J. M. Bioucas-Dias, M. A. T. Figueiredo, "A New TwIST: Two-Step Iterative Shrinkage/Thresholding 

Algorithms for Image Restoration," IEEE Transactions on Image Processing, vol. 16, no. 12, pp. 2992-

3004, Dec. 2007. 

[102] M. F. Duarte et al., "Single-Pixel Imaging via Compressive Sampling," IEEE Signal Processing 

Magazine, vol. 25, no. 2, pp. 83-91, March 2008. 

[103] R. Fergus, A. Torralba, W. T. Freeman, “Random lens imaging”, MIT CSAIL Technical Report, 

September 2006. 

[104] M. Trakimas, R. D'Angelo, S. Aeron, T. Hancock, S. Sonkusale, “A Compressed Sensing Analog-to-

Information Converter With Edge-Triggered SAR ADC Core”, IEEE Transactions on Circuits and 

Systems I: Regular Papers, pp. 1135- 1148, vol. 60, Issue: 5, 2013 

[105] M. Lustig, D.L Donoho, J.M Santos, J.M Pauly “Compressed Sensing MRI”, IEEE Signal Processing 

Magazine, 2008, vol. 25, no. 2, pp. 72-82. 

[106] M. Lustig, J.M. Santos, D.L. Donoho, and J.M. Pauly, “k-t Sparse: High frame rate dynamic MRI 

exploiting spatio-temporal sparsity,” In Proceedings of the 13th Annual Meeting ISMRM, Seattle, WA, 

2006, p. 2420. 

[107] S. Zukovic, M. Medenica, A. Draganic, I. Orovic, S. Stankovic, "A Virtual Instrument for Compressive 

Sensing of Multimedia Signals," In Proceedings of the 56th International Symposium ELMAR 2014, 

Zadar, Croatia, 2014. 

[108] M. Hong, Y. Yu, H. Wang, F. Liu, S. Crozier “Compressed sensing MRI with singular value 

decomposition-based sparsity basis”, Physics in Medicine and Biology, vol. 56 (2011), pp. 6311–6325. 

[109] D. Craven, B. McGinley, L. Kilmartin, M. Glavin, E. Jones, "Compressed Sensing for Bioelectric 

Signals: A Review," IEEE Journal of Biomedical and Health Informatics, vol. 19, no. 2, pp. 529-540, 

March 2015. 



 On Some Common Compressive Sensing Recovery Algorithms and Applications 509 

[110] Y. Liu, M. De Vos, S. Van Huffel, "Compressed Sensing of Multichannel EEG Signals: The Simultaneous 

Cosparsity and Low-Rank Optimization," IEEE Transactions on Biomedical Engineering, vol. 62, no. 8, pp. 

2055-2061, Aug. 2015. 

[111] A. M. Abdulghani, A. J. Casson, E. Rodriguez-Villegas, "Quantifying the Feasibility of Compressive Sensing in 

Portable Electroencephalography Systems," In Proceedings of the 5th International Conference on Foundations 

of Augmented Cognition. Neuroergonomics and Operational Neuroscience: Held as Part of HCI International 

2009, San Diego, CA, 2009, pp. 319-328. 

[112] S. Senay, L. F. Chaparro, M. Sun, R. J. Sclabassi, "Compressive sensing and random filtering of EEG 

signals using slepian basis," In Proceedings of the 16th European Signal Processing Conference 

(EUSIPCO 2008), Lausanne, Switzerland, 2008. 

[113] Z. Zhang, T. P. Jung, S. Makeig, B. D. Rao, "Compressed Sensing of EEG for Wireless Telemonitoring 

With Low Energy Consumption and Inexpensive Hardware," IEEE Transactions on Biomedical 

Engineering, vol. 60, no. 1, pp. 221-224, Jan. 2013. 

[114] J. K. Pant, S. Krishnan, "Reconstruction of ECG signals for compressive sensing by promoting sparsity 

on the gradient," In Proceedings of the IEEE International Conference on Acoustics, Speech and Signal 

Processing, Vancouver, BC, 2013, pp. 993-997. 

[115] L. F. Polanía, R. E. Carrillo, M. Blanco-Velasco, K. E. Barner, "Exploiting Prior Knowledge in Compressed 

Sensing Wireless ECG Systems," IEEE Journal of Biomedical and Health Informatics, vol. 19, no. 2, pp. 508-

519, March 2015. 

[116] O. Kerdjidj, K. Ghanem, A. Amira, F. Harizi, F. Chouireb, "Real ECG signal acquisition with shimmer platform 

and using of compressed sensing techniques in the offline signal reconstruction," In Proceedings of the IEEE 

International Symposium on Antennas and Propagation (APSURSI), Fajardo, 2016, pp. 1179-1180. 

[117] K. Wilhelm, Y. Massoud, "Compressive sensing based classification of intramuscular electromyographic 

signals," In Proceedings of the IEEE International Symposium on Circuits and Systems, Seoul, 2012, pp. 273-

276. 

[118] M. Brajović, I. Orović, M. Daković, S. Stanković, “Gradient-based signal reconstruction algorithm in the 

Hermite transform domain,” Electronics Letters, vol. 52, Issue 1, pp. 41-43, 2016. 

[119] M. Brajovic, I. Orovic, M. Dakovic, S. Stankovic, "On the Parameterization of Hermite Transform with 

Application to the Compression of QRS Complexes," Signal Processing, vol. 131, February 2017, Pages 113–

119. 

[120] M. Brajovic, I. Orovic, S. Stankovic, "The Optimization of the Hermite transform: Application Perspectives and 

2D Generalization," In Proceedings of the 24th Telecommunications Forum TELFOR 2016, November 2016, 

Belgrade, Serbia, 2016. 

[121] G. Teschke “Sparse Recovery and Compressive Sampling in Inverse and Ill-Posed Problems”, Lecture Notes. 

[122] J. Trzasko, A. Manduca, "Highly Undersampled Magnetic Resonance Image Reconstruction via 

Homotopic \ell _{0} -Minimization," IEEE Transactions on Medical Imaging, vol. 28, no. 1, pp. 106-

121, Jan. 2009. 

[123] P. Zhang, Z. Hu, R. C. Qiu, B. M. Sadler, "A Compressed Sensing Based Ultra-Wideband Communication 

System," In Proceedings of the IEEE International Conference on Communications, Dresden, 2009, pp. 1-5. 

[124] B. Zhang, X. Cheng, N. Zhang, Y. Cui, Y. Li, Q. Liang, “Sparse target counting and localization in sensor 

networks based on compressive sensing,” In Proceedings of the IEEE INFOCOM, 2011, pp. 2255–2263. 

[125] M. Weiss, “Passive WLAN radar network using compressed sensing,” In Proceedings of the IET 

International Conference on Radar Systems (Radar 2012), Glasgow, UK, 2012, pp. 1-6. 

[126] J. Bazerque, G. Giannakis, “Distributed spectrum sensing for cognitive radio networks by exploiting sparsity,” 

IEEE Trans. Signal Process., vol. 58, no. 3, pp. 1847–1862, Mar. 2010. 

[127] M. Brajovic, A. Draganic, I. Orovic, S. Stankovic, "FHSS signal sparsification in the Hermite transform 

domain," In Proceedings of the 24th Telecommunications Forum TELFOR 2016, November 2016, 

Belgrade, Serbia, 2016. 

[128] Y. Lu, W. Guo, X. Wang, W. Wang, "Distributed Streaming Compressive Spectrum Sensing for Wide-

Band Cognitive Radio Networks," In Proceedings of the IEEE 73rd Vehicular Technology Conference 

(VTC Spring), Yokohama, 2011, pp. 1-5. 

[129] I. Stanković, I. Orović, M. Daković, S. Stanković, “Denoising of Sparse Images in Impulsive Disturbance 

Environment,” Multimedia Tools and Applications, in print, 2017. 

[130] J. Wu, F. Liu, L. C. Jiao, X. Wang and B. Hou, "Multivariate Compressive Sensing for Image 

Reconstruction in the Wavelet Domain: Using Scale Mixture Models," IEEE Transactions on Image 

Processing, vol. 20, no. 12, pp. 3483-3494, Dec. 2011. 



510 A. DRAGANIĆ, I. OROVIĆ, S. STANKOVIĆ 

[131] J. Bobin, J. L. Starck, R. Ottensamer, "Compressed Sensing in Astronomy," IEEE Journal of Selected 

Topics in Signal Processing, vol. 2, no. 5, pp. 718-726, Oct. 2008. 

[132] J. Bobin, J.-L. Starck, “Compressed sensing in astronomy and remote sensing: a data fusion 

perspective”, In Proc. SPIE 7446, Wavelets XIII, 74460I (September 04, 2009). 

[133] R. G. Baraniuk, T. Goldstein, A. C. Sankaranarayanan, C. Studer, A. Veeraraghavan, M. B. Wakin, 

"Compressive Video Sensing: Algorithms, architectures, and applications," IEEE Signal Processing 

Magazine, vol. 34, no. 1, pp. 52-66, Jan. 2017. 

[134] I. Orovic, S. Park, S. Stankovic, "Compressive sensing in Video applications," In Proceedings of the 21st 

Telecommunications Forum TELFOR, Novembar, 2013. 

[135] L.-W. Kang, C.-S. Lu, “Distributed compressive video sensing,” In Proceedings of the IEEE 

International Conference on Acoustics, Speech and Signal Processing (ICASSP '09), 2009, pp. 1169–

1172. 

[136] X. Liao, K. Li, J. Yin, “Separable data hiding in encrypted image based on compressive sensing and 

discrete Fourier transform”, Multimedia Tools and Applications, pp. 1-15, 2016. 

[137] M. W. Fakhr, “Robust Watermarking Using Compressed Sensing Framework with Application to MP3 

Audio”, International Journal of Multimedia & Its Applications 2013. 

[138] X. Tang, Z. Ma, X. Niu, Y. Yang, "Compressive Sensing-Based Audio Semi-fragile Zero-Watermarking 

Algorithm," Chinese Journal of Electronics, vol. 24, no. 3, pp. 492-497, 07 2015. 

[139] I. Orovic, S. Stankovic, "Compressive Sampling and Image Watermarking," In Proceedings of the 55th 

International Symposium ELMAR 2013, Zadar, Croatia, Sept. 2013. 

[140] M. Orovic, T. Pejakovic, A. Draganic, S. Stankovic, "MRI watermarking in the Compressive Sensing 

context," In Proceedings of the 57th International Symposium ELMAR-2015, Zadar, Croatia, 2015. 

[141] I. Orovic, A. Draganic, S. Stankovic, "Compressive Sensing as a Watermarking Attack," In Proceedings 

of the 21st Telecommunications Forum TELFOR 2013, Novembar, 2013.