In order to improve the computational efficiency of the traditional approach to Monte Carlo simulations, we implement the semi-analytical approach (also known as the partial photon technique).^{8}^{,}^{27}^{,}^{28}^{,}^{42}^{,}^{43} Using this method, a portion of scattered intensity is collected after a photon reaches each new position within the medium. This “partial photon” intensity $Ipartial$ is calculated by multiplying the probability that the photon is scattered in the direction of the surface $F(\theta s,0)$ by the probability that it will reach surface according to the Beer-Lambert law $e\u2212(\mu s+\mu a)z$: Display Formula
$Ipartial=F(\theta s,0)\xb7e\u2212(\mu s+\mu a)z,$(34)
where $\theta s$ is the angle between $e^prop$ and the surface normal vector and $z$ is the distance to the surface. Modifying Eq. (11) under the semi-analytical approach, we find function $p(xs,ys)$ for different the different polarization channels as: Display Formula$pxx(xs,ys)=\u2211nIpartial\xb7W\xb7[1+Q(xs,ys)/I(xs,ys)],pxy(xs,ys)=\u2211nIpartial\xb7W\xb7[1\u2212Q(xs,ys)/I(xs,ys)],p++(xs,ys)=\u2211nIpartial\xb7W\xb7[1+V(xs,ys)/I(xs,ys)],p+\u2212(xs,ys)=\u2211nIpartial\xb7W\xb7[1\u2212V(xs,ys)/I(xs,ys)],$(35)
where the index of summation $n$ now represents the number of scattering events. The product of functions $p(xs,ys)\xb7pc(xs,ys)$ for the orthogonal polarization channels can be found in analogy with Eq. (36) as: Display Formula$pxy(xs,ys)\xb7pcxy(xs,ys)=\u2211nIpartial\xb7W\xb7DOCxy\xb7[1\u2212Q(xs,ys)/I(xs,ys)],p+\u2212(xs,ys)\xb7pc+\u2212(xs,ys)=\u2211nIpartial\xb7W\xb7DOC+\u2212\xb7[1\u2212V(xs,ys)/I(xs,ys)].$(36)
Scoring photons in this way enables both computational efficiency through collection of information at each scattering event and accuracy through collection of intensity in the exact backscattering direction.