Single-pixel imaging (SPI) is a unique imaging method that can capture images employing a single-pixel detector (SPD)^[1]. It illuminates the object with modulated light fields and collects the reflected or transmitted light intensities with an SPD. Compared to a pixelated detector, an SPD is cheaper and easier to fabricate especially in the invisible wavelength. Accordingly, SPI is applicable in multi-wavelength imaging, such as infrared imaging^[2], terahertz imaging^[3], and X-ray imaging^[4,5]. Additionally, the SPD demonstrates enhanced functionality, such as high detection efficiency, low dark current, and quick response. That benefits the application in remote sensing of SPI^[6,7].

Currently, illumination light fields can be generated through rotating ground glass plates^[8,9], liquid-crystal spatial light modulators (SLMs)^[10,11], digital micro-mirror devices (DMDs)^[12,13], LED arrays, and silicon-based optical phased arrays (S-OPAs)^[14,15]. With regard to the rotating ground glass and SLMs, they can achieve a light field refresh rate of no more than hundreds of hertz. DMDs are capable of achieving a refresh rate not surpassing 50 kHz and would be employed in high-speed SPI, particularly in conjunction with mechanical scanning structures^[16,17]. In 2022, Jiang et al. employed a DMD to generate illumination fields, achieving high-speed imaging of 2,000,000 frame/s (fps) for periodically repeating scenes^[12]. LED arrays can generate Hadamard light fields with a refresh rate of up to 12.5 MHz^[18]. S-OPA is integrated with a silicon chip with a multi-mode fiber (MMF). Benefiting from the electro-optic phase modulators, S-OPA ensures rapid light field refresh at frequencies exceeding 100 MHz. However, factors like significant diffraction losses in DMDs, severe heat generation and large beam divergence in LED arrays, and low power tolerance of S-OPAs restrict the application of current high-speed SPI methods in remote sensing. Despite the light sources mentioned above, researchers also reported that laser arrays can generate illumination light fields suited to SPI^[19–21]. Actually, the fiber laser array is another form of optical phased array. Combined with high-performance phase modulators, a fiber laser array can theoretically achieve a phase modulation rate exceeding 100 MHz, which corresponds to an illumination light field refresh rate of 100 MHz or higher. Additionally, the fiber laser array offers a notable advantage of high output power. Currently, the single-frequency all-fiber laser is capable of achieving power output of hundreds of watts^[22]. Moreover, the fiber laser can achieve coherent beam combining and realize power output in excess of 20 kW^[23]. This advantage would promote the application of SPI in long-distance imaging.

However, another issue must be addressed before fiber laser arrays can be employed in SPI. The configuration of sub-apertures within the laser array would have a significant impact on the normalized intensity second-order correlation function $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields. For example, the illumination light fields generated by a periodic arranged array, such as a rectangular and hexagonal array, would have many peaks in $g^{(2)}(u,v;u_0,v_0)$. This kind of illumination light field would degrade the reconstructed image due to periodic artifacts or overlapping in the reconstructed images. In 2018, Wu et al^[20]. reported a method via a sparse structured illumination source to depress the spatial periodicity of the illumination light fields. In 2018, Liu et al.^[21] proposed using a low-pixel array detector to avoid the spatial periodicity issue of illumination fields, thereby reducing the periodic artifacts in reconstructed images. In 2023, Lai et al.^[19] proposed a demonstration experiment of SPI based on the Fermat spiral laser array arrangement, which prevents the periodic artifacts in reconstructed images. Previous studies have only included demonstrative experiments and conceptual designs, without constructing an actual high-speed SPI system illuminated by a fiber laser array.

This study developed a high-quality and high-speed SPI system based on a Fermat spiral fiber laser array. Our work has thoroughly validated the feasibility of using fiber laser arrays for SPI and has pioneered the achievement of fast imaging. In comparison to rectangular and hexagonal arrays, the non-periodic configuration of the Fermat spiral ensures that the $g^{(2)}(u,v;u_0,v_0)$ of illumination light fields retain unimodal properties. This property reduces the periodic artifacts in reconstructed images, thereby enhancing the imaging quality. We have designed a high-speed data synchronous acquisition system that coordinates the synchronization of illumination light field refresh and light intensity detection at 20 kHz. Previously, we reported partial imaging results based on traditional reconstruction algorithms^[24]. In this study, combined with an untrained neural network (UNN)^[25], a $64\mathrm{pixel}\times 64$ pixel target was discernibly imaged at a sampling ratio of 4.88% with an imaging rate of 100 fps. Furthermore, due to the high-power output potential of fiber laser array, this fast SPI system can be applied in long-distance imaging.

SPI uses a series of modulated light fields to illuminate the object and collect the light intensity values reflected or transmitted from the object by an SPD. The light intensity is expressed as $$B=\iint O(u,v)P(u,v)\mathrm{d}u\mathrm{d}v,$$(1)where $P(u,v)$ represents the illumination light field intensity distribution, $O(u,v)$ denotes the object reflectivity or transmittance ratio distribution, and $B$ is the intensity value.

The image can be reconstructed using the illumination light fields and their corresponding light intensities. The classic reconstruction process is $$O_{\text{recons.}}(u,v)=\frac{⟨P(u,v)B⟩}{⟨P(u,v)⟩⟨B⟩},$$(2)where $O_{\text{recons.}}(u,v)$ is the reconstructed image. $⟨⟩$ represents the ensemble average operation.

According to the theory of statistical optics^[26], random illumination light fields can be generated through the interference of several sub-apertures within a fiber laser array. Figure 1 illustrates the generation process of the illumination light fields. The seed laser is amplified and split into several channels. Each channel is modulated separately by a high-speed electro-optic phase modulator and then injected into the collimator array. The collimator array emits laser beams into the free space. A lens is set at the emission plane in order to simulate the far-field condition, and its back focal plane is the far-field plane. The beams interfere in the back focal plane, and then the interference speckles are the random illumination light fields.

At the emission plane, each sub-aperture of the array emits a laser beam. Each beam can be regarded as a Gaussian beam truncated by a circular stop. The waist of each beam is $w_0$, and the diameter of the circular stop is $D$. The light field at the emission plane is expressed as $$\begin{gathered} E_1(x,y,z=0)=\sum _{n=1}^N{U}_n\exp .[-\frac{{(x-x_n)}^2+{(y-y_n)}^2}{w_0^2}] \\ \times \text{circ}\left[\frac{\sqrt{{(x-x_n)}^2+{(y-y_n)}^2}}{D}\right]\exp .(i{\phi }_n), \end{gathered}$$(3)where $N$ denotes the total number of the sub-apertures and subscript $n$ represents the $n$th laser beam. For the $n$th laser beam, $U_n$ represents its amplitude, $(x_n,y_n)$ is its center coordinate, and ${\phi }_n$ corresponds to the random phase introduced by the fiber system and the modulator. circ is the aperture function of the circular stop, as shown in the following equation: $$\text{circ}(x)=\{\begin{array}{cc}1& x\le 1\\ 0& x>1\end{array}.$$(4)

After passing through the lens, each beam is attached to a phase, $${\phi }_{\mathrm{lens}}=-\frac{ik(x^2+y^2)}{2f},$$(5)where $k$ means the wavenumber and $f$ is the focal length of the lens.

When the light beams propagate to the far-field plane, they interfere with each other and generate a random illumination light field. The light field can be calculated according to Fresnel–Kirchhoff diffraction theory as $$E_2(u,v,z=f)=\frac{\exp .(ikf+ik\frac{u^2+v^2}{2f})}{i\lambda f}\iint E_1(x,y,z=0)\exp .(-ik\frac{ux+vy}{f})\mathrm{d}x\mathrm{d}y,$$(6)where $\lambda$ represents the wavelength of the light and $(u,v)$ denotes the coordinate of the point in the far-field plane. The intensity distribution in the far-field plane then can be expressed as $$I_2(u,v,z=f)=E_2(u,v,z=f)·E_2^{*}(u,v,z=f).$$(7)

$I_2$ ($u,v$, $z=f$) is the illumination light field mentioned in the above equation. Once the phase modulator array renews the random phase ${\phi }_n$, there would be a refresh of the illumination light field. Thus, the refresh rate of illumination light fields corresponds to the phase modulation speed of the electro-optic phase modulators.

The principle of SPI can be further articulated as $$O_{\mathrm{recons}.}(u,v)=\iint O(u_0,v_0)\frac{⟨I(u,v)I(u_0,v_0)⟩}{⟨I(u,v)⟩⟨I(u_0,v_0)⟩}\mathrm{d}{u}_0\mathrm{d}{v}_0.$$(8)

Among the above equation, the normalized intensity second-order correlation function $g^{(2)}(u,v;u_0,v_0)$ is an important performance indicator for random illumination light fields in SPI, which is defined as $$g^{(2)}(u,v;u_0,v_0)=\frac{⟨I(u,v)I(u_0,v_0)⟩}{⟨I(u,v)⟩⟨I(u_0,v_0)⟩},$$(9)where the $I(u_i,v_i)$ denotes the intensity of $(u_i,v_i)$. $⟨⟩$ represents ensemble average operation. $g^{(2)}(u,v;u_0,v_0)$ means the correlation between $(u,v)$ and $(u_0,v_0)$ in the light field.

According to Eq. (8), an object point $(u_0,v_0)$ maps to a reconstructed image point $(u,v)$ with weight $g^{(2)}(u,v;u_0,v_0)$. If $(u_0,v_0)$ correlates strongly with a single light field point $(u,v)$ and not others, $g^{(2)}$ will have a single peak, ensuring accurate reconstruction. However, if $(u_0,v_0)$ correlates with multiple light field points $(u_1,v_1)$, $(u_2,v_2)$, and $(u_3,v_3)$, $g^{(2)}$ will have multiple peaks, leading to artifacts or overlaps in the reconstructed image. If $(u_0,v_0)$ has a strong correlation with one point $(u,v)$ and a weaker correlation with others $(u_4,v_4)$, $g^{(2)}$ will have a main peak with minor sidelobes, introducing noise in the reconstructed image at those points. To improve image quality, $g^{(2)}$ should ideally have a single dominant peak with minimal sidelobes.

$g^{(2)}(u,v;u_0,v_0)$ can be deduced referring to the work by Liu et al.^[27] as $$g^{(2)}(u,v;u_0,v_0)\propto \frac{I_{\mathrm{KF}}(u-v_0,u-v_0)-\sum _{n=1}^N{U}_n^4+{\left(\sum _{n=1}^N{U}_n^2\right)}^2}{{\left(\sum _{n=1}^N{U}_n^2\right)}^2},$$(10)$$I_{\mathrm{KF}}(u,v)=\sum _{n=1}^N\sum _{m=1}^N{U}_n^2{U}_m^2\exp .\{-i\frac{k}{f}[u(x_n-x_m)+v(y_n-y_m)]\},$$(11)where $(u,v)$ and $(u_0,v_0)$ represent the coordinates in the observation plane. $(x_{n,}{y}_n)$ and $(x_{m,}{y}_m)$ denote the center coordinates of the $n$th and $m$th laser beam in the emission plane. It can be induced that the distribution of $g^{(2)}(u,v;u_0,v_0)$ is related to the geometric configuration of the fiber laser array.

Three different arrays are selected to research the effect of the geometric configuration: a rectangle array with 36 sub-apertures, a hexagon array with 37 sub-apertures, and a Fermat spiral array with 36 sub-apertures, as shown in Figs. 2(a1)–2(a3). The Fermat spiral fiber laser array in simulations is arranged according to Ref. [19], $$\{\begin{array}{cc}{r}_n=s\sqrt{\frac{n}{\pi }},& n=1,2,3,\dots ,N\\ {\theta }_n=2\pi \beta n,& n=1,2,3,\dots ,N\end{array},$$(12)where $(r_n,{\theta }_n)$ denotes the center’s polar coordinate of the $n$th laser beam. $s$ and $\beta$ represent the radial and angular parameters deciding the compactness of the array. The angular interval $2\pi \beta$ is set to be 137.5°, and $s$ is chosen as $0.7D$. The diameter of aperture $D$ is selected as 6 mm. Additionally, the laser is simulated with a wavelength of 1064 nm and with a waist of 2.7 mm. The simulations in Sec. 2.3 also follow these parameters. It is notable that $g^{(2)}(u,v;u_0,v_0)$ is a statistical property of the illumination light fields. Infinite illumination light fields are needed to calculate the theoretical $g^{(2)}(u,v;u_0,v_0)$. However, the calculation results of $g^{(2)}(u,v;u_0,v_0)$ would approach the theoretical case when the total number of illumination light fields exceeds 1000^[27].

A total of 3000 frames illumination light fields were simulated based on the three aforementioned arrays and employed to calculate the $g^{(2)}(u,v;u_0,v_0)$ according to Eq. (9). The reference point is selected as ($u_0=0$, $v_0=0$). The 3D distributions of $g^{(2)}(u,v;u_0,v_0)$ are shown in Figs. 2(b1)–2(b3). Figures 2(c1)–2(c3) and Figs. 2(d1)–2(d3) are the distributions of $g^{(2)}(u,v;u_0,v_0)$ along the $u$-axis and $v$-axis, respectively. Figures 2(e1)–2(e3) demonstrate the vertical view for $g^{(2)}(u,v;u_0,v_0)$ of different array configurations. When the array is configured in rectangular form, the $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields has many peaks aligned in rectangular manner, as illustrated in Figs. 2(b1)–2(e1). Similarly, the $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields generated by hexagon arrays also have many hexagonal peaks, as shown in Figs. 2(b2)–2(e2). As for the Fermat spiral array, the $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields is unimodal with a prominent peak and no additional sidelobes as demonstrated in Figs. 2(b3)–2(e3). The illumination light fields generated by the Fermat spiral array satisfy the requirement for the unimodal characteristic of $g^{(2)}(u,v;u_0,v_0)$.

Furthermore, the SPI imaging simulation is performed to reveal the impact of different arrays on the imaging quality. A binary object “3-slits” and a grayscale object “drone” are selected as the targets whose resolutions are $128\mathrm{pixel}\times 128$ pixel. Differential ghost imaging (DGI)^[28] and a kind of compressive sensing algorithm based on total variation regularization prior (CS-TV) are used for image reconstruction. DGI is a linear algorithm and has less computation overhead. The CS-TV algorithm is a nonlinear iterative method capable of reconstructing target images with less sampling. A total of 5000 frames of illumination light fields are employed to acquire images with a resolution of $128\mathrm{pixel}\times 128$ pixel, corresponding to a sampling ratio of 30.52%. The sampling ratio is defined as the ratio between the number of illumination light fields and the pixel count of the image to be reconstructed. Figure 3 demonstrates that the images acquired by a rectangle array and hexagon array have periodic artifacts and overlapping deterioration. Artifacts are distributed in a rectangular and hexagonal manner for the rectangle array and hexagon array, respectively. When the array is configured along a Fermat spiral, the reconstructed images have better imaging quality with no artifacts. However, the reconstructed “3-slits” image exhibits a slight distortion, and the “drone” image shows a foggy blur when utilizing DGI as the reconstruction algorithm for the Fermat spiral array. In contrast, the nonlinear algorithm CS-TV demonstrates superior resolution and higher contrast. Based on the analysis of illumination light fields and SPI imaging, the fiber laser array is arranged along the Fermat spiral to achieve high imaging quality for subsequent work.

For Fermat spiral fiber laser arrays with different numbers of apertures $N$, we studied the properties of their illumination light fields and their SPI imaging results. A total of 3000 frames of illumination light fields are simulated according to Sec. 2.1. The $g^{(2)}(u,v;u_0,v_0)$ is calculated according to Eq. (9), and the reference point is selected as ($u_0=0$, $v_0=0$). When $N$ takes on the values of 8, 16, 24, and 32, the corresponding results are demonstrated in Fig. 4. Array configurations with $N$ apertures are shown in Figs. 4(a1)–4(a4). The 3D distribution and vertical view of $g^{(2)}(u,v;u_0,v_0)$ are illustrated in Figs. 4(b1)–4(b4) and Figs. 4(c1)–4(c4), respectively. When $N$ equals 8, the $g^{(2)}(u,v;u_0,v_0)$ has sidelobes comparable to the main lobe, as shown in Figs. 4(b1)–4(c1). As $N$ increases to 16, sidelobes still exist, while the amplitudes of sidelobes have decreased, as depicted in Figs. 4(b2)–4(c2). With the rise of $N$, the sidelobes of $g^{(2)}(u,v;u_0,v_0)$ diminish in number and amplitude, highlighting the main lobe. At $N=32$, the $g^{(2)}(u,v;u_0,v_0)$ exhibits unimodality with a main lobe and minimal sidelobes, as demonstrated in Figs. 4(b4)–4(c4).

Furthermore, the SPI imaging results for arrays with different numbers of apertures $N$ are simulated and shown in Fig. 5. A total of 5000 simulated illumination light fields were utilized for imaging a $128\mathrm{pixel}\times 128$ pixel “3-slits,” with a corresponding sampling ratio of 30.52%. DGI and CS-TV are selected as reconstruction algorithms in this part. To quantitatively evaluate the quality of the reconstructed images, the root mean square error (RMSE) is calculated as illustrated in Fig. 6. A smaller RMSE indicates that the image is more similar to the original object. As depicted in Fig. 5, reconstructed images have artifacts when $N$ equals 8. As $N$ increases, the artifacts in the images are progressively diminished, and the RMSE of the images exhibits a declining trend. When $N$ equals 32, the reconstructed images are free from artifacts and possess excellent imaging quality. As for different algorithms, DGI has fewer artifacts and lower RMSE when $N$ equals 8. With the increase of $N$, CS-TV results in a better imaging quality than DGI. Following the analysis of illumination light fields and SPI imaging, our subsequent research selected 32 as the number of apertures $N$.

The schematic of the experimental system is illustrated in Fig. 7. The seed laser is a single frequency linearly polarized fiber laser with a wavelength of 1064 nm. Following amplification via a fiber amplifier, the seed laser is then directed to a fiber splitter, where it is evenly distributed into 32 individual channels. Each channel is individually injected into a corresponding phase modulator within the modulator array and then undergoes specific phase modulation. The phase modulation signals are applied to the phase modulators via the data acquisition card 1 (DAQ 1), featuring a maximum sampling rate of 2 MS/s. Notably, equipped with ${\text{LiNbO}}_3$ electro-optic phase modulators, our modulation array can support a high phase modulation refresh rate up to 150 MHz. The modulated beams are emitted from a collimator array whose sub-apertures are configured according to the Fermat spiral mentioned in Sec. 2.2. $D$ is the aperture diameter that equals 1.1 mm, and $w_0$ is 2.92 mm. A lens with a focal length of 1.1 m is employed to achieve a far-field condition. The illumination light fields would be generated at the back focal plane of the lens. The focusing light is divided into two beams by a beam splitter. A high-speed camera, capable of recording at a maximum frame rate of 20 kfps, is set on the back focal plane to record the intensity distribution of illumination light fields. The camera resolution is $1024\times 512$, from which we extracted a $64\mathrm{pixel}\times 64$ pixel section containing the illumination light field. Meanwhile, the other beam reflected by the beam splitter illuminates the transmissive object located at the equivalent back focal plane. An SPD with a bandwidth of 260 kHz is utilized to collect the light intensity values transmitted from the object. The light intensity values are converted to digital signals by DAQ 2 with a maximum sampling rate of 1 MS/s. A PC is utilized to input data into DAQ 1 and to process data received from DAQ 2 and the camera. The communication between camera and PC is conducted via a gigabit ethernet cable. Both the DAQs and camera are synchronized by a synchronization signal generated from the DAQ 1.

According to the principle of SPI, the light intensity values are expected to match their corresponding illumination light fields. Therefore, a high-speed data synchronous acquisition system was designed as demonstrated in Fig. 8. DAQ 1 outputs the synchronization signals for the system at a frequency of 20 kHz, as shown in Fig. 8(a). Once a synchronization signal is received, the phase modulator, SPD, and high-speed camera execute their tasks simultaneously. As illustrated in Fig. 8(b), the phase modulator attaches random phase modulation to the laser beam at 20 kHz in order to refresh the illumination light fields. A high-speed camera is used to record the intensity distribution of the illumination light fields whose refresh rate is 20 kHz. A recorded illumination light field is shown in Fig. 8(c). Meanwhile, the SPD performs light intensity detection. The time-varying normalized light intensity values are demonstrated in Fig. 8(d).

In this part, 2000 frames of illumination light fields are acquired to calculate the $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields according to Eq. (9), and the reference point is also set at ($u_0=0,v_0=0$). The results are shown in Fig. 9. As mentioned in Sec. 2.2, 2000 frames are enough to capture the statistical property of the illumination light fields. An acquired illumination light field is depicted in Fig. 9(a). The 3D distribution of the $g^{(2)}(u,v;u_0,v_0)$ is illustrated in Fig. 9(b). Figures 9(c) and 9(d) show the 3D distributions of $g^{(2)}(u,v;u_0,v_0)$ from the positive direction of the $v$- and $u$-axes, respectively. We can see that $g^{(2)}(u,v;u_0,v_0)$ attains its peak value of 1.62 at the point (0, 0), with the maximum sidelobe taking the value of 1.23. Points of $g^{(2)}(u,v;u_0,v_0)$ elsewhere are approximately 1. In a conclusion, $g^{(2)}(u,v;u_0,v_0)$ of the light fields generated by the Fermat fiber laser array exhibits a unimodal property, which corresponds to the simulation results.

Furthermore, SPI imaging is performed. A transmissive “2” of $64\mathrm{pixel}\times 64$ pixel is selected as the object. The reconstruction results are illustrated in Fig. 10. Different samples, which take the values of 100, 200, 400, 600, and 2000, are employed to reconstruct the images, corresponding to the frame rates of 200, 100, 50, and 33 fps, respectively. This paper is dedicated to enhancing the refresh rate of light fields, reducing sampling time and thereby increasing imaging speed. It should be noted that the frame rate calculated here is an off-line imaging frame rate, considering only the time consumption of the sampling process. This reflects the highest temporal resolution capability supported by the current system. Despite DGI and CS-TV, UNN is introduced as a reconstruction algorithm in this part. In the Supplement 1, we provide a brief introduction to three algorithms and their runtime. In order to evaluate the imaging quality quantitatively, RMSEs are also calculated and shown below the images.

As the number of samples increases, the images from these three algorithms are increasingly clear, and the RMSE values of the images decrease. At the samples over 200, compared to the reconstructed images of DGI and CS-TV, the images from the UNN have the highest contrast, the weakest background noise, and minimal RMSE. Thus, the UNN exhibits the best imaging quality among these three algorithms. When the samples equal 100, the images reconstructed by these algorithms are indistinguishable. When the samples are 200, corresponding to a sampling ratio of 4.88%, the images reconstructed by DGI and CS-TV have heavier background noise, which reduces the recognizability. The images reconstructed by the UNN are still discernable with an RMSE value less than 0.17; hence, it is our assessment that the imaging system achieves identifiable imaging when samples are 200, which equates to a frame rate of 100 fps.

To further investigate the impact of sampling reduction on image quality, we additionally calculated the $g^{(2)}(u,v;u_0,v_0)$ of the illumination light fields at sampling volumes of 100, 200, 400, and 600. The results are shown in Fig. 11. According to the analysis in Sec. 2.2, the main lobe of $g^{(2)}(u,v;u_0,v_0)$ would correspond to the correct point in the image, and the sidelobes would contribute to the distribution of noise. Thus, the $g^{(2)}(u,v;u_0,v_0)$ of illumination light fields should ideally exhibit a dominant peak with maximum amplitude and suppressed sidelobes. With 100 samples, the main peak of the $g^{(2)}(u,v;u_0,v_0)$ for the illumination light field was 1.47, which was not prominently pronounced. Moreover, numerous sidelobes persisted, with the highest sidelobe reaching a value of 1.28, as illustrated in Figs. 11(a1) and 11(b1). There were 16 points with values exceeding 1.2, and 51 points exceeding 1.15 in the $g^{(2)}(u,v;u_0,v_0)$, as shown in Table 1. At this stage, the advantage of the Fermat spiral configuration for sidelobe suppression was not fully realized, leading to increased noise in the reconstructed image. As demonstrated in Figs. 11(a2)–11(b2), the main lobe of $g^{(2)}(u,v;u_0,v_0)$ attains a value of 1.62 at 200 samples, which is closer to the peak value at 2000 samples. The $g^{(2)}(u,v;u_0,v_0)$ distribution shows a reduction in sidelobes, with only 12 points valued above 1.2 and 40 points exceeding 1.15, as shown in Table 1. The $g^{(2)}(u,v;u_0,v_0)$ of illumination light fields exhibits a unimodal property preliminarily. With the further increase in the number of samples, there is no significant change in the height of the main lobe and the maximum sidelobe. However, the number of sidelobes with peaks exceeding 1.15 further decreases, fully leveraging the sidelobe suppression characteristic of the Fermat spiral configuration. Consequently, the $g^{(2)}(u,v;u_0,v_0)$ of illumination light fields exhibits a distinct unimodal feature. At 600 samples, $g^{(2)}(u,v;u_0,v_0)$ has 12 points over 1.2 and 21 over 1.15, slightly more than the 11 and 17 points at 2000 samples. The unimodal characteristic of the $g^{(2)}(u,v;u_0,v_0)$ is nearly equivalent to that achieved with 2000 frames. In general, the unimodal characteristic improves progressively with an increase in the number of samples. When the number of samples reaches 200, an acceptable level of unimodality is achieved.

We used the structural similarity index (SSIM) to compare the similarity of $g^{(2)}(u,v;u_0,v_0)$ at different samples to the reference $g^{(2)}(u,v;u_0,v_0)$ derived from 2000 samples. For different samples, a higher SSIM indicates that the distribution of $g^{(2)}(u,v;u_0,v_0)$ is closer to that of the reference $g^{(2)}(u,v;u_0,v_0)$. This provides a quantitative reference for assessing the unimodal characteristics of the illumination light field. The results are plotted in Fig. 12. At 100 samples, the SSIM value was only 0.4636, suggesting a poor unimodal characteristic of the $g^{(2)}(u,v;u_0,v_0)$. When the samples doubled to 200, the SSIM increased to 0.6795, indicating a better unimodal property. As the samples further increased, the SSIM continued to rise, highlighting the enhanced unimodal characteristics of the $g^{(2)}(u,v;u_0,v_0)$. At 600 samples, the SSIM reached 0.9214. This implies that when the number of samples is 600, the distribution of $g^{(2)}(u,v;u_0,v_0)$ is already very similar to that of 2000 samples, exhibiting favorable unimodal characteristics. This analysis further supports that $g^{(2)}(u,v;u_0,v_0)$ is a statistical characteristic of the illumination light fields. Limited illumination light fields fail to leverage array configuration benefit, and they inadequately control sidelobe issues. Additionally, with fewer illumination light fields, the object information obtained is limited, which further affects the quality of reconstruction. In conclusion, samples below 200 are insufficient for imaging due to the unsuppressed sidelobes of $g^{(2)}(u,v;u_0,v_0)$ and the total information acquired.

In summary, a high-quality and high-speed SPI system based on the Fermat spiral fiber laser array was developed. Via theoretical analysis, the use of the Fermat spiral array configuration reduces the periodic artifacts in reconstructed images in comparison to rectangular and hexagonal arrays, thereby enhancing the imaging quality. In the experiment, we designed and built an SPI system illuminated by the Fermat spiral fiber laser array. We developed a high-speed data synchronous acquisition system that executed illumination light field refresh and light intensity detection at 20 kHz simultaneously. In conjunction with a UNN algorithm, the system can achieve a frame rate of up to 100 fps for $64\mathrm{pixel}\times 64$ pixel images, with a sampling ratio as low as 4.88%. Under these conditions, the RMSE value of the reconstructed images is less than 0.17. Besides, we analyzed the minimum samples for imaging in the current system from the perspectives of the $g^{(2)}(u,v;u_0,v_0)$ and the total amount of information acquired. Additionally, due to the high-power output potential of fiber laser array, this fast SPI system can be applied in long-distance imaging. The next efforts can be made to acquire the illumination light fields by computation instead of recording them by the camera, which would further increase the imaging speed of SPI.