Deconvolution of 3-D Gaussian kernels
aa r X i v : . [ phy s i c s . d a t a - a n ] A ug Deconvolution of 3-D Gaussian kernels
Z. K. Silagadze
Budker Institute of Nuclear Physics and Novosibirsk State University, 630 090,Novosibirsk, Russia.
Abstract
Ulmer and Kaissl formulas for the deconvolution of one-dimensional Gaussiankernels are generalized to the three-dimensional case. The generalization isbased on the use of the scalar version of the Grad’s multivariate Hermitepolynomials which can be expressed through ordinary Hermite polynomials.
Keywords:
Deconvolution, Gaussian kernels, Multivariate Hermitepolynomials, Hermite polynomials, Laguerre polynomials.
PACS:
1. Introduction
Mathematical problem of deconvolution of Gaussian kernels was consid-ered in [1]. In one-dimensional case, such a deconvolution (inverse problem)is equivalent to the solution of the Fredholm-type integral equation φ ( x ) = ∞ Z −∞ K ( x − ξ ; σ ) ρ ( ξ ) dξ, (1)where φ ( x ) is the known function, typically a detector response in a somemeasurement process of the unknown signal ρ ( ξ ) which is blurred by the finiteresolution σ of the detector. This blurring is described by the Gaussian kernel K ( x − ξ ; σ ) = 1 σ √ π e − ( x − ξ )2 σ . (2)The integral transform (Gaussian convolution) (1) arises in various problemsof applied physics, for example in calculations of transverse profiles of photon,proton or electron beams in medical applications [1, 2, 3, 4]. Preprint submitted to Physics Letters A August 21, 2019 sually the integral equation (1) is solved by a Fourier transform method.It is well known that, due to a presence of a fast growing Gaussian functionin the deconvolution integral in this case, this method may imply ill-posed in-verse problems [1]. Therefore alternative methods of deconvolution of Gaus-sian kernels, avoiding ill-posed inverse problems, were developed [1, 2].In particular, Ulmer and Kaissl provide in [1] two different expressionsfor the inverse kernel K − ( x − ξ ; σ ) defined through the relation ρ ( x ) = ∞ Z −∞ K − ( x − ξ ; σ ) φ ( ξ ) dξ. (3)In the next section we reproduce their results by a method which can beeasily generalized to the three-dimensional case.
2. Deconvolution of 1-D Gaussian kernels
Our starting point will be an identity [5] K ( x ; σ ) = e σ d dx δ ( x ) , (4)which can be easily proved by using the Fourier integral representation ofthe Dirac’s delta function: δ ( x ) = 12 π ∞ Z −∞ e ipx dp. (5)Because of (4), we can write φ ( x ) = ∞ Z −∞ ρ ( ξ ) e σ d dξ δ ( x − ξ ) dξ = ∞ Z −∞ δ ( x − ξ ) e σ d dξ ρ ( ξ ) dξ = e σ d dx ρ ( x ) . (6)Therefore the integral transform (1) can be written in the operator form φ ( x ) = ˆ Kρ ( x ) , ˆ K = e σ d dx . (7)2ence its inverse is given by ρ ( x ) = ˆ K − φ ( x ) , ˆ K − = e − σ d dx . (8)Now we apply the same trick as in (6): ρ ( x ) = ∞ Z −∞ δ ( x − ξ ) e − σ d dξ φ ( ξ ) = ∞ Z −∞ φ ( ξ ) e − σ d dξ δ ( x − ξ ) , which is the same as ρ ( x ) = ∞ Z −∞ φ ( ξ ) e − σ d dx δ ( x − ξ ) . (9)but, because of (4), we have e − σ d dx δ ( x − ξ ) = e − σ d dx e σ d dx δ ( x − ξ ) = e − σ d dx K ( x − ξ ; σ ) . (10)Therefore ρ ( x ) = ∞ Z −∞ φ ( ξ ) e − σ d dx K ( x − ξ ; σ ) , which shows that K − ( x − ξ ; σ ) = e − σ d dx K ( x − ξ ; σ ) = ∞ X n =0 ( − σ ) n n n ! d n dx n K ( x − ξ ; σ ) . (11)Thanks to the Rodrigues formula for Hermite polynomials [6] H n ( x ) = ( − n e x d n dx n e − x , (12)derivatives of the Gaussian function in (11) can be expressed through Hermitepolynomials: d n dx n K ( x − ξ ; σ ) = σ − n H n (cid:18) x − ξσ (cid:19) K ( x − ξ ; σ ) , (13)3nd we finally obtain K − ( x − ξ ; σ ) = ∞ X n =0 ( − n n n ! H n (cid:18) x − ξσ (cid:19) K ( x − ξ ; σ ) . (14)Another version of the inverse Gaussian kernel can be obtained by the useof the following identity e − σ d dx δ ( x − ξ ) = δ ( x − ξ ) + (cid:18) e − σ d dx − e − σ d dx (cid:19) e σ d dx δ ( x − ξ ) . (15)Using this identity and the identity (4) in (9), we get ρ ( x ) = ∞ Z −∞ φ ( ξ ) (cid:20) δ ( x − ξ ) + (cid:18) e − σ d dx − e − σ d dx (cid:19) K ( x − ξ ; σ ) (cid:21) dξ = ∞ Z −∞ φ ( ξ ) " δ ( x − ξ ) + ∞ X n =1 ( − σ ) n n n ! (2 n − d n dx n K ( x − ξ ; σ ) . (16)Therefore, in light of (13), we have K − ( x − ξ ; σ ) = δ ( x − ξ ) + ∞ X n =1 ( − n n n ! (2 n − H n (cid:18) x − ξσ (cid:19) K ( x − ξ ; σ ) . (17)Equations (14) and (17) were first obtained (with some typos) in [1].
3. Deconvolution of 3-D Gaussian kernels
Generalization to the three-dimensional case is now straightforward. Allwhat is needed is to change the differential operator d dx by its three-dimensio-nal version (Laplacian) ∆ = ∇ , change x and ξ by the corresponding three-dimensional vectors r and ξ , and use a three-dimensional Gaussian kernel K ( r − ξ ; σ ) = 1 σ √ π e − ( r − ξ )2 σ . (18)Then, instead of (11), we will end up with the equation K − ( r − ξ ; σ ) = e − σ ∆ K ( r − ξ ; σ ) = ∞ X n =0 ( − σ ) n n n ! ∆ n K ( r − ξ ; σ ) . (19)4owever now the action of the powers of Laplacian ∆ = ∇ on the three-dimensional Gaussian function cannot be as simply expressed in terms of theHermite polynomials, as in the one-dimensional case. In [2] the followingapproach was suggested. Using K ( r − ξ ; σ ) = K ( x − ξ x ; σ ) K ( y − ξ y ; σ ) K ( z − ξ z ; σ ) , and e − σ ∆ = e − σ d dx e − σ d dy e − σ d dz , (20)we get from the first equation of (19) K − ( r − ξ ; σ ) = F F F K ( r − ξ ; σ ) , (21)where F = ∞ X n =0 ( − n n n ! H n (cid:18) x − ξ x σ (cid:19) , F = ∞ X n =0 ( − n n n ! H n (cid:18) y − ξ y σ (cid:19) ,F = ∞ X n =0 ( − n n n ! H n (cid:18) z − ξ z σ (cid:19) . (22)Yet this form of inverse kernel doesn’t seem convenient in applications. Es-pecially if we note that a more straightforward generalization of equations(14) and (17) is possible if we use multivariate and tensorial generalization ofHermite polynomials introduced by Grad [7, 8]. Grad’s definition generalizesthe Rodrigues formula for the univariate Hermite polynomials (12): H ( n ) i i ...i n ( r ) = ( − n e r ∇ i ∇ i · · · ∇ i n e − r . (23)Using this definition, we easily get straightforward generalization of the equa-tion (14) in the form K − ( r − ξ ; σ ) = ∞ X n =0 ( − n n n ! H n (cid:18) ( r − ξ ) σ (cid:19) K ( r − ξ ; σ ) , (24)where H n (cid:18) ( r − ξ ) σ (cid:19) = δ i i δ i i · · · δ i n − i n H (2 n ) i i ...i n − i n ( r − ξ ; σ ) (25)5s the completely contracted version of the multivariate Hermite polynomials(the so called scalar, or irreducible, Hermite polynomials [9]).Analogously, the generalization of the equation (17) looks like K − ( r − ξ ; σ ) = δ ( r − ξ ) + ∞ X n =1 ( − n n n ! (2 n − H n (cid:18) ( r − ξ ) σ (cid:19) K ( r − ξ ; σ ) . (26)
4. Calculation of scalar Hermite polynomials
Scalar Hermite polynomials can be expressed in terms of the classical(associated) Laguerre polynomials [9, 10]. Below we elaborate this connectionand show that, in fact, Scalar Hermite polynomials can be expressed throughthe ordinary Hermite polynomials.Let us introduce the following operators [11]ˆ E = 12 r , ˆ F = −
12 ∆ , ˆ H = r · ∇ + 32 . (27)Using the commutator [ ∇ i , r ] = r i r , it can be easily checked that these operators obey commutation relations ofthe sl (2 , C ) Lie algebra [11]:[ ˆ H, ˆ E ] = 2 ˆ E, [ ˆ H, ˆ F ] = − F , [ ˆ E, ˆ F ] = ˆ H. (28)It is clear from (23) and (25) that H n ( r ) = e r ∆ n e − r = e ˆ r ∆ n e − ˆ r (1) = (cid:16) e ˆ r ∆ e − ˆ r (cid:17) n (1) , (29)where ˆ r is just a multiplication operator: ˆ rg ( r ) = rg ( r ) for any function g ( r ),and an expression ˆ A (1) means that the operator ˆ A acts on the function g ( r ) =1. The trick used here is borrowed from [12] and it is useful because thelast expression in (29) allows to use the operator identity (Baker-Hausdorffformula) [13] e ˆ A ˆ B e − ˆ A = ˆ B + [ ˆ A, ˆ B ] + 12! [ ˆ A, [ ˆ A, ˆ B ]] + 13! [ ˆ A, [ ˆ A, [ ˆ A ˆ B ]]] + · · · (30)6ith the help of the commutation relations (28), we can easily check thatin the case of the last expression in (29) the corresponding Baker-Hausdorffseries (30) in fact terminates and we end up with e ˆ r ∆ e − ˆ r = e E ( − F ) e − E = − F − H + 8 ˆ E. (31)Analogously, we have e − ∆4 (4 r ) n = e − ∆4 (4ˆ r ) n e ∆4 (1) = (cid:16) e − ∆4 (4ˆ r ) e ∆4 (cid:17) n (1) , (32)and e − ∆4 (4ˆ r ) e ∆4 = e ˆ F (8 ˆ E ) e − ˆ F = 8 ˆ E − H − F . (33)Comparing (33) and (31), we see from (29) and (32) that H n ( r ) = e − ∆4 (4 r ) n = e − ∆4 (2 r ) n . (34)Note that the relation (34) is completely analogous to an alternative expres-sion for the usual Hermite polynomials [14] H n ( x ) = e − ddx (2 x ) n . (35)Although in some calculations (35) is more useful than the Rodrigues formula,it was little known for a long time (see, for example, [15, 16]).Using ∆ r n = 1 r ddr (cid:18) r dr n dr (cid:19) = n ( n + 1) r n − , it can be proved by induction that∆ m r n = (2 n + 1)!(2 n − m + 1)! r n − m ) , if n ≥ m, (36)and ∆ m r n = 0, if m > n . Then (34) gives H n ( r ) = ∞ X m =0 ( − m n − m m ! ∆ m r n = n X m =0 ( − m n − m m ! (2 n + 1)!(2 n − m + 1)! r n − m ) , which, after introduction of a new summation index i = n − m , can berewritten as follows: H n ( r ) = ( − n n X i =0 ( − i i ( n − i )! (2 n + 1)!(2 i + 1)! r i . (37)7n particular, we get for the first few scalar Hermite polynomials: H ( r ) = 2 (2 r − , H ( r ) = 4 (4 r − r + 15) , H ( r ) = 8 (8 r − r + 210 r − , H ( r ) = 16 (16 r − r + 1512 r − r + 945) . (38)In light of the Legendre duplication formula [17](2 z + 1)! = 2 z +1 π − / z ! (cid:18) z + 12 (cid:19) ! , (39)we have 4 i (2 n + 1)!( n − i )! (2 i + 1)! = 2 n n ! i ! (cid:0) n + (cid:1) !( n − i )! (cid:0) i + (cid:1) ! = (cid:18) n + n − i (cid:19) n n ! i ! , and (37) takes the form H n ( r ) = ( − n n n ! n X i =0 ( − i i ! (cid:18) n + n − i (cid:19) ( r ) i . (40)On the other hand, the associated Laguerre polynomials L kn satisfy [14] L kn ( x ) = n X i =0 ( − i i ! (cid:18) n + kn − i (cid:19) x i . (41)Comparing (40) and (41), we see that H n ( r ) = ( − n n n ! L / n ( r ) . (42)However [14] L / n ( r ) = ( − n n +1 n ! r H n +1 ( r ) , (43)and we finally get a simple expression for the scalar Hermite polynomials interms of ordinary Hermite polynomials: H n ( r ) = H n +1 ( r )2 r . (44)8 . Concluding remarks As we see, the use of Grad’s multivariate Hermite polynomials allows toperform a straightforward generalization of Ulmer and Kaissl formulas forthe inverse problem of a Gaussian convolution to the three-dimensional case.In this generalization, the scalar Hermite polynomials play the same roleas ordinary Hermite polynomials in the Ulmer and Kaissl’s one-dimensionalformulas. The scalar Hermite polynomials by themselves can be expressed interms of associated Laguerre polynomials and through them, rather surpris-ingly, through ordinary Hermite polynomials. The latter connection becomesless mysterious thanks to the new generating functions for even- and odd-Hermite polynomials [18]: F ( x, t ) = ∞ X n =0 t n n ! H n ( x ) = 1 √ t exp (cid:18) tx t (cid:19) ,G ( x, t ) = ∞ X n =0 t n n ! H n +1 ( x ) = 2 x p (1 + 4 t ) exp (cid:18) tx t (cid:19) . (45)Indeed, since(2 r ) n = 4 n ( x + y + z ) n = X n + n + n = n n ! n ! n ! n ! (2 x ) n (2 y ) n (2 z ) n , (20), (34) and (35) imply that (44) is equivalent to the curious identity X n + n + n = n n ! n ! n ! n ! H n ( x ) H n ( y ) H n ( z ) = H n +1 ( r )2 r , (46)which, however, is a simple consequence of the generating functions (45),because F ( x, t ) F ( y, t ) F ( z, t ) = G ( r, t )2 r . The task of deconvolution arises in many fields, such as geophysics, sig-nal and image processing, underwater acoustics, astrophysics, high-energyphysics, etc. A general theory is beautiful, although very abstract [19]. Formore practical-minded introduction, see [20]. Naturally many algorithms ofdeconvolution were developed. Our aim in this note was not to suggest one9ore practical algorithm and test its robustness, but indicate that the algo-rithm described in [1] can be naturally generalized for the three-dimensionalcase.We hope that this refinement of the Ulmer and Kaissl method can findsome practical applications. Besides many applications in medical radiationphysics [2], Ulmer and Kaissl method was applied to get a (one-dimensional)inverse integral kernel for diffusion in a harmonic potential of an overdampedBrownian particle [21]. The method (its three-dimensional version describedin this article) was also used to obtain a non-local generalization of theSchr¨odinger equation due to quantum-gravity effects [22]. It can be alsouseful in some problems of non-local gravity [23, 24].Finally, the relation (44), which expresses scalar Hermite polynomialsin terms of ordinary Hermite polynomials can be useful in plasma physics(in [9, 10] scalar Hermite polynomials were expressed through the Laguerrepolynomials but the connection with the ordinary Hermite polynomials wasnot noticed).
Acknowledgments
We thank Gjergji Zaimi for informing us about [18]. The work of is sup-ported by the Ministry of Education and Science of the Russian Federation.
ReferencesReferences [1] W. Ulmer and W. Kaissl, The inverse problem of a Gaussian convo-lution and its application to the finite size of the measurement cham-bers/detectors in photon and proton dosimetry, Phys. Med. Biol. ,707–727 (2003).[2] W. Ulmer, Inverse problem of linear combinations of Gaussian convo-lution kernels (deconvolution) and some applications to proton/photondosimetry and image processing, Inverse Prob. , 085002 (2010).[3] W. D. Newhauser and R. Zhang, The physics of proton therapy, Phys.Med. Biol. , R155–R209 (2015).104] W. Ulmer and A. Carabe-Fernandez, Aspects of the physical principlesof the proton therapy with inclusion of nuclear interactions, J. ProtonTher. , 311 (2017).[5] R. G. Hohlfeld, J. I. F. King, T. W. Drueding and G. V. H. Sandri,Solution of Convolution Integral Equations by the Method of DifferentialInversion, SIAM J. Appl. Math. , 154–167 (1993).[6] G. Szeg¨o, Orthogonal Polynomials, 4th edition (American MathematicalSociety, Providence, 1975).[7] H. Grad, Note on N-dimensional hermite polynomials, Commun. PureAppl. Math. (1949) 325–330.[8] B. Holmquist, The d-Variate Vector Hermite Polynomial of Order k,Linear Algebra Appl. , 155–190 (1996).[9] R. Balescu, Transport Processes in Plasmas (Elsevier Science Publishers,Amsterdam, 1988).[10] D. Pfefferl´e, E. Hirvijoki and M. Lingam, Exact collisional moments forplasma fluid theories, Phys. Plasmas , 042118 (2017).[11] H. De Bie, An Alternative Definition of the Hermite Polynomials Relatedto the Dunkl Laplacian, SIGMA , 093 (2008).[12] A. W¨unsche, Generating Functions for Products of Special Laguerre 2Dand Hermite 2D Polynomials, Appl. Math. , 61552 (2015).[13] J. D. Bjorken and S. D. Drell, Relativistic Quantum Mechanics (McGraw-Hill, New York, 1964).[14] W. W. Bell,
Special Functions for Scientists and Engineers (D. VanNostrand Company, London, 1968).[15] F. M. Fern´andez, Generating functions for Hermite polynomials of ar-bitrary order, Phys. Lett. A , 189–191 (1998)[16] A. W¨unsche, Hermite and Laguerre 2D polynomials, J. Comput. Appl.Math. , 665–678 (2001). 1117] G. B. Arfken and H. J. Weber,
Mathematical Methods for Physicists,6th edition (Elsevier Academic Press, Amsterdam, 2005).[18] F. H-Yi and Z. D.-Hui, New generating function formulae of even- andodd-Hermite polynomials obtained and applied in the context of quan-tum optics, Chinese Phys. B , 060301 (2014).[19] S. Saitoh and Y. Sawano, Theory of Reproducing Kernels and Applica-tions (Springer, Singapure, 2016).[20] J. L. Mueller and S. Siltanen,
Linear and Nonlinear Inverse Problemswith Practical Applications (SIAM, Philadelphia, 2012).[21] T. Kosugi, Inverse Integral Kernel for Diffusion in a Harmonic Potential,J. Phys. Soc. Jap. , 054002 (2014).[22] M. Maziashvili and Z. K. Silagadze, Non-local imprints of gravity onquantum theory, arXiv:1812.03993 [quant-ph].[23] G. Calcagni, L. Modesto and G. Nardelli, Non-perturbative spectrum ofnon-local gravity, Phys. Lett. B , 391-397 (2019).[24] G. Calcagni, Taming the Beast: Diffusion Method in Nonlocal Gravity,Universe , no. 9, 95 (2018); Erratum: [Universe5