This paper was awarded in the IV International Competition (1995/96) "First Step to Nobel Prize in Physics" and published in the competition proceedings (Acta Phys. Pol. A 92 Supplement, S-9 (1997)). The paper is reproduced here due to kind agreement of the Editorial Board of "Acta Physica Polonica A".

Analysis of Light Fluctuations from Photon Counting Statistics

A.V. KURASHOV

Kiev Physical and Mathematical Lyceum No. 145, Shota Rustavely 46, Kiev, Ukraine

Abstract

Reconstrution of light intensity fluctuations from photon-counting data is a typical ill-posed inverse problem for non-Hermitian transforming operator. Such problems cannot be solved with the help of eigenfunction representation which diagonalizes the operator and simplifies the regularization of a quasi-solution. We describe an analytical method to solve this inverse problem using the Poisson transform operator representation in a mixed basis which diagonalizes the Poisson transform. On the base of the proposed method we performed numerical experiments with typical intensity distributions. We showed that the reconstruction procedure is stable and may be used even when the level of statistical errors is relatively high. In conclusion we describe the results of the inversion procedure in application to the processing of experimental data on photon counts statistics.

PACS numbers: 42.50.Ar

1. Introduction

A common problem in the investigation of light fluctuations is to infer the probability density function of intensity from those of the photon counts. In the semiclassical approximation these two functions are related through the Poisson transform [1]. This problem, which is known as an inverse problem of photon counts statistics has been a subject of extensive studies. It is a typical ill-posed problem and it has the same features as a well known problem of blurred and distorted images restoration [2]. Analytical methods [3, 4] have been demonstrated to work when the corresponding functional relationship is known. Obviously, these methods are unusable for experimental data processing. Even small statistical errors cause strong variations in the resulting estimation of light intensity probability density function (PDF) when PDF of photon counts is estimated from a measured histogram. Some methods of numerical solutions of this problem have been proposed [5-9]. One procedure [5] makes use of the Taylor expansion of the light intensity PDF characteristic function whose coefficients are given by the photon counts PDF. Unfortunately, the limited number of values of the photodetected histogram makes the expression of such expansion valid only in a narrow region of characteristic function argument and prevents its application for experimental data processing. Recently, a generalization of this method with the help of Pade approximation has been proposed [8] and has demonstrated good results for low level photon rates. The problem there is the selection of a Pade approximation from the number of acceptable choices, thus additional physical constraints and consistency with experimental data is required [9]. Byrne et al. [7] gave a stable numerical procedure of the intensity PDF estimation by means of Tikhonov regularization [10] which provides the minimum norm of the solution in the Hilbert space. But as all such methods, this technique does not permit to estimate the necessary regularization level and requires additional physical considerations.

Bedard [6] photon counts statistics inverse problem solution is based on the PDF's expansion in terms of the orthogonal Laguerre series. It is a special case of the common input function orthogonal expansion which is the most useful method to solve such problems. This procedure does not permit to determine error in PDF estimation and is not optimized with respect to the Poisson transformation operator, which is its principal disadvantage.

In this work we suggest a new approach to solve the photon counts statistics inverse problem which is based on the optimized PDF's orthogonal expansions and the transform operator representation in a mixed basis. The main advantages of such approach are: (a) maximum simplicity of the solution obtained -- direct and inverse matrices in this representation are diagonal; (b) the possibility to check the regularization level and the estimation errors by simple restriction of eigenvalues taken into account in the inverse matrix; (c) the necessity to accomplish numerical calculations only on analytically defined symmetrical matrices of limited dimensions, which are given by the maximum number of registered photon counts (conventionally it does not exceed 100).

2. Probability distribution of the photoelectrons

Let us briefly consider the relationship between the statistics of the photoelectron counting and the light falling on it (the direct problem of photon counting statistics). From the experimental measurements of the photoelectron counting one can draw conclusions on the statistical properties of light which in turn have wide applications to various fields of science. The last problem is the inverse problem of photon counting statistics and it is the main problem of our consideration.

Consider a beam of light falling on a photoelectric detector. Assume that n photoelectrons are ejected in a certain time interval T. Of course, only the photoelectrons and not the photons are observable and our discussion must therefore be confined to the statistical behavior of the photoelectrons. Although it is tempting to associate the ejection of a photoelectron with the arrival of a photon, this picture becomes inadmissible by the uncertainty principle for time intervals shorter than the inverse frequency band width of light. Restriction to a ``homogeneous'' beam ensures that there are no large phase differences between different parts of the beam.

The fluctuations of the number of particles n therefore have two causes. First, there are fluctuations of the wave intensity I(t), and second, the stochastic fluctuations of particles. This twofold source of the fluctuations results in a departure from classical statistics, as was shown by Mandel [1].

Let P(n,T) denote the probability that n photoelectrons are ejected in a time interval T. It can be shown in a usual way that under the common assumptions that P(n,T) is a Poisson distribution in n:

$$P(n,T)=\exp(-\eta W)\frac{(\eta W)^n}{n!},\eqno(1)$$

where W is the integrated intensity of light beam, $W=\int_0^TI(t){\rm d}t$, $\eta$ is the quantum efficiency of the photodetector. P(n,T) is not, however a distribution that can be found experimentally because it is a stochastic function of time itself. As P(n) fluctuates stochastically the given interval T is unique and only the ensemble averages, which are equal to the time averages for an ergodic process, is observable. But the operation of averaging of P(n,T) over the time when I(t) is fluctuating will not, in general, result in some other Poisson distribution. This departure in the observed distribution from the classical form can therefore be seen to be a general consequence of the association of particles with fluctuating waves. To derive the general semi-classical expression for P(n,T) we can average Eq. (1) over the fluctuations of integrated intensity W. In this case we get

$$P(n,T)= \left\langle\exp(-\eta W)\frac{(\eta W)^n}{n!}\right... ...fty P_0 (W)\exp(-\eta W)\frac{(\eta W)^n}{n!}{\rm d}W,\eqno(2)$$

where P₀(W) is the PDF of classical integrated intensity fluctuations. A complete quantum mechanical derivation of the photoelectron counting formula gives the same expression for P(n,T), except that P₀(W) is, in general, not a positive defined function [11]. We shall not discuss further these special features. Expression (2) is known as the Mandel formula and in essence it is a Poisson transform of the function P₀(W). If the probability distribution P₀(W) is known, we can evaluate P(n,T) and so the direct problem of photon counting statistics can be solved. As we shall show in the next section, the inverse problem has a much more difficult structure.

3. Diagonal representation of Poisson transformation in mixed basis

The inverse problem of photoelectron-counting statistics consists in finding the probability density of integral intensity W from the experimental data obtained for the probability distribution P(n). Mandel's equation (2) can be expressed in the functional form as follows:

$$P(n)=\widehat LP_0(x)=\int_0^ \infty h(n,x)P_0(x){\rm d}x=\int_0^\infty\exp(-x)\frac{x^n}{n!}P_0(x){\rm d} x.\eqno(3)$$

(Without loss of generality, we assume here a detector with quantum efficiency =1 (or x=W), so that the mean of photon counts n equals that of integrated light intensity x.) The inverse problem means the evaluation of PDF P₀(x) from the known function P(n).

Effective method of regularization of any inverse problem $\widehat Ly=z$ gives the solution expansion in appropriate finite-dimensional orthogonal basis. If the mapping operator $\widehat L$ is self-conjugated it is convenient to choose its own eigenvectors as such basic functions. The main advantage of this procedure is the simplicity of quasi-solution retrieval, which is caused by diagonality of forming operator. Actually, a self-conjugated operator maps its domain onto itself, and we can use the same system of basic functions, that is the eigenfunctions of operator $\widehat L$, to represent source and outlet functions (z and y respectively). In this case decomposition coefficients of y (a_i) and z (b_i) are connected by a simple equation

$$b_i=\lambda_ia_i, \eqno(4)$$

which can be easily inverted. An advantage of this method consists not only in its simplicity, but in the possibility of inverse problem regularization by limiting of the number of eigenvalues $\lambda_i$ taken into account.

Unfortunately, this method cannot be used directly when $\widehat L^+\neq\widehat L$. In particular, a class of such problems is formed by all inverse problems in which the set of values of operator $\widehat L$ does not coincide with its domain, so that {y} and {z} sets are different. A typical problem of this class is the inverse problem of photoelectron-counting statistics. In fact, the operator $\widehat L$ in Eq. (3) maps the set of integrable positive functions P₀(x) on the set of infinite-dimensional vectors P(n), which are defined on the integer variable $n=0,\dots$ To solve such problem we propose to use mixed representation of the operator $\widehat L$ in a basis, which has the main feature of a self-conjugated operator in its own representation -- the diagonality of corresponding matrix. This representation preserves all advantages of eigenrepresentation in inverse problem regularization but the class of solvable problems is substantially longer. Indeed, Eq. (4) which relates the coefficients of the orthogonal decomposition of P₀(x) and P(n) in corresponding spaces is still valid if the basic functions are defined by the equation,

$$\widehat L\varphi_i(W)=\lambda_if_ i(n),\eqno(5)$$

where

$$\sum_0^\infty f_i(n)f_j(n)=\delta_{ij},\qquad\int_0^ \infty\varphi_i(x)\varphi_j(x){\rm d}x=\delta_{ij}. \eqno(6)$$

The existence of these functions should be proved as well as the required properties. We shall show that such basic sets do exist and they can be determined with sufficient accuracy by means of a simple numerical procedure.

At first, we have to point out an essential refinement of the transformation (3) formal definition. To determine the problem with the required mathematical accuracy, let us define the domains of P₀(x) and P(n) as finite intervals in x and n spaces, respectively, so that

$$h(n,x)= \left\{\begin{array}{ll}h(n,x),&0\le n\le N_0,~0\le x\le X,\\ 0,&n>N_0,~x>X. \end{array}\right.\eqno(7)$$

The bounds on x and n can be determined from the physical arguments. Under these restrictions the operator $\widehat L$ becomes of the Hilbert-Schmidt type [12] and therefore it is compact. But this restriction causes the output function space to be finite-dimensional. Thus the functions $\varphi_i(x)$, if they really exist, do not form a complete set for the PDF functions $P_0(x)\in L_2(0,~X)$. This feature is inherent in all regularization methods but we must constantly keep in mind that it introduces an extra systematic error.

Now it is easy to verify that Eqs. (5), (6) are satisfied if functions $\varphi_i(x)$ and vectors f_i(n) are the solutions of the following equations:

$$\int_0^x R_0(x,x')\varphi_i(x'){\rm d}x'=\lambda^2\varphi_i(x), \eqno(8)$$

$$\sum_{n'=0}^{N_0}R(n,n')f_i(n')=\lambda_i^2 f_i(n).\eqno(9)$$

Here R₀ and R are left and right iterated kernels of the operator L respectively

$$R_0(x,x')=\sum_{n=0}^{N_0} h(n,x)h(n,x'),\quad R(n,n')=\int_0^ xh(n,x)h(n',x){\rm d}x.\eqno(10)$$

Thus we can decompose P₀(x) and P(n) into the orthogonal series

$$P_0(x)= P'_0(x)+P_\perp(x)\equiv\sum_{i=0}^{N_0} a_i\varphi_i(x)+P_\perp(x),\eqno(11)$$

$$P(n)= P'_0(n)+P_\perp(n)\equiv\sum_{i=0}^{N_0}b_if_i(n)+P_\perp(n), \eqno(12)$$

where $P_\perp(x)$ and $P_\perp(n)$ are the projections of P₀(x) and P(n) on the subspaces, orthogonal to $\{\varphi_i(x)\}$ and {f_i(n)}, respectively. With the help of Eq. (5) we get the estimation $\widetilde a_i=b_i/\lambda_i$ for the coefficients a_i, which is formally equivalent to Eq. (4). The problem becomes a correct one and the solution regularization level can be simply adjusted by restriction of considered values $\lambda_i$ in the finite-dimension decomposition

$$\widetilde P_0(x)=\sum_{i=0}^N\frac{b_i}{\lambda_i}\varphi_i(x), \quad0\le N\le N_0.\eqno(13)$$

It is interesting to notice that the dimension of eigenvector space {f_i(n)} cannot be greater than N₀ and thus the dimension of the space of reconstructed functions $\widetilde P_0(x)$ is smaller than N₀. As a rule, the real single photon counting experiments are performed with the values of $N_0(10\dots100)$. Thus the numerical calculations of $\widetilde P_0(x)$ become simple because of a rather small data array. Really as it is follows from Eqs. (5), (8), (9) functions $\varphi_i(x )$ can be defined as

$$\varphi_i(x)=\frac1{\lambda_i}\sum_{n=0}^{N_0}h(n,x)f_ i(n).\eqno(14)$$

Thus there is no need to solve numerically Eq. (8) which is usually the most complicated problem in eigenfunction methods. It is sufficient to find eigenvectors of matrix R(n,n ') only and functions $\varphi_i(x)$ can be obtained from Eq. (14) with any precision in the variable x. Moreover, this procedure is performed only once and the experimental data processing can be done with standard sets of functions f_i(n) and $\varphi_i(x)$.

4. Error estimation

The question of principal importance in every method of the inverse problem solving is error estimation in the reconstructed probability density. As it follows from the solution structure, there are three main kinds of errors:

1) a systematic error P_c(x) caused by the projection of P₀(x) onto the subspace which is orthogonal to the subspace of quasi-solutions P₀(x):

$$\delta P_c(x)=P_0(x)-P'_0(x)=P_0(x)-\sum_{i=0}^{N_0}a_i\varphi_i(x),\eqno (15)$$

$$\int_0^xh(n,x)\delta P_c(x){\rm d}x=0,\quad n\le N_0;\eqno(16)$$

2) a systematic error P_N of the projection of P₀(x) onto the subspace of regularized quasi-solutions with the dimensionality N<N₀:

$$\delta P_N=P'_ 0(x)-P_N(x)=\sum_{i=N+1}^{N_0}a_i\varphi_i(x),\quad N<N_0;\eqno(17)$$

3) finally, a random error P_E which is due to the statistical fluctuations of the input data in the measured histogram $\widetilde P(n)$.

A quantitative estimation of the first type errors is rather complicated. We can use the following consideration. As it can be shown, every residual P_c(x) is mapped by the operator L on a subspace which is orthogonal to the measured functions P'(n) until n<N₀, and therefore we can determine P_c (x) only from the data for n>N₀. Because of the assumption that P(n)=0 for all n>N₀, the statistical error in this region is too large and no reasonable estimations can be obtained. It means that this error is independent of the regularization method and it can be eliminated only by a significant increase in data accuracy in the asymptotical limit. With data sample of limited size it is practically impossible. As it will be shown below, this error is sufficiently small except in extremely adverse cases, for example when P₀(x) has singularities. If P₀(x) is known, the norm of such type error can be found

$$\varepsilon_c^2=\Vert\delta P_c\Vert^2=\int _0^x[P_0-P'_0]^2... ...^x P_0^2(x){\rm d}x-\sum_{i=0}^{N_0}\vert a_i\vert^2 \eqno(18)$$

and this is the smallest error of the solution achieved for n<N₀.

The norm, N, of the finite-dimensional approximation P_N(x) error equals the remainder of the sum in P'₀(x):

$$\varepsilon_N^2=\int_0^x[P'_0(x)-P_N (x)]^2{\rm d}x=\sum_{i=N+1}^{N_0}\vert a_i\vert^2.\eqno(19)$$

Finally, to estimate the norm of random error let us assume that

$$\widetilde P(n )=P(n)+\delta P(n),\eqno(20)$$

where P(n) are statistical errors in experimental histogram. Strictly speaking, $\widetilde P(n)$ is the function of P(n). If we expand $\widetilde P(n)$ in the orthogonal series (12) and perform the inversion procedure we can find

$$\langle\varepsilon_R^2\rangle=\int_0^x\vert P_N(x)-\widetilde... ...0}^N\frac{\langle\delta b_i^2\rangle} {\lambda_i^2}.\eqno(21)$$

It has been taken into account that

$$\widetilde b_i=b_i +\delta b_i,\quad\delta b_i=\sum_{i=0}^{N... ...mbda_i}+\frac{\delta b_i}{\lambda_i}=a_i+\delta a_i.\eqno (22)$$

The resulting error estimation, $\langle\varepsilon_\Sigma^2\rangle$, with the help of Eqs. (18),(19),(21) can be found to be

$$\langle\varepsilon ^2_\Sigma\rangle=\varepsilon^2_c+\varepsi... ...=0}^N\frac{\langle \delta b_i^2\rangle}{\lambda_i^2}.\eqno(23)$$

As evident from Eq. (23), $\langle\varepsilon_\Sigma^2\rangle$ consists of three components. The first one does not depend on data accuracy and regularization level thus it determines the minimal error of the method. The last two terms determine the operational error which depends on regularization conditions. The first of them decreases and the second one increases with N so that the optimal estimation $\langle\varepsilon_\Sigma^2\rangle\rm _{min}$ of the dimension $N \rm _{opt}$ exists. If the rms errors $\langle \delta b_i^2\rangle$ and the coefficients a_i are known, the definition of $N\rm _{opt}$ should present no problems. But under the real experimental conditions these quantities are not known and they have to be estimated indirectly.

5. Computer simulation

In this section we discuss some typical models of light intensity statistics and show the results of their PDF's numerical restoration from the computer simulated histograms of photon counts.

5.1.. Bose-Einstein distribution

If polarized single mode laser light is scattered by a rough moving surface it can be predicted that the complex amplitude distribution is Gaussian, and so the distribution of intensity I(t,s) is geometrical (exponential) [13]. If it is assumed that the detector area is much smaller than the average size of speckle spot and the counting time T is much smaller than the light coherence time, the integrated intensity distribution is

$$P_0(x)=\frac1{\langle x\rangle}\exp\left(-\frac x{\langle x\rangle}\right).\eqno(24)$$

The corresponding photon counts PDF has the form of the Bose-Einstein law

$$P(n)=\frac{\langle x\rangle^n}{(1+\langle x\rangle)^{n+1}},\eqno(25)$$

where $\langle x\rangle= \langle W\rangle\langle n\rangle$. This model is of great importance not only because of its extreme simplicity, but it has many applications to practical problems of speckle optics as well. That is why this model is the most suitable for computer simulation of inverse problem solving.

We performed the computer experiments on numerical reconstruction of distribution (24) from $\widetilde P(n)$ data under various conditions. The histogram was generated with the help of random trial scheme, the size of sample M was 10⁵ or 10⁶. All experimental data were obtained for the parameters $\langle n\rangle=4,~N_0=20,~X=21.2$. The number of photon counts n was accumulated in the frequency histogram $\widetilde P(n)$. The norm of error P(n) averaged over 5 histograms was 5.12 x 10^-3 for M=10⁵ and 1.91 x 10^-3 for M=10⁶.

The typical results of PDF reconstruction are shown in Figs. 1 and 2. As can be seen from Fig. 1, essentially exact solution is obtained for very small values of N. The influence of random errors in this region is insignificant, and the resulting error is nearly equal to the projected one. Figure 2 shows the error norm as the function of the projection dimensionality N averaged over 5 realizations of histogram. The value of error for the exact data is shown too. The fluctuation error grows rapidly with N and becomes leading at N=5. It should be noted that the systematic error $\varepsilon_c^2$ in this case is extremely small $(\varepsilon_c^2<10^{-6}$ for N=20), and thus may not be considered. Such rapid divergence of the series (11) follows from the fact that the eigenfunction P₀(x) has the exponential form which matches the distribution (24). If the convergence of the series (11) is slower, the influence of systematic error $\varepsilon_c^2$ is greater, and in some cases it becomes the main factor determining the accuracy.

$\parbox{0.9\textwidth}{{\small Fig.~1.}~Reconstruction of exponentia... ...langle n\rangle=4,~N_0=20,~X=21.2,~\vert\delta P(n)\vert=4.47\times10^{- 3}$.}$

$\parbox{0.9\textwidth}{{\small Fig.~2.}~Average resulting error for ... ...2\times10^{-3},~(c)~\langle\vert\delta P(n)\vert\rangle= 2.81\times 10^{-3}$.}$

5.2.. Log-normal distribution

The next example corresponds to the physical model of coherent radiation propagating in a turbulent atmosphere. Prevailing number of experimental data indicates that the intensity fluctuation distribution in this case obeys the log-normal law [14]:

$$P_0(x)=\frac1{2\sqrt{2\pi}\sigma_\chi\langle x\rangle} \exp\... ...\rangle}+2\sigma_ \chi^2\right)^2\right],\quad x\ge0.\eqno(26)$$

It is assumed that the logarithm of field amplitude is Gaussian random variable with mean value $\langle x\rangle$ and dispersion $\sigma_\chi^2$.

We obtained the except distribution P(n) for log-normal fluctuations by direct numerical calculations of transform (3) for intensity PDF (26). Then the imitation of experimental histogram was found by means of random term P( n) which was added to exact P(n). (If the result was negative, we took $\widetilde P(n)=0$.) The norm of random value P(n) was taken 2.5 x 10^-3 which corresponded to M=10⁵. The results obtained for the inverse problem are shown in Figs. 3 and 4. It can be seen that

$\parbox{0.9\textwidth}{{\small Fig.~3.}~Reconstruction of log-normal distribution.}$

$\parbox{0.9\textwidth}{{\small Fig.~4.}~Reconstruction error for log-normal PDF for $\langle n\rangle=4$, $\vert P(n)\vert=2.5\times10^{-3}$.}$

the inversion of the received data gives almost correct results in this case too. But in this case the accurate inversion is achieved with much longer required dimensionality of the projection P_N(x).

6. Experimental results

We used the method discussed above to determine the intensity PDF for experimental data received for two partially polarized superposed light beams which were scattered by a diffuse surface. As it was shown earlier [15] this model permits to generate the variety of light intensity distributions and so it is convenient in testing other experiments. The experimental setup and the details of experimental procedure are described in [15].

The results of intensity PDF reconstruction are separated in two groups. The first deals with two-mode radiation with independent intensity fluctuations, the second corresponds to modes with partially correlated phases (non-stationary interference). Theoretical expressions for these cases can be approximately given in the form

$$P_0(x)=\frac1{p\langle n\rangle}\left\{\exp \left[-\frac{2x}... ...ft[-\frac{2x}{(1-p) \langle n\rangle}\right]\right\},\eqno(27)$$

$$P_0(x)\approx\frac1{\pi\sqrt{m ^2x^2_0-(x-x_0)^2}},\quad x_0=x_1+x_2,\eqno(28)$$

respectively. Parameters p (polarization degree), m (modulation coefficient), x_i, i=1, 2 (beam intensities) as well as the mean value $\langle n\rangle$ in the theoretical models have been chosen to fit experimental histograms with minimum error norm. Typical experimental histograms for these two types of intensity PDF are shown in Fig. 5. The results of reconstruction are shown in Figs. 6, 7. As it is evident from Fig. 6,

$\parbox{0.9\textwidth}{{\small Fig.~5.}~Experimental histograms for:... ...0.5$, and ($b$) non-stationary interference, $\langle n\rangle=10.2,~m=0.83$.}$

$\parbox{0.9\textwidth}{{\small Fig.~6.}~Reconstruction of 2-mode PDF for histogram ($a$) in Fig.~5.}$

$\parbox{0.9\textwidth}{{\small Fig.~7.}~Reconstruction of PDF for non-stationary interference (histogram ($b$) in Fig.~5).}$

the reconstruction of smooth intensity distribution is quite good no matter whether it has the monotone form or not. But it is not so for extremely singular functions. Figure 7 gives such example. In this case we deal with the model (Eq. (28)) which is mathematically poorly defined for the reconstruction procedure of any type. The most difficult problem here is that, strictly speaking, such function cannot be expanded into a convergent series because it is not square integrable. For this reason we cannot expect that any quasi-solution of the form of orthogonal expansion (11) will be close to exact one. Actually, in all our numerical experiments the reconstruction error for the PDF of the form (28) was in the range of 0.3-0.4. The main term in this result gave the systematic error $\varepsilon _c^2$. But we can see from Fig. 7 that even in this case the restored approximate PDF is somewhat similar to the initial one.

7. Conclusion

The method presented here can be used in the studies of low intensity optical fields. The results of computer simulation and experimental data processing show that the inversion procedure is easy in realization and it has the necessary stability with respect to numerical calculation errors and input data statistical fluctuations. Its advantage is also the possibility of checking the regularization level in estimation resulted. The restriction of matrix R(n,n') dimension by a finite number N₀ of input photoelectrons enables to eliminate data in the range of large values of registered counts. This permits to minimize the statistical errors and the influence of photodetector dead time effects which are of great importance for high rate photoelectron flow. It should be mentioned also that the estimations of the solution obtained are continuous functions of the argument x in spite of the fact that the initial data are discrete and the inversion is performed in finite-dimensional space.

Acknowledgment

Author wishes to thank Prof. Claude Aime (Univ. of Nice-Sophia Antipolis) for helpful discussion.

References

1. L. Mandel, in: Progress in Optics, Vol. 2, Ed. E. Wolf, North-Holland, Amsterdam 1963, p. 181

2. R.H.T. Bates, M.J. McDonnell, Image Restoration and Reconstruction, Clarendon Press, Oxford 1986

3. E. Wolf, C.L. Mehta, Phys. Rev. Lett. 13, 705 (1964)

4. B. Saleh, Photoelectron Statistics, Springer Verlag, Berlin 1978

5. J. Perina, Coherence of Light, Van Nostrand, New York 1972

6. G. Bedard, J. Opt. Soc. Am. A 57, 1201 (1967)

7. C.L. Byrne, B.M. Levine, J.C. Dainty, J. Opt. Soc. Am. A 1, 1132 (1984)

8. F. Sultani, C. Aime, H. Lantéri, Pure Appl. Opt. 4, 89 (1995)

9. F. Sultani, C. Aime, H. Lantéri, SPIE Proc. 2647, 2691 (1995)

10. A.N. Tikhonov, V.Y. Arsenin, Solutions of Ill-posed Problems, Halsted Press, Washington 1977

11. C.L. Mehta, in: Progress in Optics, Vol. 8, Ed. E. Wolf, North-Holland, Amsterdam 1970, p. 373

12. F. Riesz, B. Sz.-Nagy, Leçons d'analyse fonctionnelle, Akademiai Kiado, Budapest 1972

13. J.W. Goodman, in: Laser Speckle and Related Phenomena, Ed. J.C. Dainty, Springer Verlag, Berlin 1975, p. 8

14. Laser Beam Propagation in the Atmosphere, Ed. J.W. Strohben, Springer Verlag, Berlin 1978

15. O.I. Barchuk, V.N. Kurashov, A.I. Maschenko, A.T. Mirzaev, Ukr. J. Phys. 33, 1145 (1988)

This document was generated using the LaTeX2HTML translator Version 2K.1beta (1.47)
Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
Copyright © 1997, 1998, 1999, Ross Moore, Mathematics Department, Macquarie University, Sydney.