Light transport in homogeneous tissue with m-dependent anisotropic scattering II: Fourier transform solution and consistent relations

This paper is the second of two focusing on the analytical solutions for light transport in infinite homogeneous tissue with an azimuth-dependent (m-dependent) anisotropic scattering kernel by two approaches, Case’s singular eigenfuncions (CSEs) expansion and Fourier transform, and proving the consistence of the two solutions theoretically. In this paper, the analytical solution for the m-dependent truncated scattering kernel was derived via the Fourier transform and inversion, and expanded with the m-dependent generalized singular eigenfuncions (GSEs). Two kinds of GSEs that are defined by Ganapol in the case m=0 are extended to arbitrary azimuthal orders and proven to be consistent with CSEs both in expression forms and in intrinsic behaviors. By applying the Fourier transform inversion on the solution for the three-term recurrences, the Green’s function of radiance distributions is obtained successfully, and it conforms perfectly to the CSEs solution in the limit, which has already been discussed in our first accompanying paper. Meanwhile, as a byproduct, a series of identities about the m-dependent Chandrasekhar orthogonal polynomials were presented and will be greatly helpful for further studies.


Introduction
In our first accompanying paper we have derived Case's singular eigenfunctions (CSEs) solution of the light transport equation in the infinite homogeneous tissue with m-dependent anisotropic scattering kernel and benchmarked four computational methods for the mdependent Chandrasekhar orthogonal polynomials. In this paper we presented the corresponding Fourier transform solution and then showed the consistence of the solutions obtained from the two approaches.
The CSEs solutions are always seen as the exact solutions of integro-differential transport equations and naturally suitable for the light transport in turbid tissue [1][2][3][4][5][6][7][8][9]. Fourier transform method is a useful method and applied successfully to the transport equations in the regular space [10][11][12][13][14][15][16][17]. Ganapol [12][13][14][15][16] derived the Fourier transform solution for the infinite medium with azimuthal symmetric scattering kernel and showed the consistent theory in a more mathematically unified manner. Machida [17] used the Fourier transform and rotated reference frames technique to compute three-dimension Green's function for three-dimension linear Boltzmann equation. Their outstanding work opened up a new scope to solve the exact radiance distributions that we are seeking for the light transport in biological tissue.
Firstly, by introducing the associated Chandrasekhar polynomials of the second kind, the general analytical solutions for three-term recurrence in Fourier k-domain are extended to the azimuth-dependent or m-dependent case explicitly as the base of the Fourier transform solutions. Secondly, in m-dependent case, taking the associated Chandrasekhar polynomials of the first/second kind as the Legendre moments, the generalized singular eigenfunctions of the first/second kind (GSE1, 2) are defined, in the limit version and in the truncated version, in terms of the associated Legendre polynomials of the first kind with the same degrees and orders. And also, the implicit corresponding relations between GSEs and CSEs are pointed out and proven from the perspective of both expression forms and intrinsic behaviors. Depending on the first Legendre moment, the general analytical solution can be expressed completely and used to Fourier inversion. Lastly, in the Fourier inversion, the contributions from the poles and the branch cut are combined to the final Fourier transform solution. The poles in complex plane correspond to the discrete values in Case's method and the branch cut corresponds to the continuous spectrum, and then the final solution confirms perfectly to the CSEs solution.

M-dependent polynomials and C-D formulas
In our following derivations, many m-dependent polynomials or functions and corresponding C-D formulas will be used frequently. However, there are some aspects different from that of the case 0 m = . In this section, we give brief reviews of these polynomials and formulas.
In the real domain [20][21][22], ( )( ) In the complex plane, ( ) l P z is a Legendre orthogonal polynomial [20] and ( ) m l P z is defined as [20][21][22] ( ) ( ) and satisfies the Eq. (2) with initial terms [17] ( ) ( where z is not on the cut line −1 to 1. Similarly, we define the associated Legendre function of the second kind as [17,20] ( ) ( ) ( ) ( ) Q z should be noticed during the derivation. We have already derived these formulas, and listed here for future use. In the following formulas, L is the truncated position of the phase function, that is to say when l L > , 0 l f = . N is an integer, and N L > .
s t ϖ μ μ = is the single particle albedo, where s μ is the scattering coefficient, a μ is the absorption coefficient, and t a s μ μ μ = + is the total attenuation coefficient.

Standard Fourier transform and matrix inversion
Let's begin with the Eq. (4) in our first accompanying paper: where jl δ is the Kronecker delta,

Solution of three-term recurrence
From Eq. (4), we can also get the three-term recurrence are the coefficients to be determined. It is noted that the analytical solution form expressed in Eq. (38) may be different from that given by Machida [17], but they are identical essentially.
We assume that the first Legendre moment of  ( ) Substituting Eq. (39) into Eq. (37), to keep balance between the equations, the following relations must be satisfied Thus, the analytical solution Eq. (37) can be expressed as see that if the first Legendre moment is known, the other Legendre moments will be easily got. The first Legendre moment plays an important part during the derivation.

The first Legendre moment
The first Legendre moment is so important that Ganapol had provided several methods of solving it [12][13][14][15][16]. When 0 m = , those methods are all valid. However, when 0 m > , it is impossible to use the isotropic point source to find the first Legendre moment. So, we substitute the Eq. (48) into Eq. (30), let j m = , and use where Vol. 9, No. 9 | 1 Sep 2018 | BIOMEDICAL OPTICS EXPRESS 4038 In our first accompanying paper, we defined a holomorphic function in the complex plane excluding a cut line −1 to 1 on the real axis obvious that Eq. (52) and Eq. (53) are identical.

The generalized singular eigenfunctions of the first kind
Now, we can begin to build the generalized singular eigenfunctions of the first kind (GSE1) which are extended from the Ganapol's definition for the case 0 m = in the complex z-plane [15,16].
Let's define GSE1, called the limit version, as the expansion series of the combination of the associated Chandrasekhar polynomials of the first kind and the associated Legendre polynomials of first kind and define the corresponding truncated version Substituting Eq. (13) into Eq. (56), we get According to the continuity of ( ) m l g z at the truncated position, we have Using Eq. (22), the coefficients can be solved From Eq. (20), we get Then, as N → ∞ and v is real discrete eigenvalues, the limit version of GSE1 is From Plemelj's formula [23] we get ( ) Here P denotes the Cauchy principal value. From Eq. (6) and Plemelj's formula we get From Eq. (69) we get and Equation (85) is the evidence of the analyticity of ( ) So the consistency of eigenfurnctions form is proven. When 0 m = , our GSE1 is identical to the Ganapol's GSE1. Of course, GSE1 satisfies the normalization condition Eq. (9) in our first accompanying paper naturally.

The generalized singular eigenfunction of the second kind
Similar to GSE1, we define the limit version of the generalized singular eigenfunction of the second kind (GSE2) As the derivation in the previous section, when l L > we get (96) Similar to the derivation in previous section, as N → ∞ and z is discrete eigenvalues, we get Further, when z approaches the cut line −1 to 1 from above ( + ) and below (-), we can also get So as N → ∞ , we get the analogue of the Case's singular eigenfunctions in the continuous spectrum When 0 m = , our GSE2 is identical to the Ganapol's GSE2 [15,16]. And also, from Eq. (63) and Eq. (95), we get Considering v in Eq. (84) and Eq. (98) can be extended to any z in complex z-plane, we get

Fourier domain analytical solution
It is time to build the Fourier domain analytical solution to prepare for the Fourier transform inversion. From Eq. (27) and Eq. (48) we can build the Fourier domain analytical solution (120) Thus, from Eq. (103) and Eq. (106), we get

Fourier transform inversion
Applying Cauchy's integral theorem [24] to Fourier transform inversion is a popular and useful technique. To use the Cauchy's integral formula, we must construct a contour which includes the real integral on the real axis specific to Fourier transform inversion in the complex z-plane or k-plane. However, because Jordan's lemma [24], we will immediately conclude that as R → ∞ the integrals over R C + and R C − go to zero. As 0 ε → the integral over C ε goes to zero because k i = or 1 z = is not the pole. So we only need to consider the contributions from the poles and the branch cut.

Application of the residue theorem
Before applying the residue theorem, we need to find out all the valid poles. From Eq. It is easy to verify that the zeroes of ( ) m L z Λ are valid poles, but the zeroes of ( ) , m L z z ϕ are removable singularities. Ganapol [16] proved that is true in the case 0 m = . In the case 0 m > , the derivation is totally similar, so we will not reiterate here.

Integral over the branch cut
As shown in Fig. 1, the branch cut is [ ) , i i∞ on the imaginary axis in the upper half plane. We 1, k i y y ε ± = ± + ∈ ∞ and k + represents the integral path right side of the branch cut and k − the left side. So z i k ± ± = and the branch cut in k-plane changes to the branch cut ( ] 0,1 on the real axis in z-plane, as shown in Fig. 1 (a). So we get The integral over double sides of the branch cut is From above derivation we only need to consider the terms in Eq. (121) which include From Eq. (80) and Eq. (84), we get From Eq. (80) and Eq. (82), we get So the contributions from branch cut is

Fourier transform inversion solution
Thus, we can build the Fourier transform inversion solution in the position space from Eq.
Since now z has become real, as N → ∞ , we get Now, the last task is to confirm the norms in Eq. (142) or Eq. (143) are equal to those norms in our first accompanying paper, though it is apparently in 0 m = case. In fact, from Eq. (119), this is undoubted without much efforts due to the common zeroes.
So we can finally conclude that the Fourier transform solution is consistent with the Case's singular eigenfunctions solution. The contributions from the poles in Fourier transform inversion are equivalent to the contributions of discrete eigenvalues in Case's method, and correspondingly, the branch cut contribution is equivalent to that of the continuous spectrum.

Conclusion
Both the Case's singular eigenfunctions expansion method and Fourier transform method are applied to the light transport integro-differential equations in the infinite homogeneous tissue with m-dependent anisotropic scattering kernel, and the consistence is proven from the singular eigenfunctions and the final series solution. However, some remarks are worth to be mentioned.
The form of solutions for three-term recurrences may be the most valuable point and can be combined with other methods, such as N P method [26,27]. We can predicate that the general analytical solutions that consist of the general solution and particular solution may be applied to more approaches to seek the solutions of light transport equations in one/two/threedimension full space, half space, and slab space etc. with general scattering kernel, especially when Case's method is not suitable. These general solutions and particular solutions are usually all related to the two kinds of the associated Chandrasekhar polynomials, so the value of the polynomials of the second kind may be developed intensively.
During our derivation, the identities involving Chandrasekhar polynomials and Legendre polynomials are often confusing due to the disunity of the definitions of some polynomials, such as the associated Legendre polynomials of the second kind. It is noted that the three-term recurrences for the associated Legendre polynomials of the second kind may be different from those in other books [20][21][22] in the first term. As Ganapol [16] and other experts pointed out, some functions in the derivation, such as delta function ( )