This tutorial gives an overview of the use of the Wigner function as a tool for modeling optical field propagation. Particular emphasis is placed on the spatial propagation of stationary fields, as well as on the propagation of pulses through dispersive media. In the first case, the Wigner function gives a representation of the field that is similar to a radiance or weight distribution for all the rays in the system, since its arguments are both position and direction. In cases in which the field is paraxial and where the system is described by a simple linear relation in the ray regime, the Wigner function is constant under propagation along rays. An equivalent property holds for optical pulse propagation in dispersive media under analogous assumptions. Several properties and applications of the Wigner function in these contexts are discussed, as is its connection with other common phase-space distributions like the ambiguity function, the spectrogram, and the Husimi, P, Q, and Kirkwood–Rihaczek functions. Also discussed are modifications to the definition of the Wigner function that allow extending the property of conservation along paths to a wider range of problems, including nonparaxial field propagation and pulse propagation within general transparent dispersive media.
©2011 Optical Society of America
One of the most used mathematical tools in physics and engineering is the Fourier transformation, which converts a function of a given variable into a function of another, Fourier conjugate variable. Particularly in situations where the mathematical description of the physical process is independent of the choice of origin for the initial variable, the Fourier conjugate variable has clear physical significance. Consider, for example, a function of time representing the local air pressure for a sound or the electric field for optical radiation. Fourier transformation allows us to express this function in terms not of time but of frequency, which our senses perceive as pitch or color. When considering the spatial distribution of a monochromatic beam, on the other hand, Fourier transformation replaces the spatial variable with one associated with direction of propagation. In quantum mechanics, Fourier transformation over the position coordinates of a particle’s wave function leads to a distribution in momentum, while Fourier transformation over time gives an energy distribution.
In any of these specific contexts, the well-known mathematical properties of the Fourier transformation acquire a physical interpretation. In particular, a function and its Fourier transform cannot be simultaneously arbitrarily narrow, since the product of the widths of their squared moduli (defined as the standard deviations of their variables) must be equal to or greater than a given constant. This purely mathematical property, sometimes referred to as the space–bandwidth product theorem, has profound physical implications in the different cases mentioned earlier. Within the context of sound, it tells us that a musical note cannot be arbitrarily short and yet have a well-defined pitch. When applied to monochromatic wave beams over a transverse coordinate, it implies that a beam cannot be simultaneously arbitrarily thin and well collimated. In quantum mechanics (according to the Copenhagen interpretation), this relation gives rise to the famous Heisenberg uncertainty principle that limits the simultaneous knowledge about the position and momentum of a particle. Due to this latter interpretation, many authors refer to the general space–bandwidth relation of Fourier transform pairs as the uncertainty relation.
In physical applications, the Fourier transformation provides a change in the representation of a system, from a description in terms of one physically meaningful variable (e.g., time or position) into one in terms of another (e.g., frequency, direction, or momentum). Nevertheless, in many of these contexts, our intuition is more inclined towards a representation in terms of both Fourier conjugate variables. For example, we think of a musical piece as a sound distribution in terms of both time and frequency, represented as a musical score. Similarly, when designing and modeling macroscopic imaging or illumination optical systems, one thinks of light as a distribution over a bundle of rays, each of these specified by both its position and direction. The same is true in mechanics, where the classical model, closer than the quantum one to our everyday perception of the world, requires the specification of both the position and momentum of a particle. Notice that these three simplistic models (musical notation, ray optics, and classical mechanics) disregard the limitations imposed by the space–bandwidth product theorem or uncertainty relation.
There is a family of mathematical transformations that allow the description of a function in terms of both a variable and its Fourier conjugate. When applied to the physical situations mentioned earlier, these transformations mimic the so-called phase-space distributions used in the simplified physical models (e.g., a musical score, a radiance distribution giving the weight of all rays in an optical system, or a distribution of classical particles in terms of both their position and velocities). In this context, the term phase space refers to the extended space composed of both sets of variables (time–frequency, position–direction, position–momentum). For this reason, these mathematical transformations are sometimes referred to as phase-space quasi-probability distributions or mock-phase-space distributions , although most authors drop the cautionary qualifiers and refer to them simply as phase-space distributions, as we will do here.
One of the most widely used distributions of this type is the Wigner function, proposed by Wigner  to represent a quantum state as a distribution in terms of both position and momentum, reminiscent of a statistical–mechanical distribution of classical particles. The same mathematical transformation was independently proposed in the signal analysis community by Ville  for studying time signals, and by both Dolin  and Walther  to represent optical wave fields as ray distributions analogous to what is known in the theory of radiometry as the radiance or specific intensity  and in image science as the light field . Excellent reviews on the Wigner function and other phase-space distributions are given in the books by Cohen  and Mecklenbräuker and Hlawatsch , and the review articles by Lee , Agarwal and Wolf [11–13], Balazs and Jennings , and Hillery et al. , to name a few. For applications in optics, the reader can also consult, for example, the book edited by Testorf et al. , the book by Torre , and the review articles by Dragoman  and Sheridan et al. .
This tutorial gives a review of the Wigner function and its use within several contexts in classical optics. While many other phase-space representations have also been proposed (as discussed in Sections 3 and 6), the Wigner function presents two properties that make it unique in many applications: (i) its so-called marginal projections yield observable physical quantities, and (ii) under suitable situations, it satisfies the same conservation rules under spatial propagation or temporal evolution as the simpler physical model it mimics. For example, in the case of paraxial fields propagating through simple optical systems, the Wigner function plays the role of a distribution of weights for the rays, because its value for each ray remains constant under propagation across the system, and because its integral over all rays passing through a point leads to the local optical intensity. The main emphasis of this work is precisely on the property of conservation along paths under propagation or evolution, its usefulness in practical applications, and range of validity. Also presented here are generalizations to the standard definition of the Wigner function that allow extending this property to more general situations, including nonparaxial propagation and pulse propagation in arbitrary transparent dispersive media.
The outline of this article is as follows. Some unitary transformations (including the Fourier transformation) that are useful in the description of certain simple optical systems are described in Section 2. A very simple phase-space representation, the spectrogram, is introduced in Section 3. The definition and some properties of the Wigner function are presented in Section 4 for the case of deterministic functions, and the more general definition of this representation for the case of correlations is given in Section 5. An overview of several phase-space distributions is given in Section 6. The Wigner function’s marginal and conservation properties are discussed in Section 7 for the case of classical paraxial optical fields propagating through simple systems and in Section 9 for the case of optical pulses propagating in media presenting linear dispersion. In both these contexts (as well as in quantum optics), a method can be used to recover the coherence properties of the field based on the measurement of the marginal projections of the Wigner function at different stages. This scheme, known as phase-space tomography, is discussed briefly in Section 10. The conservation property of the Wigner function under classical propagation, however, is valid only under certain limiting cases, as discussed in Section 8. In Section 11 we discuss a procedure for defining generalizations to the Wigner function that satisfy the desired propagation and marginal properties for other systems, including nonparaxial propagation through homogeneous media and pulse propagation through more general dispersive media. Finally, some brief concluding remarks are given in Section 12.
2. The Fourier Transform and Other Unitary Linear Transformations
Before introducing the Wigner function and other phase-space distributions, let us first review some simple integral transformations that appear naturally in the approximate description of certain optical systems. All these transformations depend linearly on the function they are applied to.
2.1. The Fourier Transform
Let us start by introducing the generic definition of Fourier transformation that is used in this article. Consider a scalar, square-integrable function , where is an N-dimensional variable. The Fourier transform of this function is defined here as
- (a) time t and frequency ω, for which , , and ;
- (b) the transverse position and transverse direction of a beam at a fixed frequency ω, where , (the free-space wavenumber, given by the ratio of the frequency and the speed of light c), (the components of the position vector perpendicular to the optical axis), and p is the transverse part of a vector indicating a direction of propagation;
- (c) two quantized field quadratures, i.e., the “real” and “imaginary” parts of the complex representation of a Cartesian component of the quantized electric field, with and q and p being the two quadratures.
There are several well-known properties of the Fourier transformation (see, e.g., Ref.  for a comprehensive list). Let us mention here only a couple:
(i) Parseval–Plancherel theorem. This theorem states that the inner product of two functions f and g is the same in either representation, due to the unitarity of the Fourier transformation, i.e.,4). More general forms of the multidimensional uncertainty relation exist that account for misalignment of the axes of symmetry of the Gaussian with respect to the coordinate axes, as well as for chirping (complex Gaussian widths), but these are beyond the scope of this article. The reader can consult, e.g., Ref.  for more on the uncertainty relation.
2.2. The Fractional Fourier Transform
A generalization of the Fourier transformation that gives access to a continuum of functions joining the original function and its Fourier transform was proposed by Condon [21–24]. This transformation, referred to as the fractional Fourier transformation of degree θ, is defined here as
Notice that, for and , the fractional Fourier transform reduces, respectively, to scaled versions of the original function and its Fourier transform, i.e.,25] under the limit .) The fractional Fourier transform presents an additivity property, which takes a particularly simple form in cases where σ can be set to unity: two consecutive transformations with degrees and , respectively, are equivalent to a single transformation with degree . As will be discussed later in this article, the spatial field propagation through a quadratic gradient index waveguide or the temporal evolution of the quadratures of a quantized field are described mathematically by fractional Fourier transformation. In these contexts, the additivity property has a straightforward physical interpretation, since θ is proportional to the propagation distance or time. Due to both Eqs. (9a) and (9b) and to this additivity property, many authors use, instead of the degree θ, the order , so that corresponds to the standard Fourier transformation . This way, γ can be thought of as the power to which the Fourier transformation operator is raised before acting on a function.
It must be noted that a more general version (sometimes called anamorphic) of the fractional Fourier transformation can be written in the multidimensional case, where a different θ is used for each Cartesian component of the transformation variables . This more general definition is in turn a special case of an even more general family of transformations, described in what follows.
2.3. Linear Canonical Transformations
A more general class of unitary transformations that includes the regular and fractional Fourier transformtions is now described. These transformations, known as linear canonical transformations, are used in paraxial optics and nonrelativistic quantum mechanics to describe simple systems [27–30]. Consider a linear transformation whose kernel has constant amplitude, like those of the standard and fractional Fourier transforms:10), we will set this constant to zero. Further, it can be shown that the analog of Eq. (3) is given by
Let us now assume that v satisfies the symmetry property , so that (where 0 denotes a zero vector). It turns out to be convenient to rewrite Eq. (11) as
Note that in Eq. (17), the transformation depends on , which is a matrix in which , , and are embedded:16), it is easy to see that the remaining submatrix, , must be chosen such that
At a first glance, characterizing the linear canonical transformation in terms of the matrix instead of the three matrices , , and might seem unnecessarily convoluted. However, this parametrization is convenient, since it confers on the linear canonical transformation a group property similar to (and more general than) that of the fractional Fourier transformation: two consecutive linear canonical transformations are equivalent to a single one whose parameter is the product of the parameters of the individual transformations, i.e.,17) guarantees that the right-hand side of Eq. (23) has no phase factor in front of it.] It is easy to see from Eq. (19) that the product of two symplectic matrices is itself symplectic.
Let us now discuss four special cases that appear often in the study of optical systems, starting with the fractional Fourier transform.
2.3a. Fractional Fourier Transformation
The particular case of the fractional Fourier transformation corresponds to , where the symplectic matrix is given by
Norm-preserving scalings can be defined according to17) does not lead directly to Eq. (25). One way to show this relation is to notice that . Therefore, the fact that norm-preserving scalings are linear canonical transformations can be shown from taking the limit of a fractional Fourier transform with , and using the method of stationary phase.
It is straightforward to see that . More general, anamorphic versions of this transformation can be written where a different scaling applies in each direction. If these directions are aligned with the coordinate axes, this anamorphic scaling corresponds to a linear canonical transformation where is diagonal.
2.3c. Fresnel Transformation
The Fresnel transformation, used to model paraxial free propagation in optics or nonrelativistic free evolution in quantum mechanics, is given by
Chirping corresponds to a simple multiplication by a quadratic phase factor:
2.4. Phase-Space Shifts
A different type of linear unitary transformation corresponds to combinations of shifts in both sets of Fourier conjugate variables. Let us define an operator that causes a shift in the origin of a function as
The action of many simple systems is expressible as a combination of linear canonical transformations and phase-space shifts. This includes the general linear transformation in Eqs. (10) and (11) in which the linear coefficients and v do not vanish.
3. Simple Phase-Space Distributions: the Windowed Fourier Transform and the Spectrogram
In this and the next section, we discuss different mathematical definitions of phase-space representations, i.e., representations of a function in terms of both sets of Fourier conjugate variables. Perhaps the most intuitive representation of this type is the windowed Fourier transform. Here, the part of a function f around a prescribed location q is isolated through multiplication by a window function that differs significantly from zero only within a small neighborhood around q. This product is then Fourier-transformed, in order to give the “local frequency content” of the function at that region. The windowed Fourier transform is then a function of q, the central position of the window, and p, the Fourier variable. Mathematically, the windowed Fourier transform can be written as36), so that the corresponding expression for this transformation in terms of the Fourier transforms of f and w has a similar form: 36) as the overlap of f with a phase-space-shifted version of the window function by using the phase-space-shift operation defined in Eq. (35):
The windowed Fourier transform is linear in f and complex. For many practical applications, it is more useful to use the squared modulus of this transformation, referred to in the signal analysis community as the spectrogram:36) has no effect on the spectrogram.
In order to visualize the definition of the windowed Fourier transform and the spectrogram, consider a one-dimensional () function given by the sum of two identical Gaussians, one shifted by the action of and multiplied by a constant phase factor :Media 1) in Fig. 1 shows the real (blue) and imaginary (red) parts of for . The second row shows the window function over the same range, for varying q. The third row shows the products of each of the curves in the first row with that in the second row, i.e., the real and imaginary parts of a local detail of f isolated by the window function. The spectrogram for the corresponding value of q, given by the squared modulus of the Fourier transform of the function in the third row, is shown in the fourth row. It is clear that the spectrogram has significant contributions only when the window overlaps the regions where f is significant.
The movie frames (Media 2 and Media 3) in Figs. 2 and 3 present the spectrogram of this function over phase space for several values of , and Φ, for a rectangular window (which equals for and vanishes otherwise) and for a Gaussian window , respectively. The blue dot indicates the coordinates . Unless both and are of the order of unity or smaller, the spectrogram is composed of two separate and similar distributions, one centered at the origin (associated with ) and one at (associated with ). This example shows that the action of is indeed to shift a distribution in phase space by the amounts specified by its two arguments. The movies also show that the effect of the phase Φ between the two separate contributions on the spectrogram is negligible except when the two contributions are so close that they overlap significantly.
One drawback of the windowed-Fourier transform and the spectrogram as representations of local frequency content is that they do not have a unique definition, since they depend on the choice of a window function. In some cases, the window function is dictated by the physics of the problem, as we will discuss later. In general, though, the choice of window function limits the range of frequencies that can be described. For this reason, it can be useful to have a phase-space distribution that does not rely on the choice of an ancillary function. One such distribution is the Wigner function, described in the next section.
4. The Wigner Function
We now introduce the Wigner function as a way to represent the function f jointly in terms of both q and p, without the need for an ancillary function:5.) It is easy to show that the Wigner function can be written in terms of the Fourier transform of f in a very similar form:
The definition of the Wigner function in Eq. (41) is illustrated for the case in Media 4 (Fig. 4), for the same function as in Fig. 1, namely, the function in Eqs. (40) with and . The first and second rows show plots, respectively, of and as functions of . The value of q is varied from −12 to 12 in Media 4. In these plots, the real and imaginary parts are shown as blue and red curves, respectively. Note that, since they are plotted as functions of , these functions are stretched by a factor of 2. The third row shows the product of these two functions, and the fourth row shows the Fourier transform of this product, corresponding to the Wigner function as a function of p for the specified value of q. A comparison of the graphic representation in Fig. 4 with that for the spectrogram in Fig. 1 suggests that, in the case of the Wigner function, is somewhat analogous to a window function acting on . Note that the product has the same local frequency (i.e., the same rate of change of the phase in terms of ) for small as . Note also that, for this choice of f, the Wigner function is significant in three regions: when the contributions due to in the first two rows overlap (i.e., for ), when the contributions due to in the first two rows overlap (i.e., for ), and when the contributions due to in one row overlap with those due to in the other row (for ). This third contribution (visible in the static frame of Fig. 4) will be discussed in more detail later.
Comprehensive summaries of the properties of the Wigner function are given in Refs. [1,8–14]. Here we concentrate on a few of them that are of particular interest for the applications mentioned in what follows.
4.1. Reality and Possible Negativity
The Wigner function is always real, because it is the Fourier transform over of the Hermitian function . That is, changing the sign of the variable of integration in the integrand in Eq. (41) is equivalent to complex conjugation, so the Wigner function equals its complex conjugate. It must be noted, however, that the Wigner function is not necessarily non-negative; for example, the Wigner function plotted in the fourth row of Fig. 4 presents positive and negative oscillations when . The possible negativity of the Wigner function for some values of its arguments is sometimes seen as a conceptual problem when trying to ascribe to it a literal physical meaning as a power or probability density. However, as will be discussed later, these negative values are needed for the description of coherent effects. Note also that multiplying by a global constant phase factor does not change the resulting Wigner function, given its bilinear definition.
Among the main properties of the Wigner function are the so-called marginal relations: the integral over one set of variables of the Wigner function gives the square modulus of the function in the representation associated with the remaining variable:7): 43a) and (43b) correspond to and , respectively. This projection property was found independently by Mustard  and Lohmann , and it is sometimes referred to as the Radon–Wigner transform [33,34]. (A more general relation was found by Bertrand and Bertrand .) As will be discussed in Section 10, this marginal projection property is the basis for a technique called phase-space tomography, used in classical optics in both the position–direction  and time–frequency  phase spaces, as well as in quantum optics [35,38,39]. In fact, the formulas for calculating the Wigner function in Eqs. (41) and (42) can be generalized in terms of the fractional Fourier transform as 36], also discussed in Section 10.
From any of the marginal relations given earlier, it is easy to see that the integral of over all phase space gives the squared norm of f:
4.3. Linear Mapping under Linear Canonical Transformations and Phase-Space Shifts
The Wigner representation for the result of a linear canonical transformation on a given function is just a simple mapped version of the Wigner representation of the original function, namely,2.3. In the graphic examples, we consider a one-dimensional function () and use for simplicity. The function we will use is a Gaussian, given by
4.3a. Fractional Fourier TransformationMedia 5 (Fig. 5) for the Wigner function in Eq. (49) with .
Similarly, a scaling transformation causes a phase-space volume-preserving scaling of the Wigner function:Media 6 (Fig. 6). For anamorphic scalings (for ), each component of q and p is divided and multiplied, respectively, by the corresponding scaling factor.
4.3c. Fresnel Transformation
Fresnel transformation corresponds to a horizontal shearing of the Wigner function:Media 7 (Fig. 7).
Similarly, chirping corresponds to a vertical phase-space shearing of the Wigner function:Media 8 (Fig. 8).
4.3e. Phase-Space Shifts
A similar simple phase-space rearrangement occurs in the case of phase-space shifts: the Wigner function of is given byMedia 9 in Fig. 9. Note that using a different ordering convention in the definition of a phase-space shift would only lead to a global constant phase factor that would have no effect on the Wigner function.
4.4. Inner Product
Suppose that we have two arbitrary square-integrable functions and . The integral over all phase space of the product of their corresponding Wigner functions is always a non-negative real number, given by40,41]. A similar interpretation exists for the classical measurement of a monochromatic wave field [42–45]. In this context, the left-hand side of Eq. (55) has an intuitive interpretation, as will be discussed in Subsection 7.5.
Notice that for , we get
Given the bilinear character of this transformation, the Wigner function of the sum of two functions, , is not, in general, the sum of the corresponding two Wigner functions. Instead we get40) . The three contributions to the Wigner function are easily shown to be given by Media 4 (Fig. 4): the main contributions to the Wigner function occur not only for values of q at which a significant part of f overlaps with itself within the two top rows of that figure, but also when one significant part of f in the first row overlaps with a different significant part in the second row. At these values of q, the product shown in the third row of the figure is composed of two identical (up to an inversion and complex conjugation) contributions, roughly separated by the same distance as that between the two contributions of f. The Fourier transform of this product (i.e., the Wigner function) is then oscillatory, with frequency proportional to this distance.
The dependence of the Wigner function on the location of the second contribution and the phase difference between the two contributions is illustrated in Media 10 (Fig. 10) for the same variation in , , and Φ as in Figs. 2 and 3. It is seen there that the oscillation fringes are parallel to the straight line joining the two contributions. The fringe structure indicates the phase between the two contributions, regardless of their separation: the central fringe is positive if the two contributions are in phase, and negative if they are completely out of phase. The role of these fringes is clear from considering the marginal properties in Eq. (43a) and (43b), shown in Fig. 11(a) for the case : the interference term cancels almost entirely following integration in p, leading to a marginal distribution with two maxima at ; on the other hand, the interference term is responsible for the oscillations of the marginal in q. That is, the negative values of the Wigner function are needed to account for coherent effects like destructive interference. Figure 11(b) shows the Wigner function of a superposition of three Gaussians, one of them multiplied by a linear phase factor in order to displace it from the q axis. Notice that there is an interference term between any two contributions, and that the frequency of oscillations increases with the distance between them.
It was discussed in Subsection 4.4 that the effect of a linear canonical transformation and/or phase-space shift is a simple linear mapping of the Wigner function over phase space. These mappings are the only ones that not only preserve phase-space area, but also map any straight line in phase space onto another straight line. In other words, these mappings are the only ones that guarantee that the half-point interference contributions of the Wigner function remain as half-points following the transformation. Therefore, transformations that are not combinations of linear canonical transformations and phase-space shifts cannot correspond to mappings of the Wigner function in general. Notice that following a linear canonical transformation and/or phase-space shift, the fringes of an interference contribution remain aligned with the straight line joining the two contributions giving rise to it, and an increase/decrease in the distance between two contributions brings with it the corresponding decrease/increase in the transverse spacing of the fringes, given the constraint of phase-space area conservation.
4.6. Oscillatory Character of the Wigner Function
The interference terms described in Subsection 4.5 can cause the Wigner function to be highly oscillatory within some regions of phase space. From a computational point of view, this sometimes constitutes a major shortcoming of this representation: the Wigner function represents a function of N variables in terms of a potentially complicated, highly oscillatory function of variables. Therefore, except for a limited number of cases that admit of an analytic description, the numerical modeling of a system by using the Wigner function requires both doubling the number of variables over which to sample and potentially increasing significantly the sampling rate in each direction in order to be able to resolve the fast oscillations.
The problem with oscillations is particularly critical for functions whose Wigner functions extend over phase space regions much larger than the minimum-uncertainty-limit cell size, . One such case of particular relevance in optics and quantum mechanics is that of functions of the form , where both A and Ω are real, and A varies slowly compared to the oscillations of . In this case, one would expect the Wigner function to be localized in phase space around the N-dimensional manifold described by the equation
As described in detail by Berry and collaborators [46–48], does indeed present a positive ridge structure that mimics the manifold defined by Eq. (60), but due to the interference contributions discussed in the previous section, it also displays an intricate pattern of fast oscillations involving large values, both positive and negative, at many other regions of phase space. These regions are composed of all phase space locations that are midway between any two points along the manifold. Let us consider the case for simplicity, and use the construction proposed in Ref.  to visualize the location of these contributions. Consider a plot of the curve , shown as a black curve in Fig. 12(a). Now superpose on this plot a series of one-half-scale replicas of the same curve, where, for each replica, the anchor point for the scaling is a point along the curve, as illustrated in this figure and the accompanying movie. Let the density of these replicas be proportional to , where q is the location of the anchor point. The region occupied by these replicas indicates the region where the Wigner function is significant, as can be seen from the comparison of this picture with the Wigner function in Fig. 12(b). As explained in Ref. , the amplitude, phase, and orientation of these oscillations can be inferred from this geometrical construction: the amplitude is proportional to the density of replica intersections. The phase of the oscillations, on the other hand, is proportional to K times the area enclosed by the phase space curve and a straight line joining the anchor points of the two replicas crossing at the corresponding point, and the local direction of the fringes is approximately parallel to this line, as shown in Fig. 12. The complexity of the pattern depends on the geometry of the curve; for example, more complicated patterns result for curves including inflection points (like the one shown in Fig. 12), for which more than one pair of points along the curve can have their half-point at the same phase-space location.
The same type of construction can be used for functions associated with phase-space manifolds that are not single-valued in q. Consider the simplest such case, given by a one-dimensional Hermite–Gaussian function of the form13, this result is consistent with the replica-based construction described above. In particular, the amplitude of the oscillations is related to the density of replica intersections, with a maximum at the origin where all replicas cross.
5. Partial Coherence
5.1. Correlations and Coherent Mode Expansions
An advantage of bilinear representations like the Wigner function and the spectrogram over linear ones like the windowed Fourier transform is that they can be used to represent not only deterministic functions , but also the bilinear autocorrelation of stochastic functions at two values of their argument, i.e.,
A useful representation of F is given in terms of its eigenvalues and eigenfunctions , resulting from considering the integral Fredholm equation65) (with ) by , summing over all n, and using Eq. (67), F can be written as a diagonal sum known as a Mercer or coherent mode expansion :
5.2. Wigner Function and Spectrogram for Mixed States
The Wigner function for a partially coherent field or mixed state described by F is given by2] was precisely for mixed quantum-mechanical states, so Eq. (69) is the general definition of the Wigner function. Similarly, the spectrogram for a mixed state is given by
5.3. Overall Degree of Purity
A measure of the overall level of coherence or purity of a field is given by 52]; implies that there is only one coherent mode, that is, that F can be factored out as the product . The analog of Eq. (56) for mixed states is given by
5.4. Quasi-homogeneous Fields
As will be discussed later, there is a type of optical field known as quasi-homogeneous, for which the following factorization exists for the correlation:77). Since μ is much narrower than , is much wider than , so is significantly different from zero over a region of phase space much larger than the minimum uncertainty cell size, .
5.5. Smoothness of the Wigner Function for Highly Mixed States
As discussed in Subsection 4.6, the Wigner function of a deterministic function can be highly oscillatory and include many negative regions. However, for a significantly mixed state involving a very large number of coherent modes (i.e., for ), the Wigner function tends to be smooth and non-negative. While not a rigorous proof, we now give an argument of why this is the case, based on the following two observations: first, the Wigner function of each coherent mode of a mixed state can oscillate between positive and negative values, although overall it is “more positive than negative” given the combination of Eqs. (46) and (66):55) and (66), the following relation holds for :
To illustrate this argument, consider a Gaussian correlation of the form68) for this correlation can be found  to have the eigenvalues of the form 61), that is 71), where, according to Eqs. (63) and (51), the Wigner function of each mode is given by
The overall degree of purity can be found either from the substitution of Eq. (84) into Eq. (73) and the use of the geometric series formula, or from the substitution of Eq. (83) into Eq. (75). The result is simply83) into Eq. (75). Notice that equals the ratio between the cross-sectional phase-space area of the Wigner function for F, given in Eq. (83), to that for a fully coherent Gaussian function. That its, overall, the average relative variation of in any phase-space direction is smaller by a factor of ν than that for the Wigner function of a coherent Gaussian. Media 12 (Fig. 14) shows (a) the Wigner function for the mode of this expansion given in Eq. (86), and (b) the partial sum up to n of the superposition in Eq. (71). One can appreciate how the distribution becomes smoother as more modes are included.
It was argued in Subsection 4.6 that for “pure” or “coherent” systems described by a function , the use of the Wigner function can be problematic from a numerical point of view, due to the fast oscillations and the doubling of the number of variables. However, for “mixed” states this is not the case. First, the correlation F depends on variables, as does the Wigner function. Second, as discussed earlier, the oscillations of the Wigner function tend to wash out for significantly mixed states.
6. Other Phase-Space Representations
The Wigner function and the spectrogram are not the only representations of a function in terms of both its original variables and the Fourier conjugate ones. Extensive reviews of the many alternatives that have been proposed are given in some of the references mentioned earlier [1,8,10–14,54]. Here, I just give a brief summary of some of the standard representations and their properties.
6.1. Distributions Related to the Spectrogram
As mentioned earlier, the window function used in the definition of the spectrogram is sometimes dictated by the physics of the problem. This is the case, for example, of measurement processes that are described by an overlap integral like that in Eq. (55). Note that the definition of the spectrogram in Eqs. (39) and (38) has this same form, where the measuring device is characterized by a wave function in the space corresponding to . This “device” then makes a measurement of the process described by f, focusing around the prescribed values q and p, but respecting the limits established by the uncertainty relation. For example, consider the measurement of a function of time through absorption by resonant vibration modes of the measuring device (e.g., our hearing); in the simplest case where the mode can be modeled as a damped harmonic oscillator, the window function has the form , where Θ is the Heaviside step function and γ is the damping coefficient. For a detailed treatment of the use of the spectrogram (and other distributions) to describe realistic measurements of stochastic time signals, see Ref. .
In cases where the window function is not dictated by the physical situation, it can be chosen for mathematical convenience. For example, it is well known that the function that reaches the lower bound set by the uncertainty relation is a Gaussian, whose Fourier transform is also a Gaussian. Therefore, a Gaussian window minimizes the size of the blur in phase space. When the window is chosen as a normalized Gaussian of width a, i.e.,56]: 10.1), the lowest-order eigenstate of the system (the vacuum mode in quantum optics) is a Gaussian centered at the origin and with a specific width, referred to here as . When a Gaussian of this width is used as the window function, the Husimi function is referred to in the quantum optics literature as the Q function [57–59]: 60–62], defined mathematically such that its convolution with the Wigner function of a Gaussian window of width gives the Wigner function of F:
6.2. The Kirkwood–Rihaczek and Margenau–Hill Distributions65], is also used as a phase-space representation. Note that, unlike the spectrogram (and its special cases, the Husimi and Q functions) and the P function, the Kirkwood–Rihaczek distribution (as well as its real part) satisfies marginal relations analogous to those satisfied by the Wigner function and given in Eqs. (43a) and (43b): 44). For mixed states, the Kirkwood–Rihaczek distribution is given by
6.3. The Ambiguity Function
Another commonly used representation that is defined in a phase space is the ambiguity function [66–68]. However, this function is quite different from those described earlier, since the phase space over which it is defined is the Fourier conjugate of that for the Wigner function and other representations. The ambiguity function is the characteristic function for the Wigner function, i.e., its Fourier (or inverse Fourier) transform over all variables:
The ambiguity function satisfies a phase-space mapping under linear canonical transformations similar to that of the Wigner function. To see this, consider the substitution of Eqs. (47) into the first expression in Eq. (97):19) was used in the second step. On the other hand, phase-space shifts do not correspond to a rigid displacement of the ambiguity function, but to multiplication by a linear phase factor.
6.4. The Cohen Class of Distributions
Except for the ambiguity function, all of the representations described so far are members of the so-called Cohen class of phase-space distributions [54,69,70], which in terms of the ambiguity function can be defined as1. Note that, for all these distributions, the kernel is normalized so that , which guarantees that the integral of the distribution over all phase space is independent of the functional form of the kernel:
Since all Cohen-class representations are bilinear, they all present interference contributions like those of the Wigner function. However, these effects can be suppressed, amplified, or relocated depending on the kernel, as can be seen from Eq. (101). Since the interference contributions to the Wigner function are oscillatory, they contribute to the ambiguity function at regions away from the origin in the space. If the kernel decays away from the origin, as in the case of the spectrogram and its special cases, the Husimi and Q functions, then these interference terms are largely suppressed. For the P function, on the other hand, for which the kernel grows extremely fast away from the origin, the interference terms are greatly enhanced, causing the distribution to be singular except in cases where decays sufficiently fast as to overcome the kernel’s divergence. In the case of the Kirkwood–Rihaczek distribution, the kernel is a phase factor with unit magnitude, so the interference effects are of similar significance as those for the Wigner function, but shifted to a different region of phase space, as will be discussed later.
Notice that, by using the convolution theorem, Eq. (101) states that all members of the Cohen class are given by the convolution over phase space of the Wigner function with a blur function, i.e.,104) is the one described in Eq. (89), so . However, in the case of the P function, Φ is not well defined. The blurs for the different Cohen-class distributions are also presented in Table 1.
The inner-product property of the Wigner function in Eq. (55) is a special case of a more general relation for Cohen-class phase-space distributions, given by1,10,14].
To illustrate some of the distinctive features of these different representations, let us consider the function . As mentioned in Subsection 4.6, it would be expected that the significant values of a phase-space distribution should be concentrated around the phase-space curve . Figure 15 shows (a) the Wigner function, (b), (c) the Husimi function for two values of the width a, (d) the magnitude, and (e) the real part of the Kirkwood–Rihaczek distribution, all of which can be calculated analytically in terms of Airy functions. The Wigner function presents the characteristic interference contributions within the concavities of the curve. For the Husimi function, on the other hand, these contributions are largely suppressed, except for some small effects near the vertex of the parabola in 15(c), because the width of the blur in q is not sufficiently narrow to resolve this vertex properly. Figures 15(d) and (e) illustrate the fact that it is the real part of the Kirkwood–Rihaczek distribution (i.e., the Margenau–Hill distribution) and not its magnitude that holds some resemblance to the curve. However, even the real part of this distribution contains significant oscillatory interference contributions over large regions of phase space. These regions correspond to the intersection of all values of q occupied by and all values of p occupied by . A second example is given in Media 14–16 (Fig. 16), which show the transformation of several of these distributions under fractional Fourier transformation for a sum of two Gaussians of the form . Note that, of all the representations shown in the figure, only the Wigner and the Q function rotate rigidly in phase space with θ. While the P function would also transform in this way, it is highly singular for this case and therefore not plotted. Notice also the counterrotating interference contributions in the case of the real part of the Kirkwood–Rihaczek distribution.
Table 2 shows whether or not several of the Cohen-class distributions transform as a simple mapping under different linear unitary transformations. In the last column, the transformation represents a general unitary transformation that is not a combination of phase-space shifts and linear canonical transformations.
Figures 15 and 16 illustrate the potentially very different appearances of all the bilinear phase-space representations described so far for the case of a pure function . However, it is important to note that these differences tend to wash away in the case of mixed functions : the more incoherent or mixed the function is, the more similar all these representations become. This is particularly easy to see in the quasi-homogeneous limit, where F satisfies the factorization in Eq. (76). For such fields, the ambiguity function takes the form101), the Cohen-class distributions can be written as 79): 82). The Wigner function [already calculated in Eq. (83)], as well as the Q, P, and Kirkwood–Rihaczek functions, can be found to be given by 76), but for any function that differs significantly from zero only for values of much smaller than the scale of variation of .
As mentioned in Section 1, the main topic of this tutorial is the conservation of the Wigner function along paths associated with a simple physical model. As discussed in Subsection 4.3, there are two general classes of linear unitary transformation on a function f whose effect in phase space is a simple mapping (i.e., rearrangement of arguments) of the Wigner function: linear canonical transformations and phase-space shifts . The most general class of unitary transformation that corresponds to a mapping of the Wigner function is an arbitrary combination of these two types of transformations. Therefore, the conservation property of the standard Wigner function is restricted to systems where propagation or evolution is described by a combination of phase-space shifts and linear canonical transformations. However, as will be seen in Section 11, it is possible to retain the conservation property for other physical systems by modifying the definition of the Wigner function.
7. Wigner-Based Modeling of Beam Propagation
In optics and quantum mechanics, one of the main properties of the Wigner function is that, for a certain class of simple systems, its spatial propagation or temporal evolution follows the same rules as the simpler physical formalisms it mimics (i.e., geometrical optics or classical mechanics). In the case of paraxial classical optical propagation, this class of systems includes, but is not limited to, the so-called ABCD systems , which are characterized in both the ray and the wave models by a 4 × 4 symplectic matrix. In this section, we first give a brief review of the ray-optical and radiometric description of this type of system, and then show that, by using the Wigner function, the propagation of wave fields through these systems follows the same rules, namely, conservation along rays.
7.1. Paraxial ABCD Systems in the Ray Domain
At a plane of constant z, a ray is fully characterized by its intersection with this plane and by its direction. The intersection’s position is specified by its transverse coordinate, , while direction is specified by the optical momentum two-vector, , corresponding to the product of the local refractive index and the transverse part of a unit vector in the ray’s direction (see Fig. 17). The ray is then completely characterized by the state vector , which gives the coordinates of the ray in phase space. The paraxial approximation implies that , where θ is the angle between the ray’s tangent vector and the z axis. For a simple class of optical systems, the propagation of a paraxial ray is described mathematically by a simple relation involving linear combinations of x and p. These systems, known as linear systems or ABCD systems , are fully characterized by a 4 × 4 symplectic matrix composed of the 2 × 2 matrices , , , and , so that the initial state vector of the ray, , corresponding to a plane intersecting the optical axis at , transforms into the final state vector at , , through multiplication by this matrix:
The fact that, in this context, must be symplectic follows from a simple derivation of Eq. (111), based on the Hamiltonian formalism of ray optics [72–74]. The part of an optical system contained between the planes at and z can be described by Hamilton’s point characteristic , which equals the optical length of an extremal path (i.e., a ray) joining the points and . (This description is complete provided that there is one and only one extremal path joining these points; otherwise, a different type of characteristic function must be used .) In the paraxial regime and for a simple class of systems, sufficiently accurate results are obtained from the quadratic expansion of the point characteristic. Further simplification results from assuming that the z axis is an axis of symmetry of the system, i.e., that the system has mirror symmetry around two perpendicular planes intersecting at this axis. (This includes, as a special case, systems with rotational symmetry around the z axis.) In this case, the quadratic expansion does not include linear terms, and can be written, without loss of generality, in a form analogous to Eq. (14):15). We now arrive at Eq. (111) by using Hamilton’s equations : 21), meaning that is indeed symplectic.
To illustrate the type of system that can be modeled through a linear relation involving a symplectic ABCD matrix, consider the following three simple examples:
- (a) Propagation through a homogeneous transparent medium of index n leaves the ray’s slope unchanged but changes the ray’s transverse position in an amount proportional to both the slope (approximately equal to in the paraxial approximation) and the propagation distance , so that the matrix is given by 28).
- (b) Propagation through a perfect thin lens of focal distance f leaves the transverse position unchanged, but redirects the rays in an amount proportional to x, so that the matrix takes the form in Eq. (29):
- (c) Propagation by a distance through a gradient index (GRIN) medium with refractive index is given by the matrix 24), this matrix is, up to a scaling factor, a rotation matrix by an angle .
7.2. Radiometry and the Radiance in the Paraxial Regime
Ray optics gives a simple approximate description of the paths followed by optical radiation, but not of the amount of power they carry. For this purpose, the theory of radiometry supplements the rays with a weight, through the definition of the radiance or specific intensity, which corresponds to the amount of power per unit transverse area per unit solid angle . A closely related and perhaps more fundamental quantity, known as the basic radiance , is equal to the radiance divided by the square of the refractive index, i.e., to the amount of power per unit volume of phase space:
When losses at interfaces are negligible, the basic radiance is conserved under propagation along rays , and can therefore be regarded as a measure of the power density carried by the rays in the system. For an ABCD system, this conservation can be written as
Consider, for example, propagation through free space from an initial plane at , for which is as given in Eq. (114) with and . In this case, this propagation law reduces to121) for this matrix, i.e., 116) into Eq. (121):
7.3. Paraxial ABCD Systems in the Wave Domain
In the wave domain, the propagation of a monochromatic paraxial scalar field from an initial plane at , traveling exclusively toward larger values of z, can be described by a linear expression of the form77–79]:
The expression for the propagator in Eq. (127) was first proposed in optics by Walther [77,78]. An analogous propagator was proposed earlier by Van Vleck  to describe the evolution of solutions to Schrödinger’s equation in terms of classical particle trajectories, where instead of V one uses the classical action. This result can be considered as a form of the Huygens–Fresnel principle: the field at the initial plane can be decomposed into secondary sources, each emitting a field proportional to whose phase delay at the observation point is proportional to the optical distance from the secondary source, and whose optical intensity is proportional to the density of rays around the observation point arriving from the secondary source point [i.e., ]. Note that the original expressions by Walther and Van Vleck did not include the Maslov index correction, which was added by Gutzwiller  within the quantum context.
For situations where the ABCD formalism is valid in the ray-optical regime, sufficiently accurate results are obtained from using in Eq. (127) the quadratic expansion of the point characteristic in Eq. (112). It is straightforward to see that, in this approximation, the propagated field is proportional to a linear canonical transformation as defined in Eqs. (17) and (14) with , applied to the initial field:129) is known as Collins’ formula . For example, for propagation through free space, for which , , and is given in Eq. (114), this expression reduces to the Fresnel propagation integral: 116), is proportional to a fractional Fourier transform of degree .
7.4. Wave Propagation Through ABCD Systems in Terms of the Wigner Function
The Wigner function of the field over the transverse variables gives a representation in terms of both transverse position x and transverse momentum p:43a), the integral of the Wigner function over all directions gives the local optical intensity: 47) imply that the Wigner function is constant along rays in these systems, i.e., 132), (133), and (134) are identical to Eqs. (120), (121) and (122), with the Wigner function playing the role of the basic radiance. This implies, for example, that free propagation in the paraxial wave domain corresponds to a horizontal linear phase-space shear like that in Eq. (124) and illustrated in Fig. 7. Similarly, propagation across a thin lens imparting quadratic phase (i.e., a chirp) has the effect of a vertical linear shear on the Wigner function like that in Fig. 8, and propagation through a quadratic GRIN medium is equivalent to a phase-space rotation like that in Eq. (125), illustrated in Fig. 5. The fact that propagation through ABCD systems corresponds to a simple linear mapping of the Wigner function was first noticed by Bastiaans [70,83–86]. Given this similarity between the properties of the radiance and the Wigner function, it is tempting to state that the Wigner function constitutes a rigorous wave-optical foundation for the theory of radiometry, at least within the paraxial regime. The one conceptual objection to this association, though, is the possible negativity of the Wigner function for some bundles of rays. Under the radiometric interpretation, these rays would then be some sort of “anti-rays” that cancel some of the power transported by other rays through the system.
To illustrate these ideas, let us consider the simplest physical situation, corresponding to propagation through a region of free space surrounding the plane . The substitution in Eq. (124) of and with L replaced with in Eq. (132) gives the expressionMedia 17 and Fig. 18 for the case of a field that is independent of y, for which the phase space is two-dimensional. In this figure, the initial field is chosen to be constant for and zero otherwise. Within the paraxial regime, this field corresponds, for example, to a plane wave traveling in the z direction that is diffracted a slit of width a. The Wigner function, which can be calculated analytically, is shown in Fig. 18(b). Its marginal in p does return the desired rectangular distribution. As shown in Media 17, propagation in z is described by a horizontal linear shearing of the Wigner function, whose marginal in p gives the intensity profile at the corresponding propagation distance. This example shows the importance of the negative values of the Wigner function in order to account for destructive interference. It also gives a simple visualization of the fact that, at the propagation distance , a peak in the axial intensity is achieved, associated with the vertical alignment of some of the main positive regions of the Wigner function. This is precisely the distance at which, for example, a pinhole camera with a square pinhole of side a would achieve its best resolution. (For a circular pinhole of diameter a, the numerical factor is slightly smaller, 0.044, as found by Rayleigh.)
An alternative presentation of the same example is given in Fig. 19. Figure 19(a) shows the intensity distribution over the (scaled) physical plane. Figure 19(b) shows a superposition of lines representing rays, which are either black or white depending on whether the Wigner function for the corresponding ray is negative or positive, and whose opacity is proportional to the magnitude of the Wigner function, such that rays for which the Wigner function is close to zero are nearly invisible. These lines are superposed over a grey background (representing zero), so that the final shade of grey at each point depends on the number and weight of the positive and negative rays going through it. Media 18 shows partial superpositions of these rays for the sake of illustration. It is seen there that, once a sufficient number of rays is included, the ray superposition in Fig. 19(b) reproduces the features of the wave pattern in Fig. 19(a).
It must be noted that ABCD systems are not the most general type of system for which the rules of paraxial propagation for the Wigner function exactly mimic the radiometric model. As discussed in Section 4, the most general type of transformation that corresponds to a simple mapping of the Wigner function in phase space is a combination of linear canonical transformations and phase-space shifts. This corresponds to systems where Hamilton’s characteristic function can be approximated by a quadratic expansion like that in Eq. (112), but where this expansion can include linear terms that account for asymmetries of the system. These linear terms lead to phase-space shifts. Lateral spatial shifts can be introduced, for example, by a tilted étalon, or by a deliberate lateral shift of the optical axis in order to model a misaligned optical element. Directional shifts can be introduced by, say, a transparent wedge, or by a deliberate small rotation of the optical axis in order to model a tilted optical element.
7.5. Phase-Space Interpretation of Measurement
As discussed in , an ideal measuring device can be thought of as a source that emits a characteristic field whose amplitude and phase are such that it cancels as much of the incident field as possible. (For example, a point detector that measures the local intensity can be thought of as a point source emitting a spherical wave.) Let us consider the paraxial case. In general, let the field specific to some detector be given, at a plane z (after the detector), by the function , times an undetermined constant complex amplitude, . The power absorbed by the detector is then the difference between the power of the incident field U and that of the field past the detector, . This absorbed power is given by136) with respect to the real and imaginary parts of we find 136) leads to the expression 55) with . Therefore, according to that equation, the result of the measurement is proportional to the overlap integral of the Wigner function of the field, , and that of the detector, , at the corresponding plane. This form suggests an intuitive ray-based interpretation for the measurement process: the Wigner function of the detector indicates what rays can be accepted into the detector and in what measure. The measurement is then equal to the integral of the power of all rays from the incident field that is accepted in the detector.
Of course, this interpretation should not be taken literally, since rays are not real, independent physical objects, but rather they are just a mathematical construct. Further, the Wigner function of the field can be negative for some thin bundles of rays. However, due to the uncertainty relation, the Wigner function of the detector cannot be arbitrarily narrow, so that it cannot be focused on arbitrarily narrow bundles of rays of the field. In particular, it cannot focus on narrow bundles of rays over which the field’s Wigner function is purely negative; a sufficient amount of “positive rays” is always also accepted into the detector, so that the result of the measurement is always non-negative. Something analogous happens with possible negative regions of the Wigner function for the measuring device.
7.6. Applications of the Wigner and Ambiguity Functions in the Paraxial Regime
In the past three decades, a large number of articles have been published whose topic is the application of the Wigner and ambiguity functions to the study and visualization of optical systems within the paraxial regime. It is impossible to do justice to all of this body of work within the limited space provided by this article, and I apologize in advance to the many authors whose interesting work is not mentioned or cited here. Since the main emphasis of this article is the use of the simple conservation properties of the Wigner function, the applications singled out in this section are ones that exploit these properties.
7.6a. Imaging Systems and the Optical Transfer Function
The expression for the intensity of a field propagating through free space in terms of the Wigner function is given in Eq. (135). The corresponding expression in terms of the ambiguity function is given by139) gives 87]. As pointed out by Brenner et al. , for the field generated by a point source in object space, the ambiguity function at “radial” planar sections defined by corresponds to the OTF of the system for different amounts of defocus z, as can be seen from Eq. (140). Based on this realization, both the Wigner and ambiguity functions have been employed in the design of optical systems, e.g., to study axial resolution  or to image transparent objects with partially coherent illumination .
7.6b. Propagation-Invariant Beams and Extended Depth of Field
There is a class of solutions of the paraxial wave equation whose intensity profile is exactly invariant under propagation. These solutions do not correspond to realistic fields, since they require an infinite amount of optical power. However, they can be achieved approximately, leading to fields whose transverse intensity profile is nearly invariant under propagation within some range of distances. Perhaps the best known such solution is that of Bessel beams , given by the expression92].) This Wigner function is defined in a four-dimensional phase space, so it is difficult to visualize its horizontal shearing due to free propagation. However, it is straightforward to see that this Wigner function is invariant under such shearing, that is, , since the only occurrence of x is in a dot product with , which is perpendicular to p; so replacing x with has no effect. This invariance under horizontal shearings, though, implies that the Wigner function extends over all values of x, meaning that its integral is infinite. This is consistent with the fact that propagation-invariant fields have infinite power.
There is an interesting second class of propagation-invariant field, whose intensity profile shape is invariant under paraxial propagation, although it does not travel along a straight line but along a parabolic path. The archetypal beam of this type is the Airy beam . Let us consider first the version of this beam for a space with only one transverse variable, x, in addition to the propagation distance, z. The formula for an Airy beam is then given by20(a), is given by 15, but rotated by 90°: Media 19 (Fig. 20) that a horizontal shearing caused by propagating a distance z is equivalent to a phase-space shift. By substituting in Eq. (147), it is easy to show that this phase-space shift is given by ; the spatial part of the shift is responsible for the curvature of the beam’s path, while the directional part leads to the phase factor in Eq. (144). By using the projection properties of Airy functions , one can show that the marginal in p of this Wigner function equals the intensity distribution in Eq. (145). Note that the intriguing feature of Airy beams of following a curved trajectory becomes less mysterious under this phase-space picture: all rays composing this beam travel along straight lines, but the rays that give rise to a maximum at a given value of z are not the same ones as those that give rise to the maximum at a different propagation distance. An analogous interpretation can be given in terms of caustics . A treatment that is similar to the one presented here was published recently .
A three-dimensional version of the Airy beam can be composed as the product of the field in Eq. (144) times a similar solution in terms of y [divided by ], and the corresponding Wigner function is proportional to the products of the Wigner functions in Eq. (147) with arguments and , respectively. Approximations to these beams can be generated [96,97] by placing at the pupil of an optical system a transparent mask whose phase has a cubic dependence [as in Eq. (146)] in both Cartesian directions. Coincidentally, this cubic phase mask was proposed by Dowski and Cathey  as a way to extend the depth of field of imaging systems. With such a cubic phase, the PSF corresponds approximately to the intensity of an Airy beam, which remains nearly invariant for a range of defocus distances. Further, given this invariance and known profile, this PSF can be deconvolved computationally, roughly independently of defocus. Of course, there is a small deviation of the defocused PSFs with respect to the ideal geometrical location due to the curvature of the path followed by the maximum, but this deviation is negligibly small in imaging applications. Approximations to Bessel beams, particularly those of zeroth order, have also been used to extend depth of focus [99–101]. These are achieved, for example, by the inclusion of either annular pupils or axicons within the imaging system. Other studies of extension of the depth of field of an imaging system based on phase-space distributions are given in [102,103].
7.6c. The Talbot Effect
The ray-like propagation of the Wigner function also provides a simple, visual description of the Talbot effect. In this effect, fields that are periodic over the initial transverse plane form images of themselves at certain distances of propagation in free space, without the use of lenses or mirrors. Consider again, for simplicity, a two-dimensional field (independent of y). Let the field at the plane , , be periodic with period Λ, i.e., for any x. The Fourier transform of is then a discrete sum of weighted Dirac deltas:148) facilitates the analysis that follows.] As found by Testorf and co-workers [104–106], the Wigner function of such a field at the initial plane has the following form: 21. (Again, for a realistic field, these delta functions must be replaced by narrow finite distributions.) Notice that this spacing in p is half of that between the delta functions in in Eq. (148), due to the presence of half-point interference contributions in the Wigner function. At the discrete set of values of p where it differs from zero, the Wigner function has different periodicity in x depending on whether p is an even or odd multiple of ; it can be seen from Eq. (149) that 21, where the values of the Wigner function at even rows has periodicity , while for odd rows it oscillates between positive and negative values with period Λ.
Given its periodicity in x and discreteness in p, one can see that the Wigner funcion is also periodic under propagation in z, due to the linear shear relationMedia 20 (Fig. 21). As discussed in Refs. [104,106], this phase-space picture also helps visualize the so-called fractional Talbot effect. Consider, for example, the case when . Given the phase-space shearing, all even rows of the Wigner function, corresponding to , shift by an integer number of periods, returning to their state at . On the other hand, the odd rows, corresponding to , shift by a half-integer number of periods. Given the periodicity properties in Eq. (150), the resulting Wigner function is identical to the initial one except for a global shift in x of . This explains why, at this propagation distance, the intensity pattern is identical to that at except for a half-period shift. Extensions of self-imaging to polychromatic fields have also been studied based on the Wigner function .
Note that this phase-space picture of the Talbot effect has several features in common with that of the propagation-invariant beams discussed in Subsection 7.6b. In both cases, the Wigner function after horizontal shearing is identical to the initial one (or a shifted version of it). The only difference is that, in this case, this reproduction of the initial state happens only at a discrete set of propagation distances. Nevertheless, this similarity helps explain other common features, like the one referred to as ”self-healing.” Consider the case where a section of the initial field (a Bessel or Airy beam, or a periodic field) is perturbed somehow, e.g., through blocking by an obstacle. After some propagation, this perturbation washes away, so that the beam returns to its original unperturbed shape. As in the case of Bessel and Airy beams, local features of the field after propagation should not be thought of as being the same physical entity as the corresponding features of the initial field, even if they are identical in shape. Such misinterpretation can lead to unphysical conclusions, like the violation of special-relativistic causality.
7.6d. Sampling and Holography
The Wigner function also gives a useful visualization of why the well-known sampling theorem can be relaxed for chirped signals [108,109], e.g., for paraxial beams that have propagated away from their waist plane. Shannon’s sampling theorem states that the maximum sampling spacing over the transverse variable x that allows a complete reconstruction of the field is inversely proportional to the support of the Fourier transform of the field at the initial plane. Paraxial free propagation does not change this support, given the relation . One could therefore think that the maximum sampling spacing remains constant under propagation. Gori  showed that, for a field that is strictly limited to an interval, say, at , the maximum transverse sampling spacing needed to fully recover the field increases linearly with z. Similarly, for fields whose Fourier transform is (at least approximately) restricted to a region , the sampling frequency needed to approximately reproduce the field can be relaxed under propagation. A simple interpretation of this relaxation in terms of the Wigner function follows a derivation that is mathematically similar to that mentioned earlier for the Talbot effect, but with x and p exchanged. Consider sampling a two-dimensional field at a discrete number of equidistant points over the initial line (), i.e., , where Λ is the sampling distance. The Wigner function of this sampled field can be found to be given by [108,109]22(a) shows schematically the replicas corresponding to even , separated vertically by . To be able to reconstruct the original field, this separation must be sufficient so that the different replicas do not overlap. If the field is not chirped, this condition reduces to , where is the support of the angular spectrum. This relation can be written as , which is the standard sampling theorem. If the field is chirped, e.g., due to propagation away from its waist plane, its Wigner function is sheared, and the replicas can be stacked more densely without overlapping significantly, relaxing the minimum value of Λ.
The revised sampling theorem described earlier has been applied to digital holography [111,112]. Similarly, the Wigner function has been used to explain the principles of regular [113–116] and volume holography . Other recent applications of the Wigner function include descriptions of the moiré effect  and the detection of subsurface targets closely located beneath a randomly rough surface .
7.7. Partially Coherent Fields
The fact that paraxial propagation through ABCD systems corresponds to a simple remapping of the Wigner function makes this formalism a useful numerical tool in the case of partially coherent fields. In the frequency domain, a partially coherent field is described by its cross-spectral density, defined as the correlation of the field at two points:120]. In both the traditional and Wigner approaches, computational savings can be achieved by using fast Fourier transform algorithms or asymptotic integral estimations, whenever appropriate. For partially coherent fields involving only a reduced number of coherent modes, however, it might be more efficient to use a modal decomposition and propagate each mode independently.
While the Wigner function gives direct access to one-point field properties such as the optical intensity, its ray-like propagation properties have also been exploited for the propagation of two-point quantities like the cross-spectral density, through the inversion of Eq. (154) [121,122].
7.8. Wigner Functions that Account for Polarization in the Paraxial Regime
Wigner functions have also been proposed to model polarization effects in the paraxial regime. These functions take the form of a 2 × 2 matrix [123,124], corresponding to the substitution of the 2 × 2 cross-spectral density matrix into Eq. (154):125,126] defined a Stokes–Wigner function that takes the form of a 4-vector whose elements are a generalization of the Stokes parameters, calculated as the trace of the product of and each of the four Pauli matrices
8. The Use of the Wigner Function for More Complicated Systems and/or Nonparaxial Fields
The simple propagation properties of the Wigner function are limited to paraxial systems in which rays evolve according to linear relations. For more complicated systems involving optical elements introducing phase delays whose dependence on the transverse position are more complicated than quadratic (e.g., lenses with aberrations), propagation does not correspond to a simple rearrangement of the arguments of the Wigner function; it must be described instead either through an integral propagation kernel  involving all phase-space variables, or by a series of differential operators [10,14,41] acting on the Wigner function. This is also the case for diffraction by planar obstacles, where one can use an integral relation to model diffraction [127,128]. In the case of fairly spatially incoherent fields, a radiometric-like model can be used, where the primary effect of the obstacle is to block the rays incident on it (i.e., the value of the Wigner function for these rays is set to zero), while those missing the obstacle continue to propagate (i.e., the value of the Wigner function for these rays remains unchanged). Additionally, the edges of the obstacle behave as secondary sources emitting rays , whose Wigner function is related to that of the rays incident on the edges. Approximate paraxial propagation models based on the Wigner function have also been proposed for the case of turbulent media [130–133], as well as for nonlinear Kerr media [134–138].
Even for propagation in free space (or a homogeneous transparent medium), the Wigner function is not conserved along rectilinear paths for fields that are not paraxial. To show this, let us use the angular spectrum representation of a field traveling in free space:158) reduces to , resulting in the fact that , as expected. However, beyond the paraxial regime, the propagation of the Wigner function is not in general described by a simple rearrangement of arguments. This propagation is ruled not by a simple first-order differential equation like Eq. (134) but instead by a high-order differential equation . As mentioned in the next subsection, this difficulty is alleviated only in the case of significantly spatially incoherent fields.
8.1. Generalized Radiometry for Nonparaxial Fields
The first use of the Wigner function in classical optics was within the context of nonparaxial partially coherent fields. The objective of such work was to give a wave-optical foundation to the heuristic theory of radiative transfer [140,141] and in particular to the radiometric formalism . These efforts concentrated on deriving a wave-based radiance. This was first done, independently, by Dolin in 1964  and by Walther in 1968 . In particular, Walther proposed the use of the standard Wigner function over the transverse variables, and introduced an extra obliquity factor to account for nonparaxiality of the field:154). Several other generalized radiances, as these wave-based definitions of the radiance are referred to in the literature, have been proposed. Many articles on this topic are collected in the volume edited by Friberg , and useful summaries on the subject are given in the articles by Wolf  and Friberg . Many of these generalized radiances result from the use of other members of the Cohen class . For example, Walther proposed a second generalized radiance that turns out to be the Kirkwood–Rihaczek distribution times an obliquity factor .
There are two objections to the wave-based definition of the radiance in Eq. (159): it can be negative, and it is not rigorously conserved along rays beyond the paraxial regime. However, it is easy to show that, in the limiting case of significantly incoherent fields of the quasi-homogeneous type, these two problems disappear [145,147–148]. Consider a field whose cross-spectral density at a plane has the form159) in the form of Eq. (158) as 160), this expression becomes162) can be approximated by its second-order expansion in . It is easy to see that even-order terms in this expansion cancel, so only the linear term survives, leading to 160) but also for any field where the directional correlation differs significantly from zero only for separations of its two arguments that are much smaller than unity.
As discussed in Subsection 6.4, any definition of the generalized radiance based on a standard Cohen-class distribution would lead to the same result in this limit (including Walther’s second generalized radiance definition ). Light generated by thermal sources typically presents coherent widths that are small compared with the scale of variation of the intensity, so any Cohen-class generalized radiance can be used to justify the radiometric model for the case of thermal light.
9. Short Pulses and Propagation in Linear Dispersion Media
Another area of optics where phase-space distributions have been applied is the characterization of pulses, and the description of their propagation through dispersive media [149–153]. Most work on this topic has been limited to a one-dimensional treatment, useful for studying highly collimated fields traveling through homogeneous media or propagation modes in waveguides, particularly optical fibers, where the approximately invariant transverse structure of the field can be factorized from the more interesting longitudinal/temporal behavior. The vast majority of this work is based on the use of the phase space composed of time and frequency [149,151,154,155], often referred to as the chronocyclic phase space. However, some authors consider instead the phase space composed of the propagation distance and wavenumber [150,152,153].
As discussed in this section, within some approximations, the propagation of pulses through certain optical elements turns out to be mathematically analogous to paraxial beam propagation through ABCD systems [156,157]. Useful reviews of this analogy and its applications to the characterization of ultrashort pulses is given in the book chapter and review article by Dorrer and Walmsley [154,155]. For simplicity, most of the treatment below considers single pulses, described by a component of the electric field, , where z is the spatial propagation variable. However, for the characterization of pulsed lasers used for ultrafast applications, it is convenient to use a statistical formalism analogous to that of the theory of partial coherence. Here, each pulse in the train is treated as a member of an ensemble. Since the different pulses are generally not identical, the pulse train is characterized by the correlation function [158,159]
Within the linear domain and neglecting the effects of absorption, a dispersive medium or system accepts monochromatic solutions of the form . The temporal frequency ω and wavenumber (or spatial frequency) k are related through the dispersion relation
A pulse traveling in a dispersive medium of this type can be expressed as a superposition of monochromatic components:
9.1. Wigner Function of Time versus Frequency
Consider the Taylor expansion of k up to second order around the central frequency of the pulse, , i.e.,164).] From the substitution of Eq. (166) into Eq. (170), it is easy to show that, if the approximation in Eq. (168) is valid, the Wigner function propagates as 168), is a linear function of frequency. That is, propagation in z is described within this approximation by a simple linear shearing of the Wigner function over the phase space . One way to understand this follows from writing Eq. (166) in the following form: 130), i.e., 114), i.e., 124). The effect of the factored carrier phase is a phase-space shift.
As discussed in Refs. [154,155], devices such as quadratic temporal phase modulators introduce chirp factors and hence play a role analogous to lenses; i.e., they correspond to linear canonical transformations with matrices analogous to that in Eq. (115). The combination of dispersive media and quadratic temporal phase modulators then allows the composition of general “temporal ABCD systems” that can perform transformations such as imaging, magnification, and Fourier and fractional Fourier transformation.
9.2. Wigner Function of Position versus Wavenumber
Consider instead the parametrization of ω in terms of k, and assume that, over the range of spatial frequencies occupied by the pulse, a quadratic Taylor expansion of the frequency around (the central wavenumber) is sufficient:176) is valid. Evolution in t then corresponds to a linear shearing in the phase space , since is linear in k. Again, this is because the field can be written as a solution to a Schrödinger-like equation, although on different variables. Let us again factor out the carrier phase and express the remaining function in terms of regular and retarded times, 114), although in this case the propagation parameter is time, not distance.
10. Phase-Space Tomography
Phase-space distributions like the Wigner and ambiguity functions have been used for the recovery of the coherence properties of mixed states from measurements of the marginal projections through a set of techniques known as phase-space tomography, first proposed by Bertrand and Bertrand . This approach works particularly well for systems where phase space is two-dimensional and for which propagation or evolution corresponds either to a Fresnel-type integral (i.e., to a linear shearing of phase space) or to a fractional Fourier transform (i.e., to a rotation of phase space). The idea is that, if sufficient marginal projections are obtained through measurements at different stages, the Wigner function can be estimated through standard tomographic methods like the inverse Radon transform. Once the Wigner function is known, the correlation function can be recovered. Much of the pioneering work on this type of technique came from the same research group [36–38], who applied it to three different physical contexts, as discussed in what follows.
10.1. Phase-Space Tomography in Quantum Optics
The first experimental implementation of phase-space tomography was by Smithey et al.  (following work by Vogel and Risken ), who used this technique to recover of the quantum state of a light field. To let us discuss this technique, I give a very brief description of the use of the Wigner function in quantum optics. To make the presentation as similar as possible to that in the rest of this article, the standard quantum notation of operators and density matrices is avoided. For more standard treatments of the topic, the reader can consult the books by Mandel and Wolf , Gerry and Knight , Leonhardt , and Schleich .
In classical optics, the electric field is the “wave function” of interest, and its variables are position and/or time. Consider a Cartesian component of the electric field at a fixed point for a monochromatic mode within a cavity. In the classical picture, this component can be written as23(a). The angular coordinate of this point in this plane corresponds to the phase of the field, while the radial one is proportional to the amplitude, and therefore related to the intensity or number of photons. As time advances, this point rotates around the origin along a circular path, with frequency ω. (Note that this phase-space motion is mathematically analogous, although entirely different physically, to that for a ray propagating paraxially in a quadratic GRIN medium, where q is the transverse position, p is the transverse optical momentum, and replaces .)
In the quantum-optical description of the field, on the other hand, the quadrature q is not a simple function of time, but the variable of a wave function . Further, q and p cannot be prescribed simultaneously, since they are Fourier conjugate variables (with ). However, by using a phase-space representation like the Wigner function, we can represent the state of the field at any given time as a function of both q and p, shown in Fig. 23(b). That is, at a fixed time, instead of having a delta-like localized distribution, the field is represented by a distribution of “possible quadrature pairs” or “possible amplitudes and phases” whose extent in phase space is limited by the uncertainty relation: the narrower the distribution is in q, the wider it is in p, and vice versa. Since we have, instead of a single point, an extended distribution, the field does not have well-defined amplitude and phase: the smaller the uncertainty in the phase (related to the angular spread of the distribution with respect to the origin) of the field is, the larger the uncertainty in the number of photons (related to the radial spread), and vice versa.
As mentioned earlier, temporal evolution in the classical picture corresponds to uniform rigid rotation around the origin. The same is true in the quantum picture when using a phase-space representation like the Wigner, Q, or P functions, since temporal evolution corresponds to fractional Fourier transformation in q of , with scaling and degree . Recall from Subsection 5.2 that fractional Fourier transformation corresponds to a rigid clockwise rotation for these three representations.
Let us consider the Wigner function. Since Gaussian wave functions achieve the maximum possible joint localization in q and p, the Wigner representation of a quantum field that is closest to the phase-space representation of a classical field (a single point in phase space) corresponds to a Gaussian wave function, shifted to the desired phase-space point. This Wigner function is itself Gaussian in both q and p. If this Gaussian has circular cross-section, then the state is called a coherent state, meaning that its width in q and p does not change under temporal evolution. If the phase-space Gaussian is instead elongated in one direction, then the state is referred to as a squeezed state, and its marginal distributions will not only oscillate harmonically but also “breathe,” expanding and contracting at a frequency as time progresses. If the elongation is in the radial phase-space direction [as shown in Fig. 23(c)], then the phase of the field (related to the angular spread subtended by the distribution from the origin) is better determined than for a coherent state of the same centroid, at the cost of the number of photons being more uncertain.
By using a series of balanced homodyne measurements, one can estimate the marginal in p of the Wigner function of the quantum state corresponding to different stages of its evolution. These projections at different times are equivalent to projections in different directions of the Wigner function at a fixed time, from which the Wigner function is inferred by use of the inverse-Radon transformation. An extensive review of this method is given in the book by Leonhardt .
10.2. Phase-Space Tomography of Classical Fields
Phase-space tomography has also been applied for the spatial or temporal characterization of classical fields. In the temporal domain, it has been used for characterizing trains of ultrashort pulses by estimating the chronocyclic Wigner function [37,163]. Within the context of the spatial characterization of stationary fields, the first proposition of using phase-space tomography to recover the cross-spectral density of a paraxial optical field is due to Nugent , and the first experimental implementation of this approach was by McAlister et al. . This approach has been used for the characterization of x-ray sources [165,166] whose main variation is in one transverse coordinate. Tu and Tamura [167,168] realized that the same method can be cast in much simpler terms if, instead of the Wigner function, one employs the ambiguity function. The key to this approach is evident from the analogue of Eq. (140) for fields that depend only on one transverse variable, i.e.,169,170].
Hazak  and Gori et al.  realized that, in the more general case where the transverse cross-spectral density depends on both x and y, the information contained in the field’s intensity over a volume (a function of three variables) is not sufficient for recovering the initial cross-spectral density (a function of four variables). Within the approach of Tu and Tamura, this problem can be appreciated from the fact that Eq. (140) gives access to the ambiguity function only within a three-dimensional “slice” in which is parallel to . To address this problem, Raymer et al.  proposed introducing a fourth independent parameter in the form of the varying focal distance of a cylindrical lens inserted before the measurement of the intensity over a volume. The effect of such lenses is to give access to marginal projections of the Wigner function that correspond to the square modulus of different anamorphic fractional Fourier transforms. A system of this type has been used for beams where the field is separable in the two transverse coordinates . Recently, Rojas et al.  performed simultaneous tomographic measurements over the phase spaces of position–direction and time–frequency.
11. Generalized Wigner Functions: Achieving Exact Conservation
As discussed in the previous sections, the Wigner function allows a description of the propagation or evolution of optical wave fields in terms of a distribution that is conserved along the paths associated with simpler physical models. The condition for this description is that wave propagation in the system corresponds to a combination of linear canonical transformations and phase-space shifts. However, even in a situation as simple as free-space propagation beyond the paraxial regime, this conservation law is no longer satisfied exactly. In this section, we discuss modifications to the definition of the Wigner function for which the property of conservation along paths under propagation is preserved. The two cases that we concentrate on are nonparaxial field propagation and pulse propagation through general dispersive transparent media.
11.1. General Procedure for Constructing Conserved Wigner Functions
We now give a simple prescription for the definition of generalizations of the Wigner function that present the property of describing propagation/evolution exactly as a rearrangement of the arguments [176,177], which is valid whenever the solutions to the wave equations in question can be expressed as superpositions of exponentials of imaginary linear combinations of the spatial or spatiotemporal variables. We concentrate on two situations:
i. Nonparaxial propagation in free-space of a field without evanescent components. For the sake of illustration, let the field be monochromatic, fully coherent, scalar, and limited to two dimensions, so that it satisfies the Helmholtz equation:
ii. One-dimensional propagation of pulses through dispersive transparent media. According to Eq. (166), pulses traveling through dispersive media can be written as
Notice that, in both these situations, the field of interest is of the form189), , , , , , and , while for the pulse in Eq. (190), , , , , , and . The expression in Eq. (191) looks like a restricted Fourier synthesis, in the sense that the frequency vector g is restricted to a curve (or manifold, for higher dimensions). For the case of nonparaxial propagation, this curve is the unit circle (or the unit sphere in three dimensions) of plane-wave directions, while for pulse propagation it is the dispersion curve k versus ω. Another physical situation, not discussed here, whose solutions can be written in the form in Eq. (191) is that of free-space relativistic quantum mechanics, where the manifold described by g is a hyperboloid in the energy–momentum 4D space .
The goal is to define a generalized Wigner function for a field of the form in Eq. (191) that is conserved under propagation/evolution and whose marginal projection gives a desired bilinear property, e.g., the optical intensity. Here, ξ is a suitable variable, e.g., the angle θ with respect to the z axis of a ray passing through r in the case of nonparaxial fields, or the uniform velocity v of a light particle at position z and time t in the case of dispersive propagation. Mathematically, the conservation condition can be written as a first-order differential equation of the form192) indicates conservation along rays], and for dispersive pulse propagation [so that Eq. (192) implies the conservation of along uniform velocity trajectories]. The marginal relation, on the other hand, is given by
In order to find , the first step is to substitute Eq. (191) into the definition of the intensity (or the bilinear quantity of interest), leading to193), where 192) as long as the following relation holds: 24. Consider the curve spanned by for all values of η. Then, in order to evaluate at a given value of ξ, the integration in α must parametrize all pairs of points and that are joined by a line segment normal to . That is, and must be chosen as corresponding to all pairs of points along the curve spanned by g that are intersections with straight lines whose slopes are normal to . In other words, ξ regulates the slope of the line containing and while, for each ξ, the displacement of this line is regulated by α.
In the case of dispersion, the curve can present inflection points, in which a case there might be more than two intersections with a straight line. If this is the case, all possible pairs must be accounted for by the parametrization. On the other hand, for higher-dimensional problems where g spans not a curve but a surface or some higher-order manifold, the condition in Eq. (196) is not sufficient to find a unique change of variables, i.e., a unique definition of the generalized Wigner function with the desired properties. As discussed below, some other physically motivated criterion must be imposed in this case to remove this ambiguity.
11.2. Nonparaxial Field Propagation
For two-dimensional nonparaxial free propagation, spans a unit circle. Note that, in this case , so Eq. (196) takes the form25(a), for . The Jacobian of this change of variables equals unity, so the generalized Wigner function can be written as [178–181] 45] that it also satisfies an inner-product overlap property analogous to that in Eq. (55).
Analogous generalized Wigner functions can be defined in the three-dimensional case [182,183], where the manifold spanned by g is the surface of a unit sphere. Let this sphere be parametrized directly by the unit vector u [that is, ]. The change of variables is slightly more complicated in this case, since it is necessary to include one more integration variable in the definition of the Wigner function. The condition in Eq. (196) now becomes25(b). As in the 2D case, we choose α as the angle between and . The resulting generalized Wigner function then has the form
Wigner functions of this type have been used to describe the propagation not only of the intensity of scalar fields but also of the electric and magnetic energy densities and Poynting vector [184,185], polarization [185,186], and cross-spectral density  of nonparaxial electromagnetic fields. In the case of polarization, the Wigner functions are 3 × 3 matrices [analogous to the 2 × 2 matrix Wigner function in Eq. (155)] whose elements correspond to correlations of different Cartesian components of the field. Further, upon propagation across or reflection from interfaces between homogeneous media, these generalized Wigner functions transform following approximately a radiometric model: to leading order, they follow ray-optical laws of refraction and reflection and are multiplied by a Fresnel transmissivity or reflectivity, with corrections taking the form of derivatives of the leading term, which are less significant for more smoothly varying, less coherent fields [188,189]. Generalized Wigner functions have also been derived to model propagation through homogeneous anisotropic (e.g., birefringent and/or chiral) transparent media . For these media, g spans not one but two manifolds, one associated with each eigenpolarization of the medium. In general, these manifolds are ellipsoids, given the dependence of the refractive index on the propagation direction.
The generalized Wigner functions in Eqs. (198) and (200) are valid for fields that include components propagating in any direction. However, in the case of nonparaxial fields composed exclusively of plane waves traveling towards larger values of z (i.e., within the forward hemisphere of directions), these Wigner functions can be written in a form that makes the description of propagation formally equivalent to that for paraxial fields . To see this, consider the generalized Wigner function for two dimensions in Eq. (198). The field is assumed to be composed of only forward-propagating plane waves; i.e., the integral in Eq. (189) extends only from to . Therefore, the Wigner function in Eq. (198) vanishes for . Now perform a change in the directional parameter, so that instead of the angle θ, we specify direction by the slope . After also including a Jacobian factor , we find the expression26 for a truncated cylindrical wave, corresponding to being constant for and zero otherwise. The alternative representation in Eq. (202) allows, among other things, the definition of a generalized ambiguity function for nonparaxial fields, given by the Fourier transform over x and the inverse Fourier transform over τ of . Given the full analogy with the paraxial regime, the standard phase-space tomography techniques used in the paraxial 2D case can be extended to the nonparaxial regime .
The generalized Wigner functions in Eqs. (198) and (200) are formally equivalent to wave-based definitions of the radiance proposed by Ovchinnikov and Tatarskii  and Pedersen et al. [193–195]. The derivation followed by those authors, however, was quite different: first, the Wigner function over all three spatial variables is calculated, yielding a distribution over a six-dimensional phase space. This distribution turns out to be confined to the interior of a sphere in the Fourier variables, and to diverge at the surface of this sphere. Then, a radial projection of this distribution over the 3D Fourier variable is performed, with an appropriate weight factor, to eliminate the redundant variable corresponding to the magnitude of the direction vector. The resulting distribution presents the desired properties of conservation along rectilinear rays, as well as a directional marginal related to the field’s flux. However, to compute this distribution, one must first know the field in all space. A very different derivation for what is essentially the same distribution (up to an obliquity factor) was also proposed by Littlejohn and Winston . Wolf et al.  also proposed a Fourier-subspace radial projection for the case of two spatial dimensions, but by expressing the field as a weighted superposition of plane waves, they found the expression in Eq. (198). This result was then extended to fields in 3D space and with any state of coherence, leading to a family of generalized Wigner functions [179,182,183] whose marginals yield different local properties like the intensity, flux, and energy density. Sheppard and Larkin [180,181] found a more direct derivation in the two-dimensional case through a change of the angular variables to centroid and difference, consistent with the construction discussed here. The generalized Wigner functions in both two and three dimensions can also be calculated in terms of a series of derivatives acting on the standard Wigner function .
11.3. Pulse Propagation Through Transparent Dispersive Media
In the case of propagation through dispersive media, the generalized Wigner function behaves as a statistical–mechanical distribution of independent “light particles” traveling at different velocities v, since it is defined to obey the propagation and marginal properties [176,177]196). Since in this case and , this constraint takes the form
11.3a. Lorentz Model Dispersion
Consider first a nonabsorbing approximation to the Lorentz model of dispersion for a medium presenting a single resonance at , given by176], the change of variables that results from imposing Eq. (207) can be written as 207) for k given in Eq. (208). The corresponding Wigner function then takes the form Media 25 (Fig. 27) for the case where the pulse’s spectral amplitude is constant within the range and vanishes outside of it, for a medium where . The marginal projection of this distribution in v, equal to the spatial intensity profile, is also shown.
11.3b. Waveguide Dispersion
As a second example, consider the case of a perfect metallic hollow waveguide, whose dispersion is given by207) for this dispersion relation is given by Media 24 in Fig. 28 for a pulse with a spectrum constrained to the range , within which is constant.
In the examples we just presented, we plotted the Wigner function over the phase space of z versus v for different propagation times t. However, there are other ways to present these Wigner functions. For example, if the temporal intensity profiles are of more interest than the spatial ones, it might be more useful to use a variant of these Wigner functions in terms of t and a “slowness” variable , as was done in Ref. . This way, propagation in z corresponds to a horizontal linear shearing in the phase space . Alternatively, if it were preferable to express these Wigner functions in the chronocyclic phase space, one can substitute and include a Jacobian factor for normalization, so that the conserved chronocyclic Wigner function is given by199–201].
11.4. Reduction to the Standard Wigner Function for Parabolic Manifolds
Notice that, for the simple case in which g spans a parabola with any location and orientation, the generalized Wigner function reduces to a form of the standard Wigner function. For a general parabolic curve, this vector can be parametrized as196) can be written as
The reduction of the generalized Wigner function to the standard one for parabolic curves or manifolds explains the usefulness of the standard Wigner function in the case of paraxial free propagation and narrow-bandwidth pulses. In these cases, only a small segment of the manifold spanned by g is of interest, which can then be accurately approximated by a parabola: for paraxial optics we can approximate , while for narrow-band pulse propagation we can use either of the quadratic approximations in Eq. (168) or (176). The same is true in the case of quantum mechanics where, following the Klein–Gordon or Dirac equations, the manifold for free evolution of a particle is a (two-sheeted) hyperboloid of “radius” equal to the particle’s mass in the energy–momentum space. In the nonrelativistic limit, only the region near the (positive energy) vertex of this hyperboloid is of importance, so it can be approximated by a paraboloid. The standard Wigner function is therefore useful in the nonrelativistic limit.
12. Concluding Remarks
An overview of the use of the Wigner function in optics was presented, with an emphasis on its application to model field propagation. As discussed briefly in Section 6, there are many other phase-space distributions that are commonly used in physics and engineering, each of them with a unique set of properties that make it well suited for some purposes. However, there are two properties that, together, confer the Wigner function a particular connection to physical intuition. The first is that the Wigner function for an optical field after spatial propagation or temporal evolution corresponds to a simple mapping of arguments of the Wigner function for the initial field. This property is valid in cases where propagation or evolution correspond mathematically to an arbitrary combination of linear canonical transformations and phase-space shifts. Such cases include propagation of paraxial stationary fields through optical systems that are linear in the ray domain, propagation of time-dependent fields through transparent media with linear dispersion, and the temporal evolution of quantum states of light. This mapping of the Wigner function can be interpreted as conservation along paths associated with a simpler physical model, like ray optics, free particle dynamics, or classical field evolution. The second key property is that, at any point of propagation over such systems, the marginal projection of the Wigner function over the Fourier variables is associated with an observable, physically meaningful local property of the field. Given these two properties, the Wigner function’s behavior mimics that of the phase-space distributions used in the simpler physical models for these situations (e.g., the basic radiance from the theory of radiometry). However, the Wigner function incorporates the fundamental limitations imposed by the space–bandwidth product theorem, and this causes it to be potentially negative within some regions of phase space, and to have strict restrictions in its functional form.
By modifying the definition of the Wigner function according to a simple geometrically motivated prescription, the exact satisfaction of its two key properties can be extended at least for certain physical situations that are not described by linear canonical transformations and phase-space shifts. These situations include all cases where the wave equation accepts a complete set of solutions corresponding to purely imaginary exponentials of linear combinations of the variables. This is the case of nonparaxial field propagation (if evanescent waves are excluded) and of pulse propagation through homogeneous media with arbitrary dispersion (if absorption is neglected). This extension brings to mind the following question: what is the most general class of system and field for which a generalized Wigner function can be defined that satisfies exactly both the marginal relation to the intensity and the property of conservation along paths? This general question was not addressed in this article. Let us consider the case of propagation in general transparent inhomogeneous linear media. There are inhomogeneous refractive index distributions that allow a rigorous description of the propagation of intensity as a conserved weight distribution for rays, even beyond the paraxial regime. One trivial example of such case is that of inhomogeneous and anisotropic media that are mappings of homogeneous media through what is known as transformation optics , since the Wigner function would then be the result of the corresponding mapping of the Wigner function for a homogeneous medium. On the other hand, there are other distributions for which this type of representation is certainly impossible. Take the case of a gradient index waveguide that is symmetric with respect to the z axis but that is anharmonic, i.e., that present nonlinear mode dispersion. In the ray picture, propagation in z would correspond to a swirling of phase space, since different rays oscillate around the axis with different periods. In the wave domain, on the other hand, fields that are concentrated spatially and directionally far from the axis of the waveguide are known to return to their initial state after some propagation distance, due to the phenomenon known as a wave packet revival . Therefore, there seems to be an irreconcilable difference in the behavior of the rays and the waves.
In many cases, phase-space distributions like the Wigner function are valuable more for conceptual reasons than for computational ones. They provide an intuitive framework that leads to useful ideas and interpretations, even if the methods resulting from them can then be expressed in a more direct way that bypasses the phase-space picture. Computationally, the increased number of variables and complicated functional form of the Wigner function or its generalizations discussed here can be important obstacles towards numerical implementation, particularly in the case of coherent fields or pure states. For partially coherent fields or mixed states, on the other hand, the value of Wigner functions can be not only conceptual but also computational, because using standard propagation formulas can become computationally prohibitive in this case. If used properly, the Wigner function can allow significant computational savings in these cases, as discussed in Ref. .
This article is dedicated to the memory of Luis Edgar Vicent (1972–2011), friend, colleague and coauthor of Ref. .
References and Notes
1. N. L. Balazs and B. K. Jennings, “Wigner’s function and other distribution functions in mock phase spaces,” Phys. Rep. 104, 347‒ 391 (1984). [CrossRef]
2. E. P. Wigner, “On the quantum correction for thermodynamic equilibrium,” Phys. Rev. 40, 749‒ 759 (1932). [CrossRef]
3. J. Ville, “Thèorie et applications de la notion de signal analytique,” Cables Transm. 2A, 6174 (1948).
4. L. S. Dolin, “Beam description of weakly inhomogeneous wave fields,” Izv. Vyssh. Uchebn. Zaved. Radiofiz. 7, 559‒ 563 (1964).
5. A. Walther, “Radiometry and coherence,” J. Opt. Soc. Am. 58, 1256‒ 1259 (1968). [CrossRef]
6. R. W. Boyd, Radiometry and the Detection of Optical Radiation (Wiley, 1983), pp. 13‒ 27).
7. P. Moon and D. E. Spencer, The Photic Field (MIT Press, 1981).
8. L. Cohen, Time–Frequency Analysis (Prentice Hall, 1995).
9. W. Mecklenbräuker and F. Hlawatsch, The Wigner Distribution: Theory and Applications in Signal Processing (Elsevier, 1997).
10. H. W. Lee, “Theory and application of the quantum phase-space distribution functions,” Phys. Rep. 259, 147‒ 211 (1995). [CrossRef]
11. G. S. Agarwal and E. Wolf, “Calculus for functions of noncommuting operators and general phase-space methods in quantum mechanics. I. Mapping theorems and ordering of functions of noncommuting operators,” Phys. Rev. D 2, 2161‒ 2186 (1970). [CrossRef]
12. G. S. Agarwal and E. Wolf, “Calculus for functions of noncommuting operators and general phase-space methods in quantum mechanics. II. Quantum mechanics in phase space,” Phys. Rev. D 2, 2187‒ 2205 (1970). [CrossRef]
13. G. S. Agarwal and E. Wolf, “Calculus for functions of noncommuting operators and general phase-space methods in quantum mechanics. III. A generalized Wick theorem and multitime mapping,” Phys. Rev. D 2, 2206‒ 2225 (1970). [CrossRef]
14. M. Hillery, R. F. O’Connell, M. O. Scully, and E. P. Wigner, “Distribution functions in physics: Fundamentals,” 106 121‒ 167 (1984).
15. M. E. Testorf, B. M. Hennelly, and J. Ojeda-Castañeda, Phase-Space Optics: Fundamentals and Applications (McGraw-Hill, 2009).
16. A. Torre, Linear Ray and Wave Optics in Phase Space: Bridging Ray and Wave Optics via the Wigner Phase-Space Picture (Elsevier, 2005).
17. D. Dragoman, “The Wigner distribution function in optics and optoelectronics,” in Progress in Optics XXXVII, E. Wolf (Elsevier, 1997), pp. 1‒ 56.
18. J. T. Sheridan, W. T. Rhodes, and B. M. Hennelly, “Wigner Optics,” Proc. SPIE 5827, 627‒ 638 (2005).
19. A. Papoulis, The Fourier Integral and its Applications (McGraw-Hill, 1962).
20. G. Folland and A. Sitaram, “The uncertainty principle: a mathematical survey,” J. Fourier Anal. Appl. 3, 207‒ 238 (1997). [CrossRef]
21. E. U. Condon, “Immersion of the Fourier transform in a continuous group of functional transformations,” Proc. Natl. Acad. Sci. 23, 158‒ 163 (1937). [CrossRef]
22. V. Namias, “The fractional order Fourier transform and its application to quantum mechanics,” IMA J. Appl. Math. 25, 241‒ 265 (1980). [CrossRef]
23. A. W. Lohmann, D. Mendlovic, and Z. Zalevsky, E. Wolf, “Fractional transformations in optics,” in Progress in Optics XXXVIII1998), pp. 263‒ 242.
24. H. M. Ozaktas, Z. Zalevsky, and M. A. Kutay, The Fractional Fourier Transform with Applications in Optics and Signal Processing (John Wiley & Sons, 2001).
25. A. Erdélyi, Asymptotic Expansions (Dover, 1956).
26. D. Mendlovic, Y. Bitran, R. G. Dorsch, C. Ferreira, J. Garcia, and H. M. Ozaktaz, “Anamorphic fractional Fourier transform: optical implementation and applications,” Appl. Opt. 34, 7451‒ 7456 (1995). [CrossRef] [PubMed]
27. S. A. Collins, “Lens-system diffraction integral written in terms of matrix optics,” J. Opt. Soc. Am. 60, 1168‒ 1177 (1970). [CrossRef]
28. M. Moshinsky and C. Quesne, “Linear canonical transformations and their unitary representations,” J. Math. Phys. 12(8), 1772‒ 1783 (1971). [CrossRef]
29. K. B. Wolf, Integral Transforms in Science and Engineering (Plenum Press, 1979), Ch. 9, 10.
30. T. Alieva and M. J. Bastiaans, “Properties of the linear canonical integral transformation,” J. Opt. Soc. Am. A 24, 3658‒ 3665 (2007). [CrossRef]
31. D. Mustard, “The fractional Fourier transform and the Wigner distribution,” J. Aust. Math. Soc. B-Appl. Math. 38, 209‒ 219 (1996) Published earlier as Applied Mathematics Preprint AM89/6 School of Mathematics, UNSW, Sydney, Australia (1989). [CrossRef]
32. A. W. Lohmann, “Image rotation, Wigner rotation, and the fractional Fourier transform,” J. Opt. Soc. Am. A 10, 2181‒ 2186 (1993). [CrossRef]
33. A. W. Lohmann and B. H. Soffer, “Relationships between the Radon–Wigner and fractional Fourier transforms,” J. Opt. Soc. Am. A 11, 1798‒ 1801 (1994). [CrossRef]
34. See Chapter 4 by G. Saavedra and W. Furlan in Ref. , pp. 107–164.
35. J. Bertrand and P. Bertrand, “A tomographic approach to Wigner’s function,” Found. Phys. 17, 397‒ 405 (1987). [CrossRef]
36. D. F. McAlister, M. Beck, L. Clarke, A. Mayer, and M. G. Raymer, “Optical phase retrieval by phase space tomography and fractional-order Fourier transforms,” Opt. Lett. 20, 1181‒ 1183 (1995). [CrossRef] [PubMed]
37. M. Beck, M. G. Raymer, I. A. Walmsley, and V. Wong, “Chronocyclic tomography for measuring the amplitude and phase structure of optical pulses,” Opt. Lett. 18, 2041‒ 2043 (1993). [CrossRef] [PubMed]
38. D. T. Smithey, M. Beck, M. G. Raymer, and A. Faridani, “Measurement of the Wigner distribution and the density matrix of a light mode using optical homodyne tomography: application to squeezed states and the vacuum,” Phys. Rev. Lett. 70, 1244‒ 1247 (1993). [CrossRef] [PubMed]
39. U. Leonhardt, Measuring the Quantum State of Light (Cambridge U. Press, 1997).
40. C. Cohen-Tannoudji, B. Diu, and F. Laloe, Quantum Mechanics, Vol. 1 (Wiley, 1977), pp. 214‒ 227.
41. J. E. Moyal, “Quantum mechanics as a statistical theory,” Math. Proc. Camb. Phil. Soc. 45, 99‒ 124 (1949). [CrossRef]
42. R. G. Littlejohn and R. Winston, “Generalized radiance and measurement,” J. Opt. Soc. Am. A 12, 2736‒ 2743 (1995). [CrossRef]
45. M. A. Alonso, “Measurement of Helmholtz wave fields,” J. Opt. Soc. Am. A 17, 1256‒ 1264 (2000). [CrossRef]
46. M. V. Berry, “Semi-classical mechanics in phase space: a study of Wigners function,” Philos. Trans. R. Soc. Lond. 287, 237‒ 271 (1977). [CrossRef]
47. M. V. Berry and N. L. Balazs, “Evolution of semiclassical quantum states in phase space,” J. Phys. A Math. Phys. 12, 625‒ 642 (1979). [CrossRef]
48. M. V. Berry, “Quantum scars of classical closed orbits in phase space,” Proc. R. Soc. Lond. Ser. A 423, 219‒ 231 (1989). [CrossRef]
49. M. A. Alonso and G. W. Forbes, “Phase-space distributions for high-frequency fields,” J. Opt. Soc. Am. A 17, 2288‒ 2300 (2000). [CrossRef]
50. L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge University Press, 1995), Sec. 4.7.
51. M. J. Bastiaans, “Uncertainty principle for partially coherent light,” J. Opt. Soc. Am. 73, 251‒ 255 (1983). [CrossRef]
52. A. Starikov, “Effective number of degrees of freedom of partially coherent sources,” J. Opt. Soc. Am. 72, 1538‒ 1544 (1982). [CrossRef]
53. See Ref. , p. 261.
54. L. Cohen, “Time–frequency distributions—a review,” Proc. IEEE 77, 941‒ 981 (1989). [CrossRef]
55. J. H. Eberly and K. Wódkiewicz, “The time-dependent physical spectrum of light,” J. Opt. Soc. Am. 67, 1252‒ 1261 (1977). [CrossRef]
56. K. Husimi, “Some formal properties of the density matrix,” Proc. Phys. Math. Soc. Jpn. 22, 264‒ 314 (1940).
57. Y. Kano, “A new phase-space distribution function in the statistical theory of the electromagnetic field,” J. Math. Phys. 6, 1913‒ 1915 (1965). [CrossRef]
58. C. L. Mehta and E. C. G. Sudarshan, “Relation between Quantum and Semiclassical Description of Optical Coherence,” Phys. Rev. 138, B274‒ B280 (1965). [CrossRef]
59. R. J. Glauber, “Optical coherence and photon statistics,” in Quantum Optics and Electronics, C. Dewitt, A. Blandin, and C. Cohen-Tannoudji (Gordon and Breach, 1965), p. 65).
60. R. J. Glauber, “Coherent and incoherent states of the radiation field,” Phys. Rev. 131, 2766‒ 2788 (1963). [CrossRef]
61. J. R. Klauder, “Continuous representation theory. I. Postulates of continuous representation theory,” J. Math. Phys. 4, 1055‒ 1058 (1963). [CrossRef]
62. E. C. G. Sudarshan, “Equivalence of semiclassical and quantum mechanical descriptions of statistical light beams,” Phys. Rev. Lett. 10, 277‒ 279 (1963). [CrossRef]
63. J. G. Kirkwood, “Quantum statistics of almost classical ensembles,” Phys. Rev. 44, 31‒ 37 (1933). [CrossRef]
64. A. W. Rihaczek, “Signal energy distribution in time and frequency,” IEEE Trans. Info. Theory 14, 369‒ 374 (1968). [CrossRef]
65. H. Margenau and R. N. Hill, “Correlation between measurements in quantum theory,” Prog. Theoret. Phys. 26, 722‒ 738 (1961). [CrossRef]
66. P. M. Woodward, Probability and Information Theory with Applications to Radar (Pergamon, 1953).
67. A. Papoulis, “Ambiguity function in Fourier optics,” J. Opt. Soc. Am. 64, 779‒ 788 (1974). [CrossRef]
68. L.-P. Guigay, Ambiguity Function in Optical Imaging, Chapter 2 of Ref. , pp. 45–62.
69. L. Cohen, “Generalized phase-space distributions,” J. Math. Phys. 7, 781‒ 786 (1966). [CrossRef]
70. M. J. Bastiaans, Wigner Distribution in Optics, Chapter 1 of Ref. , pp. 1–44.
71. R. K. Luneburg, Mathematical Theory of Optics (University of California Press, 1966), pp. 246‒ 257.
72. See Ref. , pp. 103–110.
73. M. Born and E. Wolf, Principles of Optics, 7th Ed. (Cambridge University Press, 1999), pp. 142‒ 144.
74. H. A. Buchdahl, Hamiltonian Optics (Dover, 1993), pp. 7‒ 12.
75. See Ref. , pp. 13–25.
76. See Ref. , pp. 75–80.
77. A. Walther, “Lenses, wave optics and eikonal functions,” J. Opt. Soc. Am. 59, 1325‒ 1333 (1969). [CrossRef]
78. A. Walther, The Ray and Wave Theory of Lenses (Cambridge University Press, 1995), pp. 169‒ 187.
79. M. A. Alonso and G. W. Forbes, “Semigeometrical estimation of Green’s functions and wave propagators in optics,” J. Opt. Soc. Am. A 14, 1076‒ 1086 (1997). [CrossRef]
80. J. H. Van Vleck, “The correspondence principle in the statistical interpretation of quantum mechanics,” Proc. Natl. Acad. Sci. USA 14, 178‒ 188 (1928). [CrossRef]
81. M. C. Gutzwiller, “Periodic orbits and classical quantization conditions,” J. Math. Phys. 12, 343‒ 358 (1971). [CrossRef]
82. D. Mendlovic and H. M. Ozaktas, “Fractional Fourier transforms and their optical implementation: I,” J. Opt. Soc. Am. A 10, 1875‒ 1881 (1993). [CrossRef]
83. M. J. Bastiaans, “Wigner distribution function applied to optical signals and systems,” Opt. Commun. 25, 26‒ 30 (1978). [CrossRef]
84. M. J. Bastiaans, “Transport equations for the Wigner distribution function,” J. Mod. Opt. 26, 1265‒ 1272 (1979).
85. M. J. Bastiaans, “The Wigner distribution function and Hamilton’s characteristics of a geometric–optical system,” Opt. Commun. 30, 321‒ 326 (1979). [CrossRef]
86. M. J. Bastiaans, “Wigner distribution function and its application to first-order optics,” J. Opt. Soc. Am. 69, 1710‒ 1716 (1979). [CrossRef]
87. J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, 1988), pp. 101‒ 136.
88. K.-H. Brenner, A. W. Lohmann, and J. Ojeda-Castañeda, “The ambiguity function as a polar display of the OTF,” Opt. Commun. 44, 323‒ 326 (1983). [CrossRef]
90. S. B. Mehta and C. J. R. Sheppard, “Using the phase-space imager to analyze partially coherent imaging systems: bright-field, phase contrast, differential interference contrast, differential phase contrast, and spiral phase contrast,” J. Mod. Opt. 57, 718‒ 739 (2010). [CrossRef]
92. W. P. Schleich, J. P. Dahl, and S. Varró, “Wigner function for a free particle in two dimensions: a tale of interference,” Opt. Commun. 283, 786‒ 789 (2010). [CrossRef]
93. M. V. Berry and N. L. Balazs, “Nonspreading wave packets,” Am. J. Phys. 47, 264‒ 267 (1979). [CrossRef]
94. M. V. Berry and F. J. Wright, “Phase-space projection identities for diffraction catastrophes,” J. Phys. A. Math. Gen. 13, 149‒ 160 (1980). [CrossRef]
95. R.-P. Chen, H.-P. Zheng, and C.-Q. Dai, “Wigner distribution function of an Airy beam,” J. Opt. Soc. Am. A 28, 1307‒ 1311 (2011). [CrossRef]
97. G. Siviloglou, J. Broky, A. Dogariu, and D. N. Christodoulides, “Observation of accelerating Airy beams,” Phys. Rev. Lett. 99, 213901-1-213901-4 (2007). [CrossRef]
99. W. T. Welford, “Use of annular apertures to increase focal depth,” J. Opt. Soc. Am. 50, 749 (1960). [CrossRef]
100. C. J. R. Sheppard, D. K. Hamilton, and I. J. Cox, “Optical microscopy with extended depth of field,” Proc. R. Soc. Lond. A 387, 171‒ 186 (1983). [CrossRef]
103. Q. Yang, L. Liu, J. Sun, Y. Zhu, and W. Lu, “Analysis of optical systems with extended depth of field using the Wigner distribution function,” Appl. Opt. 45, 8586‒ 8595 (2006). [CrossRef] [PubMed]
104. M. Testorf and J. Ojeda-Castañeda, “Fractional Talbot effect: analysis in phase space,” J. Opt. Soc. Am. A 13, 119‒ 125 (1996). [CrossRef]
105. M. Testorf, “Designing Talbot array illuminators with phase-space optics,” J. Opt. Soc. Am. A 23, 187‒ 192 (2006). [CrossRef]
106. M. E. Testorf, Self-imaging in Phase SpaceChapter 9 of Ref. , pp. 279–307.
107. J. Lancis, E. E. Sicre, A. Pons, and G. Saavedra, “Achromatic white-light self-imaging phenomenon: an approach using the Wigner distribution function,” J. Mod. Opt. 42, 425‒ 434 (1995). [CrossRef]
108. A. Stern and B. Javidi, “Sampling in the light of Wigner distribution,” J. Opt. Soc. Am. A 21, 360‒ 366 (2004) errata, 21, 1602–1612 (2004). [CrossRef]
109. B. M. Hennelly, J. J. Healy, and J. T. Sheridan, Sampling and Phase Space, Chapter 10 of Ref. , pp. 309–336.
110. F. Gori, “Fresnel transform and sampling theorem,” Opt. Commun. 39, 293‒ 297 (1981). [CrossRef]
111. A. Stern and B. Javidi, “Improved-resolution digital holography using the generalized sampling theorem for locally band-limited fields,” J. Opt. Soc. Am. A 23, 1227‒ 1235 (2006). [CrossRef]
112. A. Stern and B. Javidi, “Space-bandwidth conditions for efficient phase-shifting digital holographic microscopy,” J. Opt. Soc. Am. A 25, 736‒ 741 (2008). [CrossRef]
113. K. B. Wolf and A. L. Rivera, “Holographic information in the Wigner function,” Opt. Commun. 144, 36‒ 42 (1997). [CrossRef]
114. A. W. Lohmann, M. E. Testorf, and J. Ojeda-Castañeda, “Holography and the Wigner function,” in The Art and Science of Holography: A Tribute to Emmett Leith and Yuri Denisyuk, H. J. Caulfield (SPIE Press, 2004), pp. 127‒ 144.
115. A. W. Lohmann, M. E. Testorf, J. Ojeda-Castañeda, and A. W. Lohmann, “The space–bandwidth product, applied to spatial filtering and to holography,” in Selected Papers on Phase-Space Optics (SPIE Press, 2006), pp. 11‒ 32.
118. M. Testorf, “Analysis of the moiré effect by use of the Wigner distribution function,” J. Opt. Soc. Am. A 17, 2536‒ 2542 (2000). [CrossRef]
119. N. Morelle, M. E. Testorf, N. Thirion, and M. Saillard, “Electromagnetic probing for target detection: rejection of surface clutter based on the Wigner distribution,” J. Opt. Soc. Am. A 26, 1178‒ 1186 (2009). [CrossRef]
120. L. E. Vicent and M. A. Alonso, “Generalized radiometry as a tool for the propagation of partially coherent fields,” Opt. Commun. 207, 101‒ 112 (2002). [CrossRef]
123. K. Duan and B. Lü, “Wigner-distribution-function matrix and its application to partially coherent vectorial nonparaxial beams,” J. Opt. Soc. Am. B 22, 1585‒ 1593 (2005). [CrossRef]
124. R. Castañeda and J. Garcia-Sucerquia, “Electromagnetic spatial coherence wavelets,” J. Opt. Soc. Am. A 23, 81‒ 90 (2006). [CrossRef]
125. A. Luis, “Spatial-angular Mueller matrices,” Opt. Commun. 263, 141‒ 146 (2006). [CrossRef]
126. A. Luis, “Ray picture of polarization and coherence in a Young interferometer,” J. Opt. Soc. Am. A 23, 2855‒ 2860 (2006). [CrossRef]
127. R. Castañeda, J. Carrasquilla, and J. Herrera, “Radiometric analysis of diffraction of quasi-homogeneous optical fields,” Opt. Commun. 273, 8‒ 20 (2007). [CrossRef]
129. M. A. Alonso, “Diffraction of paraxial partially coherent fields by planar obstacles in the Wigner representation,” J. Opt. Soc. Am. A 26, 1588‒ 1597 (2009). [CrossRef]
130. A. Wax and J. E. Thomas, “Measurement of smoothed Wigner phase-space distributions for small-angle scattering in a turbid medium,” J. Opt. Soc. Am. A 15, 1896‒ 1908 (1998). [CrossRef]
131. H. T. Yura, L. Thrane, and P. E. Andersen, “Closed-form solution for the Wigner phase-space distribution function for diffuse reflection and small-angle scattering in a random medium,” J. Opt. Soc. Am. A 17, 2464‒ 2474 (2000). [CrossRef]
132. C.-C. Cheng and M. G. Raymer, “Propagation of transverse optical coherence in random multiple-scattering media,” Phys. Rev. A 62, 023811 (2000). [CrossRef]
133. A. C. Fannjiang, “White-noise and geometrical optics limits of Wigner–Moyal equation for wave beams in turbulent media,” Commun. Math. Phys. 254, 289‒ 322 (2005). [CrossRef]
136. B. Hall, M. Lisak, D. Anderson, R. Fedele, and V. E. Semenov, “Statistical theory for incoherent light propagation in nonlinear media,” Phys. Rev. E 65, 035602R (2002). [CrossRef]
137. M. Lisak, L. Helczynski, and D. Anderson, “Relation between different formalisms describing partially incoherent wave propagation in nonlinear optical media,” Opt. Commun 220, 321‒ 323 (2003). [CrossRef]
138. H. Gao, L. Tian, B. Zhang, and G. Barbastathis, “Iterative nonlinear beam propagation using Hamiltonian ray tracing and Wigner distribution function,” Opt. Lett. 35, 4148‒ 4150 (2010). [CrossRef] [PubMed]
139. S. Abe and J. T. Sheridan, “Wigner optics in the metaxial regime,” Optik 114, 139‒ 141 (). [CrossRef]
140. Yu. A. Kravtsov and L. A. Apresyan, “Radiative transfer: new aspects of the old theory,” in Progress in OpticsVol. XXXVI, E. Wolf (North Holland, 1996), pp. 179‒ 244.
141. L. A. Apresyan and Yu. A. Kravtsov, Radiation Transfer: Statistical and Wave Aspects (Gordon and Breach, 1996).
142. A. T. Friberg, Selected Papers on Coherence and Radiometry, Milestone Series Vol. MS69 (SPIE Optical Engineering Press, 1993).
143. E. Wolf, “Coherence and radiometry,” J. Opt. Soc. Am. 68, 6‒ 17 (1978). [CrossRef]
144. A. T. Friberg, “Effects of coherence in radiometry,” Opt. Eng. 21, 927‒ 936 (1982).
145. G. S. Agarwal, J. T. Foley, and E. Wolf, “The radiance and phase-space representations of the cross-spectral density operator,” Opt. Commun. 62, 67‒ 72 (1987). [CrossRef]
146. A. Walther, “Radiometry and coherence,” J. Opt. Soc. Am. 63, 1622‒ 1623 (1973). [CrossRef]
147. A. T. Friberg, “On the generalized radiance associated with radiation from a quasihomogeneous planar source,” Opt. Acta 28, 261‒ 277 (1981). [CrossRef]
148. J. T. Foley and E. Wolf, “Radiometry as a short-wavelength limit of statistical wave theory with globally incoherent sources,” Opt. Commun. 55, 236‒ 241 (1985). [CrossRef]
149. J. Azaña, “Time–frequency (Wigner) analysis of linear and nonlinear pulse propagation in optical fibers,” EURASIP J. Appl. Signal Process. 2005, 1554‒ 1565 (2005). [CrossRef]
150. P. Loughlin and L. Cohen, “A Wigner approximation method for wave propagation,” J. Acoust. Soc. Am. 118, 1268‒ 1271 (2005). [CrossRef]
151. J. Ojeda-Castañeda, J. Lancis, C. M. Gómez-Sarabia, V. Torres-Company, and P. Andrés, “Ambiguity function analysis of pulse train propagation: applications to temporal Lau filtering,” J. Opt. Soc. Am. A 24, 2268‒ 2273 (2007). [CrossRef]
152. P. Loughlin and L. Cohen, “Approximate wave function from approximate non-representable Wigner distributions,” J. Mod. Opt. 55, 3379‒ 3387 (2008). [CrossRef]
153. L. Cohen, P. Loughlin, and G. Okopal, “Exact and approximate moments of a propagating pulse,” J. Mod. Opt. 55, 3349‒ 3358 (2008). [CrossRef]
154. C. Dorrer and I. A. Walmsley, Phase space in ultrafast optics, Chapter 11 of Ref. , pp. 337–383.
155. I. A. Walmsley and C. Dorrer, “Characterization of ultrashort electromagnetic pulses,” Adv. Opt. Phot. 1, 308‒ 437 (2009). [CrossRef]
156. B. H. Kolner, “Space–time duality and the theory of temporal imaging,” IEEE J. Quantum Electron. 30, 1951‒ 1963 (1994). [CrossRef]
157. B. E. A. Saleh and M. C. Teich, Fundamentals of Photonics, 2nd ed. (Wiley, 2007), p. 188.
158. H. Lajunen, J. Tervo, J. Turunen, P. Vahimaa, and F. Wyrowski, “Spectral coherence properties of temporally modulated stationary light sources,” Opt. Express 11, 1894‒ 1899 (2003). [CrossRef] [PubMed]
159. B. J. Davis, “Observable coherence theory for statistically periodic fields,” Phys. Rev. A 76, 043843 (2007). [CrossRef]
160. K. Vogel and H. Risken, “Determination of quasiprobability distributions in terms of probability distributions for the rotated quadrature phase,” Phys. Rev. A 40, 2847‒ 2849 (1989). [CrossRef] [PubMed]
161. C. C. Gerry and P. L. Knight, Introductory Quantum Optics (Cambridge University Press, 2005), pp. 56‒ 71.
162. W. P. Schleich, Quantum Optics in Phase Space (Wiley-VCH, 2001).
165. C. Q. Tran, A. G. Peele, A. Roberts, K. A. Nugent, D. Paterson, and I. McNulty, “X-ray imaging: a generalized approach using phase-space tomography,” J. Opt. Soc. Am. A 22, 1691‒ 1700 (2005). [CrossRef]
166. C. Q. Tran, A. P. Mancuso, B. B. Dhal, K. A. Nugent, A. G. Peele, Z. Cai, and D. Paterson, “phase space reconstruction of focused x-ray fields,” J. Opt. Soc. Am. A 23, 1779‒ 1786 (2006). [CrossRef]
167. J. Tu and S. Tamura, “Wave field determination using tomography of the ambiguity function,” Phys. Rev. E 55, 1946‒ 1949 (1997). [CrossRef]
168. J. Tu and S. Tamura, “Analytic relation for recovering the mutual intensity by means of intensity information,” J. Opt. Soc. Am. A 15, 202‒ 206 (1998). [CrossRef]
170. “Partially coherent ambiguity functions for depth-variant point spread function design,” presentation during PIERS, Marrakesh, March, 2011.
172. F. Gori, M. Santarsiero, and G. Guattari, “Coherence and the spatial distribution of intensity,” J. Opt. Soc. Am. A 10, 673‒ 679 (1993). [CrossRef]
174. A. Cámara, T. Alieva, J. A. Rodrigo, and M. L. Calvo, “Phase space tomography reconstruction of the Wigner distribution for optical beams separable in Cartesian coordinates,” J. Opt. Soc. Am. A 26, 1301‒ 1306 (2009). [CrossRef]
176. J. C. Petruccelli and M. A. Alonso, “Phase space distributions tailored for dispersive media,” J. Opt. Soc. Am. A 27, 1194‒ 1201 (2010). [CrossRef]
177. Generalized Wigner Functions, Ph.D. Thesis, University of Rochester (Rochester, NY, 2010).
178. K. B. Wolf, M. A. Alonso, and G. W. Forbes, “Wigner functions for Helmholtz wave fields,” J. Opt. Soc. Am. A 16, 2476‒ 2487 (1999). [CrossRef]
179. M. A. Alonso, “Radiometry and wide-angle wave fields. I. Coherent fields in two dimensions,” J. Opt. Soc. Am. A 18, 902‒ 909 (2001). [CrossRef]
180. C. J. R. Sheppard and K. G. Larkin, “Wigner function for nonparaxial wave fields,” J. Opt. Soc. Am. A 18, 2486‒ 2490 (2001). [CrossRef]
182. M. A. Alonso, “Radiometry and wide-angle wave fields. II. Coherent fields in three dimensions,” J. Opt. Soc. Am. A 18, 910‒ 918 (2001). [CrossRef]
183. M. A. Alonso, “Radiometry and wide-angle wave fields. III. Partial coherence,” J. Opt. Soc. Am. A. 18, 2502‒ 2511 (2001). [CrossRef]
185. M. A. Alonso, “Wigner functions for nonparaxial, arbitrarily polarized electromagnetic wave fields in free-space,” J. Opt. Soc. Am. A. 21, 2233‒ 2243 (2004). [CrossRef]
186. J. C. Petruccelli, N, J. Moore, and M. A. Alonso, “Two methods for modeling the propagation of the coherence and polarization properties of nonparaxial fields,” Opt. Commun. 283, 4457‒ 4466 (2010). [CrossRef]
187. J. C. Petruccelli and M. A. Alonso, “Ray-based propagation of the cross-spectral density,” J. Opt. Soc. Am. A 25, 1395‒ 1405 (2008). [CrossRef]
188. J. C. Petruccelli and M. A. Alonso, “Propagation of partially coherent fields through planar dielectric boundaries using angle-impact Wigner functions I. Two dimensions,” J. Opt. Soc. Am. A 24, 2590‒ 2603 (2007). [CrossRef]
189. J. C. Petruccelli and M. A. Alonso, “Propagation of nonparaxial partially coherent fields across interfaces using generalized radiometry,” J. Opt. Soc. Am. A 26, 2012‒ 2022 (2009). [CrossRef]
190. J. C. Petruccelli and M. A. Alonso, “Generalized radiometry model for the propagation of light within anisotropic and chiral media,” J. Opt. Soc. Am. A 28, 791‒ 800 (2011). [CrossRef]
191. S. Cho and M. A. Alonso, “Ambiguity function and phase-space tomography for nonparaxial fields,” J. Opt. Soc. Am. A 28, 897‒ 902 (2011). [CrossRef]
192. G. I. Ovchinnikov and V. I. Tatarskii, “On the problem of the relationship between coherence theory and the radiation-transfer equation,” Radiophys. Quantum Electron. 15, 1087‒ 1089 (1972). [CrossRef]
193. H. M. Pedersen, J. H. Eberly, L. Mandel, and E. Wolf, “Exact geometrical description of free space radiative energy transfer for scalar wavefields,” in Coherence and Quantum Optics VI (Plenum, 1990), pp. 883‒ 887.
194. H. M. Pedersen, “Exact theory of free-space radiative energy transfer,” J. Opt. Soc. Am. A 8, 176‒ 185 (1991) errata, 8, 1518 (1991). [CrossRef]
195. H. M. Pedersen, “Geometric theory of fields radiated from three-dimensional, quasi-homogeneous sources,” J. Opt. Soc. Am. A 9, 1626‒ 1632 (1992). [CrossRef]
196. R. G. Littlejohn and R. Winston, “Corrections to classical radiometry,” J. Opt. Soc. Am. A 10, 2024‒ 2037 (1993). [CrossRef]
197. S. Cho, J. C. Petruccelli, and M. A. Alonso, “Wigner functions for paraxial and nonparaxial fields,” J. Mod. Opt. 56, 1843‒ 1852 (2009). [CrossRef]
198. M. A. Alonso, T. Setälä, and A. T. Friberg, “Optimal pulses for arbitrary dispersive media,” J. Eur. Opt. Soc. R.P. 6, 1100 (2011).
199. A. Papandreou-Suppappola, F. Hlawatsch, and G. Boudreaux-Bartels, “Quadratic time–frequency representations with scale covariance and generalized time-shift covariance: A unified framework for the affine, hyperbolic, and power classes,” Digital Signal Processing 8, 3‒ 48 (1998). [CrossRef]
200. F. Hlawatsch, A. Papandreou-Suppappola, and G. Boudreaux-Bartels, “The power classes-quadratic time–frequency representations with scale covariance and dispersive time-shift covariance,” IEEE Trans. Signal Process. 47, 3067‒ 3083 (1999). [CrossRef]
201. A. Papandreou-Suppappola, R. L. Murray, B.-G. Iem, and G. F. Boudreaux-Bartels, “Group delay shift covariant quadratic time–frequency representations,” IEEE Trans. Signal Process. 49, 2549‒ 2564 (2001). [CrossRef]
203. R. W. Robinett, “Quantum wave packet revivals,” Phys. Rep. 392, 1‒ 119 (2004). [CrossRef]
Miguel A. Alonso received the degree of Engineer in Physics in 1990 from the Universidad Autonoma Metropolitana in Mexico City and a Ph.D. in Optics in 1996 from The Institute of Optics, University of Rochester. After postdoctoral work at Macquarie University in Sydney, Australia, and a faculty position at the National Autonomous University of Mexico (UNAM) in Cuernavaca, Mexico, he joined the faculty of The Institute of Optics in 2003, where he is currently an Associate Professor. His main research interest is in finding new mathematical models to describe the propagation of waves through a variety of systems. He is a Deputy Editor of Optics Express, Chair of Spotlight on Optics, and a fellow of the Optical Society of America.