**3. Results**

## *3.1. Numerical Simulation*

Numerical simulations were used to validate the proposed method. The simulation is depicted in Figure 2. The radius of the annular array was set to *Ra* = 40 mm, and there were *N* = 576 elements in the array. The acoustic scattering layer consisted of 120 acoustic scatterers with a diameter of *d* = 0.8 mm. The inner and outside radius of the scattering layer were *R*1 = 10 mm and *R*2 = 15 mm, respectively. The sound velocity and the density of these scatterers were 5200 m/s and 7,870 kg/m3, respectively. The sound velocity of the surrounding medium was 1500 m/s, the density of which was 1,000 kg/m3, which is close to that of soft tissues or water. Therefore, the scatterers had a severe impedance mismatch with the background media, which resulted in strong scattering. The imaging targets were five optical absorbers with a diameter of *d* = 0.8 mm in the ROI within the scattering layer. The coordinates of these targets were (−2,−2), (−2,2), (2,−2), (2,2), and (0,0). Under the illumination of the laser pulse, these optical absorbers in the scattering layer emitted ultrasonic pulses with a central frequency of 2.0 MHz and a −3 dB bandwidth of 1.26–2.68 MHz to the surrounding medium. The PA signals were subsequently detected by the annular array enclosing the region of interest. The PA signals awee sampled with a frequency of 40 Mhz. The imaging targets had the same sound velocity and the same density as the scatterers.

**Figure 2.** The schematic of the scenario of the simulation. Five identical optical absorbers with a diameter of *d* = 0.8 mm are placed in the ROI within a scattering layer. The inner and outside radius of the scattering layer are *R*1 = 10 mm and *R*2 = 15 mm, respectively. The ultrasound signal is detected by an annular array with a radius of *R*a = 40 mm.

Figure 3 depicts the results of the simulation with the scattering layer. Figure 3a illustrates the typical PA signals detected by one transducer. The measured broadband signal was partly overwhelmed by the noise and the amplitude was unstable. The PA signals detected by all transducers are shown in Figure 3b to visualize the effect of the scattering layer on the received signals, where the gray dotted line corresponds to the results in Figure 3a. The coherence of the signals was broken by the randomness of the scattering. The presence of a scattering layer inevitably causes a drop in the image quality.For the sake of comparison, Figure 3c illustrates the imaging results reconstructed by the classic delay-and-sum (DAS) beamforming method. As shown, for the DAS beamforming method, the image reconstructed in the case with a scattering layer (Figure 3c) suffered from low resolution, and the background showed strong noise. As a result, the contrast between the targets and the background decreased distinctly. Figure 3d depicts the image recovered by the time reversal method combined with a correlation matrix filter. The image shows good contrast, as the background noise has been significantly reduced after filtering.

**Figure 3.** Simulative results. (**a**) The detected photoacoustic (PA) signals by a single transducer. (**b**) The detected PA signals by all transducers. (**c**) The image reconstructed by the delay and sum (DAS) beamforming method. (**d**) The image reconstructed by the proposed matrix filtering method.

We used the signal-to-noise ratio (SNR) and the full width at half maximum (FWHM) of the targets' image to quantify the image quality. The SNR is defined as the ratio of the mean signal intensity in the signal area and the mean noise intensity in the background area. We defined the area inside the circles with a diameter of 0.8 mm (equal to the actual size of the targets) at the locations of the five targets as the signal area and the area outside as the background area. The averaged FWHM along the *X* and *Z* axes of all targets was calculated to quantify the degree of the distortion. For the results in Figure 3c, the SNR of the image was 17.2 dB, and the averaged FWHM along the *X* and *Z* axes were 1.42 mm and 2.6 mm, respectively. The difference between the FWHM along the X and Z axes were caused by the non-uniform distribution of the scatterers in the scattering layer. The SNR of the image in Figure 3d increased to 21.0 dB due to the matrix filter. In addition, the FWHM along the *X* and *Z* axes were 0.87 mm and 1.63 mm, respectively, which are both closer to the actual size of the targets. The results demonstrate that the matrix filtering method can more accurate depict the shape and position of the targets in the ROI.
