## 2 The Methods of PAI and LUI
### Forward and Inverse Problems of PAT and LUT--Reconstruction Algorithms
#### 2.2.1 Algorithms for PAT
For illustration, we follow [75], and take a discrete representation of the integral in (2) in spherical coordinates, expressing the surface element as \(dS^{\prime}(\mathbf{r}^{\prime},t)=\left|\mathbf{r}-\mathbf{r}^{\prime}\right|^{2}\cos \phi d\phi d\theta\)  
\[\int_{S^{\prime}(\mathbf{r},t)}\frac{p_{0}(\mathbf{r}^{\prime})}{\left|\mathbf{r}-\mathbf{r} ^{\prime}\right|}dS^{\prime}(\mathbf{r}^{\prime},t)=\int_{\phi}\int_{\theta}p_{0 }(\mathbf{r}^{\prime})|\mathbf{r}-\mathbf{r}^{\prime}|\cos\phi d\phi d\theta\approx\Delta \phi\Delta\theta\sum_{l}p_{0}\big{(}\mathbf{r}^{\prime}_{l}\big{)}|\mathbf{r}-\mathbf{r}^{ \prime}|\cos\phi_{l}, \tag{9}\]  
Here, \(p_{0}\big{(}\mathbf{r}^{\prime}_{l}\big{)}\) denotes the values of initial pressure in the points of spatial grid, \(\phi\) and \(\theta\) are the spherical coordinates of the surface element \(dS^{\prime}\), and \(\Delta\phi\) and \(\Delta\theta\) are its angular sizes. Combining Equations (2) and (9), discretizing in time, a system of linear equations can be written:  
\[p\big{(}\mathbf{r}_{i},t_{j}\big{)}\approx\sum_{l=1}a_{l}^{ij}p_{0}\big{(}\mathbf{r}^{ \prime}_{l}\big{)}\rightarrow\mathbf{p}=\mathcal{A}\mathbf{P}_{0}, \tag{10}\]  
where the second equation is a matrix form of the first equation, \(p\big{(}\mathbf{r}_{i},t_{j}\big{)}\) is the pressure at time instant \(t_{j}\) recorded by the sensor located at \(\mathbf{r}_{i}\), \(\mathbf{p}=\big{\{}p\big{(}\mathbf{r}_{i},t_{j}\big{)}\big{\}}\). Here, \(\mathcal{A}=\left\{a_{l}^{ij}\right\}\) is the model matrix, which is defined by the geometry of the setup and the properties of the medium. The vector \(\mathbf{P}_{0}\) denotes the initial pressure in grid points, and is to be found, as it gives the distribution of optical absorption. A minimization problem  
\[\mathbf{P}_{0}=\text{arg}\underset{\mathcal{P}}{\text{min}}\|\mathbf{p}-\mathcal{A} \mathcal{P}\|^{2}+\mathcal{R}, \tag{11}\]  
where \(\mathcal{R}\) is a regularization term, has to be solved to find vector \(\mathbf{P}_{0}\). The numerical methods, such as LSQR, Moore-Penrose pseudoinverse [76], and others, can be used for this purpose. Summing up the benefits of different reconstruction algorithms, when a small number of receivers is used (\(N\lesssim 100\)), the images reconstructed by the time-domain algorithms have less artifacts than those reconstructed by Fourier domain methods. If large number of receivers is used, then the situation is reversed [73]. Fourier domain algorithms are generally fast, due to use of fast Fourier transform. However, they are applicable for real-time imaging only when planar detector geometry is used, and cannot account for non-uniform distribution of speed of sound. Model-based algorithms can be modified to account for different a priori factors and provide more accurate images, but requirements for memory and computational capacities are high. In our studies, we used back-projection algorithm for PAT, as it can be applied to different detector surfaces and can be parallelized for real-time operation. Another advantage is that obtained images are easily interpreted, that enables modification to account for non-uniform speed of sound distribution in a set of important cases.