## 2 Materials and Methods
SSIM is used to measure the similarities between two given images in three dimensions: brightness, contrast and structure. The specific form of SSIM index is:  
\[\text{SSIM}(x,y)=\frac{\left(2\mu_{x}\mu_{y}+C_{1}\right)\left(2\sigma_{xy}+C_ {2}\right)}{\left(\mu_{x}^{2}+\mu_{y}^{2}+C_{1}\right)\left(\sigma_{x}^{2}+ \sigma_{y}^{2}+C_{2}\right)} \tag{15}\]  
\(\mu_{x}\), \(\mu_{y}\), \(\sigma_{x}^{2}\), \(\sigma_{y}^{2}\) and \(\sigma_{x}\sigma_{y}\) are the average value, variance and covariance of two non-negative image signals \(x\) and \(y\), respectively. \(C_{1}\) and \(C_{2}\) are constants to maintain stability with \(C_{1}=\left(0.01L\right)^{2}\) and \(C_{2}=\left(0.03L\right)^{2}\), where \(L\) is the dynamic range of the pixel values. Simultaneously, the MSSIM index is used to evaluate the overall image quality from local image features:  
\[\text{MSSIM}(X,Y)=\frac{1}{MN}\sum_{i=1}^{M}\sum_{j=1}^{N}\text{SSIM}(x_{i},y_ {j}) \tag{16}\]  
where \(X\) and \(Y\) are the reference and the distorted images; \(x_{j}\) and \(y_{j}\) are the image contents at the \(i\)th row and \(j\)th column local window; \(M\), \(N\) is the number of local windows of the image. Here, we use a symmetric Gaussian filter matrix of \(11*11\) with a standard deviation of 1.5 as a weighting window to filter the image, which is equivalent to sliding the window on the distorted image, so the MSSIM can be obtained by averaging the matrix mapped by the local SSIM index. PSNR is used to measure the quality of the degraded image:  
\[\text{PSNR}=10\log_{10}\left(\frac{L^{2}}{\text{MSE}}\right) \tag{17}\]  
\[\text{MSE}=\frac{1}{N}\sum_{i=1}^{N}(x_{i}-y_{i})^{2} \tag{18}\]  
where MSE is mean squared error, \(N\) is the total number of pixels in the image. Generally, the larger the SSIM, MSSIM and PSNR, the higher the quality of the degraded image and the more similar the two images are. Scattering: photon packets collide after the free path \(l_{R}\) with the scattering angle \(\theta\) and the azimuth angle \(\gamma\). The transmission direction after collision converted into the global coordinate system is expressed as \(\left(u_{x}^{\prime},u_{y}^{\prime},u_{z}^{\prime}\right)\), denoted by the original direction \(\left(u_{x},u_{y},u_{z}\right)\):  
\[ \begin{bmatrix}u_{x}^{\prime}\\ u_{y}^{\prime}\\ u_{z}^{\prime}\end{bmatrix} = \begin{bmatrix}u_{x}u_{z}/\sqrt{1-u_{z}^{2}}&-u_{y }/\sqrt{1-u_{z}^{2}}&u_{x}\\ u_{y}u_{z}/\sqrt{1-u_{z}^{2}}&u_{x}/\sqrt{1-u_{z}^{2}}&u_{y}\\ -\sqrt{1-u_{z}^{2}}&0&u_{z}\end{bmatrix}\begin{bmatrix}\sin\theta\cos\gamma\\ \sin\theta\sin\gamma\\ \cos\theta\end{bmatrix}  \tag{19}\]  
Scattering is also accompanied by the loss of energy: \(E_{s}=\omega E_{pre}\) where the albedo \(\omega=b/c\), \(E_{s}\) is the energy after collision while \(E_{pre}\) is the energy before it. The bidirectional transmission of photons is a cyclic process with the termination conditions:  
1. The photon packet returns to the detector within the field of view;
2. The photon packet reaches the target plane but does not hit the target area;
3. The photon packet energy is less than the threshold, namely, it is lost during the transmission process;
4. The number of photon packet scattering exceeds the preset maximum number. The initialized photon packet participates in the calculation of the modules in turn, once the termination conditions are met, the tracking is interrupted and the next photon packet is emitted.