Fraunhofer diffraction effects on total power for a planckian source.

An algorithm for computing diffraction effects on total power in the case of Fraunhofer diffraction by a circular lens or aperture is derived. The result for Fraunhofer diffraction of monochromatic radiation is well known, and this work reports the result for radiation from a Planckian source. The result obtained is valid at all temperatures.

Key words: diffraction; Fraunhofer; Planckian; power; radiometry.

Accepted: August 28, 2001

Available online: http://www.nist.gov/jres

1. Introduction

Fraunhofer diffraction by a circular lens or aperture is a ubiquitous phenomenon in optics in general and radiometry in particular. Figure 1 illustrates two practical situations in which Fraunhofer diffraction occurs. In the first example, diffraction limits the ability of a lens or other focusing optic to focus light. According to geometrical optics, it is possible to focus rays incident on a lens to converge at a focal point. In practice, even with the focal point at the center of a finite-sized circular detector, some of the light incident on the lens fails to be collected, giving rise to a "diffraction loss."

The second example can occur when using a source, such as a blackbody cavity source viewed through a small "pinhole" aperture, to calibrate optical systems intended for observing celestial objects or for use in remote-sensing applications. In this case, the angle subtended by the cavity opening at the aperture is larger than the angle subtended by the detector pupil on the dark side of the aperture (1). Because the former angle is larger than the latter one, geometrical optics suggests that the total power detected depends only on the blackbody temperature and a geometrical factor related to the pinhole aperture, detector optics, and relative separation, because the detector pupil is overfilled. However, diffraction leads to losses in the total power reaching the detector.

All of the above diffraction losses have been a subject of considerable interest, and they have been considered by Blevin (2), Boivin (3), Steel, De, and Bell (4), and Shirley (5). The formula for the relative diffraction loss in spectral power is already well known, and it is noted below. However, a formula for the diffraction loss in the total power for the case of a Planckian source such as a star or blackbody appears to have been given only in the high-temperature limit. This article reports a formula for the diffraction loss in the power for such a source at all temperatures. The formula is useful for predicting the loss in cases such as the examples discussed, and it can be used as an independent check for numerically calculated diffraction losses.

For the radiation present at a given wavelength [lambda] the diffraction loss can be described in terms of a unit-less parameter, v = 2[pi][PSI]R/[lambda]. R is the radius of the lens or aperture. The angle [PSI] is defined in either of two ways: either 2[PSI] is the full angle subtended by the detection pupil at the focusing optic, or 2 [PSI] is the full angle subtended by the blackbody cavity opening at the pinhole aperture. For a source in thermal equilibrium at temperature T, the effects of temperature on the distribution of spectral power enter relevant equations through the ratio, [c.sub.2]/([lambda]T) Av. Here, [c.sub.2] = hc/k = 1.438 7752(25) X [10.sup.-2] m K is the second radiation constant, where the number in parenthesis is the one-standard-deviation uncertainty in the last two digits. This implicitly defines another unit-less parameter, A = [c.sub.2]/(2[pi][psi]RT). The relative effects of diffraction on the total power only depend on this one parameter, A.

As a final note, if the detection pupil subtends a considerable angle, the diffraction losses can be approximately given as the weighted average of losses arising, in the case of an infinitesimal detector, for values of [psi] varying between [psi] minus the half-angle subtended by the detector and [psi] plus the same half-angle. Analogous weighting can be used, in the context of the first example, to account for the finite angular diameter of a source (smaller than 2[psi]). Such averaging is discussed elsewhere [5] but is not discussed further in this work, which henceforth assumes a single effective value of [psi].

2. Derivation of Formula

In the Fraunhofer case, the ratio of the true spectral power to the spectral power predicted by geometrical optics is

F(v) = 1 - [J.sup.2.sub.0](v) - [J.sup.2.sub.1](v). (1)

Here, [J.sub.m](v) is a cylindrical Bessel function. It is necessary to incorporate this wavelength-dependent factor when evaluating the ratio of the total power detected to the total power predicted by geometrical optics. For a source whose spectral output obeys the Planck distribution, this unit-less ratio is given by

[F.sup.-](A) = [[integral].sup.[infinity].sub.0] [dvv.sup.3] [[exp(Av) - 1].sup.-1] [1 - [J.sup.2.sub.0](v) - [J.sup.2.sub.1](v)] / [F.sup.-](A) = [[integral].sup.[infinity].sub.0] [dvv.sup.3] [[exp(Av) - 1].sup.-1] (2)

The value of the denominator is familiar, being given by

[[integral].sup.[infinity].sub.0]dv[v.sup.3][[exp(Av) - 1].sup.-1] = [summation over (n=1/[infinity])][[integral].sup.[infinity].sub.0]dv[v.sup.3] exp(-nAv)

= 6[summation over (n=1/[infinity])][(nA).sup.-4] = 6[A.sup.-4][zeta](4). (4)

Here, [zeta](n) is the Riemann zeta function. When evaluating the numerator, two techniques have been found to be helpful: one for the "low-temperature" case (defined as A > 4) and one for the "high-temperature" case (defined as A < 4). By such definitions, "high" and "low" temperature cases arise depending on the average diffraction loss for a source at temperature T. This depends on T and the geometry through A. In either case, one first makes use of the relation,

[J.sub.0][2v sin([theta]/2)] = [J.sup.2.sub.0](v)+2[summation over (m=1/[infinity])][J.sup.2.sub.m](v)cos(m[theta]), (4)

from which follows,

[J.sup.2.sub.0](v) + [J.sup.2.sub.1](v) = 1/2[pi][[integral].sup.2[pi].sub.0]d[theta](1+cos[theta])[J.sub.0][2v sin([theta]/2)]. (5)

Low-Temperature Case

To evaluate the numerator in the low-temperature case, series expansion of the latter Bessel function in Eq. (5) yields

[J.sup.2.sub.0]+[J.sup.2.sub.1](v)=1/[pi][summation over (s=0/[infinity])][(-1).sup.s][V.sup.2s]/[(s!).sup.2][[integral].sup.2 [pi].sub0]d[theta][1 - [sin.sup.2]([theta]/2)][sin.sup.2s]([theta]/2)

=2[summation over (s=0/[infinity])][(-1).sup.s][v.sup.2s]/[(s!).sup.2] [(2s -1)!!/(2s)!! - (2s+1)!!/(2s+2)!!]

=2[summation over (s=0/[infinity])][(-1).sup.s][v.sup.2s]/[(s!).sup.2][(2s - 1)!!/(2s + 2)!!] (6)

This permits one to rewrite Eq. (2) as follows:

Q(A) [equivalent to] 6[A.sup.-4][zeta](4)[1 - F(A)]

=[summation over (n=1/[infinity])][[integral].sup.[infinity].sub.0]dv[v.sup.3]exp(-nAv )[[J.sup.2.sub.0](v)+[J.sup.2.sub.1](v)]

=2[summation over (s=1/[infinity])][(-1).sup.s](2s+3)/[(s!).sup.2][A.sup.2 s+4][(2s-1)!!/(2s+2)!!]

=-[summation over (s=0/[infinity])][(2[pi]).sup.2s+4](2s+3)!(2s-1)!![B.sub.2s+4]/[(s!). sup.2](2s+4)!(2s+2)!![A.sup.2s+4)], (7)

In this expression, [B.sub.i] is a Bernoulli number. The sum converges for all A > 2, and one obtains an error in F(A) estimated to be as small as [10.sup.-14] for A > 4 if one sums up to s = 28. (6)

High-Temperature Case

To evaluate the numerator in the high-temperature case, the relation in Eq. (2) is rewritten to render

Q(A)=6[A.sup.-4][zeta](4)[1 - F(A)]

=[summation over (n=1/[infinity])][[integral].sup.[infinity].sub.0]dv[v.sup.3]exp(-nAv )[[J.sup.2.sub.0](v)+[J.sup.2.sub.1](v)]

=[summation over (n=1/[infinity])][Q.sub.1](nA), (8)

with

[Q.sub.1](A)=[[integral].sup.[infinity].sub.0]dv[v.sup.3]exp(-Av)[[J. sup.2.sub.0](v)+[J.sup.2.sub.1](v)]. (9)

Use of Eq. (5) and the relation,

[[integral].sup.[infinity].sub.0]dvexp(-[alpha]v)[J.sub.0]([beta]v)=1 /[square root of ([[alpha].sup.2]+[[beta].sup.2])] (10)

one has

[Q.sub.1](A)=[(-d/dA).sup.3]1/2[pi][[integral].sup.2[pi].sub.0]d[thet a]1+cos[theta]/[square root of ([A.sup.2]+4[sin.sup.2]([theta]/2))]. (11)

Making the abbreviation, w=4/([A.sup.2]+4), one has

[Q.sub.1](A)=[(-d/dA).sup.3][4/[pi][square root of ([A.sup.2]+4)][[integral].sup.[pi]/2.sub.0]d[phi][cos.sup.2][phi]/[sq uare root of (1 - w[cos.sup.2][phi])]

=[(-d/dA).sup.3]{4/[pi][square root of ([A.sup.2]+4)][K(w)-E(w)/w]}

=[(-d/dA).sup.3]{[square root of ([A.sup.2]+4)]/[pi][K(w)-E(w)]}, (12)

where E(w) and K(w) are respectively complete elliptic integrals of the first and second kind, defined according to the convention, (7)

E(w)=[[integral].sup.[pi]/2.sub.0]d[phi][square root of (1-w[sin.sup.2][phi])],

K(w)=[[integral].sup.[pi]/2.sub.0]d[phi]/[square root of (1-w[sin.sup.2][phi]. (13)

(In a different convention, w is replaced by [w.sup.2] in the integrand but nowhere else.) From the properties of elliptic integrals [8], one may obtain

[Q.sub.1](A)=4/[pi][A.sup.3][4(2+[A.sup.2])E(w)-[A.sup.2]K(w)/[([A.su p.2]+4).sup.3/2]], (14)

and

[Q.sub.1](A) = -d/dA[4E(w)/[pi][A.sup.2][square root of ([A.sup.2]+4)]. (15)

The latter result is easily obtained, because [Q.sub.1](A) has already been written as the third derivative of an expression with respect to A.

The closed-from result in Eq. (14) may be directly substituted in Eq. (8), after which one may sum over n. Because the summand can be integrated according to Eq. (15), however, the Euler-Maclaurin formula may be used to obtain a more easily evaluated, asymptotic expression for Q(A):

Q(A) = 4[zeta](3)/[pi][A.sup.3] + [summation over (N-1/n=1)][Q.sub.1](nA) + 1/A [[integral].sup.[infinity].sub.NA] dA'[Q.sub.1](A') + 1/2[Q.sub.1](NA)

-[[[summation over (m-1/s=1)] [B.sub.2s]/(2s)![A.sup.2s-1][(d/dA').sup.2s-1][Q.sub.1](A')]\.sub.A'= NA] + [R.sub.m](N,A), (16)

where [B.sub.i] is again a Bernoulli number. Here a new function,

[Q.sub.1](A) = [Q.sub.1](A) - 4/[pi][A.sup.3], (17)

has been introduced as the summand, so that the largest part of [Q.sub.1](A) is extracted from the sum prior to summation.

The parameters m and N are both positive integers. [R.sub.m](N,A), the remainder (i.e., error) term, usually decreases initially with increasing m, but then rapidly diverges. It is therefore good to choose a suitable combination of N and m, for a given value of A, to have an efficient yet accurate result. Increasing N should usually improve results, and one obtains an error in F(A) estimated to be as small as [10.sup.-14] for A<4 if one uses N = 60 and m = 5.

Explicitly, Q(A) is given by

Q(A) = 4[zeta](3)/[pi][A.sup.3] + [summation over (N-1/n=1)][Q.sub.1](nA)

+2/[pi][[alpha].sup.3]{[8(N+1)+2[[alpha].sup.2](N+2)]E([omega])-[[alp ha].sup.2]K([omega])/[(4+[[alpha].sup.2]).sup.3/2]-N-1}

+1/[pi][[alpha].sup.3](-1/N+1/3[N.sup.3]-1/3[N.sup.5]+3/5[N.sup.7]- ...)

+1/[pi][[alpha].sup.3] [summation over (m-1/s=1)] [B.sub.2s]/(2s)E([omega])-[[alpha].sup.2][k.sub.s ]([alpha])K([omega])][f.sub.s]/[N.sup.2s-1][(4+[[alpha].sup.2]).sup.( 4s+1)/2]

+ [R.sub.m](N,A). (18)

One should include only the first m - 1 terms in the parenthetic expression preceding the second sum. Two abbreviations, [alpha] = NA and [omega] = 4/(4 + [[alpha].sup.2]), have been introduced. In the second sum, the first several constants are given by

[B.sub.2]/2!=1/12, [B.sub.4]/4!=-1/720, [B.sub.6]/6!=1/30240, [B.sub.8]/8!=-1/1209600,

[f.sub.1]=4, [f.sub.2]=24, [f.sub.3]=48, [f.sub.4]=2880. (19)

Likewise, the first several e- and k-functions are given by

[e.sub.1]([alpha]) = 96 + 68[[alpha].sup.2] + 19[[alpha].sup.4],

[k.sub.1]([alpha]) = 12 + 7[[alpha].sup.2],

[e.sub.2]([alpha]) = 5120 + 6112[[alpha].sup.2] + 2880[[alpha].sup.4] + 630[[alpha].sup.6] + 117[[alpha].sup.8],

[k.sub.2]([alpha]) = 640 + 624[[alpha].sup.2] + 216[[alpha].sup.4] + 57[[alpha].sup.6],

[e.sub.3]([alpha]) = 1720320 + 2908160[[alpha].sup.2] + 2093376[[alpha].sup.4] + 829552[[alpha].sup.6] + 196764[[alpha].sup.8] + 19917[[alpha].sup.10] + 3726[[alpha].sup.12],

[k.sub.3]([alpha]) = 215040 + 316480[[alpha].sup.2] + 193072[[alpha].sup.4] + 62892[[alpha].sup.6] + 8265[[alpha].sup.8] + 2046[[alpha].sup.10],

[e.sub.4]([alpha]) = 33030144 + 72310784[[alpha].sup.2] + 70087936[[alpha].sup.4] + 39429888[[alpha].sup.6] + 14163680[[alpha].sup.8] 3350128[[alpha].sup.10] + 572853[[alpha].sup.12] + 21062[[alpha].sup.14] + 6121[[alpha].sup.16],

[k.sub.4]([alpha]) = 4128768 + 8135680[[alpha].sup.2] + 6993408[[alpha].sup.4] + 3421440[[alpha].sup.6] + 1037920[[alpha].sup.8] + 219516[[alpha].sup.10] + 8820[[alpha].sup.12] + 3601[[alpha].sup.14]. (20)

3. Evaluation of Formula

In the high-temperature limit, this gives a result of the form

F(A) = 1 - [4A[zeta](3)/[pi] + O([A.sup.3] ln A]/6[zeta](4), (21)

and evaluating the successive terms in Eq. (21) may be both difficult and unfruitful, because the series appears to be both very complicated and slow to converge. Note that the leading terms in Eq. (21) are consistent with other works cited, such as the work of Blevin. In the low-temperature limit, one has

F(A) = 10[[pi].sup.2]/21[A.sup.2] - [[pi].sup.4]/4[A.sup.4]+ .... (22)

Figure 2 shows F(A) for a range of values of A bridging the low- and high-temperature regions. The solid curve indicates "exact" results obtained using Eq. (7) for A > 4 and Eq. (18) for A < 4. The dashed curves indicate approximate results for low- and high-temperature limits, obtained using the respective formulas in Eqs. (21-22). Likewise, Table 1 shows sample values of F(A) over a similar range. For 0 < A < 0.2 and A > 8, inclusion of terms shown in Eqs. (21-22) yields an error in F(A) smaller than 0.001. While this immediate discussion helps provide a sense of the behavior of F(A), accurate values should be found using Eq. (7) and Eq. (18).

Acknowledgment

The author is grateful to A. W. Smith and M. A. Edwards for critical readings of this manuscript.

About the author: Eric L. Shirley is a physicist in the Optical Technology Division of the NIST Physics Laboratory. The National Institute of Standards and Technology is an agency of the Technology Administration, U.S. Department of Commerce.

4. References

(1.)In practice, such optics may be those being calibrated or may be a reference detector used to calibrate the blackbody.

(2.) W. R. Blevin, Diffraction Losses in Radiometry and Photometry, Metrologia 6, 39-44 (1970).

(3.) L. P. Boivin, Diffraction Losses Associated with Tungsten Lamps in Absolute Radiometry, Appl. Opt. 14, 197-200 (1975); Diffraction corrections in radiometry: comparison of two different methods of calculation, ibid. 14, 2002 (1975).

(4.) W. H. Steel, M. De, and J. A. Bell, Diffraction corrections in radiometry, J. Opt. Soc. Am. 62, 1099-1103 (1972).

(5.) E. L. Shirley, Revised formulas for diffraction effects with point and extended sources, Appl. Opt. 37, 6581-6590 (1998).

(6.) This error, as well as the error in the high-temperature case, is estimated by comparing results of the two techniques at and around A = 4. It is surmised that this error arises because of computer round-off. Because the current relative standard uncertainty of [c.sub.2] is 1.7 X [10.sup.-6], the present error is acceptable.

(7.) M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables, John Wiley and Sons, New York (1972) p. 590.

(8.) See, for instance, P. F. Byrd and M. D. Friedman, Handbook of Elliptic Integrals for Engineers and Physicists, Springer-Verlag, Berlin (1954).

[Figure 2 omitted]

Table 1

F(A) at sample values of A

 A    F(A)

0.2  0.9526
0.4  0.9041
0.6  0.8543
0.8  0.8033
1.0  0.7514
1.2  0.6994
1.4  0.6481
1.6  0.5984
1.8  0.5508
2.0  0.5060
2.2  0.4643
2.4  0.4259
2.6  0.3907
2.8  0.3587
3.0  0.3297
3.2  0.3034
3.4  0.2797
3.6  0.2582
3.8  0.2389
4.0  0.2214
4.2  0.2055
4.4  0.1912
4.6  0.1781
4.8  0.1663
5.0  0.1555
5.2  0.1457
5.4  0.1367
5.6  0.1285
5.8  0.1210
6.0  0.1141
6.2  0.1077
6.4  0.1018
6.6  0.0964
6.8  0.0914
7.0  0.0867
7.2  0.0824
7.4  0.0784
7.6  0.0746
7.8  0.0712
8.0  0.0679