a method of compressed sensing imaging includes acquiring a sparse digital image b , said image comprising a plurality of intensities corresponding to an i - dimensional grid of points , initializing points , y ), wherein x is an element of a first expanded image x defined by b = rφ − 1 x , wherein r is a fourier transform matrix , φ is a wavelet transform matrix , y is a point in ∂ ) 2 ) 1 / 2 , ∇ i is a forward finite difference operator for a i th coordinate , and k is an iteration counter ; calculating a first auxiliary variable s from x - τ 1 + φ ⁢ ⁢ r * - b ) ) , wherein τ 1 , α are predetermined positive scalar constants , the sum is over all points n in x , and l * is an adjoint of operator l =; calculating a second auxiliary variable t n from y n + τ 2 l n φ − 1 x , wherein τ 2 is a predetermined positive scalar constant ; updating x from sign ) max |− τ 1 β }, wherein β is a predetermined positive scalar constant ; and updating y n from min ⁢ { 1 τ 2 ,  t n ( k )  2 } ⁢ t n ( k )  t n ( k )  2 .