Contents

Page updated: May 19, 2021
Author: Curtis Mobley
View PDF

Beam and Point Spread Functions

This page deﬁnes two equivalent quantities that describe light propagation in absorbing and scattering media and that are fundamental to image analysis.

The Beam Spread Function

Consider a collimated source emitting spectral power $P$ (units of $W n m^{- 1}$ ) in direction $𝜃 = 0$ as shown in Fig. 1. As the beam passes through the medium, scattering will spread out the beam as illustrated by the green arrows in the ﬁgure, and absorption will reduce the beam power. The combined eﬀects of scattering and absorption give some spectral plane irradiance $E (r, 𝜃)$ (units of $W m^{- 2} n m^{- 1}$ ) on the surface of a sphere of radius $r$ at an oﬀ-axis angle $𝜃$ relative to the direction of the emitted beam. The irradiance sensor used to measure $E (r, 𝜃)$ has a cosine response for angles relative to the normal to the detector surface. It is assumed here that the emitted light is unpolarized and that the medium is isotropic, so that there is no azimuthal dependence of the detected irradiance. The green arrows in the ﬁgure represent one possible path for light between the source and the detector.

Figure 1: Geometry for deﬁnition of the Beam Spread Function.

The beam spread function (BSF) is then deﬁned as the detected irradiance normalized by the emitted power:

B S F (r, 𝜃) \equiv \frac{E (r, 𝜃)}{P} [\frac{W m^{- 2} n m^{- 1}}{W n m^{- 1}} = m^{- 2}] .

(1)

Recall that the Volume Scattering Function (VSF) describes a single scattering event. Two such scatterings are shown for the green arrows in Fig. 1. The BSF on the other hand describes the cumulative eﬀects on the emitted beam of all of the scattering and absorption events between the source and the detector. The BSF thus depends on the distance between the source and the detector, whereas the VSF depends only on the optical properties of the medium and is independent of the locations of the source and detector.

The Point Spread Function

Now suppose that there is a source at the location of the detector in Fig. 1 and that this source is emitting spectral intensity $I (γ)$ (units of $W s r^{- 1} n m^{- 1}$ ) with an angular pattern given by

I (γ) \equiv \{\begin{matrix} \frac{P}{π} cos (γ) 0 \leq γ \leq π ∕ 2 [W s r^{- 1} n m^{- 1}] \\ 0 π ∕ 2 < γ \leq π . \end{matrix}

For this intensity pattern the total emitted power is

\int_{4 π} I (γ) d Ω = 2 π \int_{0}^{π ∕ 2} \frac{P}{π} cos (γ) sin (γ) d γ = P .

The short green arrows in Fig. 2 illustrate the emission of light by this cosine source. As illustrated by the long green arrows in this ﬁgure, the emitted intensity will give rise to a radiance $L (r, 𝜃)$ , where $r$ is the distance from the source and $𝜃$ is direction measured from the $γ = 0$ axis of the emitted intensity. The light emitted by this source is then detected by a well collimated radiance sensor that can scan past the source as shown in Fig. 2.

Figure 2: Geometry for deﬁnition of the Point Spread Function.

The Point Spread Function (PSF) is then deﬁned as the detected radiance normalized by the maximum of the emitted intensity:

P S F (r, 𝜃) \equiv \frac{L (r, 𝜃)}{P ∕ π} [\frac{W m^{- 2} s r^{- 1} n m^{- 1}}{W s r^{- 1} n m^{- 1}} = m^{- 2}] .

The PSF can be visualized as the “glow” of light around distant street light seen through a foggy atmosphere as seen in Fig. 3. Although a street light is not a cosine-emitting point source, the angular pattern of the glow seen around the light gives a qualitative feeling for the PSF, which is just the pattern of the radiance distribution centered on the source.

Figure 3: The glow of light around a distance street light seen through fog is essentially the PSF for the given distance and atmospheric properties.

Figure 4 shows PSFs computed by Monte Carlo simulations for a homogeneous water body with a Petzold average particle phase function and a single-scattering albedo of $ω_{o} = 0.8$ . The left panel shows the PSF for the ﬁrst 15 deg of $𝜃$ and for nondimensional optical distances between the source and detector of $τ$ = 0.1, 1, 5, 10, and 30. These curves show that the magnitude of the PSF decreases as $τ$ increases because of absorption. The colored dots are the centers of the angular bins used to tally the light rays in the Monte Carlo simulations.

The right panel of the ﬁgure shows the PSF values normalized to 1 at $𝜃 = 0$ . These curves show that the shape of the PSF starts out very highly peaked near $𝜃 = 0$ for small $τ$ and eventually becomes relatively ﬂat in $𝜃$ as the optical distance increases. This is progression of PSF shapes can be understood as follows.

As $τ$ goes to zero, there is almost no chance for light to scatter. As the radiometer seen in Fig. 2 scans past the point source, it sees either almost nothing, or it sees the point source at $𝜃 = 0$ . The PSF then approaches a Dirac delta function in $𝜃$ as $τ \to 0$ . As the optical distance increases, there is more and more chance for scattering until many photons have been scattered once. The PSF then begins to look similar to the scattering phase function, which describes the redistribution of radiance by single scattering. The red curve lying close to the curve for $τ = 5$ is the normalized Petzold phase function used in these simulations. (The PSF would not in general ever have exactly the same shape as the phase function. Even if, at some distance, most light rays have been scattered once, others will not yet be scattered, and others will have been scattered more than once. Thus a PSF never describes just single scattering.)

Finally, as $τ$ becomes very large, all rays have been scattered many times, and the resulting radiance distribution, hence the PSF, approaches the shape of the asymptotic radiance distribution $L_{\infty}$ . The black curve near the $τ = 30$ PSF shows $L_{\infty}$ as computed by HydroLight for the IOPs of these simulations. The $τ = 30$ curve is still noticeably diﬀerent from $L_{\infty}$ . It is computationally expensive to trace enough rays to large optical distances to reproduce the asymptotic distribution. A computation for $τ = 50$ required emitting $1 0^{9}$ initial rays, of which about 0.05% reached $τ = 50$ for $𝜃 < 15$ deg, but the resulting PSF (not shown) is much closer to $L_{\infty}$ . Similarly, statistical noise can become large for large $𝜃$ angles if the angular bins are small. However, it is the small $𝜃$ values that are most important for image analysis (after all, you are usually looking generally toward an object, not away from it), so noise at large angles is seldom a problem in practical applications.

Figure 4: Example point spread functions. Left panel: PSF as a function of optical distance for

𝜃 = 0

to 15 deg. Right panel: the same curves normalized to 1 at

𝜃 = 0

. In the right panel, the red line near the

τ = 5

curve is the normalized Petzold phase function used in the Monte Carlo simulations. The black line at the top is the asymptotic radiance distribution.

Equivalence of the BSF and PSF

In Maxwell’s Equations for the propagation of electromagnetic energy, if you change time to minus time, nothing changes except that the direction of propagation is reversed. This means that if light is propagating from the source to the receiver along the path shown in Fig. 1, then light emitted from the receiver location in the same direction as the detected beam will propagate back along the same path to the source, as shown by the green arrows in Fig. 2. This is known as the principle of electromagnetic reciprocity, or Helmholtz reciprocity, or the “if I can see you, you can see me” theorem. Figures 1 and 2 were drawn to highlight the symmetry between the BSF and the PSF.

The BSF and PSF are superﬁcially diﬀerent: a collimated source and a cosine detector vs a cosine source and a collimated detector. However, reciprocity suggests that the BSF and PSF contain equivalent information. Indeed, the BSF and PSF are numerically the same. Figure 5 shows PSF and BSF measurements from two very similar water bodies obtained by Mertens and Replogle (1977). The instruments were mounted on a 20 m long underwater frame, which allowed a maximum range $r$ between source and detector of approximately 19 m. Although the measurements in the ﬁgure were made on diﬀerent days and the water IOPs were slightly diﬀerent, the closeness of the PSF and BSF over several orders of magnitude suggests that the BSF and PSF are numerically equal, as they stated without proof. Gordon (1994b) started with the radiative transfer equation and used reciprocity to show that the BSF and PSF are indeed numerically equal. It is thus customary refer to just the PSF, even if the geometry of a problem corresponds to that of the BSF. (This is the case in the derivation of the Lidar equation, for example.)

Figure 5: The BSF and PSF at 488 nm measured on two diﬀerent days in very clear Bahamas water. The range was approximately 19 m. The BSF was measured when the water beam attenuation length was

1 ∕ c (488) = 6.8 m

. The PSF was measured on the following day when

1 ∕ c (488) = 7.0 m

. Curves extracted from Figs. 4 and 6 of Mertens and Replogle (1977).

Importance of the PSF in Image Analysis

The PSF plays a fundamental role in the prediction of how an object appears when seen through an absorbing and scattering medium such as water. Suppose you are looking in a particular direction at point $(x_{o}, y_{o})$ in Fig. 6 from a distance $z$ away from the shark. Most of the light detected in this viewing direction probably comes from point $(x_{o}, y_{o})$ in the image. However, every other point in the scene contributes at least a small amount of light to the detected signal in the viewing direction $(x_{o}, y_{o})$ because of scattering from those other directions into the viewing direction. How much light those other points contribute to a given direction is determined by the value of their PSFs at the angular distances away from the viewing direction.

Figure 6: Illustration of how the PSF describes the contributions of all points

(x, y)

in an image to the signal at a particular point

(x_{o}, y_{o})

. The white circles surrounding points

(x_{1}, y_{1})

and

(x_{2}, y_{2})

represent the PSFs of those points.

Let $S (x, y, 0)$ represent the “bright-dark” pattern of the shark image when seen at distance $z = 0$ , or as seen through a vacuum. Thus for a gray-scale digital image and 8 bit resolution, $S = 0$ for a black pixel and $S = 255$ for white pixel. When viewing a particular point $(x_{o}, y_{o})$ in the image through a vacuum at some distance $z$ from the shark, the detected radiance comes only from point $(x_{o}, y_{o})$ and has the value $S (x_{o}, y_{o}, 0)$ . However, when viewing point $(x_{o}, y_{o})$ in the through a scattering medium, every point in the scene contributes to the radiance seen at $(x_{o}, y_{o})$ . The concentric white circles in the ﬁgure represent contours of $P S F (z, 𝜃)$ . Point $(x_{1}, y_{1})$ thus contributes a radiance of $S (x_{1}, y_{1}, 0) P S F (z, 𝜃_{1})$ , where angle $𝜃_{1} = {tan}^{- 1} (r_{1} ∕ z)$ . Likewise point $(x_{2}, y_{2})$ contributes radiance $S (x_{2}, y_{2}, 0) P S F (z, 𝜃_{2})$ , where angle $𝜃_{2} = {tan}^{- 1} (r_{2} ∕ z)$ . The total radiance at point $(x_{o}, y_{o})$ is given by summing up the contributions from all points in the image:

S (x_{o}, y_{o}, z) = \int_{a l l x} \int_{a l l y} S (x, y, 0) P S F (x_{o} - x, y_{o} - y, z) d x d y .

(2)

The connection between points in the image at distance $z$ away and the distance and angular variables used for PSFs is

r = \sqrt{{(x_{o} - x)}^{2} + {(y_{o} - y)}^{2}}, 𝜃 = {tan}^{- 1} (\frac{r}{z})

so that

P S F (x_{o} - x, y_{o} - y, z) = P S F (z, 𝜃) .

Thus, given the inherent appearance of an object (the appearance at distance $z = 0$ or through a vacuum) and the PSF of the medium, the appearance of the object as seen through the medium can be computed. In other words, the PSF completely characterizes the eﬀect of the environment on an image.

Note that in a vacuum and with a perfect optical system, the PSF reduces to a 2D Dirac delta function: $P S F (x_{o} - x, y_{o} - y) = δ (x_{o} - x) δ (y_{o} - y)$ , in which case

S (x_{o}, y_{o}, z) = \int_{a l l x} \int_{a l l y} S (x, y, 0) δ (x_{o} - x) δ (y_{o} - y) d x d y = S (x_{o}, y_{o}, 0) .

Thus the image at distance $z$ is the same as at distance 0.

Integrals of the form

h (x_{o}) = \int f (x) g (x_{o} - x) d x \equiv f \otimes g

(3)

are called convolution integrals, and $h$ is called the convolution of $f$ and $g$ . Equation (2) is a two-dimensional convolution showing that the image seen from a distance $z$ is the convolution of the image seen from zero distance and the PSF of the medium at distance $z$ . This is an extremely important result and is the foundation of image analysis through absorbing and scattering media. As will be seen in the Photometry and Visibility Chapter (under development), the Fourier transform of Eq. (2) has very useful properties. Rather than compute the double integral of Eq. (2) for every point $(x_{o}, y_{o})$ of the scene (which is computationally ineﬃcient), the Fourier transforms of $S (x, y, 0)$ and $P S F (x, y, z)$ can be computed (computationally eﬃcient using Fast Fourier Transforms), multiplied, and then the inverse transform of the product gives the ﬁnal image. This topic will be developed in the Level 2 pages of the Photometry and Visibility chapter. Meanwhile see Mobley (2020), Goodman (1996), or the Mertens and Replogle paper cited above. Hou et al. (2007) used this equation as the starting point for a rigorous analysis of the visibility of a Secchi disk.

Because of the importance of the PSF in image analysis, considerable eﬀort has been expended to develop models of the PSF as a function of water IOPs. A number of such models are reviewed in Hou et al. (2008). The page on The Lidar Equation shows the fundamental role of the BSF in lidar remote sensing.

Comment on terminology

he term “convolution” is frequently misused. An integral of the form $\int f (x) g (x) d x$ is not a convolution. It can be called the inner product of $f$ and $g$ , or the integral of $f$ weighted by $g$ (or of $g$ weighted by $f$ ), but it returns a number. The convolution of $f$ and $g$ returns a functional $h$ as shown in Eq. (3). At best, you could call $\int f (x) g (x) d x$ the value of the convolution at $x_{o} = 0$ , but that is missing the central idea that a convolution shows how two functions “overlap” so that, in the case of an image, each point in an image contributes to all other points.]