Optimization of light fields in ghost imaging using dictionary learning

O p t i m i z a t i o n o f l ig ht fi el ds in gh os t i m ag in g u s i n g d ic ti on ar y l ea rn in g C H E N Y U H U , 1 , 2 Z H I S H E N T O N G , 1 , 2 Z H E N TA O L I U , 1 Z E N G F E N G H U A N G , 3 J I A N W A N G 3 , * A N D S H E N S H E N G H A N 1 , 2 1 Key Laborat or y f or Quantum Optics and Center f or Cold Atom Physics of CAS, Shanghai Institute of Optics and Fine Mechanics, Chinese Academy of Sciences, Shanghai 201800, China 2 Center of Materials Science and Optoelectr onics Engineering, U niv ersity of Chinese Academy of Sciences, Bei jing 100049, China 3 Sc hool of Data Science, F udan Univ er sity, Shanghai 200433, China * jian_wang@fudan.edu.cn Abstract : Ghost imaging (GI) is a no v el imaging technique based on the second-order cor relation of light fields. Due to limited number of samplings in practice, traditional GI methods often reconstruct objects with unsatisfactory quality . T o impro v e the imaging results, many reconstruction methods ha v e been dev eloped, y et the reconstruction quality is still fundamentall y restricted by the modulated light fields. In this paper , w e propose to impro v e the imaging quality of GI by optimizing the light fields, which is realized via matr ix optimization f or a learned dictionary incor porating the sparsity prior of objects. A closed-f orm solution of the sampling matrix, which enables successive sampling, is der iv ed. Through simulation and e xperimental results, it is sho wn that the proposed scheme leads to better imaging quality compared to the state-of-the-art optimization methods f or light fields, especially at a low sampling rate. © 2022 Optical Society of Amer ica under the ter ms of the OSA Open Access Publishing Agreement 1. Introduction As a nov el technique for optical imaging, ghost imaging (GI) was initially implemented with quantum-entangled photons two decades ago [1, 2]. In recent y ears, o wing to its realization with thermal light and other ne w sources [3 – 8], GI has gained new attention and dev eloped applications in various imaging areas, such as remote sensing [9, 10], imaging through scatter ing media [11, 12], spectral imaging [13, 14], photon-limited imaging [15, 16] and X -ra y imaging [17, 18]. Different from conv entional imaging techniq ues that are based on the first-order cor relation of light fields, GI extracts information of an object by calculating the second-order correlation betw een the light fields of the reference and the object arms [4, 5]. Theoretically , calculation of the second-order cor relation requires infinite number of samplings of the light fields at both arms. In practice, how e v er , the number of samplings is alwa y s finite, which often leads to reconstructed imag es of degraded signal-to-noise ratio (SNR) [19 – 21]. T o address this issue, much effor t has been made in designing more effectiv e reconstruction methods for GI. On the one hand, approaches impro ving the second-order correlation ha v e been proposed, which can increase the SNR of reconstructed images with theoretical guarantees [21 – 24]. On the other hand, b y e xploiting sparsity of the objects’ images in transf orm bases (e.g., wa v elets [25]), methods built upon the compressed sensing (CS) theor y [26 – 28] ha v e also been dev eloped [29, 30]. In g eneral, the CS based methods ha v e superior per f ormance ov er those relying on the second-order correlation, especially for imaging smooth objects [31]. While improving the reconstruction methods has greatly promoted the practical applications of GI, there has been increasing e vidence that the reconstruction quality of GI may be fundamentally restricted by the sampling efficiency [32, 33], i.e., how well inf ormation of objects is acquired in the samplings. T o enable a satisf actory reconstruction from limited number of samplings in GI, a natural w a y is to enhance the sampling efficiency . In fact, this can be realized b y optimizing the light fields of GI; see [32 – 35] and the references therein. Considering an orthogonal sparsifying basis, X u et al. [35] optimized the sampling matrix in order for their product, so called the equiv alent sampling matrix , to ha v e the minimum mutual coherence, which results in much refinement of the imaging quality . Though or thogonal basis is widely suitable f or sparse representation of natural images, f or images from a specific categor y , dictionar y learning [36, 37] usuall y produces much sparser representation coefficients, suggesting room f or fur ther impro v ements of the reconstruction quality . Motivated b y this, in this paper we propose to optimize the light fields of GI for a sparsifying basis obtained via dictionar y learning. By minimizing the mutual coherence of the equivalent sampling matr ix, the proposed scheme enhances the sampling efficiency and thus achie v es an impro v ed recons truction quality . In comparison with the state-of-the-art optimization methods f or light fields in GI, the super iority of our scheme is confirmed via both simulations and experiments. The main advantag es of the proposed scheme is summar ized as f ollo w s: • Inspired from some pre vious researc hes in CS [38, 39], we f ormulate the problem of minimizing the mutual coherence of the equivalent sampling matrix in GI as a Frobenius- norm minimization problem, which yields a closed-form solution that depends on the sparsifying basis only . T o the best of our kno w ledg e, the suggested solution of the light fields is the first closed-f orm result in the GI optimization field. • The proposed scheme enables successiv e samplings. In GI, successiv e samplings means that when more samplings are av ailable (or needed), one can simply augment new ro w s to the cur rentl y optimized sampling matr ix in order to f orm a ne w one, without the need to per f orm additional optimization ov er the entire matr ix. Such f eature can br ing great con v enience to the practical applications of GI and was not addressed in pre vious works. It is w or th mentioning that matrix optimization based on dictionary lear ning has also been studied in the CS literature, see, e.g., [38 – 40]. Ho we v er , the optimizations in [38, 39] w ere carr ied out ov er sampling matr ices of fix ed sizes, which does not allow successive samplings. Although Duar te ’ s method [40] dealt with the matr ix optimization problem of alterable sampling size, it is also not compatible to GI because of the demanding quantization accuracy . Moreov er, those methods all fail to cope with the non-negativ e nature of sampling matrices in GI. 2. The Proposed Scheme The detection process in GI can be approximatel y formulated as [30] y = Φx + n , (1) where y ∈ R M stands for the signal measured by the detector in the object arm, Φ ∈ R M × N is the sampling matr ix consisting of the light-field intensity distribution recorded by the detector in the ref erence arm, x ∈ R N signifies the object’ s information to be retr ie ved, and n denotes the detection noise. Let Ψ be the sparsifying basis obtained via dictionar y learning, in which x can be sparsely represented as x = Ψz , where z is the sparse coefficient vector . Also, consider the equiv alent sensing matr ix D : = ΦΨ . Then, (1) can be rewritten as y = Dz + n , (2) Evidences from the CS theory ha v e rev ealed that a matr ix D w ell preserving inf or mation of the sparse vector z guarantees a faithful reconstruction [26 – 28]. As a pow er ful measure of inf or mation preservation, the mutual coherence µ ( D ) characterizes how incoherent each column pairs in D are [26, 41], namely , µ ( D ) = max 1 ≤ i < j ≤ K    d i , d j      d i   2   d j   2 (3) with d i being the i -th column of D , K the number of columns in D and k · k 2 the ` 2 -norm. For its simplicity and ease of computation, the mutual coherence µ ( D ) has been widel y used to descr ibe the per f ormance guarantees of CS reconstruction algorithms. For e xample, e xact reco v ery of sparse signals via or thogonal matching pursuit (OMP) [42] is ensured by µ ( D ) < 1 2 k − 1 [43], where k is the sparsity le v el of input signals. In this w ork, with the aim of enhancing the sampling efficiency , we emplo y µ ( D ) as the objective function to be minimized in our optimization sc heme. In par ticular , our proposed scheme consists of the follo wing tw o main steps: • Firs tly , an o ver -complete dictionar y Ψ is learned from a collection of images, under the constraint that its first column has identical entries N − 1 / 2 , while each of the other columns has entries summing to zero. Specifically , giv en X = [ x ( 1 ) , x ( 2 ) , · · · , x ( K ) ] ∈ R N × L , in which each column is a reshaped v ector of the training image sample, the sparsifying dictionary Ψ ∈ R N × K is lear ned by sol ving the f ollo wing problem: min Ψ , Z k X − ΨZ k 2 F subject to Ψ 11 = · · · = Ψ N 1 = N − 1 / 2 , k z i k 0 ≤ T 0 , i = 1 , · · · , L , (4) where k · k F and k · k 0 are the Frobenius- and ` 0 -norm, respectivel y , Z =  z 1 , z 2 , · · · , z L  ∈ R K × L represents the sparse coefficient matr ix of training images, and T 0 denotes the predetermined sparsity le vel of v ectors z i . In this work, w e shall emplo y K -SVD as a representativ e method to per f or m the dictionary learning task, which results in simultaneous sparse representation of input images in the learned dictionar y Ψ . Readers are referred to [37] f or more details of the K -SVD method. • Secondly , the sampling matr ix Φ is optimized by minimizing the mutual coherence of the equiv alent sampling matr ix D . Put f or mall y , min Φ µ ( D ) subject to Φ i j ≥ 0 and D = ΦΨ . (5) The non-negativ e constraint Φ i j ≥ 0 is imposed due to the fact that the intensity of light fields is alwa ys non-negativ e. W e now proceed to solv e the optimization problem in (5). Without loss of generality , assume that matr ix D has ` 2 -normalized columns, that is, k d i k 2 = 1 f or i = 1 , · · · , K . Then, µ ( D ) = max 1 ≤ i < j ≤ K | h d i , d j i | . (6) T o optimize µ ( D ) , it suffices to minimize the off-diagonal entries of the Gram matrix D > D , each of which cor responds to the coherence between two different columns in D (i.e., ( D > D ) i j = |  d i , d j  | , i , j ). In particular, w e w ould like the Gram matr ix to be as close to the identity matrix as possible, namely , Ψ > Φ > ΦΨ ≈ I . Since replacing the identity matr ix with Ψ > Ψ yields a sampling matrix robust to the sparse representation er ror of images [44], w e propose to optimize Φ via min Φ   Ψ > Φ > ΦΨ − Ψ > Ψ   2 F . (7) By multiplying Ψ and Ψ > on the left- and r ight-hand sides of both ter ms inside the Frobenius norm, respectiv ely , one has min Φ   ΨΨ > Φ > ΦΨΨ > − ΨΨ > ΨΨ >   2 F . (8) After substituting ΨΨ > with its eig env alue decomposition VΛV > , and also denoting W : = ΛV > Φ > , (8) can be re written as min W   VWW > V > − VΛ 2 V >   2 F , (9) or equiv alently , min W      Λ 2 − M Õ i = 1 w i w > i      2 F where W = [ w 1 , · · · , w M ] . (10) Denoting Λ = [ r 1 , · · · , r N ] , (10) fur ther becomes min W      N Õ j = 1 r j r > j − M Õ i = 1 w i w > i      2 F . (11) Clearl y problem (11) has the solution c W = Λ > 1 , where Λ 1 is the matr ix consisting of the first M columns of Λ , which is obtained by setting w k = r k , k = 1 , · · · , M . R ecalling that W : = ΛV > Φ > , the optimized sampling matrix Φ can be simply calculated as b Φ = c W T  Λ − 1  T V T = Λ 1  Λ − 1  T V T = h I M × M 0 i       V T 1 V T 2       = V T 1 , (12) where matr ix V 1 consists of the first M columns of V . When more samplings become av ailable, interestingl y , it suffices to update b Φ b y augmenting more row s of V > to the previous one, thereb y enabling a successiv e sampling. As af orementioned, the feature of successiv e sampling is of vital impor tance to the practical applications of GI. Further more, due to the fact that the intensity of light fields is alwa ys non-negativ e, additional treatments are needed to make sure that elements of the sampling matrix are non-negativ e (NN). T o the end, we propose a NN lifting, which adds a constant matr ix to the optimized matr ix b Φ in (12) as b D =  b Φ + c 1 M × N  Ψ , (13) where 1 M × N is an M -by - N matrix with entr ies being ones and c : = ( − min i , j b Φ i j if min i , j b Φ i j < 0 , 0 if min i , j b Φ i j ≥ 0 . (14) As af orementioned, the first column of Ψ has identical entr ies and other columns hav e entr ies summing to zero. Thus, b D = b ΦΨ + c N − 1 / 2  1 M × 1 , 0 , . . . , 0 | {z } M ×( K − 1 )  . (15) It can be noticed that after the NN lifting, b D and b ΦΨ differ onl y in the firs t column. Ne v er theless, the mutual coherence µ ( b D ) is not much affected b y the NN lifting, as confir med b y our e xtensiv e empirical test. 3. Results 3.1. Simulations T o ev aluate the effectiveness of the proposed optimization scheme, both simulations and e xper iments are per f ormed. MNIST handwritten digits of size 28 × 28 pixels [45] are chosen to (a) (b) Fig. 1. (Left) A subset of atoms in the lear ned dictionar y Ψ ; (r ight) a subset of the optimized lgiht-field intensity distributions T able 1. A summary of test methods Method Sampling matrix Dictionary Φ i j ≥ 0 Proposed Eq. (12) K -S VD NN-lifting Gaussian Random Gaussian K -S VD NN-lifting Duarte [40] Matrix optimization K -S VD NN-lifting X u [35] Matrix optimization DCT Zero-f orcing NGI [24] Random Gaussian None NN-lifting be the imaging objects, and the dictionar y Ψ is learned based on 20,000 digits randomly selected from the training set. Moreo ver , the optimized sampling matr ix b Φ is obtained from (12), f ollo w ed by the NN lifting. A subset of atoms in the lear ned dictionar y Ψ and the optimized light-field intensity distr ibutions b Φ are sho wn as Fig. 1(a) and Fig. 1(b), respectiv ely . For comparativ e pur poses, our simulation includes other four methods: 1) Gaussian method, 2) Duar te ’ s method [40], 3) Xu ’ s method [35] and 4) nor malized GI (NGI) method [24]. T able 1 giv es a br ief summar y of the methods under test. In the Gaussian method, the sampling matr ices are random Gaussian matrices, whose entr ies are drawn independently from the standard Gaussian distribution ( Φ i j ∼ N ( 0 , 1 ) ). T o meet the NN constraint Φ i j ≥ 0 of GI, the matr ices Φ ’ s generated from the Gaussian and Duar te ’s methods are also inflicted with the NN lifting. F or the Gaussian, Duar te ’ s and our proposed methods, the images are retriev ed in two steps. F irstl y , the sparse coefficient v ector ˆ z of image under the learned dictionary Ψ is obtained by solving the ` 0 -minimization problem: ˆ z = arg min z k b Dz − y k 2 2 subject to k z k 0 ≤ T 0 (16) Gaussian Duarte Proposed Xu Ground Truth SR= 0.10 SR= 0.20 SR= 0.51 NGI Gaussian Duarte Proposed Xu Ground Truth SR = 0.15 SR = 0.38 SR = 0.64 NGI (a) 0 0.2 0.4 0.6 0.8 1 Sampling rate 0 5 10 15 20 25 30 PSNR (dB) Proposed Gaussian Xu Duarte NGI (b) 0 0.2 0.4 0.6 0.8 1 Sampling rate 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 SSIM Proposed Gaussian Xu Duarte NGI (c) Fig. 2. Simulation results with sampling matrix of high accuracy . Fig. 2(a) show s the reconstructed imag es reconstructed via different methods under different SR. Fig. 2(b) and 2(c) illustrate the PSNR and SSIM of reconstructed imag es via different methods as a function of SR, respectivel y . via the OMP algor ithm [42]. Secondly , the object’ s imag e is reconstructed as b x = Ψ b z . (17) For X u’ s method, the Discrete Cosine T ransform (DCT) basis is chosen as the or thogonal basis. And the images in methods 3) and 4) are reconstructed via approaches proposed in their corresponding ref erences. In our simulation, w e first adopt matr ices with entr ies of double-type in MA TLAB. Fig. 2 sho ws the simulation results of different methods. In Fig. 2(a), the reconstructed images at the different sampling ratios (SR’ s) (i.e., SR = 0 . 10 , 0 . 20 , and 0 . 51 ) are display ed, where the SR is computed b y dividing the number of samplings by the number of imag e pix els. Fig. 2(b) and 2(c) depict the peak signal-to-noise ratio (PSNR) and structural similar ity (SSIM) inde x of the reconstructed images as functions of the SR, respectiv ely . Given the ref erence image X and the reconstructed image Y , the PSNR and SSIM are defined as f ollow s, MSE ( X , Y ) = 1 mn m Õ i = 1 n Õ j = 1 [ X ( i , j ) − Y ( i , j )] 2 , (18a) PSNR ( X , Y ) = 10 log 10  B 2 MSE ( X , Y )  , (18b) SSIM ( X , Y ) = ( 2 µ X µ Y + c 1 ) ( 2 σ XY + c 2 )  µ 2 X + µ 2 Y + c 1   σ 2 X + σ 2 Y + c 2  , (18c) where the pix el size of imag e is m × n , B denotes the dynamic range of image pixels, which takes the v alue 255 in this paper , ( µ X , σ X ) and ( µ Y , σ Y ) are the means and variances of X and Y , respectiv ely , σ XY is the cov ar iance of X and Y , c 1 = ( 0 . 01 B ) 2 and c 2 = ( 0 . 03 B ) 2 . These tw o metrics measure the difference and similar ity betw een the reconstructed imag es and the original ones, respectiv ely . For each SR under test, the PSNR and SSIM are av eraged o ver 500 reconstructed digit images to plot the cur v es. By compar ing them, the reconstruction quality of different methods are compared empirically . From Fig. 2, it can be obser v ed that the reconstruction quality of the proposed scheme is unif or ml y better than the other methods under test. In par ticular , it achie ves 2 dB to 4 dB gain of PSNR and up to 10 % higher SSIM ov er the Gaussian method, owing to the optimized light fields. Compared to Xu ’ s method [35], the Gaussian method ha v e a notable adv antage in the lo w SR region, which gradually conv erg es as the SR approaches one. The performance gap is mainly attr ibuted to the utilization of dictionar y learning that can better incor porate the sparsity prior of imag es. Among all tes t methods, the PSNR and SSIM of the NGI method lie in the lo west lev el. This is mainly due to the image noise, which often happens to the cor relation-based reconstruction methods of GI, especially when the SR is low . It can also be observed from Fig. 2 that Duar te ’s method [40] per f or ms comparably with the proposed method in the lo w SR region, but deteriorating dramatically when the SR increases. Such phenomenon seems unreasonable at first glance, but can be inter preted from the condition number perspective. T o be specific, the sampling matrix Φ of Duar te ’ s method has larg er condition number as the SR increases (see detailed e xplanations in Footnote 7 of [40]). Thus, when Φ multiplies with the representation er ror e = x − Ψz , it could significantl y amplify this er ror , and ev entually degrade the reconstruction quality . Indeed, in this case one would need to reconstruct the sparse vector z from the samplings y = Φx + n = ΦΨz + Φe + n = Dz + ( Φe + n ) , which can be difficult since the larg ely amplified er ror Φe essentiall y becomes par t of noise f or the reconstruction. In practice, detectors measure the intensity signals with quantization, which means that Φ is actually a quantized sampling matr ix. Thus, we also simulate the case where the sampling matrices are quantized to 8 -bit of precision and plot the results in Fig. 3. Similarl y , Fig. 3(b) and 3(c) show the curves of PSNR and SSIM with av eraged values o ver 500 reconstruction trails, respectiv ely . W e observe that the o verall beha vior is similar to that of Fig. 2(b) and 2(c) ex cept that both the PSNR and SSIM curves of Duarte’ s method [40] fluctuate in the lo w est le vel f or the whole SR region. Accordingl y , Duar te ’ s method [40] also f ails to retr ie ve the imag es in Fig. 3(a). This is mainly because Duar te ’ s method [40] is demanding in the quantization accuracy . When larg e quantization er rors are introduced, sparse coefficients of the test images may not be correctly calculated by the reconstruction algorithm. The PSNR and SSIM curves of the NGI method lie in a lo w le v el similar to that of Duar te ’ s method, and the images are onl y v aguely reconstructed, as shown in Fig.3(a). W e also obser v e that the per f ormance of Xu ’s method [35] becomes worse in the quantized case, although the images can still be retr ie v ed. Overall, our proposed method per f or ms the best f or both the high accurate as well as the quantized scenar ios. Gaussian Duarte Proposed Xu Ground Truth SR= 0.10 SR= 0.20 SR= 0.51 NGI Gaussian Duarte Proposed Xu Ground Truth SR = 0.15 SR = 0.38 SR = 0.64 NGI (a) 0 0.2 0.4 0.6 0.8 1 Sampling rate 0 5 10 15 20 25 30 PSNR (dB) Proposed Gaussian Xu Duarte NGI (b) 0 0.2 0.4 0.6 0.8 1 Sampling rate 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 SSIM Proposed Gaussian Xu Duarte NGI (c) Fig. 3. Simulation results with 8 -bit quantized sampling matrix. Fig. 3(a) sho ws the reconstructed images via different methods under different SR. F ig. 3(b) and 3(c) sho w PSNR and SSIM of the reconstructed images via different methods as a function of SR, respectiv ely . 3.2. Experimental Results W e also experimentally compare the proposed method with the Gaussian method, Duar te ’s method [40] and NGI method [24]. The schematic diag ram of experimental setup is shown as Fig. 4. The light-field patterns are first display ed on the digital micro-mir ror device (DMD) after preloaded via the computer . Ne xt, light from a light-emitting diode (LED) source is modulated b y a Kohler illumination system to be ev enly incident on DMD. The light reflected b y DMD is then projected onto the imaging object by a lens system. Finall y , the whole light reflected from the object is collected b y the lens and measured by the detector. In our e xper iments, the light-field patterns are displa y ed at a rate of 10Hz to av oid frame dropping of the detector, so that the sampling procedure lasts f or one minute or so. After the samplings, the subsequent imag e  !"#"$#%& '()*"&+ ,**-./01#/%0 234"$# 56! 5"07+8 5"07+9 !:! ;%.<-#"& Fig. 4. Schematic diagram of experimental setup. Light-field patter ns are generated b y DMD and projected onto the object, afterwards collected b y the detector . T able 2. Running time of test methods Method Matrix Optimization (sec.) Reconstruction (sec.) Proposed 0.365 0 . 037 to 0 . 150 Gaussian – 0 . 037 to 0 . 158 Duarte [40] 0 . 345 0 . 039 to 0 . 157 X u [35] 90 . 18 (f or 100 iterations) 0 . 028 to 0 . 360 NGI [24] – 0 . 001 to 0 . 007 retriev al steps f or each method are the same as those in the simulation test. The reconstruction was car ried out on an industrial computer with 32 GB RAM and Intel(R) Core(TM)-I7 2600 CPU @ 3 . 4 GHz, and the consuming time of matr ix optimization and reconstruction f or different methods is specified in T able 3.2. The reconstruction time f or each method is giv en as a time period, since it varies according to the sampling rate. The compar ison of reconstructed images b y different methods is sho wn as Fig. 5(a), where the ground tr uth is obtained b y pixel-wise detection and ser v e as a reference imag e. By compar ing the reconstructed images with the ground tr uth, again, the PSNR and SSIM are calculated and plotted in Fig. 5(b) and 5(c), respectivel y . Ov erall, the experimental results demonstrate that the reconstructed quality of the proposed optimization scheme is superior to that of other methods under test, which w ell matches our simulation results. Gaussian Duarte Proposed Ground Truth SR = 0.15 SR = 0.32 SR = 0.64 NGI (a) 0 0.2 0.4 0.6 0.8 Sampling rate 6 8 10 12 14 16 18 20 22 PSNR (dB) Proposed Gaussian Duarte NGI (b) 0 0.2 0.4 0.6 0.8 Sampling rate 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 SSIM Proposed Gaussian Duarte NGI (c) Fig. 5. Exper imental results. F ig. 5(a) show s the reconstructed images via different methods under different SR. Fig. 5(b) and 5(c) show PSNR and SSIM of the reconstructed images via different methods as a function of SR, respectivel y . 3.3. Discussions W e would like to point out some interes ting points that ar ise from the simulation and experimental results. • Firs tly , the superior ity of the proposed method is mainl y owing to tw o factors: i) optimization of sampling matr ix and ii) dictionar y learning. Indeed, the proposed optimization scheme outperforms the Gaussian method, ev en though they share the same spasifying basis obtained b y dictionar y lear ning. This is because our method essentially perf orms a “global” optimization of light fields that incor porates the image statis tics captured in the dictionar y learning process, thereb y enhancing the sampling efficiency . Besides, the improv ement of the proposed method ov er Xu ’ s method can be attr ibuted to the use of both dictionar y learning and our optimized sampling matrix. • Secondly , the PSNR and SSIM cur v es of the proposed method tend to wards flat after the SR reaches a critical value. This in turn implies that the inherent inf or mation of the imaging object acquired at this v er y SR value already suffices to produce a satisfactory reconstruction. The critical SR can thus be utilized to e valuate the capability of information acquisition and also allow s the compar ison of different approaches. • Thirdl y , w e would like to point out a practical limitation of the proposed scheme in handling images of large size due to the use of dictionary learning. Specifically , while dictionar y learning in our scheme can br ing in some per f ormance gain, it is usually demanding in the requirements of storage and computational cost. Thus the patch size used in dictionary learning should not be large, which, how ev er , poses a limitation to the imag e size that w e can handle. Ne vertheless, efficient dictionar y learning methods dealing with imag es of larg er scales ha v e recently been proposed [46 – 48], in which the handled image size can go be y ond 64 × 64 pixels. T o demonstrate the effectiveness of the proposed method f or imag es of larg er size, w e car ry out simulations ov er the LFWcrop database [49], which consists of more than 13 , 000 imag es of 64 × 64 pix els. The dictionary is trained offline using the algor ithm in [47] with 12 , 000 images, which takes about 20 hours in our industrial computer . The results of the proposed method and the Gaussian method, which in v ol v e dictionary learning, are sho wn as Fig. 6 f or comparison. • Finall y , w e mention that if one wish to deal with images of e v en larg er size such that e xisting dictionar y learning methods f ail to handle or cannot lear n the images offline, then the proposed light-field optimization scheme can still be applied by using explicit dictionaries (e.g., Cropped W av elets [47]) to incorporate the sparse pr ior of imag es. 4. Conclusion In this paper , an optimization scheme of light fields has been proposed to improv e the imaging quality of GI. The k e y idea is to minimize the mutual coherence of the equivalent sampling matrix in order to enhance the sampling efficiency . A closed-f or m solution of the sampling matrix has been der iv ed, which enables successiv e sampling. Simulation and e xperimental results ha v e sho wn that the proposed scheme is v ery effectiv e in impro ving the reconstruction quality of images, compared to the state-of-the-art methods f or GI. The proposed scheme can thus be used to imaging specific targ ets with higher quality . W e w ould also like to point out a technical limitation in our scheme. Recall that we ha v e employ ed a NN lifting to cope with the constraint Φ i j ≥ 0 . This operation, how ev er , ma y sev erely influence the incoherence of the equiv alent sampling matr ix in the wors t case, though such situation rarely happens as confir med b y our empir ical test. Deriving analytical results can better address the non-negativ e issue while also w ould require a bit more effort, and our future work will be directed tow ards this a venue. Funding N ational Ke y Researc h and Dev elopment Program of China (2017YFB0503303, 2017YFB0503300); National Natural Science Foundation of China (NSFC) (11627811). References 1. T . B. Pittman, Y . H. Shih, D. V . Strekalo v , and A. V . Sergienk o, “Optical imaging by means of 2-photon quantum entanglement, ” Phy s. Re v. A At. Mol. & Opt. Phy s. 52 , R3429–R3432 (1995). 2. D. Strekalov , A. Sergienko, D. Kly shko, and Y . Shih, “Observation of two-photon “ghost ” inter f erence and diffraction, ” Phy s. Re v. Lett. 74 , 3600 (1995). 3. R. S. Bennink, S. J. Bentley , and R. W . Boyd, ““two-photon ” coincidence imaging with a classical source, ” Ph ys. Re v. Lett. 89 , 113601 (2002). Gaussian Proposed Ground Truth SR = 0.15 SR = 0.32 SR = 0.64 (a) 0 0.2 0.4 0.6 0.8 1 Sampling rate 15 17 19 21 23 25 PSNR (dB) Proposed Gaussian (b) 0 0.2 0.4 0.6 0.8 1 Sampling rate 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 SSIM Proposed Gaussian (c) Fig. 6. Simulation results of LFWcrop face images. Fig. 6(a) show s the reconstructed imag es via different methods under different SR. Fig. 6(b) and 6(c) show PSNR and SSIM of the reconstructed images via different methods as a function of SR, respectivel y . 4. A. Gatti, E. Brambilla, M. Bache, and L. A. Lugiato, “Ghost imaging with ther mal light: compar ing entanglement and classicalcorrelation, ” Ph ys. Re v. Lett. 93 , 093602 (2004). 5. J. Cheng and S. Han, “Incoherent coincidence imaging and its applicability in x-ray diffraction, ” Phy s. Re v. Lett. 92 , 093903 (2004). 6. D. Zhang, Y .-H. Zhai, L.-A. W u, and X.-H. Chen, “Correlated two-photon imaging with true ther mal light, ” Opt. Lett. 30 , 2354–2356 (2005). 7. R. I. Khakimov , B. Henson, D. Shin, S. Hodgman, R. Dall, K. Baldwin, and A. T r uscott, “Ghost imaging with atoms,” Nature 540 , 100 (2016). 8. S. Li, F . Cropp, K. Kabra, T . Lane, G. W etzstein, P . Musumeci, and D. Ratner , “Electron ghost imaging, ” Phy s. Re v. Lett. 121 , 114801 (2018). 9. M. Malik, O. S. Magaña-Loaiza, and R. W . Boyd, “Quantum-secured imaging, ” Appl. Phys. Lett. 101 , 241103 (2012). 10. C. Zhao, W . Gong, M. Chen, E. Li, H. W ang, W . Xu, and S. Han, “Ghost imaging lidar via sparsity constraints, ” Appl. Ph ys. Lett. 101 , 141123 (2012). 11. W . Gong and S. Han, “Cor related imaging in scattering media, ” Opt. Lett. 36 , 394–6 (2011). 12. M. Bina, D. Magatti, M. Molteni, A. Gatti, L. A. Lugiato, and F . Ferr i, “Backscattering differential ghost imaging in turbid media, ” Phy s. Re v. Lett. 110 , 083901 (2013). 13. Y . W ang, J. Suo, J. F an, and Q. Dai, “Hyperspectral computational ghost imaging via temporal multiple xing,” IEEE Photonics T echnol. Lett. 28 , 288–291 (2016). 14. Z. Liu, S. T an, J. Wu, E. Li, X. Shen, and S. Han, “Spectral camera based on ghost imaging via sparsity constraints, ” Sci. R eports 6 , 25718 (2016). 15. P . A. Mor ris, R. S. Aspden, J. E. Bell, R. W . Bo yd, and M. J. Padgett, “Imaging with a small number of photons, ” Nat. Commun. 6 , 5913 (2015). 16. X. Liu, J. Shi, X. W u, and G. Zeng, “Fas t first-photon ghost imaging,” Sci. Reports 8 , 5012 (2018). 17. D. P elliccia, A. Rack, M. Scheel, V . Cantelli, and D. M. Paganin, “Experimental x-ray ghost imaging, ” Phy s. Re v. Lett. 117 , 113902 (2016). 18. H. Y u, R. Lu, S. Han, H. Xie, G. Du, T . Xiao, and D. Zhu, “Fourier -transform ghost imaging with hard x ray s, ” Ph ys. Re v. Lett. 117 , 113901 (2016). 19. X. Shen, Y . Bai, T . Qin, and S. Han, “Experimental inv estigation of quality of lensless ghos t imaging with pseudo-thermal light, ” Chin. Phy s. Lett. 25 , 3968 (2008). 20. B. I. Erkmen and J. H. Shapiro, “Signal-to-noise ratio of gaussian-state ghost imaging,” Phy s. Re v. A 79 , 023833 (2009). 21. F . Ferr i, D. Magatti, L. Lugiato, and A. Gatti, “Differential ghost imaging,” Phy s. Re v. Lett. 104 , 253603 (2010). 22. W . Gong and S. Han, “ A method to impro v e the visibility of ghost images obtained by thermal light, ” Phy s. Lett. A 374 , 1005–1008 (2010). 23. G. Brida, M. Genov ese, and I. R. Berchera, “Exper imental realization of sub-shot-noise quantum imaging,” Nat. Photonics 4 , 227 (2010). 24. B. Sun, S. S. W elsh, M. P . Edgar , J. H. Shapiro, and M. J. Padgett, “Normalized ghost imaging, ” Opt. Express 20 , 16892–16901 (2012). 25. M. Antonini, M. Barlaud, P . Mathieu, and I. Daubechies, “Image coding using wa v elet transf or m. ” IEEE Trans. Image Process. 1 , 205–220 (1992). 26. D. L. Donoho, “Compressed sensing, ” IEEE Trans. Inform. Theor y 52 , 1289–1306 (2006). 27. E. J. Candès and T . T ao, “Decoding by linear programming, ” IEEE Trans. Inform. Theor y 51 , 4203–4215 (2005). 28. E. J. Candès and T . T ao, “Near -optimal signal recov ery from random projections: U niv ersal encoding strategies?” IEEE T rans. Inf orm. Theor y 52 , 5406–5425 (2006). 29. O. Katz, Y . Bromberg, and Y . Silberberg, “Compressiv e ghost imaging, ” Appl. Phys. Lett. 95 , 131110 (2009). 30. S. Han, H. Y u, X. Shen, H. Liu, W . Gong, and Z. Liu, “ A revie w of ghost imaging via sparsity constraints, ” Appl. Sci. 8 , 1379 (2018). 31. W . Gong, Z. Bo, E. Li, and S. Han, “Exper imental inv estig ation of the quality of ghost imaging via sparsity constraints, ” Appl. Opt. 52 , 3510–3515 (2013). 32. M. Chen, E. Li, and S. Han, “ Application of multi-cor relation-scale measurement matr ices in ghost imaging via sparsity constraints, ” Appl. Opt. 53 , 2924–2928 (2014). 33. S. M. Khamoushi, Y . Nosrati, and S. H. T av assoli, “Sinusoidal ghost imaging, ” Opt. Lett. 40 , 3452–3455 (2015). 34. E. Li, M. Chen, W . Gong, H. Y u, and S. Han, “Mutual information of ghost imaging sys tems,” Acta Opt. Sinica 33 , 93–98 (2013). 35. X. Xu, E. Li, X. Shen, and S. Han, “Optimization of speckle patterns in ghost imaging via sparse constraints by mutual coherence minimization,” Chin. Opt. Lett. 13 , 071101 (2015). 36. B. A. Olshausen and D. J. Field, “Natural image statistics and efficient coding,” Netw ork: Comput. Neural Sys t. 7 , 333–339 (1996). 37. M. Aharon, M. Elad, and A. Bruckstein, “K -svd: An algor ithm for designing ov ercomplete dictionaries f or sparse representation, ” IEEE T rans. Signal Process. 54 , 4311–4322 (2006). 38. M. Elad, “Optimized projections for compressed sensing, ” IEEE T rans. Signal Process. 55 , 5695–5702 (2007). 39. V . Abolghasemi, S. Ferdo wsi, and S. Sanei, “ A gradient-based alter nating minimization approach f or optimization of the measurement matr ix in compressive sensing, ” Signal Process. 92 , 999–1009 (2012). 40. J. M. Duar te-Carva jalino and G. Sapiro, “Lear ning to sense sparse signals: Simultaneous sensing matr ix and sparsifying dictionary optimization, ” IEEE T rans. Image Process. 18 , 1395–1408 (2009). 41. D. L. Donoho and M. Elad, “Optimall y sparse representation in g eneral (nonorthogonal) dictionar ies via l 1 minimization, ” Proc. Natl. Acad. Sci. 100 , 2197–2202 (2003). 42. Y . C. Pati, R. Rezaiifar, and P . S. Kr ishnaprasad, “Or thogonal matching pursuit: recursive function appro ximation with applications to wa v elet decomposition, ” in Proceedings of 27th Asilomar Confer ence on Signals, Sys tems and Computers, (IEEE, 1993), pp. 40–44. 43. J. A. Tropp, “Greed is good: Algor ithmic results for sparse appro ximation, ” IEEE Trans. Inform. Theor y 50 , 2231–2242 (2004). 44. N. Cleju, “Optimized projections for compressed sensing via rank -constrained nearest correlation matr ix, ” Appl. Comput. Harmon. Analy sis 36 , 495–507 (2014). 45. L. Deng, “ The mnist database of handwritten digit imag es f or machine lear ning research [best of the web], ” IEEE Signal Process. Mag. 29 , 141–142 (2012). 46. L. Le Magoarou and R. Gr ibon val, “Chasing butter flies: In search of efficient dictionar ies, ” in 2015 IEEE International Confer ence on Acous tics, Speec h and Signal Processing (ICASSP), (IEEE, 2015), pp. 3287–3291. 47. J. Sulam, B. Ophir, M. Zibulev sky , and M. Elad, “T rainlets: Dictionar y lear ning in high dimensions,” IEEE T rans. Signal Process. 64 , 3180–3193 (2016). 48. C. F . Dantas, M. N. Da Costa, and R. da Rocha Lopes, “Learning dictionar ies as a sum of kroneck er products, ” IEEE Signal Process. Lett. 24 , 559–563 (2017). 49. C. Sanderson and B. C. Lov ell, “Multi-region probabilistic histograms for robust and scalable identity inf erence, ” in International Conf erence on Biometrics, (Spr ing er, 2009), pp. 199–208.