In order to improve spatial resolution, conventional optical telescopes were usually designed with large apertures, leading to large volume, heavy weight and high power consumption [1, 2]. With the rapid development of photonic integrated circuit (PIC) technology and interferometric imaging technology, the concept known as the segmented planar imaging detector for electro-optical reconnaissance (SPIDER) has been proposed, which is of great importance for its potential applications in astronomy and space science [3]. For photonic integrated interferometric imaging (PIII), the large optical components, devices and mounting structures in conventional system are replaced by a densely packed interferometer array realized in an extremely thin layer of lenslets and PIC components, which offers substantial reductions in system size, weight, and cost. Using on-chip components, the high quality image is reconstructed by collecting and processing the interferometric information of Fourier domain at the output end [4–8]. At present, PIII technology is at the primary research stage and has not been transformed into production. Relevant research institutions only carried out some demonstrative experiments and theoretical investigations [9]. Key technologies concerning PIII need to be addressed, which include the principle of optical interferometric imaging and signal transmission mechanism, spatial frequency coverage and image reconstruction of the scene.

This paper focuses on spatial frequency coverage and image reconstruction algorithm of PIII. The former aspect determines the direct restored image, and the latter one influences the reconstructed image quality. Section 2 discusses the basic theory and concept of PIII. Section 3 describes a novel lenslet array arrangement and lenslet pairing approach that effectively improve spatial frequency coverage of the scene. Meanwhile, it introduces an image reconstruction approach based on compressive sensing (CS). Section 4 shows image reconstruction from a simulated PIII observation. Section 5 has the discussion. The conclusion is given in Section 6.

PIII system is to put a linear array of lenslets onto a PIC card, and the PIC cards are mounted as radial spokes on a disc. PIII system consists of three working modules: signal receiving module, signal transmission module and signal processing module. The signal receiving module mainly includes the lenslet array and its components. The signal transmission module consists of optical waveguides, arrayed waveguide gratings (AWGs), phase shifters, and multi-mode interferometers (MMIs). The signal processing module contains balanced detectors and digital signal processing. PICs transform the target intensity distribution into Fourier domain information. By scanning the phase delay on one of the arms in the MMI, the balanced detectors receive interferometric information, the phase and the amplitude of the complex visibility will be calculated. Then the image is restored by processing Fourier transform of the complex visibility information [10, 11]. Figure 1 shows the traditional structure of system and the working principle of a PIC.

Figure 1. The traditional structure of system and the working principle of a photonic integrated circuit (PIC) [12,13].

PIII is a new implementation approach of optical interferometric technology. According to the Van Cittert–Zernike theory [14], the measured complex visibility y between the paired lenslets (x₁, y₁) and (x₂, y₂) can be expressed as

(1)$$y(u,v)={\displaystyle \int {\displaystyle \int x(\xi ,\eta )\rho \left(\xi ,\eta \right)\exp \left[-j2\pi (u\xi +v\eta )d\xi d\eta \right]}},$$

where ρ is the lenslet coupling efficiency projected onto the object plane. (ξ, η) is the coordinate of image intensity x. (u, v) is the spatial frequency corresponding to the distance vector between two lenslets (x₁, y₁) and (x₂, y₂):

(2)$$u=\frac{\Delta x}{\lambda z}=\frac{x_2-x_1}{\lambda z},$$

(3)$$u=\frac{\Delta y}{\lambda z}=\frac{y_2-y_1}{\lambda z},$$

where Δx, Δy are the distances of the chosen pair of lenslets in the x and y directions. λ is the current working wavelength of the system, and z is the scene distance.

For the PIII system, the paired beam passing through the phase shifters would interfere in MMI, and the interferometric information would be measured by the balanced detectors. As shown in Fig. 1, the optical signals input into the 90° optical hybrids can be expressed as E_S and E_R respectively, and the output signals of the 90° optical hybrids can be expressed as [15]

(4)$$E_{out1}=\frac{1}{2}\left(E_S+E_R\right),$$

(5)$$E_{out2}=\frac{1}{2}\left(E_S-E_R\right),$$

(6)$$E_{out3}=\frac{1}{2}\left(E_S+j{E}_R\right),$$

(7)$$E_{out4}=\frac{1}{2}\left(E_S-j{E}_R\right).$$

The output electrical signals I and Q of balanced detectors can be derived as follows:

(8)$$I=I_1-I_2=E_{out1}{E}_{out1}^{*}-E_{out2}{E}_{out2}^{*}=\mathrm{Re}(E_S{E}_R^{*}),$$

(9)$$Q=Q_1-Q_2=E_{out3}{E}_{out3}^{*}-E_{out4}{E}_{out4}^{*}=\mathrm{Im}(E_S{E}_R^{*}),$$

where E* is the conjugate of E. From equations above, the phase φ and amplitude E of the measured complex visibility can be expressed as

(10)$$E={\left[I^2+Q^2\right]}^{\frac{1}{2}},$$

(11)$$\phi =\arctan \left(\frac{Q}{1}\right).$$

The measured complex visibility y can be obtained by using the output information from the system as follows:

(12)$$y(u,v)=I+jQ=E\exp (-j\phi ).$$

According to the Van Cittert-Zernike theorem, the complex visibility can be found by a two-dimensional Fourier transformation of the intensity distribution across the source [14]. Therefore, the restored image can be expressed as

(13)$$x(\xi ,\eta )=F^{-1}\left[y\left(u,v\right)\right].$$

3.1. Structure and Spatial Frequency Coverage

We designed a novel arrangement of lenslet array as shown in Fig. 2. This lenslet array arrangement is composed of thirty-seven long interference arms. Between two adjacent long interference arms, pairs of short interference arms are evenly distributed. Each long interference arm consists of twenty-nine lenslets (red color) and short interference arm contains fifteen lenslets (yellow color). The interval between two adjacent lenslets is equal to the diameter d of the lenslet. Considering the system is composed of P₁ long interference arms, the angle between adjacent long interference arms is α₁ = 2π / P₁. Considering the system is composed of P₂ short interference arms, the angle between adjacent short interference arms is α₂ = 2π/(P₁ + P₂).

Figure 2. The schematic diagram of the new system: (a) the new structure of system and (b) the novel arrangement of lenslet array.

The spatial frequency coverage is influenced by the lenslet array arrangement and pairing approach of the lenslets in each radial spoke. As shown in Fig. 3, 1 × 2 beam splitters are designed to follow each lenslet except the middle and last lenslets. Two 1 × 3 beam splitters are designed to follow the middle and last lenslets respectively. Erbium-doped fibre amplifiers (EDFAs) are designed to amplify the signal intensity output from each 1 × 3 beam splitter to maintain equal signal intensity corresponding to the 1 × 2 beam splitter. The PIII system takes advantage of AWGs to divide broadband light into multiple narrow spectral channels to improve spatial frequency coverage. The paired narrow spectral beam passing through the phase shifters would be interfered in MMI, and the interferometric information would be measured by the balanced detectors. By setting up 1 × 2 and 1 × 3 beam splitters following the lenslets, all the lenslets can be paired twice in each radial spoke. For example, assuming each radial spoke contains 5 lenslets, the pairing approach is (1, 5), (2, 4), (3, 3), (1, 4), (2, 3), (5, 5).

Figure 3. The schematic diagram of the optimized photonic integrated device. (a) The working principle of a PIC, (b) the 1 × 2 beam splitter, (c) the 1 × 3 beam splitter, and (d) 5 lenslets pairing.

Based on setting up 1 × 2 and 1 × 3 beam splitters following the lenslets in Fig. 3, the novel pairing details of lenslets can be described as follows. Assuming each interference arm contains N lenslets, the lenslets were paired firstly following (1, N − 1), (2, N − 2), (3, N − 3)... ((N − 1) / 2, (N + 1) / 2), (N, N) as shown in Fig. 4(a). After that, the lenslets were secondly matched as (1, N), (2, N − 1), (3, N − 2)... ((N − 1) / 2, (N + 3) / 2), ((N + 1) / 2, (N + 1) / 2), as shown in Fig. 4(b). Zero spatial frequency sampling with the image average brightness information is obtained by pairing lenslets N and (N + 1) / 2 with itself [16]. Through this novel lenslet pairing approach, the length of the baseline is 0, d, 2d, 3d, …, (N − 1)d. The spatial frequency coverage is approximately twice than that of traditional lenslet pairing [17]. Therefore, this novel lenslet pairing approach can effectively improve the spatial frequency coverage of the scene. Consequently, low frequency and high frequency information of the scene can be achieved simultaneously, which is beneficial for the image restoration of both outline and detail of the scene.

Figure 4. The Schematic diagram of the novel lenslet pairing method: (a) the first pairing and (b) the second pairing.

3.2. Compressive Sensing Theory

CS is a novel signal sampling theory based on the sparsity and compression of an observed signal. CS is a robust framework for signal reconstruction. The theoretical framework of CS motivates sparse regularisation for solving inverse problems, through which, the quality of image reconstruction can be improved effectively [18].

3.2.1. Measurement Model

According to Van Cittert-Zernike theorem, the measurement equation takes the form of a Fourier transform. By inverse Fourier transform of the complex visibility measured by the system, intensity distribution of the scene will be achieved. Due to the imaging character of the PIII system, it is not possible to sample the entire (u, v) domain, leading to an ill-posed inverse problem of reconstructing image x from measured values y. The measurement equation for the PIII system can be represented as

(14)$$y=Ax+n,$$

where x ϵ R^T×1 is intensity of image. T is the total number of image pixels. y ϵ C^M×1 is the measured complex visibility. M represents the number of measurements. A ϵ C^M×T is the linear transformation matrix from the spatial domain to the frequency domain. n represents independently and identically distributed noise. Noise is not considered in this manuscript. A is the linear transformation matrix from the spatial domain to the frequency domain:

(15)$$A=SF,$$

where S ϵ R^M×T is the binary matrix for spatial frequency sampling of the object. F ϵ C^T×T is the discrete Fourier transform. Then, we have the following linear model:

(16)$$y=SFx+n,$$

where y can be calculated from the equation (12). S ϵ R^M×T can be obtained through the spatial frequency coverage of the scene. F ϵ C^T×T can be gained through the discrete Fourier transform basis.

3.2.2. Reconstruction Method

Fourier frequencies are not measured adequately through the PIII system, fitting the data alone does not uniquely define the image sought. Therefore, the image is reconstructed by using the total variation augmented Lagrangian alternating direction algorithm. The basic framework for image reconstruction is shown in Fig. 5.

Figure 5. The basic framework for image reconstruction.

Sparse regularization contributes to sparse solutions in terms of solving ill-posed inverse problems. To be specific, the ill-posed inverse problem can be solved by sparse regularization through imposing a priori constraints to select a unique image among all those that are consistent with the data. Moreover, this inverse problem of reconstructing image x from measured values y can be solved by building a total variation regularization mode, which can be expressed as

(17)$$\underset{x}{\min }{\displaystyle \sum _i‖D_{i}x‖subjecttoAx=y,}$$

where D_ix is gradient of x at pixel i. The slack variable w_i is introduced. The model of image reconstruction problem can be expressed as

(18)$$\underset{w_i,x}{\min }{\displaystyle \sum _i‖w_i‖subjecttoAx=y,}{D}_{i}x=w_{i}foralli.$$

The constrained model (18) can be converted to an unconstrained objective function by the augmented Lagrangian method, which can be expressed as

(19)$${\tilde{\mathcal{l}}}_A(w_i,x)={\displaystyle \sum _i\left(‖w_i‖-v_i^T\left(D_{i}x-w_i\right)+\frac{{\beta }_i}{2}{‖D_{i}x-w_i‖}_2^2\right)-k^T(Ax-y)+\frac{\mu }{2}}{‖Ax-y‖}_2^2,$$

where ν_i, κ are Lagrangian multipliers. β_i, μ are penalty parameters. As shown in Fig. 5, the image reconstruction problem can be decomposed into two simple sub-problems by the alternating direction algorithm, the sub-problem w_i and the sub-problem x. Through the iterative method, the sub-problem w_i is solved firstly, then the sub-problem x is done. The sub-problem is to minimize the augmented Lagrangian function, which can be represented as

(20)$$\underset{w_i,x}{\min }{\tilde{\mathcal{l}}}_A\left(w_i,x\right)={\displaystyle \sum _i\left(‖w_i‖-v_i^T\left(D_{i}x-w_i\right)+\frac{{\beta }_i}{2}{‖D_{i}x-w_i‖}_2^2\right)}-k^T(Ax-y)+\frac{\mu }{2}{‖Ax-y‖}_2^2.$$

Supposing x^k denotes the minimizer of the sub-problem x at the k-th iteration, w_i^k+1 can be obtained by minimizing the objective function with respect to w_i on each edge as follows:

(21)$$\underset{w_i}{\min }{\displaystyle \sum _i\left(‖w_i‖-v_i^T\left(D_i{x}^k-w_i\right)+\frac{{\beta }_i}{2}{‖D_i{x}^k-w_i‖}_2^2\right).}$$

Having w_i^k+1 at hand, x^k+1 can be achieved by minimizing the objective function with respect to x as follows:

(22)$$\underset{x}{\min }{\tilde{\mathcal{l}}}_A\left(w_i^{k+1},x\right)={\displaystyle \sum _i\left(‖w_i^{k+1}‖-v_i^T\left(D_{i}x-w_i^{k+1}\right)+\frac{{\beta }_i}{2}{‖D_{i}x-w_i^{k+1}‖}_2^2\right)}-k^T(Ax-y)+\frac{\mu }{2}{‖Ax-y‖}_2^2.$$

Lagrangian multipliers are updated through the well known formulas:

(23)$$v_i^{k+1}=v_i^k-{\beta }_i^k(D_i{x}^{k+1}-w_i^{k+1}),$$

(24)$${\kappa }^{k+1}={\kappa }^k-{\mu }^k(A{x}^{k+1}-y).$$

New penalty parameters μ and β_i are chosen. The most effective way to choose new penalty parameters is to try different values and compare the corresponding image quality under different parameter values as shown in Fig. 6. The selected penalty parameters have better image quality.

Figure 6. Image quality under different parameter values: (a) image quality under different μ parameter values and (b) image quality under different β_i parameter values.

Variables w_i and x are solved alternatively through updating Lagrangian multipliers and choosing penalty parameters to ensure that the constraint condition is met in the final solution [19].

The PIII system parameters used for the simulations are listed in Table 1. The sampling pattern in the uv-plane is shown in Fig. 7. Since the acquisition of spatial frequency information is an under-sampling process, this results in the loss of the partial image content. Therefore, the PSNR, the NMSE and SSIM are used as the evaluation parameters of image reconstruction. For a test image Y and its original image X, the PSNR, NMSE and SSIM metrics are defined as [20]

TABLE 1. System parameters used for the simulations.

Parameter	Symbol of Parameters	Parameter Values
Spectral Coverage	λ	380–760 nm
Lenslet Diameter	d	1 mm
Longest Baseline	Bmax	28 mm
Number of Long PIC Spokes	P1	37
Number of Lenslets per Long PIC Spokes	N1	29
Number of Short PIC Spokes	P2	74
Number of Lenslets per Short PIC Spokes	N2	15
Scene Distance	z	100 km
Number of Spectral Bins	m	10

Figure 7. The sampling pattern of PIII system in the uv-plane.

(25)$$PSNR=10{\log }_{10}\left[\frac{MA{X}^2}{MSE}\right],$$

(26)$$NMSE=\frac{{\displaystyle \sum _{i=1}^W{\displaystyle \sum _{j=1}^H{\left[X\left(i,j\right)-Y\left(i,j\right)\right]}^2}}}{{\displaystyle \sum _{i=1}^W{\displaystyle \sum _{j=1}^H{\left[X\left(i,j\right)\right]}^2}}},$$

(27)$$SSIM=\frac{(2{\mu }_X{\mu }_Y+C_1)(2{\sigma }_{XY}+C_2)}{({\mu }_X{}^2+{\mu }_Y{}^2+C_1)({\sigma }_X{}^2+{\sigma }_Y{}^2+C_2)},$$

where MAX is the maximum possible pixel value of the image. MSE is the mean square error of the image. μ_X and μ_Y are the averages of X and Y. σ²_X and σ²_Y are the variances of X and Y. σ_XY is the covariance of X and Y. C₁ and C₂ are two variables used to stabilize division involving a weak denominator. The larger PSNR and SSIM are, the better image quality is acquired. On the contrary, the smaller NMSE is, the better image quality is achieved.

Traditional lenslet array arrangement (as shown in Fig. 1) consists of radial spokes with equal length on a disc. The traditional lenslet pairing follows symmetrical scheme: (1, N), (2, N−1), (3, N−2)…(N/2, N/2+1). For comparison, the imaging simulation experiments were carried out for the PIII system in traditional and novel lenslet array arrangements and lenslet pairing method. All the images acquired through different approaches are shown in Figs. 8–10. The graphs of spatial frequency coverage are shown in Fig. 11. The calculated PSNR, NMSE and SSIM are listed in Table 2.

TABLE 2. The calculated PSNR, NMSE and SSIM.

The Image	PSNR(/dB)	NMSE	SSIM
USAF1951(b)	10.37	0.45	0.47
USAF1951(c)	10.81	0.41	0.53
USAF1951(d)	19.64	0.05	0.96
USAF1951(e)	24.09	0.01	0.98
The Saturn(b)	24.53	0.40	0.67
The Saturn(c)	24.67	0.39	0.68
The Saturn(d)	33.45	0.05	0.96
The Saturn(e)	34.42	0.04	0.97
The Nebula(b)	30.55	0.50	0.68
The Nebula(c)	30.74	0.48	0.69
The Nebula(d)	37.84	0.09	0.94
The Nebula(e)	38.06	0.08	0.95

Figure 8. USAF1951: (a) the original image, (b) direct restored image with traditional lenslet array arrangement and symmetrical lenslet pairing method, (c) direct restored image with novel lenslet array arrangement and lenslet pairing method presented in this paper, (d) compressive sensing (CS) reconstruction image with traditional lenslet array arrangement and symmetrical lenslet pairing method, and (e) CS reconstruction image with novel lenslet array arrangement and lenslet pairing method.

Figure 9. The Saturn: (a) the original image, (b) direct restored image with traditional lenslet array arrangement and symmetrical lenslet pairing method, (c) direct restored image with novel lenslet array arrangement and lenslet pairing method presented in this paper, (d) compressive sensing (CS) reconstruction image with traditional lenslet array arrangement and symmetrical lenslet pairing method, and (e) CS reconstruction image with novel lenslet array arrangement and lenslet pairing method.

Figure 10. The Nebula: (a) the original image, (b) direct restored image with traditional lenslet array arrangement and symmetrical lenslet pairing method, (c) direct restored image with novel lenslet array arrangement and lenslet pairing method presented in this paper, (d) compressive sensing (CS) reconstruction image with traditional lenslet array arrangement and symmetrical lenslet pairing method, and (e) CS reconstruction image with novel lenslet array arrangement and lenslet pairing method.

Figure 11. The graphs of spatial frequency coverage: (a) the graph of spatial frequency coverage constructed by a conventional system, (b) the graph of spatial coverage constructed by our system, and (c) the graph of increased spatial frequency coverage.

As listed in Table 2, the PSNR of direct restored image using the novel lenslet array arrangement and pairing method increases approximately 0.25 dB than that of the direct restored image using the traditional lenslet array arrangement and symmetrical pairing method. The NMSE of the direct restored image using the novel lenslet array arrangement and pairing method proposed by this paper decreases approximately 0.02 more than that of the direct restored image using the traditional lenslet array arrangement and pairing method. The SSIM of the direct restored image using the novel lenslet array arrangement and pairing method increases approximately 0.02 than that of direct restored image by using traditional lenslet array arrangement and symmetrical pairing method.

Image reconstruction based on CS theory is carried out and uses the novel lenslet array arrangement and pairing method. The PSNR of CS reconstructed image increases by approximately 10 dB from that of the direct restored image. Meanwhile, the NMSE calculated from the direct restored image is approximately 0.38 higher than that of the reconstructed image based on CS. The SSIM of the CS reconstructed image increases approximately 0.33 than that of the direct restored image. The PSNR of the CS reconstructed image using the novel lenslet array arrangement and pairing method increases approximately 1.88 dB from that of CS reconstructed image using the traditional lenslet array arrangement and symmetrical pairing method. Meanwhile, the NMSE calculated from the CS reconstructed image using traditional lenslet array arrangement and symmetrical pairing method increases approximately 0.02 than that of CS reconstructed image using the novel lenslet array arrangement and the pairing method. The SSIM of the CS reconstructed image using the novel lenslet array arrangement and pairing method increases approximately 0.01 than that of the CS reconstructed image using traditional lenslet array arrangement and symmetrical pairing method.

Figures 8–10 are the images acquired through different approaches. According to Figs. 8–10, the image restoration and reconstruction simulation results show that the novel lenslet array arrangement and pairing method presented in this paper could effectively improve the spatial frequency coverage and have influence on direct restored image quality. Meanwhile, image reconstruction based on CS theory could improve the image quality of interferometric imaging.

The novel lenslet array arrangement adds lenslet arrays in more directions, which are different from the directions of the existing radial lenslets spokes. The novel lenslet pairing approach adds more spatial frequency sampling points. As shown in Fig. 11, the traditional lenslet array arrangement and lenslet pairing method can collect interferometric information in 37 directions (radial spokes) and obtains 5180 measured points, while this novel design can collect interferometric information in 111 directions and obtains 20720 measured points. The obtained spatial frequency points based on the novel lenslet array arrangement and lenslet pairing approach are approximately four times than that of traditional design. Therefore, this novel design improves spatial frequency coverage of the system. The novel lenslet array arrangement and lenslet pairing method designed in this paper can improve the image quality because spatial frequency coverage has the direct influence on image quality.

This paper proposed the novel approach for improving spatial frequency coverage and image reconstruction. As shown in Fig. 1, the traditional lenslet array arrangement is composed of radial spokes with equal length on a disc. However, the novel lenslet array arrangement consists of different length spokes on a disc, as shown in Fig. 2. Compared with the traditional lenslet array arrangement, the obtained spatial frequency points based on the novel lenslet array arrangement are approximately twice than that of traditional design, which can effectively improve the spatial frequency coverage of the scene.

The traditional lenslet pairing follows symmetrical scheme: (1, N), (2, N−1), (3, N−2)…(N/2, N/2+1), the length of the baseline obtained through the traditional lenslet pairing approach is d, 3d, 5d, …, (N−1)d. Nevertheless, the length of the baseline obtained by the novel lenslet pairing approach is 0, d, 2d, 3d, …, (N−1)d, and the number of the baseline is approximately twice than that of traditional lenslet pairing. Therefore, the novel lenslet pairing approach can acquire the more spatial frequency information, which efficiently improve the image quality of the PIII system.

The traditional reconstruction algorithms mainly include CLEAN algorithm and maximum entropy method (MEM) algorithm. CLEAN algorithm can be seen as a greedy pursuit algorithm with an over-complete dictionary. Therefore, CLEAN algorithm seems unstable. The core idea of MEM algorithm is to obtain the least informative image that is consistent with the data [21, 22]. CS algorithm is applied to the PIII system in this paper, and it has significant impacts in other applications such as medical imaging, synthetic aperture radar imaging, remote sensing imaging, face recognition, speech recognition. Compared with traditional imaging methods, CS provides a framework for recovering image from the fewer measurements and has great advantages in flexibility [23, 24].

In this paper, a spatial frequency coverage method and image reconstruction algorithm for the PIII system were proposed. The novel lenslet array arrangement and lenslet pairing method positively influence the spatial frequency coverage. The acquired dense spatial frequencies efficiently improve the quality of the direct restored image. Based on CS theory, the image reconstruction model for the PIII system was established. By using this image reconstruction model, interferometric image with better quality could be acquired. The research content of this paper can provide technological support for optical interferometric imaging.

Although the image quality of the PIII system is improved based on the research content of this paper, there is still potential space for improving the quality of the image. In the future, the arrangement of lenslet array and the approach of lenslet pairing, which can increase the spatial frequency coverage, will still continue to be researched.