1 Introduction

There are at least three logically independent descriptions of quantum mechanics: the Hilbert space formalism [31], the path-integral method [49], and the phase-space approach such as given by the Wigner function [24, 32, 52, 69, 75, 84, 109, 111, 125]. The phase-space formulation of quantum mechanics was initiated by Wigner in his ground-breaking work [123] from 1932, in which the Wigner function of a spinless non-relativistic quantum particle was introduced as a quasi-probability distribution. The Wigner function can be used to express quantum-mechanical expectation values as classical phase-space averages. More than a decade later, Groenewold [63] and Moyal [99] formulated quantum mechanics as a statistical theory on a classical phase space by mapping a quantum state to its Wigner function and they interpreted this correspondence as the inverse of the Weyl quantization [119,120,121].

Coherent states have become a natural way to extend phase spaces to more general physical systems [5, 8,9,10,11, 13, 21, 54, 100]. In this regard, a new focus on phase-space representations for coupled, finite-dimensional quantum systems (as spin systems) [53, 76,77,78,79,80,81,82, 87, 106, 108, 115] and their tomographic reconstructions [81, 85, 86, 107] has emerged recently. A spherical phase-space representation of a single, finite-dimensional quantum system has been used to naturally recover the infinite-dimensional phase space in the large-spin limit [78, 81]. These spherical phase spaces have been defined in terms of quantum-mechanical expectation values of rotated parity operators [77, 78, 81, 87, 108, 115] (as discussed below) in analogy with displaced reflection operators in flat phase spaces. But in the current work, we exclusively focus on the (usual) infinite-dimensional case which has Heisenberg–Weyl symmetries [21, 54, 90, 100]. This case has been playing a crucial role in characterizing the quantum theory of light [59] via coherent states and displacement operators [3, 4, 22, 23] and has also been widely used in the context of time–frequency analysis and pseudo-differential operators [16,17,18, 28, 29, 43, 44, 62]. Many particular phase spaces have been unified under the concept of the so-called Cohen class [28, 29, 44] (see Definition 2), i.e., all functions which are related to the Wigner function via a convolution with a distribution (which is also known as the Cohen kernel).

Phase-space distribution functions are mostly described by one of the following three forms: (a) convolved derivatives of the Wigner function [43, 44], (b) integral transformations of a pure state (i.e., a rapidly decaying, complex-valued function) [16,17,18, 28, 29, 43, 44, 123], or (c) as integral transformations of quantum-mechanical expectation values [3, 4, 22, 23]. Also, Wigner functions (and the corresponding Weyl quantization) are usually described by integral transformations. But the seminal work of Grossmann [44, 64] (refer also to [102]) allowed for a direct interpretation of the Wigner function as a quantum-mechanical expectation value of a displaced parity operator \(\mathrm {\Pi }\) (which reflects coordinates \(\mathrm {\Pi }\, \psi (x) = \psi (-x)\) of a quantum state \(\psi \)). In particular, Grossmann [64] showed that the Weyl quantization of the delta distribution determines the parity operator \(\mathrm {\Pi }\). This approach has been widely adopted [15, 25, 33, 36, 52, 91, 103, 104].

However, parity operators similar to the one by Grossmann and Royer [44, 64, 102] have still been lacking for general phase-space distribution functions. (Note that such a form appeared implicitly for s-parametrized distribution functions in [22, 98].) In the current work, we generalize the previously discussed parity operator \(\mathrm {\Pi }\) [44, 64, 102] for the Wigner function by introducing a family of parity operators \(\mathrm {\Pi }_\theta \) (refer to Definition 3) which is parametrized by a function or distribution \(\theta \). This enables us to specify general phase-space distribution functions in the form of quantum-mechanical expectation values (refer to Definition 4) as

$$\begin{aligned} F_\rho (\mathrm {\Omega },\theta ) := (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_\theta \mathcal {D}^\dagger (\mathrm {\Omega }) ]. \end{aligned}$$

We will refer to the above operator \(\mathrm {\Pi }_\theta \) as a parity operator following the lead of Grossmann and Royer [64, 102] and given its resemblance and close analogy to the reflection operator \(\mathrm {\Pi }\) discussed in prior work [15, 77, 78, 81, 87, 108, 115]. Here, \(\mathcal {D}(\mathrm {\Omega })\) denotes the displacement operator and \(\Omega \) describes suitable phase-space coordinates (see Sect. 3.1). (Recall that \(\hbar = h / (2\pi )\) is defined as the Planck constant h divided by \(2\pi \).) The quantum-mechanical expectation values in the preceding equation give rise to a rich family of phase-space distribution functions \(F_\rho (\mathrm {\Omega },\theta )\) which represent arbitrary (mixed) quantum states as given by their density operator \(\rho \). In particular, this family of phase-space representations contains all elements from the (above-mentioned) Cohen class and naturally includes the pivotal Husimi Q and Born–Jordan distribution functions.

We would like to emphasize that our approach to phase-space representations averts the use of integral transformations, Fourier transforms, or convolutions as these are subsumed in the parity operator \(\mathrm {\Pi }_\theta \) which is independent of the phase-space coordinate \(\Omega \). Although our definition also relies on an integral transformation given by a Fourier transform, it is only applied once and is completely absorbed into the definition of a parity operator, thereby avoiding redundant applications of Fourier transforms. This leads to significant advantages as compared to earlier approaches:

  • conceptual advantages (see also [76, 81, 98, 108, 115]):

    • The phase-space distribution function is given as a quantum-mechanical expectation value. This form nicely fits with the experimental reconstruction of quantum states [7, 14, 46, 68, 81, 93, 107].

    • All the complexity from integral transformations (etc.) is condensed into the parity operator \(\mathrm {\Pi }_\theta \).

    • The dependence on the distribution \(\theta \) and the particular phase space is separated from the displacement \(\mathcal {D}(\mathrm {\Omega })\).

  • computational advantages:

    • The repeated and expensive computation of integral transformations (etc.) in earlier approaches is avoided as \(\mathrm {\Pi }_\theta \) has to be determined only once. Also, the effect of the displacement \(\mathcal {D}(\mathrm {\Omega })\) is relatively easy to calculate.

In this regard, the current work can also be seen as a continuation of [81] where the parity operator approach has been emphasized, but mostly for finite-dimensional quantum systems. Moreover, we connect results from quantum optics [22, 23, 59, 88], quantum harmonic analysis [30, 37,38,39,40, 43, 44, 74, 118], and group-theoretical approaches [21, 54, 90, 100]. It is also our aim to narrow the gap between different communities where phase-space methods have been successfully applied.

On the other hand, a major contribution of our work is the analysis of existence properties of generalized phase-space distributions and their parity operators. While the Wigner function has been known to exist for the general class of tempered distributions (a class of generalized functions that includes the pivotal \(L^2\) space), we further illuminate which classes of Cohen kernels yield well-defined generalized phase-space distribution functions. Such existence questions are fully absorbed into the parity operators and precise conditions are used to guarantee their mathematical existence.

Similar to the parity operator \(\mathrm {\Pi }\) (which is the Weyl quantization of the delta distribution), we show that its generalizations \(\mathrm {\Pi }_\theta \) are Weyl quantizations of the corresponding Cohen kernel \(\theta \) (refer to Sect. 4.3 for the precise definition of the Weyl quantization used in this work). We discuss how these general results reduce to well-known special cases, and discuss properties of phase-space distributions in relation to their parity operators \(\mathrm {\Pi }_\theta \). In particular, we consider the class of s-parametrized distribution functions [22, 23, 59, 98], which include the Wigner, Glauber P, and Husimi Q functions, as well as the \(\tau \)-parametrized family which has been proposed in the context of time–frequency analysis and pseudo-differential operators [16,17,18, 43]. We derive spectral decompositions of parity operators for all of these phase-space families, including the Born–Jordan distribution. Relations of the form \(\mathrm {\Pi }_{\theta }=A_\theta \circ \mathrm {\Pi }\) motivate the name “parity operator” as they are in fact compositions of the usual parity operator \(\mathrm {\Pi }\) followed by some operator \(A_\theta \) that usually corresponds to a geometric or physical operation (which commutes with \(\mathrm {\Pi }\)). In particular, \(A_\theta \) is a squeezing operator for the \(\tau \)-parametrized family and corresponds to photon loss for the s-parametrized family (assuming \(s<0\)). This structure of the parity operators \(\mathrm {\Pi }_{\theta }\) connects phase spaces to elementary geometric and physical operations (such as reflection, squeezing operators, photon loss), and these concepts are central to applications: The squeezing operator models a nonlinear optical process which generates non-classical states of light in quantum optics [60, 88, 94]. These squeezed states of light have been widely used in precision interferometry [61, 96, 110, 124] or for enhancing the performance of imaging [92, 116]; also, the gravitational-wave detector GEO600 has been operating with squeezed light since 2010 [1, 65].

The Born–Jordan distribution and its parity operator constitute a most peculiar instance among the phase-space approaches. This distribution function has convenient properties, e.g., it satisfies the marginal conditions and therefore allows for a probabilistic interpretation [43]. The Born–Jordan distribution is, however, difficult to compute. But most importantly, the Born–Jordan distribution and its corresponding quantization scheme have a fundamental importance in quantum mechanics. In particular, there have been several attempts in the literature to find the “right” quantization rule for observables using either algebraic or analytical techniques. In a recent paper [42], one of us has analyzed the Heisenberg and Schrödinger pictures of quantum mechanics, and it is shown that the equivalence of both theories requires that one must use the Born–Jordan quantization rule (as proposed by Born and Jordan [20])

$$\begin{aligned} \text {(BJ)} \quad x^{m}p^{\ell } \mapsto \frac{1}{m{+}1} \sum _{k=0}^{m} \hat{x}^{k} \hat{p}^{\ell } \hat{x}^{m-k}, \end{aligned}$$

instead of the Weyl rule

$$\begin{aligned} \text {(Weyl) } \quad x^{m} p^{\ell } \mapsto \frac{1}{2^{m}} \sum _{k=0}^{m} \left( {\begin{array}{c}m\\ k\end{array}}\right) \hat{x}^{k} \hat{p}^{\ell } \hat{x}^{m-k} \end{aligned}$$

for monomial observables. The Born–Jordan and Weyl rules yield the same result only if \(m<2\) or \(\ell <2\); for instance, in both cases the quantization of the product xp is \(\frac{1}{2}(\hat{x}\hat{p} + \hat{p}\hat{x})\). It is, however, easy to find physical examples which give different results. Consider, for instance, the square of the z component of the angular momentum: It is given by

$$\begin{aligned} \ell _{z}^{2}=x^{2}p_{y}^{2}+y^{2}p_{x}^{2}-2xp_{x}yp_{y} \end{aligned}$$

and its Weyl quantization is easily seen to be

$$\begin{aligned} \text {Op}_{\text {Weyl}}(\ell _{z}^{2}) =\hat{x}_{x}^{2} \hat{p}_{y}^{2} + \hat{x}_{y}^{2} \hat{p}_{x}^{2} - \tfrac{1}{2}(\hat{x}_{x} \hat{p}_{x} {+} \hat{p}_{x} \hat{x}_{x} ) (\hat{x}_{y} \hat{p}_{y} {+} \hat{p}_{y} \hat{x}_{y}) \end{aligned}$$
(1)

while its Born–Jordan quantization is the different expression

$$\begin{aligned} \text {Op}_{\text {BJ}}(\ell _{z}^{2}) =\hat{x}_{x}^{2} \hat{p}_{y}^{2} + \hat{x}_{y}^{2} \hat{p}_{x}^{2} - \tfrac{1}{2}(\hat{x}_{x} \hat{p}_{x} {+} \hat{p}_{x} \hat{x}_{x} ) ( \hat{x}_{y} \hat{p}_{y} {+} \hat{p}_{y} \hat{x}_{y} ) - \tfrac{1}{6}\hbar ^{2}. \end{aligned}$$
(2)

(Recall that the operators \(\hat{x}_\eta \) and \(\hat{p}_\kappa \) satisfy the canonical commutation relations \([\hat{x}_\eta ,\hat{p}_\kappa ] = i\hbar \delta _{\eta \kappa }\) using the spatial coordinates \(\eta ,\kappa \in \{x,y,z\}\) and the Kronecker delta \(\delta _{\eta \kappa }\).) One of us has shown in [45] that the use of (2) instead of (1) solves the so-called angular momentum dilemma [34, 35].

To a general observable a(xp), the Weyl rule associates the operator

$$\begin{aligned} \text {Op}_{\text {Weyl}}(a) =(2\pi \hbar )^{-1}\int \mathcal {F}_{\sigma }a(x,p) \mathcal {D}(x,p) \, \textrm{d}x \, \textrm{d}p \end{aligned}$$

where \(\mathcal {F}_{\sigma }a\) is the symplectic Fourier transform of a and \(\mathcal {D}(x,p)\) the displacement operator (see Sect. 3.1); in the Born–Jordan case, this expression is replaced with

$$\begin{aligned} \text {Op}_{\text {BJ}}(a) =(2\pi \hbar )^{-1}\int \mathcal {F}_{\sigma }a(x,p) K_{\text {BJ}}(x,p) \mathcal {D}(x,p) \, \textrm{d}x \, \textrm{d}p \end{aligned}$$

where the filter function \(K_{\text {BJ}}(x,p)\) is given by

$$\begin{aligned} K_{\text {BJ}}(x,p)=\textrm{sinc}( \tfrac{ px }{2\hbar } ) = \frac{\sin (\tfrac{px}{2\hbar })}{px/(2 \hbar )}. \end{aligned}$$

We obtain significant, new results for the case of Born–Jordan distributions and therefore substantially advance on previous characterizations. In particular, we derive its parity operator \(\mathrm {\Pi }_{\text {BJ}}\) in the form of a weighted average of geometric transformations

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} = \tfrac{1}{4 \pi \hbar } \int \textrm{sinc}( \tfrac{ px }{2\hbar } ) \, \mathcal {D}(x,p) \, \textrm{d}x \, \textrm{d}p = \left[ \, \tfrac{1}{4} \int _{-\infty }^{\infty } \, \textrm{sech}(\tfrac{\xi }{2}) S(\xi ) \,\textrm{d} \xi \,\right] \; \mathrm {\Pi }, \end{aligned}$$
(3)

where \(\mathcal {D}(x,p)\) is the displacement operator and \(S(\xi )\) is the squeezing operator (see Eq. (46)) with a real squeezing parameter \(\xi \). We have used the sinus cardinalis \(\textrm{sinc}(x):= \textrm{sin}(x)/x\) and the hyperbolic secant \( \textrm{sech}(x):=1/\textrm{cosh}(x)\) functions. The parity operator \(\mathrm {\Pi }_{\text {BJ}}\) in Eq. (3) decomposes into a product \(\mathrm {\Pi }_{\theta }=A_\theta \circ \mathrm {\Pi }\) containing the usual reflection operator \(\mathrm {\Pi }\). This is another example of the above-discussed motivation for our terminology of parity operators. We prove in Proposition 2 that \(\mathrm {\Pi }_{\text {BJ}}\) is a bounded operator on the Hilbert space of square-integrable functions and therefore gives rise to well-defined phase-space distribution functions of arbitrary quantum states. We derive a generalized spectral decomposition of this parity operator based on a continuous family of generalized eigenvectors that satisfy the following generalized eigenvalue equation for every real E (see Theorem 5):

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \, | \psi ^E_{\pm } \rangle = \pm \pi /2 \, \textrm{sech}(\pi E) \, | \psi ^E_{\pm } \rangle . \end{aligned}$$

Facilitating a more efficient computation of the Born–Jordan distribution, we finally derive explicit matrix representations in the so-called Fock or number-state basis, which constitutes a natural representation for bosonic quantum systems such as in quantum optics [60, 88, 94]. In this case, the parity operator \(\mathrm {\Pi }_{\text {BJ}}\) of the Born–Jordan distribution is not diagonal in the Fock basis—as compared to the diagonal parity operators of s-parametrized phase spaces (cf. [81]) that enable the experimental reconstruction of distribution functions from photon-count statistics [7, 14, 46, 93] in quantum optics. We calculate the matrix elements \([\mathrm {\Pi }_{\text {BJ}}]_{mn}\) in the Fock or number-state basis and provide a convenient formula for a direct recursion, for which we conjecture that the matrix elements are completely determined by eight rational initial values. This recursion scheme has significant computational advantages for calculating Born–Jordan distribution functions as compared to previous approaches and allows for an efficient implementation. In particular, large matrix representations of the parity operator \(\mathrm {\Pi }_{\text {BJ}}\) can be well approximated using rank-9 matrices. We finally illustrate our results for simple quantum states by calculating their Born–Jordan distributions and comparing them to other phase-space representations. Let us summarize the main results of the current work:

  • quantum-mechanical expectation values of the parity operators \(\mathrm {\Pi }_\theta \) from Definition 3 define distribution functions (see Definition 4) and form the Cohen class (Theorem 1);

  • existence properties of parity operators and generalized phase-space functions are clarified in Sect. 4. We refer in particular to the crucial Lemma 2;

  • the parity operators \(\mathrm {\Pi }_\theta \) are Weyl quantizations of the corresponding Cohen convolution kernels \(\theta \) (Sect. 4.3);

  • parity operators for important distribution functions are summarized in Sect. 4.4 along with their operator norms (Theorem 2) and generalized spectral decompositions in Sect. 5.2;

  • the Born–Jordan parity operator is a weighted average of displacements (Theorem 3) or, equivalently, a weighted average of squeezing operators (Theorem 4), and it is bounded (Proposition 2);

  • the Born–Jordan parity operator admits a generalized spectral decomposition (Theorem 5);

  • its matrix representation is calculated in the number-state basis in Theorem 6; and an efficient, recursion-based computation scheme is proposed in Conjecture 1.

Our work has significant implications: General (infinite-dimensional) phase-space functions can now be conveniently and effectively described as natural expectation values. We provide a much more comprehensive understanding of Born–Jordan phase spaces and means for effectively computing the corresponding phase-space functions. Working in a rigorous mathematical framework, we also facilitate future discussions of phase spaces by connecting different communities in physics and mathematics.

We start by recalling precise definitions of distribution functions and quantum states for infinite-dimensional Hilbert spaces in Sect. 2. In Sect. 3, we discuss phase-space translations of quantum states using coherent states, recall one known formulation of translated parity operators, and relate a general class of phase spaces to Wigner distribution functions and their properties. We note that an experienced reader can skip most of the introductory Sects. 2 and 3 and jump directly to our results. These preparations will, however, guide our study of phase-space representations of quantum states as expectation values of displaced parity operators in Sect. 4. We present and discuss our results for the case of the Born–Jordan distribution and its parity operator in Sect. 5. Formulas for the matrix elements of the Born–Jordan parity operator are derived in Sect. 6. Explicit examples for simple quantum systems are discussed and visualized in Sect. 7, before we conclude. A larger part of the proofs are relegated to “Appendices.”

2 Distributions and Quantum States

All of our discussion and results in this work will strongly rely on precise notions of distributions and related descriptions of quantum states in infinite-dimensional Hilbert spaces. Although most (or all) of this material is quite standard and well known [44, 66, 73, 101], we find it prudent to shortly summarize this background material in order to fix our notation and keep our presentation self-contained. This will also help to clarify differences and connections between divergent concepts and notations used in the literature. We hope this will also contribute to narrowing the gap between different physics communities that are interested in this topic.

2.1 Schwartz Space and Fourier Transforms

We will now summarize function spaces that are central to this work, refer also to [44, Ch. 1.1.3] and to [56, 112]. The set of all smooth, complex-valued functions on \(\mathbb R^n\) that decrease faster (together with all of their partial derivatives) than the reciprocal of any polynomial is called the Schwartz space and is usually denoted by \(\mathcal {S}(\mathbb {R}^n)\), refer to [101, Ch. V.3] or [73, Ch. 6]. More precisely, a function \(\psi :\mathbb R^n\rightarrow \mathbb C\) is called fast decreasing if the absolute values \(| x^{\beta } \partial _x^{\alpha } \psi (x) |\) are bounded for each multi-index of natural numbers \(\alpha :=(\alpha _1, \dots , \alpha _n)\) and \(\beta :=(\beta _1, \dots , \beta _n)\), where by definition \(x^{\beta }:=x_1^{\beta _1} \cdots x_n^{\beta _n}\) and \(\partial _x^{\alpha }:=\partial _{x_1}^{\alpha _1} \cdots \partial _{x_n}^{\alpha _n}\), refer to [44, Ch. 1.1.3]. This gives rise to a family of seminorms \(\Vert \psi \Vert _{\alpha ,\beta }:=\sup _{x\in \mathbb R^n}|x^\beta \partial _x^\alpha \psi (x)|\) which turn \(\mathcal S(\mathbb R^n)\) into a topological space which is even a Fréchet space [101, Thm. V.9].

The topological dual space \(\mathcal S'(\mathbb R^n)\) of \(\mathcal S(\mathbb R^n)\) is often referred to as the space of tempered distributions, and we will denote the distributional pairing for \(\phi \in \mathcal S'(\mathbb R^n)\) and \(\psi \in \mathcal S(\mathbb R^n)\) as \(\langle \phi , \psi \rangle :=\phi (\psi )\in \mathbb C\). In Sect. 2, we will consistently use the symbol \(\phi \) to denote distributions and \(\psi ,\psi '\) to denote Schwartz or square-integrable functions. Also, note that \(\mathcal {S}(\mathbb {R}^n)\) is dense in \( L^2 (\mathbb {R}^n)\),Footnote 1 and that tempered distributions naturally include the usual function spaces \(\mathcal {S}(\mathbb {R}^n)\subset L^2 (\mathbb {R}^n) \subset \mathcal S'(\mathbb R^n)\) via distributional pairings in the form of an integral \( \langle \phi , \psi \rangle = \int _{\mathbb R^n}\phi ^*(x)\psi (x)\,\textrm{d}x\), where \(\phi ^*(x)\) is the complex conjugate of \(\phi (x) \in L^2(\mathbb {R}^n)\) or \(\phi (x) \in \mathcal S(\mathbb {R}^n)\). This inclusion is usually referred to as a rigged Hilbert space [26, 57] or the Gelfand triple.

Remarkably, every tempered distribution is the derivative of some polynomially bounded continuous function, that is, given \(\phi \in \mathcal S'(\mathbb R^n)\) there exists \(g:\mathbb R^n\rightarrow \mathbb C\) continuous such that \(|g(x)|\le C(1{+}x^2)^m\) for some \(C,m\ge 0\) and all \(x\in \mathbb R^n\), and a multi-index \(\alpha \) such that \( \langle \phi ,\psi \rangle {=}(-1)^{|\alpha |}\int _{\mathbb R^n}g^*(x) (\partial _x^\alpha \psi )(x)\,\textrm{d}x\) for all \(\psi \in \mathcal S(\mathbb R^n)\)—for short one can write \(\phi =\partial _x^\alpha g\) [101, Thm. V.10].

In particular, one can construct tempered distributions by considering smooth functions \(\phi \) that (together with all of their partial derivatives) grow slower than certain polynomials. More precisely, a smooth map \(\phi :\mathbb R^n\rightarrow \mathbb C\) is said to be slowly increasing or of slow growth if for every \(\alpha =(\alpha _1, \dots , \alpha _n)\) there exist constants C, m, and A such that \(|\partial _x^\alpha \phi (x)|\le C\Vert x\Vert ^m\) for all \(\Vert x\Vert >A\), where \(\Vert x \Vert \) is the Euclidean norm in \(\mathbb R^n\), refer to [73, Ch. 6.2]. A standard example of such functions are polynomials. In particular, every slowly increasing function \(\phi (x)\) generates a tempered distribution \(\langle \phi , \psi \rangle = \int _{\mathbb R^n} \phi ^*(x) \psi (x)\,\textrm{d}x\) for all \(\psi \in \mathcal S(\mathbb R^n)\), and therefore, such functions are usually denoted as \(\phi (x) \in \mathcal {S}'(\mathbb {R}^n)\) (refer to [73, Ch. 6.2]).

Example 1

This motivates the delta distribution \(\langle \delta _b, \psi \rangle := \psi (b)\) which is in its integral representation commonly written as \( \int _{\mathbb R^n} \delta (x{-}b) \psi (x) \,\textrm{d}x = \psi (b) \). We emphasize that the notation \(\delta (x)\) is, however, only formal, cf. [101, Eq. (V.3)]. Moreover, this tempered distribution is generated by the second derivative of the polynomially bounded continuous function \(g(x):=x{-}b\) for \(x\ge b\) and zero otherwise, i.e., \(\langle \delta _b,\psi \rangle =\int _\mathbb Rg(x)\psi ''(x)\,\textrm{d}x\) for all \(\psi \in \mathcal S(\mathbb R)\) [101, Ch. V, Ex. 8]. This generating function is not unique as, for example, one also has \(\delta _b=({d^2}/{\textrm{d}x^2})|x{-}b|/2\).

For the rest of our work, we will restrict the general space \(\mathbb R^n\) to the case of \(\mathbb R\) which is most relevant for the applications we highlight. This simplifies our notation, even though many statements could be generalized.

Recall that for all \(a\in \mathcal S(\mathbb R^2)\) the symplectic Fourier transform \([\mathcal {F}_{\sigma }a](x,p)\) (see App. B in [44]) is related to the usual Fourier transformFootnote 2

$$\begin{aligned}{}[\mathcal {F}a](x,p):= (2 \pi \hbar )^{-1}\int e^{- \tfrac{i}{\hbar } ( x'x + p'p )} a(x',p') \, \textrm{d}x'\, \textrm{d}p', \end{aligned}$$

up to a coordinate transformation \([\mathcal {F}_{\sigma }a](x,p) = [\mathcal {F}a](p,-x)\) where

$$\begin{aligned} {[}\mathcal {F}_{\sigma } a](x,p):= (2 \pi \hbar )^{-1}\int e^{-\tfrac{i}{\hbar } (x'p - xp')} a(x',p') \, \textrm{d}x'\, \textrm{d}p'. \end{aligned}$$
(4)

Note that the square \([\mathcal {F}_{\sigma } \mathcal {F}_{\sigma } a](x,p) = a(x,p)\) is equal to the identity, and that the Fourier transform of every function in \(\mathcal {S}(\mathbb {R}^n)\) is also in \(\mathcal {S}(\mathbb {R}^n)\), cf. [101, Ch. IX.1]. The fact that \(\mathcal {F}_\sigma \) is Hermitian, i.e., \(\langle \mathcal {F}_\sigma \phi ,\psi \rangle _{L^2}=\langle \phi ,\mathcal {F}_\sigma \psi \rangle _{L^2}\) for all \(\phi ,\psi \in \mathcal S(\mathbb R^2)\) (see Sect. 2.2) motivates us to define the symplectic Fourier transform of tempered distributions via the distributional pairing \(\langle \mathcal {F}_{\sigma } \phi , \psi \rangle := \langle \phi , \mathcal {F}_{\sigma } \psi \rangle = \phi (\mathcal {F}_{\sigma } \psi )\) for \(\phi \in S'(\mathbb R^2)\) and \(\psi \in \mathcal S(\mathbb R^2)\). Thus, this is the extension of \(\mathcal {F}_\sigma \) with respect to the distributional pairing in our sense, cf. also “Appendix A.” In particular, the symplectic Fourier transform generalizes to phase-space distribution functions a(xp) without further adjustment and all the properties of \(\mathcal {F}_\sigma \) on \(\mathcal S(\mathbb R^2)\) transfer to \(\mathcal S'(\mathbb R^2)\).

Let us come back to our previous example: the delta distribution can be identified formally via the brackets \( \langle \delta _0, \mathcal {F}_{\sigma } \psi \rangle = [\mathcal {F}_{\sigma } \psi ](0) =(2\pi \hbar )^{-1} \langle 1, \psi \rangle \) as the Fourier transform \(\delta (x) = (2 \pi \hbar )^{-1} \mathcal {F}_{\sigma }[1] \) of the constant function, refer to [73, Ch. 6.4].

2.2 Quantum States and Expectation Values

Let us denote the abstract state vector of a quantum system by \(| \psi \rangle \) which is an element of an abstract, infinite-dimensional, separable complex Hilbert space (here and henceforth denoted by) \(\mathcal H\). The Hilbert space \(\mathcal H\) is known as the state space and it is equipped with a scalar product \(\langle \,\cdot \,|\,\cdot \, \rangle \) [66]. Considering projectors \(\mathcal {P}_{\psi }:=| \psi \rangle \langle \psi |\) defined via the open scalar products \(\mathcal {P}_{\psi } = \langle \psi | \,\cdot \, \rangle \, | \psi \rangle \), an orthonormal basis of \(\mathcal H\) is given by \(\{ | \phi _n \rangle , n\in \mathbb {N} \} \) if \( \langle \phi _n | \phi _m \rangle = \delta _{nm}\) for all \(m,n\in \mathbb N\) and \(\sum _{n = 0}^\infty \mathcal {P}_{\phi _n} = \mathbbm {1}\) in the strong operator topology. For a broader introduction to this topic we refer to [66].

Depending on the given quantum system, explicit representations of the state space can be obtained by specifying its Hilbert space [58]. In the case of bosonic systems, the Fock (or number-state) representation is widely used. There a quantum state \(| \psi \rangle \) is an element of the Hilbert space \(\ell ^2\) of square-summable sequences of complex numbers [66], and it is characterized by its expansion \( | \psi \rangle = \sum _{n=0}^{\infty } \psi _n | n \rangle \) into the Fock basis \(\{ | n \rangle , n=0,1, \dots \}\) of number states using the expansion coefficients \(\psi _n = \langle n | \psi \rangle \in \mathbb {C}\), refer to, e.g., [23] and [66, Ch. 11]. The scalar product \(\langle \psi | \psi ' \rangle \) then corresponds to the usual scalar product of vectors, i.e., to the absolutely convergent sum \(\sum _{n=0}^{\infty } (\psi _n)^* \psi '_n =: \langle \psi | \psi ' \rangle _{\ell ^2}\). The corresponding norm of vectors is then given by \(\Vert \psi \Vert _{\ell ^2} = \Vert ( | \psi \rangle ) \Vert _{\ell ^2} = [\langle \psi | \psi \rangle _{\ell ^2}]^{1/2}\).

For a quantum state \(| \psi \rangle \), the coordinate representation \(\psi (x)\in \mathcal {S}(\mathbb {R})\) and its Fourier transform (or momentum representation) \(\psi (p)\in \mathcal {S}(\mathbb {R})\) are given by complex, square-integrable, and smooth functions that are also fast decreasing. The quantum state \(| \psi \rangle = \int _{\mathbb {R}} \psi (x) | x \rangle \, \textrm{d}x\) of \(\psi (x) = \langle x | \psi \rangle \) is then defined via coordinate eigenstatesFootnote 3\( | x \rangle \). The coordinate representation of a coordinate eigenstate is given by the distribution \(\delta (x'{-}x) \in \mathcal {S}'(\mathbb {R})\), refer to [58, 66]. The scalar product \(\langle \psi | \psi ' \rangle \) is then fixed by the usual \(L^2\) scalar product, i.e., by the convergent integral \(\int _{\mathbb {R}}\psi ^*(x) \psi '(x) \, \textrm{d}x =: \langle \psi | \psi ' \rangle _{L^2}\). This integral induces the norm of square-integrable functions via \(\Vert \psi (x)\Vert _{L^2} = [\langle \psi | \psi \rangle _{L^2}]^{1/2}\).

The two examples given above are particular representations of the state space, which are convenient for specific physical systems; however, these representations are well known to be equivalent via

$$\begin{aligned} \mathcal {H} \simeq \ell ^2 \simeq L^2(\mathbb {R},\textrm{d}x) \simeq L^2(\mathbb {R},\textrm{d}p), \end{aligned}$$
(5)

refer to Theorem 2 in [58]. In particular, any coordinate representation \(\psi (x)\in {L^2}(\mathbb {R})\) of a quantum state \(| \psi \rangle \) can be expanded in the number-state basis into \(\psi (x) = \sum _{n=1}^{\infty }\psi _n\, \psi ^{\textrm{Fock}}_{n}(x)\) via \(\psi _n = \int _{\mathbb {R}} [\psi ^{\textrm{Fock}}_{n}(x)]^*(x) \psi (x) \, \textrm{d}x\) where \(\psi ^{\textrm{Fock}}_{n}(x) \in \mathcal {S}(\mathbb {R})\) are eigenfunctions of the quantum harmonic oscillator. For any \(\psi (x),\psi '(x)\in {L^2}(\mathbb {R})\), the \(L^2\) scalar product is equal to the \(\ell ^2\) scalar product

$$\begin{aligned} \int _{\mathbb {R}}\psi ^*(x) \psi '(x) \, \textrm{d}x = \sum _{n,m=1}^{\infty }\psi ^*_n \psi '_m \int _{\mathbb {R}} [\psi ^{\textrm{Fock}}_{n}(x)]^* \psi ^{\textrm{Fock}}_{m}(x) \, \textrm{d}x = \sum _{n=1}^{\infty }\psi ^*_n \psi '_n, \end{aligned}$$
(6)

and it is independent of the chosen orthonormal basis as any two orthonormal bases are related via a unitary transformation. The Plancherel formula \(\int _{\mathbb {R}}\psi ^*(x) \psi '(x) \, \textrm{d}x =\int _{\mathbb {R}}\psi ^*(p) \psi '(p) \, \textrm{d}p\) yields the equivalence \(L^2(\mathbb {R},\textrm{d}x) \simeq L^2(\mathbb {R},\textrm{d}p)\).

Motivated by the invariance of the scalar product under the choice of representation, in the following we will consistently use the notation \(\langle \,\cdot \,|\,\cdot \, \rangle \) for scalar products in Hilbert space, without specifying the type of representation. However, in order to avoid confusion with different types of operator or Euclidean norms, we will use in the following the explicit norms \(\Vert \psi (x)\Vert _{L^2}\) and \(\Vert \, |\psi \rangle \Vert _{\ell ^2}\), despite their equivalence.

Finally, let us summarize some of the main concepts on operators on infinite-dimensional Hilbert spaces, refer to [97, Ch. 15 & 16] for a comprehensive introduction. We start with the set of bounded linear operators \(\mathcal B(\mathcal H)\), that is, the collection of all \(A:\mathcal H\rightarrow \mathcal H\) linear for which the operator norm

$$\begin{aligned} \Vert A\Vert _\textrm{sup}:= \sup _{\Vert \,|\psi \rangle \Vert _{\mathcal H}=1} \Vert A | \psi \rangle \Vert _{\mathcal H} \end{aligned}$$
(7)

is finite. As discussed previously—after translating A into an equivalent operator on \(L^2\) or \(\ell ^2\)—this formalism encompasses \(\Vert A\Vert _\textrm{sup}= \sup _{\Vert \,|\psi \rangle \Vert _{\ell ^2}=1} \Vert A | \psi \rangle \Vert _{\ell ^2}\) for the Hilbert space \(\ell ^2\) (number-state representation), as well as \(\Vert A\Vert _\textrm{sup}=\sup _{\Vert \psi (x) \Vert _{L^2} =1 } \Vert A \psi (x)\Vert _{L^2}\) for square-integrable functions \(\psi (x)\) (coordinate representation). Next one looks at the set of all compact operatorsFootnote 4\(\mathcal K(\mathcal H)\subseteq \mathcal B(\mathcal H)\) as every compact A can be written as \(\sum _{n\in N}s_n(A)|f_n\rangle \langle g_n|\) where \(\{|f_n\rangle \}_{n\in N}\), \(\{|g_n\rangle \}_{n\in N}\), \(N\subseteq \mathbb N\) are orthonormal systems and \(s_n(A)>0\) are the unique singular values of A [97, Prop. 16.3]. This lets one define the trace class \( \mathcal B^1(\mathcal H):=\{A\in \mathcal K(H):\sum _{n\in N} s_n(A)<\infty \} \) and the Hilbert–Schmidt operators \( \mathcal B^2(\mathcal H):=\{A\in \mathcal K(H):\sum _{n\in N} s_n(A)^2<\infty \} \).

For all \(A\in \mathcal B^1(\mathcal H)\), one then defines the trace via the absolutely convergent sum \(\textrm{Tr}{(A)}:=\sum _{n=1}^\infty \langle \psi _n|A\psi _n \rangle \) where the right-hand side is independent of the chosen orthonormal basis \(\lbrace |\psi _n\rangle \,,\,n\in \mathbb N\rbrace \) of \(\mathcal H\). The name “trace class” is due to the fact that it is the largest subset of \(\mathcal B(\mathcal H)\) where the trace can be reasonably defined [97, Prop. 16.18]. Equipped with the trace norm \(\Vert A\Vert _1:=\sum _{n\in N} s_n(A)\), the trace class is a Banach space, and the Hilbert–Schmidt operators even form a Hilbert space under the (well-defined) inner product \(\langle A,B\rangle _\textrm{HS}:=\textrm{Tr}{(A^\dagger B)}\); here, \(A^{\dagger }\) is the adjoint of A (which is in finite dimensions given by the complex conjugated and transposed matrix). Trace-class operators \(A\in \mathcal B^1(\mathcal H)\) have the important property that their products with bounded operators \(B\in \mathcal B(\mathcal H)\) are also in the trace class, i.e., \(A B,BA\in \mathcal B^1(\mathcal H)\). With this, one finds that the trace is linear and continuous with respect to the trace norm, and one has the following important trace inequality: \(|\textrm{Tr}{(AB)}|\le \Vert A\Vert _1\Vert B\Vert _\textrm{sup}\) for all \(A\in \mathcal B^1(\mathcal H)\), \(B\in \mathcal B(\mathcal H)\) [97, Lemma 16.23].

Thus, one defines a density operator or state \(\rho \in \mathcal B^1(\mathcal H)\) to be positive semi-definiteFootnote 5 with \(\textrm{Tr}{(\rho )}=1\). It therefore admits a spectral decomposition [97, Prop. 16.2], i.e., there exists an orthonormal system \(\lbrace |\psi _n\rangle ,n\in \mathbb {N} \rbrace \) in \(\mathcal H\) such that

$$\begin{aligned} \rho = \sum _{n = 1}^\infty \, p_n | \psi _n \rangle \langle \psi _n |. \end{aligned}$$
(8)

The probabilities \(\lbrace p_n, n\in \mathbb {N} \rbrace \) satisfy \(p_1\ge p_2\ge \ldots \ge 0\) and \(\sum _{n = 1}^\infty p_n=1\). As expectation values of observables are computed via the trace \( \langle O \rangle _\rho = \textrm{Tr}{(\rho \,O)} = \sum _{n=1}^\infty p_n \langle \psi _n | O \psi _n \rangle \) where \(O\in \mathcal B(\mathcal H)\) is self-adjoint, as a simple consequence of the trace inequality stated earlier one finds:

Lemma 1

The expectation value of an observable O in a mixed quantum state is upper bounded by the operator norm \( | \textrm{Tr}{(\rho \,O)} | \le \Vert O\Vert _\textrm{sup}\) for arbitrary density operators \(\rho \).

3 Coherent States, Phase Spaces, and Parity Operators

We continue to fix our notation by discussing an abstract definition of phase spaces that relies on displaced parity operators. This usually appears concretely in terms of coherent states [21, 54, 90, 100], for which we consider two equivalent but equally important parametrizations of the phase space using the coordinates \(\alpha \) or (xp) (see below). This definition of phase spaces can be also related to convolutions of Wigner functions which is usually known as the Cohen class [28, 29, 44]. We also recall important postulates for Wigner functions as given by Stratonovich [21, 113], and these will be considered later in the context of general phase spaces.

3.1 Phase-Space Translations of Quantum States

We will now recall a definition of the phase space for quantum-mechanical systems via coherent states, refer to [21, 50, 54, 72, 90, 100]. We consider a quantum system which has a specific dynamical symmetry group given by a Lie group G. The Lie group G acts on the Hilbert space \(\mathcal {H}\) using an irreducible unitary representation \(\mathcal {D}\) of G. By choosing a fixed reference state as an element \(| 0 \rangle \in \mathcal {H}\) of the Hilbert space, one can define a set of coherent states as \( | g \rangle := \mathcal {D}(g) | 0 \rangle \) where \(g\in G\). Considering the subgroup \(H\subseteq G\) of elements \(h \in H\) that act on the reference state only by multiplication \(\mathcal {D}(h)| 0 \rangle := e^{i\phi }| 0 \rangle \) with a phase factor \(e^{i\phi }\), any element \(g \in G\) can be decomposed into \(g = \mathrm {\Omega }h\) with \(\mathrm {\Omega }\in G /H\). The phase space is then identified with the set of coherent states \(| \mathrm {\Omega }\rangle := \mathcal {D}(\mathrm {\Omega }) | 0 \rangle \). In the following, we will consider the Heisenberg–Weyl groupFootnote 6\(H_3\), for which the phase space \(\mathrm {\Omega }\in H_3 /\mathbb R\) is a plane.

Next, we introduce the corresponding displacement operators that generate translations of the plane. These operators are also known as Heisenberg–Weyl operators [44] or, in the physics literature, simply as Weyl operators [6, 41, 70]. In particular, for harmonic oscillator systems, the phase space \(\mathrm {\Omega }\equiv \alpha \in \mathbb {C}\) is usually parametrized by the complex eigenvalues \(\alpha \) of the annihilation operator \(\hat{a}\) and Glauber coherent states can be represented explicitly [23] in the so-called Fock (or number-state) basis as

$$\begin{aligned} | \alpha \rangle = e^{-|\alpha |^2/2} \sum _{n=0}^{\infty } \tfrac{\alpha ^n}{\sqrt{n!}} | n \rangle = e^{\alpha \hat{a}^\dagger -\alpha ^* \hat{a}} | 0 \rangle =: \mathcal {D}(\alpha ) | 0 \rangle . \end{aligned}$$
(9)

Here, the second equality specifies the displacement operator \(\mathcal {D}(\alpha )\) as a power series of the usual bosonic annihilation \(\hat{a}\) and creation \(\hat{a}^\dagger \) operators, which satisfy the commutation relation \([\hat{a},\hat{a}^\dagger ] = 1\), refer to Eq. (2.11) in [23]. In particular, the number-state representation of displacements is given by [23]

$$\begin{aligned}{}[\mathcal {D}(\alpha )]_{m n} := \langle m | \mathcal {D}(\alpha ) | n \rangle = (\tfrac{n!}{m!})^{1/2} \alpha ^{m-n} e^{-|\alpha |^2/2} L_n^{(m-n)}(|\alpha |^2), \end{aligned}$$
(10)

forFootnote 7\(m\ge n\) where \(L_n^{(m-n)}(x)\) are generalized Laguerre polynomials. This is the usual formulation for bosonic systems (e.g., in quantum optics) [88], where the optical phase space is the complex plane and the phase-space integration measure is given by \(\textrm{d}\mathrm {\Omega }=2\hbar \,\textrm{d}^2\alpha =2 \hbar \, \textrm{d}\Re (\alpha ) \,\textrm{d}\Im (\alpha )\) (where one often sets \(h=2\pi \hbar =1\), cf., [21,22,23]). The real and imaginary parts of \(\alpha \) are denoted by \(\Re (\alpha )\) and \(\Im (\alpha )\), respectively. The annihilation operator admits a simple decomposition

$$\begin{aligned} \hat{a}= 2 \hbar \, \int _{-\infty }^{\infty } \int _{-\infty }^{\infty } \alpha \, | \alpha \rangle \langle \alpha | \, \textrm{d}\Re (\alpha ) \,\textrm{d}\Im (\alpha ) \end{aligned}$$

with respect to its eigenvectors, see, e.g., [23, Eqs. (2.21)–(2.27)].

Let us now consider the coordinate representation \(\psi (x) \in \mathcal {S}(\mathbb {R})\) of a quantum state. The phase space is parametrized by \(\mathrm {\Omega }\equiv (x,p) \equiv z \in \mathbb {R}^2\) and the integration measure is \(\textrm{d}\mathrm {\Omega }= \textrm{d}z= \textrm{d}x \, \textrm{d}p\). The displacement operator acts via (see also [119,120,121])

$$\begin{aligned} \mathcal {D}(x_0,p_0) \psi (x) := e^{\tfrac{i}{\hbar }(p_{0}x-\tfrac{1}{2}p_{0}x_{0})} \psi (x{-}x_{0}) = e^{-\tfrac{i}{\hbar } (x_0 \hat{p} - p_0 \hat{x}) } \psi (x), \end{aligned}$$
(11)

where \(x,x_0,p_0 \in \mathbb {R}\). The right-hand side of Eq. (11) specifies the displacement operator as a power series of the usual operators \(\hat{x}\) and \(\hat{p}\), which satisfy the commutation relation \([\hat{x},\hat{p}] = i\hbar \), refer to [44, Sec. 1.2.2., Def. 2].

The most common representations of these two unbounded operators are \(\hat{x} \psi (x)= x \psi (x)\) and \(\hat{p} \psi (x)= -i \hbar {\partial \psi (x)}/{\partial x}\). Displacements of tempered distributions \(\phi (x) \in \mathcal {S}'(\mathbb {R})\) are understood via the distribution pairings \( (\mathcal {D}(\mathrm {\Omega })\phi )( \psi ):= \phi ( \mathcal {D}(-\Omega ) \psi )\) where \(-\Omega =(-x_0,-p_0)\). This definition guarantees thatFootnote 8

$$\begin{aligned} \mathcal D(\Omega )[\langle \phi ,\cdot \,\rangle ](\psi )=\langle \mathcal D(\Omega )\phi ,\psi \rangle \end{aligned}$$

as integrals from Sect. 2.1 (cf. Example 3(2), “Appendix A”) for all \(\phi :\mathbb R\rightarrow \mathbb C\) such that \(\langle \phi ,\cdot \,\rangle \in \mathcal S'(\mathbb R)\), and all \(\psi \in \mathcal S(\mathbb R)\), \(\Omega \in \mathbb R^2\). In particular it does not matter whether \(\mathcal D(\Omega )\) acts on a function \(\phi :\mathbb R\rightarrow \mathbb C\) or on the induced functional \(\psi \mapsto \langle \phi ,\psi \rangle \).

The two (above mentioned) physically motivated examples are particular representations of the displacement operator for the Heisenberg–Weyl group in different Hilbert spaces that rely on different parametrizations of the phase space. Let us now highlight the equivalence of these two representations. In particular, we obtain the formulas \(\hat{a}_{\lambda }=(\lambda \hat{x} + i \lambda ^{-1} \hat{p} )/\sqrt{2 \hbar }\) and \( \hat{a}^\dagger _{\lambda }=(\lambda \hat{x} - i \lambda ^{-1} \hat{p} )/\sqrt{2 \hbar }\) for any nonzero real conversion factor \(\lambda \) with physical dimension \(\sqrt{\hbar }/[x]\) (where [x] denotes the physical dimension of x), refer to Eqs. (2.1–2.2) in [23]. In the context of quantum optics, the operators \(\hat{x}\) and \(\hat{p}\) are the so-called optical quadratures [88]. The operators \(\hat{a}_{\lambda }\) and \( \hat{a}^\dagger _{\lambda }\) are now defined on the Hilbert space \(L^2(\mathbb R)\), whereas \(\hat{a}\) and \(\hat{a}^\dagger \) act on elements of the Hilbert space \(\ell ^2\). For any \(\lambda \ne 0\), they reproduce the commutator \([\hat{a}_{\lambda }, \hat{a}^\dagger _{\lambda }]={\text {id}}_{L^2}\), i.e., \([\hat{a}_{\lambda }, \hat{a}^\dagger _{\lambda }] \, \psi (x) = \psi (x)\) for all \(\psi (x)\in L^2(\mathbb R)\), and they correspond to raising and lowering operators of the quantum harmonic oscillatorFootnote 9 eigenfunctions \(\psi ^{\textrm{Fock}}_{n}(x)\), refer to [66]. Substituting now \(\hat{x}=\sqrt{\hbar /2} \lambda ^{-1} (\hat{a}_{\lambda }{+} \hat{a}^\dagger _{\lambda })\) and \(\hat{p}=-i \sqrt{\hbar /2} \lambda (\hat{a}_{\lambda }{-} \hat{a}^\dagger _{\lambda })\) into the exponent on the right-hand side of (11) yields

$$\begin{aligned} -\tfrac{i}{\hbar } (x_0 \hat{p} - p_0 \hat{x}) = \hat{a}^\dagger _{\lambda }(x_0 \lambda + i \lambda ^{-1} p_0 )/\sqrt{2 \hbar } - \hat{a}_{\lambda }(x_0 \lambda - i \lambda ^{-1} p_0 )/\sqrt{2 \hbar }. \end{aligned}$$

This then confirms the equivalence

$$\begin{aligned} \mathcal {D}(x_0,p_0) \psi (x) = e^{-\tfrac{i}{\hbar } (x_0 \hat{p} - p_0 \hat{x}) } \psi (x) = e^{ \hat{a}^\dagger _{\lambda }\alpha - \hat{a}_{\lambda }\alpha ^* } \psi (x) = \mathcal {D}(\alpha ) \psi (x), \end{aligned}$$
(12)

where the phase-space coordinate \(\alpha \) is defined by \(\alpha := (x_0 \lambda + i \lambda ^{-1} p_0 )/\sqrt{2 \hbar }\). Note that the corresponding phase-space element is then \(\textrm{d}\mathrm {\Omega }= 2 \hbar \, \textrm{d}\Re (\alpha ) \,\textrm{d}\Im (\alpha ) = \textrm{d} x \, \textrm{d}p \) which is independent of the choice of \(\lambda \). Let us also recall two properties of the displacement operator [119,120,121] (see, e.g., [44, p. 7]):

$$\begin{aligned} \mathcal {D}(x_0,p_0) \mathcal {D}(x_1,p_1)&= e^{\tfrac{i}{\hbar }(p_0 x_1 - p_1 x_0)} \mathcal {D}(x_1,p_1) \mathcal {D}(x_0,p_0) \end{aligned}$$
(13)
$$\begin{aligned} \mathcal {D}(x_0{+}x_1,p_0{+}p_1)&= e^{-\tfrac{i}{2\hbar }(p_0 x_1 - p_1 x_0)} \mathcal {D}(x_0,p_0) \mathcal {D}(x_1,p_1). \end{aligned}$$
(14)

In the following, we will use both phase-space coordinates \(\alpha \) and (xp) interchangeably. The displacement operator is obtained in both parametrizations, and they are equivalent via (12). Motivated by the group definition, we will also use the parametrization \(\mathrm {\Omega }\) for the phase space via \(\mathcal {D}(\mathrm {\Omega })\), where \(\mathrm {\Omega }\) corresponds to any representation of the group, including the ones given by the coordinates \(\alpha \) and (xp).

3.2 Phase-Space Reflections and the Grossmann–Royer Operator

Recall that the parity operator \(\mathrm {\Pi }\) reflects wave functions via \(\mathrm {\Pi }\psi (x):= \psi (-x)\) and \(\mathrm {\Pi }\psi (p):= \psi (-p)\) for coordinate-momentum representations [15, 44, 64, 91, 102], and \(\mathrm {\Pi }| \mathrm {\Omega }\rangle := | {-} \mathrm {\Omega }\rangle \) for phase-space coordinates of coherent states [15, 22, 91, 102]. This parity operator is obtained as a phase-space average

$$\begin{aligned} \mathrm {\Pi }:= (4 \pi \hbar )^{-1}\int \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }\end{aligned}$$
(15)

of the displacement operator from (11). One finds for all \(\psi \in \mathcal S(\mathbb R)\), \(x\in \mathbb R\) that

$$\begin{aligned} {[}\mathrm {\Pi }\psi ](x)&=(4 \pi \hbar )^{-1}\int [\mathcal {D}(\mathrm {\Omega })\psi ](x) \, \textrm{d}\mathrm {\Omega }\\&=\tfrac{1}{2}\cdot (2\pi \hbar )^{-1}\int e^{-\frac{i}{\hbar }(xp'-x'p)}[\mathcal {D}(\mathrm {\Omega })\psi ](x) \, \textrm{d}\mathrm {\Omega }\;\Big |_{x'=p'=0}\\&=\tfrac{1}{2} \big \{ \mathcal {F}_{\sigma } [\mathcal {D}\psi (x)](\mathrm {\Omega }')\big \} \big |_{\mathrm {\Omega }'=0}, \end{aligned}$$

or \(\mathrm {\Pi }=\tfrac{1}{2} \{ [\mathcal {F}_{\sigma } \mathcal {D}](\mathrm {\Omega }')\} |_{\mathrm {\Omega }'=0}\) for short. Thus, the parity operator equals evaluating the symplectic Fourier transform of the displacement operator at the phase-space point \(\mathrm {\Omega }' = 0\). This is related to the Grossmann–Royer operator

$$\begin{aligned} \tfrac{1}{2} [\mathcal {F}_{\sigma } \mathcal {D}](-\mathrm {\Omega }) = \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }\mathcal {D}^\dagger (\mathrm {\Omega }), \end{aligned}$$
(16)

which is the parity operator transformed by the displacement operator [15, 44, 64, 91, 102]. Here, we use in both (15) and (16) an abbreviated notation for formal integral transformations of the displacement operator.

Remark 1

This abbreviation in Eq. (16) is justified as the existence of the corresponding integral \( (\mathrm {\Pi }\phi )(\psi )= (4 \pi \hbar )^{-1} \int \phi ( \mathcal {D}^\dagger (\mathrm {\Omega })\psi )\, \textrm{d}\mathrm {\Omega }=(4 \pi \hbar )^{-1} \int \phi ( \mathcal {D}(\mathrm {\Omega })\psi )\, \textrm{d}\mathrm {\Omega }\) is guaranteed by, e.g., [44, Sec. 1.3., Prop. 8] for all \(\phi \in \mathcal {S}'(\mathbb {R})\). In the following, we will use this abbreviated notation for formal integral transformations of the displacement operator, i.e., by dropping \(\phi \). However, we might need to restrict the domain of more general parity operators to ensure the existence of the respective integrals.

3.3 Wigner Function and the Cohen Class

The Wigner function \(W_{\psi }(x,p)\) of a pure quantum state \(| \psi \rangle \) was originally defined by Wigner in 1932 [123], and it is (in modern terms) the integral transformation of a pure state \(\psi \in L^2(\mathbb {R})\), i.e.,

$$\begin{aligned}{} & {} W_{\psi }(x,p) = (2 \pi \hbar )^{-1}\int e^{- \tfrac{i}{\hbar } py} \psi ^*(x{-}\tfrac{1}{2}y) \psi (x{+}\tfrac{1}{2}y) \, \textrm{d}y \\{} & {} \quad = (\pi \hbar )^{-1} \langle \psi , \, \mathcal {D}(x,p) \mathrm {\Pi }\mathcal {D}^\dagger (x,p) \psi \rangle = (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, (|\psi \rangle \langle \psi |) \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }\mathcal {D}^\dagger (\mathrm {\Omega }) ]. \end{aligned}$$

The second and third equalities specify the Wigner function using the Grossmann–Royer operator [44, 64] from (16), refer to [44, Sec. 2.1.1., Def. 12]. We use this latter form to extend the definition of the Wigner function to mixed quantum states as in [4, 15, 22, 102].

Definition 1

The Wigner function \(W_\rho (\mathrm {\Omega })\in L^2(\mathbb {R}^2)\) of an infinite-dimensional density operator (or quantum state) \(\rho =\sum _n p_n | \psi _n \rangle \langle \psi _n |\in \mathcal B^1(L^2(\mathbb R))\) is proportional to the quantum-mechanical expectation value

$$\begin{aligned} W_\rho (\mathrm {\Omega }) := (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }\mathcal {D}^\dagger (\mathrm {\Omega }) ] = \sum _n p_n W_{\psi _n}(\mathrm {\Omega }) \end{aligned}$$
(17)

of the Grossmann–Royer operator from (16), which is the parity operator \(\mathrm {\Pi }\) transformed by the displacement operator \(\mathcal {D}(\mathrm {\Omega })\), refer also to [4, 15, 22, 44, 91, 102].

The square-integrable cross-Wigner transform \( W_{\psi ,\psi '}(\mathrm {\Omega }) \in L^2(\mathbb {R}^2)\) of two functions \(\psi ,\psi ' \in L^2(\mathbb {R})\) used in time–frequency analysis [43, 44] is obtained via the finite-rank operator \(A = | \psi \rangle \langle \psi '|\) in the form \(W_{\psi ,\psi '}(\mathrm {\Omega }):=W_A (\mathrm {\Omega })\). Furthermore, as \(\lbrace \mathcal D(\mathrm {\Omega }):\mathrm {\Omega }\in \mathbb C\rbrace \) forms a subgroup of the unitary group, the range of \(W_\rho \) is a subset of the \(\rho \)-numerical range of \((\pi \hbar )^{-1}\mathrm {\Pi }\) [47].

The Wigner representation is in general a bijective, linear mapping between the set of density operators (or, more generally, the trace-class operators) and the phase-space distribution functions \(W_\rho \) that satisfy the so-called Stratonovich postulates [21, 113]:

$$\begin{aligned}&\text {Postulate (i):}\; \rho \mapsto W_\rho \text { is linear and injective}, \\&\text {Postulate (ii):}\; W_{\rho ^\dagger }=W_\rho ^* \;\text {(reality)}, \\&\text {Postulate (iiia):}\; \textrm{Tr}( \rho ) = \int W_\rho \, \textrm{d} \mathrm {\Omega }\;\text {(normalization)}, \\&\text {Postulate (iiib):}\; \textrm{Tr}( A^\dagger \rho ) = \int a^* \, W_{\rho } \, \textrm{d} \mathrm {\Omega }\;\text {(traciality)}, \\&\text {Postulate (iv):}\; W_{\mathcal {D}(\mathrm {\Omega }')\rho }(\mathrm {\Omega })=W_\rho (\mathrm {\Omega }{-} \mathrm {\Omega }') \;\text {(covariance)}. \end{aligned}$$

The not necessarily boundedFootnote 10 operator A is the Weyl quantization of the phase-space function (or distribution) \(a(\mathrm {\Omega }) \in \mathcal {S}'(\mathbb {R}^2)\), refer to Sect. 4.3. Based on these postulates, the Wigner function was defined for phase-spaces of quantum systems with different dynamical symmetry groups via coherent states [21, 54, 79, 81, 100, 115].

Before finally presenting the definition of the Cohen class for density operators following [44, Sec. 8.1., Def. 93] or [29], let us first recall the concept of convolutions. Given Schwartz functions \(a,\phi \in \mathcal S(\mathbb R^2)\), one defines their convolution via

$$\begin{aligned} \phi *a:= 2\pi \hbar \, \mathcal {F}_\sigma [(\mathcal {F}_\sigma \phi )\, (\mathcal {F}_\sigma a)] \end{aligned}$$
(18)

which is again in \(\mathcal S(\mathbb R^2)\). In principle, this formula extends to general functions, although convergence may become an issue. These extensions are used in Theorem 1 as well as Sect. 4.3. Now, Eq. (18) as well as the fact that

$$\begin{aligned} (\phi *a)(\mathrm {\Omega })=\int \phi (\mathrm {\Omega }')a(\mathrm {\Omega }{-}\mathrm {\Omega }')\,\mathrm d\mathrm {\Omega }'=\langle \phi ^*,[\mathcal {T}(\Omega ) a]^\vee \rangle \end{aligned}$$

are, for example, shown in [101, Thm. IX.3], where \(a^\vee (\Omega ):= a(-\Omega )\) and \(\mathcal {T}(\Omega )\) is the operator which translates a function by \(\Omega \) [i.e., \(\mathcal {T}(\Omega ) a(\Omega '):= a(\Omega '{-}\Omega )\)]. With this in mind, one arrives at an extension of the convolution to tempered distributions [62, Eq. (4.37) ff.]: Given \(\theta \in \mathcal S'(\mathbb R^2)\), \(a\in \mathcal S(\mathbb R^2)\) set

$$\begin{aligned} (\theta *a)(\mathrm {\Omega }):=\theta \big ([\mathcal {T}(\Omega ) a]^\vee \big ) \end{aligned}$$
(19)

for all \(\mathrm {\Omega }\in \mathbb R^2\). This definition extends in a natural way to general linear functionals \(\theta :D_\theta \rightarrow \mathbb C\) on some subspace \(D_\theta \subseteq (\mathbb R^2\rightarrow \mathbb C)\), and general functions \(a:\mathbb R^2\rightarrow \mathbb C\) as long as \([\mathcal {T}(\Omega ) a]^\vee \in D_\theta \) for all \(\Omega \in \mathbb R^2\).

Defining the convolution via Eq. (19) is consistent with the distributional pairing in the sense that \(\langle \phi ^*|*a\equiv \phi *a\), if \(\langle \phi ^*| (\psi ):=\langle \phi ^*,\psi \rangle \) on \(\mathcal S(\mathbb R^2)\). Moreover, one readily verifies the identity \(\langle (\theta *a)^*,\psi \rangle =\theta (a^\vee *\psi )\) for all \(\theta \in \mathcal S'(\mathbb R^2)\), \(a,\psi \in \mathcal S(\mathbb R^2)\). This shows that Eq. (19) is equivalent to other extensions of convolutions commonly found in the literature, e.g., [101, p. 324]. Be aware that \(\theta *a\) is always a function of slow growth, that is, \(\langle (\theta *a)^*,\cdot \rangle \in \mathcal S'(\mathbb R^2)\) for all \(\theta \in \mathcal S'(\mathbb R^2)\), \(a\in \mathcal S(\mathbb R^2)\) [101, Thm. IX.4].

Definition 2

The Cohen class is the set of all linear mappings from density operators to phase-space distributions that are related to the Wigner function via a convolution. More precisely, a linear map \(F:\mathcal B^1(L^2(\mathbb R))\rightarrow (\mathbb R^2\rightarrow \mathbb C)\), \(\rho \mapsto F_\rho \) maps to the phase-space distributions if \(\langle F_\rho ,\cdot \,\rangle \in \mathcal S'(\mathbb R^2)\) for all \(\rho \in \mathcal B^1(L^2(\mathbb R))\). Then, \(F\) belongs to the Cohen class if there existsFootnote 11\(\theta \in \mathcal S'(\mathbb R^2)\) (called “Cohen kernel”) such that

$$\begin{aligned} F_\rho (\mathrm {\Omega })=[\theta *W_\rho ](\mathrm {\Omega }). \end{aligned}$$

This is a generalization of the definition commonly found in the literature [44, Def. 93]: There one restricts the domain of \(F\) from the full trace class to only rank-one operators \(\rho =|\phi \rangle \langle \psi |\) for some \(\phi ,\psi \in L^2(\mathbb R)\) or even \(\in \mathcal S(\mathbb R)\). As a simple example [44, p. 90], the Wigner function is in the Cohen class: To see this, choose \(\theta =\delta \) in the above definition: \( [\delta *W_\rho ](\mathrm {\Omega })=\delta ([\mathcal T(\mathrm {\Omega })W_\rho ]^\vee )=W_\rho (\mathrm {\Omega }).\)

Remark 2

Given some \(\theta \in \mathcal S'(\mathbb R^2)\) associated with an element \(F\) of the Cohen class, one formally obtains \( \mathcal {F}_\sigma [F_\rho ]=\mathcal {F}_\sigma [\theta *W_\rho ]=K_\theta \mathcal {F}_\sigma [W_\rho ] \) if the symplectic Fourier transform of \(\theta \) is generated by a function \(K_\theta :\mathbb R^2\rightarrow \mathbb C\) via the usual distributional pairing (we will call this “admissible” later, cf. Sect. 4.1). The reason we make this observation is that this object always exists: It is a product of two classical functions where \(\mathcal {F}_\sigma [W_\rho ]\) is a bounded and square-integrable function, i.e., \(|\mathcal {F}_\sigma [W_\rho ](\mathrm {\Omega })| = |\textrm{Tr}[\mathcal D(\mathrm {\Omega }) \rho ]| \le \Vert \mathcal D(\mathrm {\Omega })\Vert _\textrm{sup}\Vert \rho \Vert _1 =1\) due to unitarity of \(\mathcal D(\mathrm {\Omega })\), and \(W_\rho \in L^2(\mathbb {R}^2)\) [44, Proposition 68] so the same holds true for its Fourier transform. Thus—while the expression \(\theta *W_\rho \) may be ill defined for certain \(\theta \in \mathcal S'(\mathbb R^2)\), \(\rho \in \mathcal B^1(L^2(\mathbb R))\)—going to the Fourier domain yields a well-defined object which can be studied rather easily.

4 Theory of Parity Operators and Their Relation to Quantization

4.1 Phase-Space Distribution Functions via Parity Operators

We propose a definition for phase-space distributions and the Cohen class based on parity operators, the explicit form of which will be calculated in Sect. 4.4. A similar form has already appeared in quantum optics for the so-called s-parametrized distribution functions, see, e.g., [22, 98]. In particular, an explicit form of a parity operator that requires no integral transformation appeared in (6.22) of [23], including its eigenvalue decomposition which was later re-derived in the context of measurement probabilities in [98], refer also to [91, 102]. Apart from those results, mappings between density operators and their phase-space distribution functions have been established only in terms of integral transformations of expectation values, as in [3, 4, 22].

For a convolution kernel \(\theta \in S'(\mathbb R^2)\), we introduce the corresponding filter kernel

$$\begin{aligned} K_\theta := 2\pi \hbar \, \mathcal {F}_\sigma (\theta ) \end{aligned}$$
(20)

where \(\mathcal {F}_\sigma \) denotes the symplectic Fourier transform (see Sect. 2.1). Henceforth, we say \(\theta \in \mathcal S'(\mathbb R^2)\) is admissible if its filter kernel is generated by a function via the usual integral form of the distributional pairing \(\langle \phi , \psi \rangle =\phi (\psi )\in \mathbb C\) for \(\phi \in \mathcal S'(\mathbb R)\) and \(\psi \in \mathcal S(\mathbb R)\) (see Sect. 2.1): More precisely, \(\theta \) is admissible if there exists a function \(K_\theta \) from \(\mathbb R^2\) to \(\mathbb C\) such that \(2\pi \hbar \, \mathcal {F}_\sigma (\theta ) (\psi ) =\langle K_\theta ^*, \psi \rangle \) for all \(\psi \in \mathcal S(\mathbb R)\). In this case, we call \(K_\theta \) the filter function associated with \(\theta \).

Most importantly, if the convolution kernel is admissible and itself is generated by a function, i.e., if we consider \(\langle \theta ^*,\cdot \,\rangle \in S'(\mathbb R^2)\) admissible, then Eq. (20) simplifies to

$$\begin{aligned} K_\theta (\Omega ) = 2\pi \hbar \, [\mathcal {F}_\sigma \theta ^\vee ](\Omega )= 2\pi \hbar \, [\mathcal {F}_\sigma \theta ](-\Omega ) \end{aligned}$$
(21)

for all \(\Omega \in \mathbb R^2\). As before \(\theta ^\vee (\Omega ) = \theta (-\Omega )\). The technical condition of \(\theta \) being admissible is always satisfied in practice (cf. Tables 2 and 3). The advantage of only considering admissible kernels is that the definition of the (generalized) parity operator makes for an obvious generalization of the parity operator from Sect. 3.2. For an even more general definition, we refer to Remark 12 in “Appendix A.”

Definition 3

Given any admissible convolution kernel \(\theta \in S'(\mathbb R^2)\) with associated filter function \(K_\theta \), we define a parity operator \(\mathrm {\Pi }_\theta \) on \(\mathcal S(\mathbb R)\) via

$$\begin{aligned} \mathrm {\Pi }_\theta := (4 \pi \hbar )^{-1}\int K_\theta (\mathrm {\Omega }) \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }, \end{aligned}$$
(22)

that is, \([\mathrm {\Pi }_\theta \psi ](x):=(4\pi \hbar )^{-1}\int K_\theta (\Omega )[\mathcal D(\Omega )\psi ](x)\,\mathrm d\mathrm {\Omega }\) for all \(\psi \in \mathcal S(\mathbb R)\), \(x\in \mathbb R\). This extends to a parity operator on the tempered distributions \(\langle \Pi _\theta |:D_\theta \rightarrow \mathcal S'(\mathbb R)\) via

$$\begin{aligned} \langle \Pi _\theta | := (4 \pi \hbar )^{-1}\int K_\theta ^*(\mathrm {\Omega }) \mathcal {D}^\dagger (\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }\end{aligned}$$
(23)

(where the notation \(\langle \Pi _\theta |\) is replaced below with \(\Pi _\theta \)) with domain

$$\begin{aligned} D_\theta := \{ \phi \in \mathcal {S}'(\mathbb {R}) \,\text { s.t. }\, \int K_\theta ^*(\mathrm {\Omega })\phi [ \mathcal {D}(-\mathrm {\Omega })(\cdot )]\, \textrm{d}\mathrm {\Omega }\in \mathcal {S}'(\mathbb {R}) \}. \end{aligned}$$
(24)

We remark that the operator (22) has already appeared in Eq. (33) of [12] for the special case \(K_\theta (0)=1\). The latter, however, does not avoid potential domain problems, cf. Example 2.

The derivation of the extension (23) of \(\Pi _\theta \) to tempered distributions is detailed in “Appendix A.” Displacements of tempered distributions \(\phi \in \mathcal {S}'(\mathbb {R})\) are understood via the distributional pairing \( \langle \mathcal {D}(\mathrm {\Omega })\phi , \psi \rangle = \langle \phi , \mathcal {D}^\dagger (\mathrm {\Omega }) \psi \rangle \) and (23) gives rise to a well-defined linear operator \(\langle \Pi _\theta |\) from \(D_\theta \) to \(\mathcal S'(\mathbb R)\) acting on \(\psi \in \mathcal S(\mathbb R)\) via

$$\begin{aligned} ( \langle \Pi _\theta |\phi )(\psi )=(4\pi \hbar )^{-1}\int K_\theta ^*(\mathrm {\Omega }) \langle \phi ,\, e^{-\tfrac{i}{\hbar }(p_{0}x+\tfrac{1}{2}p_{0}x_{0})} \psi (x{+}x_{0}) \rangle \, \textrm{d}\mathrm {\Omega }. \end{aligned}$$
(25)

The definition of \(\Pi _\theta \) is independent of the object it acts on (see “Appendix A”): \(\langle \Pi _\theta |\langle \phi ,\cdot \,\rangle =\langle \Pi _\theta \phi ,\cdot \,\rangle \) for all \(\phi \in \mathcal S(\mathbb R)\) where \(\langle \phi ,\cdot \,\rangle \) denotes the functional \(\psi \mapsto \langle \phi ,\psi \rangle \in \mathcal S'(\mathbb R)\). All filter functions used in practice (refer to Tables 2 and 3) obey \(K_\theta ^*(x_0,p_0)=K_\theta (x_0,-p_0)\) for all \(x_0,p_0\in \mathbb R\). In this case, \(\langle \Pi _\theta |\) is not only compatible with the inner product on \(L^2(\mathbb R)\), but also with the embedding \(\mathcal S(\mathbb R)\hookrightarrow \mathcal S'(\mathbb R)\) usually employed in mathematical physics (see Lemma 3 in “Appendix A”). This motivates us to henceforth write \(\Pi _\theta \) in the case of both (22) and (23) (instead of \(\langle \Pi _\theta |\)).

While our definition above is pleasantly intuitive, we have to explicitly consider the domain of the parity operator. For a general (admissible) kernel \(\theta \), one needs to restrict the domain \(D_\theta \subseteq \mathcal {S}'(\mathbb {R})\) of \(\mathrm {\Pi }_\theta \) to tempered distributions for which the integral in Eq. (22) exists, as done in Eq. (24) and already hinted at in Remark 1.

Example 2

Domain considerations are illustrated using the standard ordering with \(K_\theta (\Omega _0)=\exp [{i p_0 x_0/(2 \hbar ) }]\) (see Table 2). Given any \(\phi ,\psi \in \mathcal S(\mathbb R)\), we have

$$\begin{aligned} \langle \phi ,\Pi _\theta \psi \rangle&=(\Pi _\theta \langle \psi ,\cdot \rangle )(\phi )^* =(4\pi \hbar )^{-1}\iiint \phi ^*(x{+}x_0)e^{\tfrac{i}{\hbar }p_0(x+x_0)}\psi (x)\,\textrm{d}x_0\,\textrm{d}x\,\textrm{d}p_0\nonumber \\&=(8\pi \hbar )^{-1/2}\Big (\int [\mathcal {F}\phi ](p_0)\,\textrm{d}p_0\Big )^*\Big (\int \psi (x)\,\textrm{d}x\Big )\nonumber \\&=\sqrt{\tfrac{\pi \hbar }{2}}\mathcal {F}[\mathcal {F}\phi ]^*(\tilde{p})\big |_{\tilde{p}=0}[\mathcal {F}\psi ](\tilde{x})\big |_{\tilde{x}=0} =\sqrt{\tfrac{\pi \hbar }{2}}\phi ^*(0)[\mathcal {F}\psi ](0). \end{aligned}$$
(26)

This reproduces known properties as in Eq. (5.39) of [30] (cf. Remark 3); however, we emphasize that, although Eq. (26) exists for all functions \(\phi ,\psi \) as long as \([\mathcal {F}\psi ](0)\) exists, this expression is only equal to \(\langle \phi ,\Pi _\theta \psi \rangle \) if in addition \(\phi \) and \(\mathcal {F}\phi \) are both in \(L^1\) (else the Fourier inversion formula used in the last step cannot be applied). In other words, a function \(\phi :\mathbb R\rightarrow \mathbb R\) is in the domain \(D_\theta \) of \(\Pi _\theta \) if and only if its Fourier transform exists and is in \(L^1(\mathbb R)\) if and only if (26) (resp. Eq. (5.39) of [30]) equals \(\langle \phi ,\Pi _\theta \psi \rangle \) for all suitable \(\psi \). In particular, \(D_\theta \) contains all Schwartz functions confirming that \(\Pi _\theta \) is densely defined. However, the functional \(\langle \phi ,\cdot \,\rangle \in \mathcal S'(\mathbb R)\) fails to be in \(D_\theta \) for most functions \(\phi :\mathbb R\rightarrow \mathbb C\) of slow growth including nonzero constant ones such as \(\phi := 1 \in \mathcal {S}'(\mathbb {R})\). In particular, \(\Pi _\theta \) does not extend to a well-defined operator on \(L^2(\mathbb R)\) as not all square-integrable functions will be contained in \(D_\theta \).

Following this line of thought, we investigate the well-definedness and boundedness of \(\Pi _\theta \) on the Hilbert space \(L^2(\mathbb R)\). As in Example 2, we observe that \(\mathcal S(\mathbb R)\subseteq D_\theta \) for all filter functions \(K_\theta \) which is particularly relevant for applications. This follows by interpreting \(\Pi _\theta \) as a Weyl quantization (cf. Sect. 4.3) whereby \(\theta \mapsto \Pi _\theta \) is specified as a map from \(\mathcal S'(\mathbb R^2)\) to the linear maps between \(\mathcal S(\mathbb R)\) and \(\mathcal S'(\mathbb R)\) (cf. Chapter 6.3 in [43] or Lemma 14.3.1 in [62]). Consequently, every parity operator has a well-defined matrix representation in the number-state basis (which is a subset of \(\mathcal S(\mathbb R)\), cf. Sect. 2.2). The following stronger statement is shown in “Appendix B.1”:

Lemma 2

Given any convolution kernel \(\theta \in \mathcal S'(\mathbb R^2)\), the following are equivalent:

(i,a) \(\Pi _\theta :L^2(\mathbb R)\rightarrow L^2(\mathbb R)\) is a well-defined linear operator, that is, the mapping \(x\mapsto \frac{1}{2}\theta (\mathcal {F}_\sigma [ \mathcal D\psi (x) ])\) (cf. Remark 12, “Appendix A”) is in \(L^2(\mathbb R)\) for all \(\psi \in L^2(\mathbb R)\).

(i,b) \([\theta *W_{\psi }](0,0)\) exists for all \(\psi \in L^2(\mathbb R)\), i.e., \([\theta *W_{\psi }](0,0)<\infty \).

Also, the following statements are equivalent:

(ii,a) \(\sup _{\psi ,\phi \in L^2(\mathbb R),\Vert \psi \Vert =\Vert \phi \Vert =1}| [\theta *W_{\phi ,\psi }](0,0)|<\infty \).

(ii,b) \((\phi ,\psi )\mapsto \theta *W_{\phi ,\psi }\) is weakly continuous on \(L^2(\mathbb R)\) in the sense that there exists \(C>0\) such that \(| [\theta *W_{\phi ,\psi }](0,0)|\le C\Vert \phi \Vert \Vert \psi \Vert \) for all \(\phi ,\psi \in L^2(\mathbb R)\).

(ii,c) \(\Pi _\theta \in \mathcal B(L^2(\mathbb R))\).

Moreover, if \(\theta \) is admissible, then (i,a), (i,b) and (ii,a), (ii,b), (ii,c) are all equivalent.

Recall from Sect. 3.3, \(W_{\phi ,\psi }\) is the usual cross-Wigner transform given by

$$\begin{aligned} W_{\phi ,\psi }(x,p)&=(\pi \hbar )^{-1}\langle \psi ,\mathcal D(x,p)\Pi \mathcal D^\dagger (x,p)\phi \rangle \\&=(2\pi \hbar )^{-1}\int _{-\infty }^\infty e^{-\tfrac{i}{\hbar }py}\psi ^*(x{-}\tfrac{y}{2}) \phi (x{+}\tfrac{y}{2})\,dy. \end{aligned}$$

Let us highlight that condition (ii,b) in Lemma 2 is a known sufficient condition from time–frequency analysis to ensure that a tempered distribution \(\theta \) is an element of the Cohen class, cf. Theorem 4.5.1 in [62]. Now, the almost magical result of Lemma 2 is that \(\Pi _\theta \) being well defined on \(L^2(\mathbb R)\) automatically implies boundedness as long as \(\theta \) is admissible. This can also be attributed to the folklore that unbounded operators “cannot be written down explicitly”: As the operator \(\Pi _\theta \) for admissible kernels is defined via an explicit integral, one gets the boundedness of \(\Pi _\theta \) “for free.” Indeed, the proof that all five statements from the above lemma are equivalent breaks down if one considers not only admissible but arbitrary kernels.

We define a general class of phase-space distribution functions \(F_\rho (\mathrm {\Omega },\theta )\) via the (formal) expression \( (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_\theta \mathcal {D}^\dagger (\mathrm {\Omega }) ]\). For general \(\theta \), however, this only makes sense if all displaced quantum states \(\mathcal D^\dagger (\mathrm {\Omega })\rho \mathcal D(\mathrm {\Omega })\) are supported on \(D_\theta \). We avoid these technicalities by restricting the definition to those filter functions which give rise to operators \(\Pi _\theta \) that are bounded on \(L^2(\mathbb R)\) and thereby allow for general \(\rho \).

Definition 4

Given any \(\theta \in \mathcal S'(\mathbb R^2)\) such that \(\Pi _\theta \in \mathcal B(L^2(\mathbb R))\), we define a linear mapping \(F_\rho (\cdot ,\theta )\) on the density operators \(\rho \in \mathcal B^1(L^2(\mathbb R))\) in the form of the quantum-mechanical expectation value

$$\begin{aligned} F_\rho (\mathrm {\Omega },\theta ) := (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_\theta \mathcal {D}^\dagger (\mathrm {\Omega }) ]. \end{aligned}$$
(27)

While our definition considers the practically most important case of bounded parity operators, in “Appendix C,” we give a detailed account of the extension of \(F_\rho (\mathrm {\Omega },\theta )\) to arbitrary \(\theta \in \mathcal S'(\mathbb R^2)\) whereby the associated parity operators may be unbounded. This is of importance for, e.g., the standard and antistandard orderings as shown in Example 2. The prototypical case where these extensions may not apply due to \(\theta \notin \mathcal S'(\mathbb R^2)\) is the case of the Glauber P function which is well known to be singular except for classical thermal states. However, most other convolution kernels appearing in practice are induced by a tempered distribution and thus fit into the framework of either Definition 4 or its extension in “Appendix C.”

Either way Definition 4 has many conceptual and computational advantages as we have detailed in the introduction. To further clarify the scope of said definition we now—similar to the proof of Lemma 2—relate the distribution functions \(F_\rho (\mathrm {\Omega },\theta )\) from Eq. (27) to the Cohen class (see Definition 2 and [44, Ch. 8]) by considering the filter function associated with any admissible kernel.

Theorem 1

Given any \(\theta \in S'(\mathbb R^2)\) such that \(\Pi _\theta \in \mathcal B(L^2(\mathbb R))\), the corresponding phase-space distribution function \(F_\rho (\mathrm {\Omega },\theta ) \in \mathcal {S}'(\mathbb {R}^2)\) as defined in Eq. (27) is an element of the Cohen class. In particular, \(F_\rho (\mathrm {\Omega },\theta )\) is related to the Wigner function \(W_\rho (\mathrm {\Omega })\) via the convolution

$$\begin{aligned} F_\rho (\mathrm {\Omega },\theta ) = [\theta *W_\rho ](\mathrm {\Omega }). \end{aligned}$$
(28)

If the convolution kernel \(\theta \in \mathcal S'(\mathbb R^2)\) is additionally admissible—meaning it is the reflected symplectic Fourier transform \(\theta =(2\pi \hbar )^{-1}\langle \mathcal {F}_\sigma K_\theta ^*|\) of its filter function \(K_\theta \)—then in analogy to (16) one finds

$$\begin{aligned} \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_\theta \mathcal {D}^\dagger (\mathrm {\Omega }) =\tfrac{1}{2} \mathcal {F}_{\sigma } [ K_\theta (\cdot ) \mathcal {D}(\cdot ) ](-\mathrm {\Omega }). \end{aligned}$$
(29)

The proof of Theorem 1 is given in “Appendix B.2.” The construction of a particular class of phase-space distribution functions was detailed in [4], where the term “filter function” also appeared in the context of mapping operators. However, these filter functions were restricted to nonzero, analytic functions. Definition 4 extends these cases to the Cohen class via Theorem 1 which allows for more general phase spaces. For example, the filter function of the Born–Jordan distribution has zeros (see Theorem 3), and is therefore not covered by [4]. Most of the well-known distribution functions are elements of the Cohen class. We calculate important special cases in Sect. 4.4. The Born–Jordan distribution and its parity operator are detailed in Sect. 5.

Our approach to define phase-space distribution functions using displaced parity operators also nicely fits with the characteristic [22, 29, 88] or ambiguity [44, Sec. 7.1.2, Prop. 5] function \(\chi (\Omega ) \in L^2(\mathbb {R}^2)\) of a quantum state that is defined as the expectation value \(\chi (\Omega ):=\textrm{Tr}[\rho \mathcal {D}(\mathrm {\Omega })]=[\mathcal {F}_{\sigma } W_\rho ](\Omega )\) or, equivalently, as the symplectic Fourier transform of the Wigner function \(W_\rho (\Omega )\). By multiplying the characteristic function \(\chi (\Omega )\) with a suitable filter function \(K_\theta (\Omega )\) and applying the symplectic Fourier transform, one obtains the Cohen class of phase-space distribution functions.

Remark 3

Definitions 3 and 4 for the parity operator and the phase-space function can be compared to prior work where special cases or similar parity operators have implicitly appeared and where similar restrictions on their existence must be observed. For example, the integral definition [28] of phase-space functions

$$\begin{aligned}&F_{|\phi \rangle \langle \phi |}(x,p,\theta )\nonumber \\&\quad =(4\pi ^2\hbar ^2)^{-1}\iiint \phi ^*(x'{-}\tfrac{y}{2}) \phi (x'{+}\tfrac{y}{2})K_\theta (-y,p')e^{-\tfrac{i}{\hbar }(xp'+ yp -x'p' )}\,\textrm{d}x'\,\textrm{d}y\,\textrm{d}p' \end{aligned}$$
(30)
$$\begin{aligned}&\quad =(\pi \hbar )^{-1}\langle \phi ,\mathcal D(x,p)\Pi _\theta \mathcal D^\dagger (x,p)\phi \rangle . \end{aligned}$$
(31)

as givenFootnote 12 in Eq. (5.2) of [30] translates into the definition (31) with the parity operator. Both Eqs. (30) and (31) need to respect domain restrictions as discussed in Example 2 and neither equation is well defined for tempered distributions in \(S'(\mathbb R)\) or square-integrable functions in \(L^2(\mathbb R)\) that are not contained in the domain \(D_\theta \).

4.2 Common Properties of Phase-Space Distribution Functions

We now detail important properties of \(F_\rho (\mathrm {\Omega },\theta )\) and their relation to properties of \(K_\theta (\mathrm {\Omega })\) and \(\mathrm {\Pi }_{\theta }\). These properties will guide our discussion of parity operators and this allows us to compare the Born–Jordan distribution to other phase spaces. Table 1 provides a summary of these properties, and the proofs are deferred to “Appendix D.” Recall that we are dealing exclusively with convolution kernels \(\theta \in \mathcal S'(\mathbb R^2)\) which give rise to bounded operators \(\Pi _\theta \) so the induced phase-space distribution \(F_\rho (\mathrm {\Omega },\theta )\) is well defined everywhere.

Property 1

Boundedness of phase-space functions \(F_\rho (\mathrm {\Omega },\theta )\): The phase-space distribution function \(F_\rho (\mathrm {\Omega },\theta )\) is bounded in its absolute value, i.e., \( \pi \hbar \, | F_\rho (\mathrm {\Omega },\theta ) | \le \Vert \mathrm {\Pi }_\theta \Vert _\textrm{sup}\) for all quantum states \(\rho \), refer to Lemma 1. In particular, then \(F_\rho (\Omega ,\theta )\in \mathcal S'(\mathbb R^2)\). Moreover, one finds that square-integrable filter functions give rise to bounded parity operators due to \( \Vert \mathrm {\Pi }_\theta \Vert _\textrm{sup} \le \Vert K_\theta \Vert _{L^2}/\sqrt{ 8 \pi \hbar }\). The proof of Property 1 in “Appendix D” implies the even stronger statement that \(\mathrm {\Pi }_{\theta }\) is a Hilbert–Schmidt operator if and only if \(K_\theta \) is square integrable.

Property 2

Square integrability: The phase-space distribution function \(F_\rho (\mathrm {\Omega },\theta )\) is square integrable [i.e., \(F_\rho (\mathrm {\Omega },\theta ) \in L^2(\mathbb {R}^2)\)] for all \(\rho \in \mathcal {B}^1(L^2(\mathbb R))\) if the absolute value of the filter function is bounded [i.e., \(K_\theta (\mathrm {\Omega }) \in L^\infty (\mathbb {R}^2)\)]. In particular, this implies \( F_\rho (\Omega ,\theta )\in \mathcal S'(\mathbb R^2)\).

Property 3

Postulate (iv): The phase-space distribution function \(F_\rho (\mathrm {\Omega },\theta )\) satisfies, by definition, the covariance property. In particular, a displaced density operator \(\rho ':=\mathcal {D}(\mathrm {\Omega }') \rho \mathcal {D}^\dagger (\mathrm {\Omega }')\) is mapped to the inversely displaced distribution function \(F_{\rho '}(\mathrm {\Omega },\theta ) = F_\rho (\mathrm {\Omega }{-} \mathrm {\Omega }',\theta )\).

Property 4

Rotational covariance: Let us denote a rotated density operator \(\rho ^{\phi } = U_\phi \rho U^\dagger _\phi \), where the phase-space rotation operator is given by \(U_\phi := \exp {(-i \phi \hat{a}^\dagger \hat{a})}\) in terms of creation and annihilation operators. The phase-space distribution function is covariant under phase-space rotations,Footnote 13 i.e., \(F_{\rho ^{\phi }}(\mathrm {\Omega },\theta )=F_\rho (\mathrm {\Omega }^{-\phi },\theta )\), if the filter function \(K_\theta (\mathrm {\Omega })\) (or equivalently the parity operator \(\mathrm {\Pi }_{\theta }\)) is invariant under rotations. Here, \(\mathrm {\Omega }^{-\phi }\) is the inversely rotated phase-space coordinate, e.g., \(\alpha ^{-\phi } = \exp {(i \phi )} \alpha \). As a consequence of this symmetry, the corresponding parity operators are diagonal in the number-state representation, i.e., \(\langle n | \mathrm {\Pi }_\theta | m \rangle \propto \delta _{nm}\).

Table 1 Properties of phase-space distribution functions from Definition 4

Property 5

Postulate (ii): The phase-space distribution function \(F_\rho (\mathrm {\Omega },\theta )\) is real if \(\mathrm {\Pi }_\theta \) is self-adjoint. This condition translates to the symmetry \(K^*_\theta (-\mathrm {\Omega })=K_\theta (\mathrm {\Omega })\) of the filter function.

Property 6

Postulate (iiia): The trace of a trace-class operator \( \textrm{Tr}\,[ \rho ]\) is mapped to the phase-space integral \(\int F_\rho (\mathrm {\Omega },\theta ) \, \textrm{d} \mathrm {\Omega }\) if the corresponding filter function satisfies \(K_\theta (0) = 1\). Note that this property also implies that the trace exists, i.e., \(\textrm{Tr}(\mathrm {\Pi }_{\theta }) = K_\theta (0)/2\), in some particular basis, even though \(\mathrm {\Pi }_{\theta }\) might not be of trace class.

Property 7

Marginals: An even more restrictive subclass of the Cohen class satisfies the marginal properties \( \int F_{|\psi \rangle \langle \psi |}(x,p,\theta ) \,\textrm{d} x = |\psi (p)|^2 \) and \( \int F_{|\psi \rangle \langle \psi |}(x,p,\theta ) \,\textrm{d} p = |\psi (x)|^2 \) if and only if \([K_\theta (x,p)]|_{p=0}=1\) and \([K_\theta (x,p)]|_{x=0}=1\). This follows, e.g., directly from Proposition 14 in Sec. 7.2.2 of [43].

4.3 Relation to Quantization

The Weyl quantization of a tempered distribution \(\mathfrak {a}\in \mathcal S'(\mathbb R^2)\) is obtained from the Grossmann–Royer operator in Eq. (16) (cf. [43, Sec. 6.3., Def. 7 and Prop. 9]), i.e., 

$$\begin{aligned} \text {Op}_{\text {Weyl}}(\mathfrak {a})=(\pi \hbar )^{-1}\mathfrak {a}(\mathcal {D}\Pi \mathcal {D}^\dagger ). \end{aligned}$$
(32)

More precisely, \(\text {Op}_{\text {Weyl}}(\mathfrak {a}):\mathcal S(\mathbb R)\rightarrow \mathcal S'(\mathbb R)\) is the well-defined linear map

$$\begin{aligned} {[}\text {Op}_{\text {Weyl}}(\mathfrak {a})\psi ](x)=(\pi \hbar )^{-1}\mathfrak {a}[(\mathcal {D}\Pi \mathcal {D}^\dagger \psi )(x)] \end{aligned}$$
(33)

for all \(\psi \in \mathcal S(\mathbb R)\), \(x\in \mathbb R\), where the argument of \(\mathfrak {a}\) is the Schwartz function \(\mathrm {\Omega }\mapsto (\mathcal D(\mathrm {\Omega })\Pi \mathcal D^\dagger (\mathrm {\Omega })\psi )(x)\) on \(\mathbb R^2\) [43, Sec. 6.3., Prop. 13]. If \(\mathfrak {a}\) is generated by a phase-space function \(a:\mathbb R^2\rightarrow \mathbb C\), i.e., \(\mathfrak a\equiv \langle a^*,\cdot \rangle \), then

$$\begin{aligned} \text {Op}_{\text {Weyl}}(\mathfrak a) = (\pi \hbar )^{-1} \int a(\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }\mathcal {D}^\dagger (\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }= (2 \pi \hbar )^{-1} \int a_\sigma (\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }, \end{aligned}$$
(34)

where the symplectic Fourier transform \(a_\sigma (\mathrm {\Omega }) = [\mathcal {F}_{\sigma } a(\cdot )](\mathrm {\Omega })\) is used for the second equality. Thus, \(\text {Op}_{\text {Weyl}}\) is similar to the generalized parity operator in the sense that it maps a function (or tempered distribution) to a linear operator which acts on real-valued functions. This is not by chance as these two objects are very much related to each other: recall that quantizations associated with the Cohen class \(\text {Op}_{\theta }(\mathfrak a)\) are essentially Weyl quantizations of convolved phase-space functions up to coordinate reflection, i.e., \(\text {Op}_{\theta }(\mathfrak a):= \text {Op}_{\text {Weyl}}(\theta ^\vee *a)\) where \(\theta ^\vee (\psi ):=\theta (\psi ^\vee )\) for all \(\psi \in \mathcal S(\mathbb R^2)\). If \(\theta \) is an admissible kernel in the sense of Sect. 4.1, then formally

$$\begin{aligned} \text {Op}_{\theta }(\mathfrak a)&= (\pi \hbar )^{-1} \int (\theta ^\vee *a) (\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }\mathcal {D}^\dagger (\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }\end{aligned}$$
(35)
$$\begin{aligned}&=(2\pi \hbar )^{-1} \int a_\sigma (\mathrm {\Omega }) K_\theta (\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }, \end{aligned}$$
(36)

cf. [43, Sec. 7.2.4, Prop. 17]. The symplectic Fourier transform \( [\mathcal {F}_{\sigma } (\theta ^\vee *a) ](\mathrm {\Omega }) = a_\sigma (\mathrm {\Omega }) K_\theta (\mathrm {\Omega })\) (as functionals on \(\mathcal S(\mathbb R^2)\) so in particular \( \langle \theta ^\vee *a,\cdot \rangle \in \mathcal S'(\mathbb R^2) \)) from Theorem 1 is used for the second equality, refer to §7.2.4 in [43].

Proposition 1

Let \(\theta ,\mathfrak {a}\in \mathcal S'(\mathbb R^2)\) be given such that \(\theta \) is admissible and the parity operator \(\Pi _\theta \) from Definition 3 is in \(\mathcal B(L^2(\mathbb R))\). If \(\mathfrak a\) is generated by a phase-space function \(a:\mathbb R^2\rightarrow \mathbb C\) (i.e., \(\mathfrak a\equiv \langle a^*,\cdot \rangle \)) and if \( \langle \theta ^\vee *a,\cdot \rangle \in \mathcal S'(\mathbb R^2) \), then

$$\begin{aligned} \text {Op}_{\theta }(\mathfrak a) = (\pi \hbar )^{-1} \int a(\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_{\theta } \mathcal {D}^\dagger (\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }\end{aligned}$$

in analogy to (34) as quadratic forms on \(\mathcal S(\mathbb R)\).

Proof

The Plancherel formula \(\int a(\mathrm {\Omega }) b(\mathrm {\Omega })\, \textrm{d} \mathrm {\Omega }= \int a_\sigma (\mathrm {\Omega }) b_\sigma (-\mathrm {\Omega })\, \textrm{d} \mathrm {\Omega }\) implies that

$$\begin{aligned} \text {Op}_{\theta }(\mathfrak a)&=(2\pi \hbar )^{-1} \int a_\sigma (\mathrm {\Omega }) K_\theta (\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }\\&=(2\pi \hbar )^{-1} \int \mathcal {F}_\sigma [a_\sigma (\cdot )](\mathrm {\Omega }) \, \mathcal {F}_\sigma [ K_\theta (\cdot ) \, \mathcal {D}(\cdot )](-\mathrm {\Omega }) \, \textrm{d} \mathrm {\Omega }, \end{aligned}$$

where the equality \(\mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_\theta \mathcal {D}^\dagger (\mathrm {\Omega }) =\tfrac{1}{2} \mathcal {F}_{\sigma } [ K(\cdot ) \mathcal {D}(\cdot ) ](-\mathrm {\Omega })\) follows from (29) (Theorem 1) and \(\mathcal {F}_\sigma [a_\sigma (\cdot )](\mathrm {\Omega }) = a(\mathrm {\Omega })\) is applied. \(\square \)

This result motivates the following extension of Eq. (32), special cases of which have already appeared in [12, Eq. (36)] as well as [55, (3.24)]:

Definition 5

The quantization \(\text {Op}_{\theta }(\mathfrak a)\) of any \(\mathfrak a \in \mathcal {S}'(\mathbb {R}^2)\) is defined to be

$$\begin{aligned} \text {Op}_{\theta }(\mathfrak {a})=(\pi \hbar )^{-1}\mathfrak {a}(\mathcal {D}\Pi _\theta \mathcal {D}^\dagger )\ \end{aligned}$$

in the sense of Eq. (33).

One can consider single Fourier components \(e^{i (p_0 x - x_0 p)/\hbar }=:f_{\mathrm {\Omega }_0}(\mathrm {\Omega })\), for which the Weyl quantization yields the displacement operator \(\text {Op}_{\text {Weyl}}( f_{\mathrm {\Omega }_0}) = \mathcal {D}(\mathrm {\Omega }_0)\) from Sect. 3.1, refer to Proposition 51 in [44] or Proposition 11 in Sec. 6.3.2 of [43]. Let us now consider the \(\theta \)-type quantization of a single Fourier component, which results in the displacement operator being multiplied by the corresponding filter function via (35)–(36). Substituting \(a_\sigma (\mathrm {\Omega })=\mathcal {F}_\sigma [f_{\mathrm {\Omega }_0}](\mathrm {\Omega })\) into (35), one obtains

$$\begin{aligned} \text {Op}_{\theta }( f_{\mathrm {\Omega }_0}) = K_{\theta }(\mathrm {\Omega }_0) \mathcal {D}(\mathrm {\Omega }_0) \quad \text {and} \quad \mathrm {\Pi }_{\theta } = (4 \pi \hbar )^{-1}\int \text {Op}_{\theta }( f_{\mathrm {\Omega }_0}) \, \textrm{d}\mathrm {\Omega }_0. \end{aligned}$$

The second equality follows from (22), and it specifies the parity operator as a phase-space average of quantizations of single Fourier components.

But it is even more instructive to consider the case of the delta distribution \(\delta ^{(2)} \), the Weyl quantization of which yields the Grossmann–Royer parity operator \(\text {Op}_{\text {Weyl}}( \delta ^{(2)} ) = \pi \hbar \, \mathrm {\Pi }\) (as obtained in [64]). Applying (35), the Cohen quantization of the delta distribution yields the parity operator from (22). In particular, the operator \(\mathrm {\Pi }_{\theta }\) from Definition 3 is a \(\theta \)-type quantization of the delta distribution as

$$\begin{aligned} \mathrm {\Pi }_{\theta } = ( \pi \hbar )^{-1} \, \text {Op}_{\theta }( \delta ^{(2)}) = (\pi \hbar )^{-1} \, \text {Op}_{\text {Weyl}}( \theta ^\vee ), \end{aligned}$$
(37)

or equivalently, the Weyl quantization of the Cohen kernel, up to coordinate reflection. Since \(\mathrm {\Pi }_{\theta }\) is the Weyl quantization of the tempered distribution \(\theta ^\vee \in \mathcal {S}'(\mathbb {R}^2)\), one can adapt results contained in [38] to precisely state conditions on \(\theta \), for which bounded operators \(\mathrm {\Pi }_{\theta }\) are obtained via their Weyl quantizations, refer also to Property 1. For example, square-integrable \(\theta \in L^2(\mathbb {R}^2)\) result in Hilbert–Schmidt operators \(\mathrm {\Pi }_{\theta }\), absolutely integrable \(\theta \in L^1(\mathbb {R}^2)\) result in compact operators \(\mathrm {\Pi }_{\theta }\), and Schwartz functions \(\theta \in \mathcal {S}(\mathbb {R}^2)\) result in trace-class operators \(\mathrm {\Pi }_{\theta }\), refer to [38].

Table 2 Common operator orderings, their defining filter functions \(K_{\tau s}(\mathrm {\Omega }_0)\), and the corresponding single Fourier component quantizations \(\text {Op}_{\tau {}s}( f_{\mathrm {\Omega }_0})\) as displacement operators with \(e^{i (p_0 x - x_0 p)/\hbar }=:f_{\mathrm {\Omega }_0}(\mathrm {\Omega })\), refer to, e.g., [3, 4, 22, 43, 102]. Coordinates \(\mathrm {\Omega }_0 \simeq \alpha _0 \simeq (x_0,p_0)\) with subindex 0 are used for clarity

We consider now a class of explicit quantization schemes along the lines of [3, 22, 23, 43] which are motivated by different \((\tau ,s)\)-orderings of non-commuting operators \(\hat{x}\) and \(\hat{p}\) or \(\hat{a}\) and \(\hat{a}^\dagger \) (\(-1 \le s \le 1\) and \(0 \le \tau \le 1\)). This class is obtained via the \((\tau ,s)\)-parametrized filter function (where the relation \(\alpha = (\lambda x + i p/\lambda )/\sqrt{2 \hbar }\) from Sect. 3.1 is used)

$$\begin{aligned} K_{\tau s}(\mathrm {\Omega }) := \exp {[ \tfrac{2\tau -1}{4} (\alpha ^2 {-} (\alpha ^*)^2) + \tfrac{s}{2} |\alpha |^2]} =\exp {[ \tfrac{i(2\tau -1) px}{2 \hbar } + \tfrac{s( \lambda ^2 x^2 + \lambda ^{-2} p^2)}{4 \hbar } ]},\nonumber \\ \end{aligned}$$
(38)

which admits the symmetries \(K_{\tau s}(\mathrm {\Omega }) = K_{\tau s}(-\mathrm {\Omega })\) and \( K_{\tau s}^*(x,p)=K_{\tau s}(x,-p) \). The corresponding \((\tau ,s)\)-parametrized quantizations of a single Fourier component are given by the operators \(\text {Op}_{\tau {}s}( f_{\mathrm {\Omega }_0}):= K_{\tau {}s}(\mathrm {\Omega }) \mathcal {D}(\mathrm {\Omega })\), which are central in ordered expansions into non-commuting operators. Also, note that for \(s \le 0\), the resulting parity operators are bounded as is readily verified; hence, the corresponding distribution functions are in the Cohen class with \(K_{\tau s}(\mathrm {\Omega }) \in \mathcal {S}'(\mathbb {R}^2)\) due to Theorem 1. Important, well-known special cases are summarized in Table 2, refer also to [3, 4, 22, 43].

4.4 Explicit Form of Parity Operators

Expectation values of displaced parity operators

$$\begin{aligned} \mathrm {\Pi }_{\tau s} = (4 \pi \hbar )^{-1}\int K_{\tau s}(\mathrm {\Omega }) \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }= (4 \pi \hbar )^{-1}\int \text {Op}_{\tau s}( f_{\mathrm {\Omega }_0}) \, \textrm{d}\mathrm {\Omega }\end{aligned}$$
(39)

are obtained via the kernel function in (38) and recover well-known phase-space distribution functionsFootnote 14 for particular cases of \(\tau \) or s, which are motivated by the ordering schemes \(\text {Op}_{\tau s}( f_{\mathrm {\Omega }_0})\) from Table 2. Important special cases of these distribution functions and their corresponding filter functions and Cohen kernels are summarized in Table 3.

In particular, the parameters \(\tau =1/2\) and \(s=0\) identify the Wigner function with \(K_{1/2, 0}(\mathrm {\Omega }) \equiv 1\) and (39) reduces to (15). Note that the corresponding Cohen kernel \(\theta \) from Theorem 1 is the two-dimensional delta distribution \(\delta ^{(2)}(\mathrm {\Omega })\) and that convolving with \(\delta ^{(2)}(\mathrm {\Omega })\) is the identity operation, i.e., \(\delta ^{(2)} *W_\rho = W_\rho \) [see (28)].

Table 3 Well-known phase-space distribution functions and their corresponding Cohen kernels recovered for particular values of \(\tau \) or s via expectation values of displaced parity operators from (39)

The filter function \(K_{\tau s}\) from (38) for a fixed parameter of \(\tau =1/2\) results in the Gaussian \(K_{s}(\mathrm {\Omega }):= K_{1/2, s}(\mathrm {\Omega }) = \exp {[ \tfrac{s}{2} |\alpha |^2]}\). The corresponding parity operators are diagonal in the number-state representation (refer to Property 4), and they can be specified for \(-1 \le s < 1\) in terms of number-state projectors [22, 98, 102] as

$$\begin{aligned} \mathrm {\Pi }_s := \mathrm {\Pi }_{1/2,s} = (4 \pi \hbar )^{-1} \int e^{s |\alpha |^2/2} \, \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }= \sum _{n=0}^\infty (-1)^n \frac{(1{+}s)^n}{(1{-}s)^{n+1}} | n \rangle \langle n | , \end{aligned}$$
(40)

where the second equality specifies \(\mathrm {\Pi }_s\) in the form of a spectral decomposition. We further discuss the representation of \(\mathrm {\Pi }_s\) in terms of creation and annihilation operators in Remark 9 (Sect. 5.3). This form has implicitly appeared in, e.g., [22, 98, 102]. We provide a more concise proof in “Appendix E.” Equation (40) readily implies \(\Vert \mathrm {\Pi }_s\Vert _\textrm{sup}=(1{-}s)^{-1}\) for \(s\le 0\), and for \(s<0\) one even finds that \(\mathrm {\Pi }_s\) are trace-class operators due to

$$\begin{aligned} \Vert \mathrm {\Pi }_s\Vert _1=\sum _{n=0}^\infty \frac{(1{+}s)^n}{(1{-}s)^{n+1}}=\frac{1}{(1{-}s)-(1{+}s)}=(2|s|)^{-1}. \end{aligned}$$

Note that for \(s>0\) the corresponding filter functions lie outside of our framework as then \(K_s\not \in \mathcal S'(\mathbb R^2)\) due to its superexponential growth. While one can still formally write down their distribution functions, one runs into convergence problems resulting in singularities. However, their symplectic Fourier transform always exists and it is related to the Wigner function via \(K_s(\mathrm {\Omega }) \mathcal {F}_\sigma [W_\rho ](\mathrm {\Omega })\) by multiplying with the filter function \(K_{s}(\mathrm {\Omega })\) (cf. Remark 2). This class of s-parametrized phase-space representations has gained widespread applications in quantum optics and beyond [32, 60, 94, 111, 125], and they correspond to Gaussian convolved Wigner functions

$$\begin{aligned} F_\rho (\mathrm {\Omega }, s) = F_{|0\rangle }(\mathrm {\Omega }, s{+}1) *W_\rho (\mathrm {\Omega }), \end{aligned}$$

for \(s<0\) such as the Husimi Q function for \(s=-1\). Note that the Cohen kernel \(\theta _{s} \) via Theorem 1 corresponds to the vacuum state \( F_{|0\rangle }(\mathrm {\Omega }, s{+}1) \) of a quantum harmonic oscillator [22]. Gaussian deconvolutions of the Wigner function are formally obtained for \(s>0\), which includes the Glauber P function for \(s=1\) [22]. Due to the rotational symmetry of its filter function \(K_{s}(\mathrm {\Omega })\), the s-parametrized distribution functions are covariant under phase-space rotations, refer to Property 4.

Another important special case is obtained for the fixed parameter \(s=0\) which results in Shubin’s \(\tau \)-distribution, refer to [16,17,18, 43]. Its filter function from (38) reduces to the chirp function

$$\begin{aligned} K_{\tau {}0}(\mathrm {\Omega }) =: K_{\tau }(\mathrm {\Omega }) =\exp {[ i (2\tau {-}1) px /(2 \hbar ) ]} \end{aligned}$$

while relying on the parametrization with x and p. The resulting distribution functions \(F_{\rho }(\mathrm {\Omega },\tau )\) are in the Cohen class due to Theorem 1, and they are square integrable following Property 2 as the absolute value of \( K_{\tau }(\mathrm {\Omega })\) is bounded. We calculate the explicit action of the corresponding parity operator \(\mathrm {\Pi }_{\tau }\).

Theorem 2

The action of the \(\tau \)-parametrized parity operator \(\mathrm {\Pi }_\tau :=\mathrm {\Pi }_{\tau {}0}\) on some coordinate representation \(\psi (x) \in L^2(\mathbb {R})\) is explicitly given for any \(\tau \ne 1\) by

$$\begin{aligned} \mathrm {\Pi }_\tau \psi (x) = \frac{ 1}{2 |\tau {-}1|} \psi (\tfrac{\tau x }{ \tau {-} 1}), \end{aligned}$$
(41)

which for the special case \(\tau = 1/2\) reduces (as expected) to the usual parity operator \(\mathrm {\Pi }\). It follows that \(\mathrm {\Pi }_\tau \) is bounded for every \(0< \tau < 1\) (or in general for every real \(\tau \) that is not equal to 0 or 1) and its operator norm is given by \(\Vert \mathrm {\Pi }_\tau \Vert _\textrm{sup}= 1/\sqrt{4 (\tau {-} \tau ^2)}\).

Proof

By (39), the parity operator \(\mathrm {\Pi }_{\tau }\) acts on the coordinate representation \(\psi (x)\) via

$$\begin{aligned} \mathrm {\Pi }_{\tau } \psi (x) = (4 \pi \hbar )^{-1} \int K_{\tau {}0}(\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \psi (x) \, \textrm{d}\mathrm {\Omega }. \end{aligned}$$

This integral can be evaluated using the explicit form of \(K_{\tau {}0}(\mathrm {\Omega }) \) from (38), and the action of \(\mathcal {D}\) on coordinate representations \(\psi (x)\) from (11) yields

$$\begin{aligned} \int K_{\tau {}0}(\mathrm {\Omega }) \, \mathcal {D}(\mathrm {\Omega }) \psi (x) \, \textrm{d}\mathrm {\Omega }&= \int e^{\tfrac{i}{2 \hbar } (2\tau -1) p_0 x_0 } e^{\tfrac{i}{\hbar }(p_{0}x-\tfrac{1}{2}p_{0}x_{0})} \psi (x{-}x_{0}) \, \textrm{d}x_0 \, \textrm{d}p_0 \\&= \int \left[ \int e^{\tfrac{i}{\hbar } p_{0}(x+(\tau -1)x_{0})}\, \textrm{d}p_0 \right] \psi (x{-}x_{0}) \, \textrm{d}x_0 \\&= \frac{1}{|\tau {-}1|} \int \left[ \int e^{\tfrac{i}{\hbar } p_{0}y}\, \textrm{d}p_0 \right] \psi (\tfrac{\tau x {-} y}{ \tau {-} 1}) \, \textrm{d}y, \end{aligned}$$

where the change of variables \(y=x+(\tau {-}1)x_{0}\) with \(x_{0}=(y{-}x)/(\tau {-}1)\) and \(\textrm{d}y= |\tau {-}1|\,\textrm{d}x_0 \) was used. Therefore, the right-hand side is

$$\begin{aligned} \frac{ 2\pi \hbar }{|\tau {-}1|} \int \delta (y)\, \psi (\tfrac{\tau x {-} y}{ \tau {-} 1}) \, \textrm{d}y = \frac{ 2\pi \hbar }{|\tau {-}1|} \psi (\tfrac{\tau x }{ \tau {-} 1}). \end{aligned}$$

Now, let \(\tau \in (0,1)\). Recalling that the operator norm \(\Vert \mathrm {\Pi }_\tau \Vert _{\textrm{sup}}\) is calculated via \(\textrm{sup}_{\Vert \phi (x) \Vert _{L^2}=1} \Vert \mathrm {\Pi }_\tau \phi (x) \Vert _{L^2}\), for an arbitrary square-integrable \(\phi (x)\) with \(L^2\) norm \(\Vert \phi (x) \Vert _{L^2}=1\) one obtains

$$\begin{aligned} \Vert \mathrm {\Pi }_\tau \phi \Vert ^2_{L^2}&= \langle \mathrm {\Pi }_\tau \phi | \mathrm {\Pi }_\tau \phi \rangle = (2 |\tau {-}1|)^{-2} \int _{\mathbb {R}} \phi ^*(\tfrac{\tau x }{ \tau - 1}) \phi (\tfrac{\tau x }{ \tau - 1}) \, \textrm{d}x\\&=(2 |\tau {-}1|)^{-2} \tfrac{ | \tau - 1 |}{| \tau |} \Vert \phi (x) \Vert _{L^2} =\tfrac{ 1}{ 4 |\tau -1| | \tau |} = \tfrac{1 }{ 4 (\tau - \tau ^2)} \end{aligned}$$

by applying a change of variables. This results in \(\Vert \mathrm {\Pi }_\tau \Vert _{\textrm{sup}} = [ 4 (\tau {-} \tau ^2)]^{-1/2}\). \(\square \)

This parity operator is bounded for every \(0< \tau < 1\), and its expectation value gives rise to well-defined distribution functions (Property 1) which are also integrable as \(K_{\tau }(0)=1\) (Property 6). Note that this family of distribution functions \(F_\rho (\mathrm {\Omega },\tau ,0)\) for \(\tau \ne 1/2\) does not satisfy Property 5, i.e., self-adjoint operators \(\rho \) are mapped to complex functions. In particular, the symmetry \(K^*_{\tau }(\mathrm {\Omega }) = K_{1-\tau }(\mathrm {\Omega })\) implies that

$$\begin{aligned} F_{\rho ^\dagger }(\mathrm {\Omega },\tau ) = F^*_\rho (\mathrm {\Omega }, 1{-} \tau ) . \end{aligned}$$

In the following, we will rely on this \(\tau \)-parametrized family to construct and analyze the parity operator of the Born–Jordan distribution.

5 The Born–Jordan Distribution

5.1 Parity Operator Description of the Born–Jordan Distribution

The Born–Jordan distribution \(F_\rho (\mathrm {\Omega }, \text {BJ}) \) is an element of the Cohen class [17, 28, 43] and is obtained by averaging over the \(\tau \)-distributions \(F_{\rho }(\mathrm {\Omega },\tau ) \in L^2(\mathbb {R}^2)\):

$$\begin{aligned} F_\rho (\mathrm {\Omega }, \text {BJ}) := \int _0^1 F_\rho (\mathrm {\Omega },\tau ) \, \textrm{d}\tau . \end{aligned}$$
(42)

As in Definition 4, this distribution function is also obtained via the expectation value of a parity operator.

Theorem 3

The Born–Jordan distribution \(F_\rho (\mathrm {\Omega }, \text {BJ})\) of a density operator \(\rho \in \mathcal B^1(L^2(\mathbb R))\) is an element of the Cohen class, and it is obtained as the expectation value

$$\begin{aligned} F_\rho (\mathrm {\Omega },\text {BJ}) = (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) \mathrm {\Pi }_{\text {BJ}} \mathcal {D}^\dagger (\mathrm {\Omega }) ] \end{aligned}$$
(43)

of the (displaced) parity operator \(\mathrm {\Pi }_{\text {BJ}}\) that is defined by the relation

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} := (4 \pi \hbar )^{-1} \int K_{\text {BJ}}(\mathrm {\Omega }) \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }, \end{aligned}$$
(44)

where \(K_{\text {BJ}}(\mathrm {\Omega }) = \textrm{sinc}(a)=\sin (a)/a\) is the cardinal sine function with the argument \(a=(2\hbar )^{-1} \, px = i [(\alpha ^*)^2{-}\alpha ^2]/4\). Here, one applies the substitution \(\alpha = (\lambda x + i p/\lambda )/\sqrt{2 \hbar }\) from Sect. 3.1 and the expression for a is independent of \(\lambda \).

Proof

Combining Eqs. (42) and (27), the Born–Jordan distribution is the expectation value

$$\begin{aligned} F_\rho (\mathrm {\Omega }, \text {BJ}) = (\pi \hbar )^{-1} \, \textrm{Tr}\,[ \, \rho \, \mathcal {D}(\mathrm {\Omega }) (\int _0^1 \mathrm {\Pi }_{\tau {}0} \, \textrm{d}\tau ) \mathcal {D}^\dagger (\mathrm {\Omega }) ], \end{aligned}$$

and the corresponding parity operator can be expanded as

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}}= (4 \pi \hbar )^{-1} \int \left[ \int _0^1 K_{\tau {}0}(\mathrm {\Omega }) \, \textrm{d}\tau \right] \mathcal {D}(\mathrm {\Omega }) \, \textrm{d}\mathrm {\Omega }. \end{aligned}$$

Using the explicit form of \(K_{\tau {}0}(\mathrm {\Omega })\) from (38), the evaluation of the integral

$$\begin{aligned} \int _0^1 K_{\tau {}0}(\mathrm {\Omega }) \, \textrm{d}\tau = \int _0^1 \exp {[ \tfrac{i}{2 \hbar } (2\tau {-}1) px ]} \, \textrm{d}\tau = \textrm{sinc}[(2\hbar )^{-1} \, px] \end{aligned}$$

over \(\tau \) concludes the proof. \(\square \)

This confirms that the Born–Jordan distribution \(F_\rho (\mathrm {\Omega }, \text {BJ}) \in L^2(\mathbb {R}^2)\) is square integrable following Property 2 as the absolute value of its filter function is bounded, i.e., \(| \textrm{sinc}[(2\hbar )^{-1} \, px] \,| \le 1\) for all \((x,p)\in \mathbb R^2\). The filter function \(K_{\text {BJ}}\) satisfies \(K_{\text {BJ}}(x,0) =K_{\text {BJ}}(0,p)=1\), and the Born–Jordan distribution therefore gives rise to the correct marginals as quantum-mechanical probabilities (Property 7). In particular, integrating over the Born–Jordan distribution reproduces the quantum-mechanical probability densities, i.e., \( \int F_\rho (x,p,\text {BJ}) \,\textrm{d} x = |\psi (p)|^2 \) and \( \int F_\rho (x,p,\text {BJ}) \,\textrm{d} p = |\psi (x)|^2 \).

Most importantly the operator \(\mathrm {\Pi }_{\text {BJ}}\) is bounded, meaning Born–Jordan distributions are well defined and bounded for all quantum states, refer to Property 1. Also, the largest (generalized) eigenvalue of \(\mathrm {\Pi }_{\text {BJ}}\) is exactly \(\pi /2\) as shown in Theorem 5.

Proposition 2

The Born–Jordan parity operator \(\mathrm {\Pi }_{\text {BJ}}\) is bounded and an upper bound of its operator norm is given by \(\Vert \mathrm {\Pi }_{\text {BJ}}\Vert _\textrm{sup}\le \pi /2\).

Proof

Using the \(\mathrm {\Pi }_\tau \)-representation of \(\mathrm {\Pi }_{\text {BJ}}\), we compute

$$\begin{aligned}&\Vert \mathrm {\Pi }_{\text {BJ}}\psi (x)\Vert _{L^2}= \big \Vert \int _0^1 \mathrm {\Pi }_\tau \psi (x) \, \textrm{d}\tau \big \Vert _{L^2}\le \int _0^1 \Vert \mathrm {\Pi }_\tau \psi (x)\Vert _{L^2} \, \textrm{d}\tau \\&\le \int _0^1 \frac{1}{\sqrt{4(\tau {-}\tau ^2)}} \, \textrm{d}\tau \,\Vert \psi (x)\Vert _{L^2}=\frac{\pi }{2}\Vert \psi (x)\Vert _{L^2} \end{aligned}$$

for arbitrary \(\psi (x)\in L^2(\mathbb R)\). In the second-to-last step, we used Theorem 2. \(\square \)

It is well known that the Born–Jordan distribution is related to the Wigner function via a convolution with the Cohen kernel \(\theta _{\text {BJ}}\), refer to [43, 44]. However, calculating this kernel, or the corresponding parity operator directly might prove difficult. In the following, we establish a more convenient representation of the Born–Jordan parity operator which is an “average” of \(\mathrm {\Pi }_\tau \) from Theorem 2 via the formal integral transformation

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} = \int _0^1 \mathrm {\Pi }_\tau \, \textrm{d}\tau , \end{aligned}$$
(45)

which—as in Sect. 3.2—is interpreted as \(\mathrm {\Pi }_{\text {BJ}} \psi (x) = \int _0^1 \mathrm {\Pi }_\tau \psi (x) \, \textrm{d}\tau \) for all \(\psi (x)\in L^2(\mathbb R)\). Recall that the parity operator \(\mathrm {\Pi }_\tau \) is well defined and bounded for every \(0< \tau < 1\).

Remark 4

Note that evaluating \(\mathrm {\Pi }_{\text {BJ}}\psi (x)\) at \(x=0\) for some \(\psi (x)\in L^2(\mathbb R)\) with \(\psi (0)\ne 0\) leads to a divergent integral in (45). This comes from the singularity at \(\tau =1\) in (41). However, we will later see that this is harmless as it only happens on a set of measure zero (so one can define \(\mathrm {\Pi }_{\text {BJ}}\psi (x)|_{x=0}\) to be 0 or \(\psi (0)\) or arbitrary) and, more importantly, that \(\psi (x)\in L^2(\mathbb R)\) implies \(\mathrm {\Pi }_{\text {BJ}}\psi (x)\in L^2(\mathbb R)\) (Proposition 2).

Following [27] and Chapter 2.3 in [88], the squeezing operator is defined to be

$$\begin{aligned} S(\xi ) := \exp [ \tfrac{i \xi }{2} ( \hat{p} \hat{x} {+} \hat{x}\hat{p}) ]= \exp [ \tfrac{i \xi }{2} ( \hat{a}^2 {-} (\hat{a}^\dagger )^2) ] \end{aligned}$$
(46)

(where \(\xi \in \mathbb R\)), and it acts on a coordinate representation via \(S(\xi ) \psi (x) = e^{\xi /2} \psi (e^{\xi } x)\).

Theorem 4

The Born–Jordan parity operator

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} = [\, \tfrac{1}{4} \int _{-\infty }^{\infty } \, \textrm{sech}(\xi /2) S(\xi ) \,\textrm{d} \xi \,] \; \mathrm {\Pi }\end{aligned}$$
(47)

is a composition of the reflection operator \(\mathrm {\Pi }\) followed by a squeezing operator (and the two operations commute), and this expression is integrated with respect to a well-behaved weight function \(\textrm{sech}(\xi /2) =2/(e^{\xi /2}{+}e^{-\xi /2})\). Note that the function \(\textrm{sech}(\xi /2) \in \mathcal {S}(\mathbb {R}) \) is fast decreasing and invariant under the Fourier transform (e.g., as Hermite polynomials).

Proof

The explicit action of \(\mathrm {\Pi }_{\text {BJ}}\) on a coordinate representation \(\psi (x) \in L^2(\mathbb {R})\) is given by (see Theorem 2)

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \psi (x) = \int _0^1 \mathrm {\Pi }_\tau \psi (x) \, \textrm{d}\tau = \int _0^1 \frac{ 1}{2 |\tau {-}1|} \psi (\tfrac{\tau x }{ \tau - 1}) \, \textrm{d}\tau . \end{aligned}$$

Applying a change of variables \(e^{\xi } = \tau /(1{-}\tau )\) with \(\xi \in \mathbb {R}\) yields the substitutions \(\tau = 1/(1 {+} e^{-\xi })\), \({1}/(2 |\tau {-}1|) = (1 {+} e^\xi )/2\), and \(\textrm{d} \tau = e^\xi /(1{+}e^\xi )^2\,\textrm{d} \xi \). One obtains

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \psi (x) = \int _{-\infty }^{\infty } \frac{e^{\xi }}{2(1{+}e^{\xi })} \psi ( - e^{\xi } x) \,\textrm{d} \xi . \end{aligned}$$

Let us recognize that \(\psi ( - e^{\xi } x)= e^{- \xi /2} S(\xi )\, \mathrm {\Pi }\, \psi (x)\) is the composition of a coordinate reflection and a squeezing of the pure state \(\psi (x)\); also, the two operations commute. This results in the explicit action

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \psi (x) = \int _{-\infty }^{\infty } \frac{e^{\xi /2} }{2(1{+}e^{\xi })} S(\xi )\, \mathrm {\Pi }\, \psi (x) \,\textrm{d} \xi , \end{aligned}$$

where \( e^{\xi /2} /[ 2(1{+}e^{\xi }) ] = [ 2(e^{-\xi /2} {+}e^{\xi /2}) ]^{-1}\) concludes the proof. \(\square \)

The expression for the parity operator in Theorem 4 is very instructive when compared to Theorem 3, and this confirms that the parity operator \(\mathrm {\Pi }_{\text {BJ}}\) decomposes into the usual parity operator \(\mathrm {\Pi }\) followed by a geometric transformation, refer also to Sect. 5.3. In the case of the Born–Jordan parity operator this geometric transformation is an average of squeezing operators.

Remark 5

The Born–Jordan distribution is covariant under squeezing, which means that the squeezed density operator \(\rho ' = S(\xi ) \rho S^\dagger (\xi )\) is mapped to the inversely squeezed phase-space representation \(F_{\rho '}(x,p, \text {BJ}) = F_{\rho }(e^{-\xi } x, e^{\xi } p, \text {BJ}) \).

The form of \(\mathrm {\Pi }_{\text {BJ}}\) given in (47) allows for an alternative proof of Proposition 2:

Remark 6

Recalling \( \Vert \mathrm {\Pi }_{\text {BJ}}\Vert _\textrm{sup}= \sup _{\Vert \phi (x)\Vert _{L^2}=1} \Vert \mathrm {\Pi }_{\text {BJ}} \phi (x) \Vert _{L^2}, \) the norm of the function \(\mathrm {\Pi }_{\text {BJ}} \phi (x)\) for any \(\phi (x)\) with \(\Vert \phi (x)\Vert _{L^2}=1\) can be expressed as

$$\begin{aligned} \Vert \mathrm {\Pi }_{\text {BJ}} \phi (x) \Vert _{L^2}^2&= \langle \phi | \mathrm {\Pi }_{\text {BJ}}^\dagger \mathrm {\Pi }_{\text {BJ}} | \phi \rangle \\&=\tfrac{1}{16} \iint \textrm{sech}(\xi /2) \textrm{sech}(\xi '/2) \langle \phi | S(\xi '{-}\xi ) | \phi \rangle \,\textrm{d} \xi \, \textrm{d} \xi '\\&\le \tfrac{1}{16} \iint \textrm{sech}(\xi /2) \textrm{sech}(\xi '/2) | \langle \phi | S(\xi '{-}\xi ) | \phi \rangle | \,\textrm{d} \xi \, \textrm{d} \xi ', \end{aligned}$$

and it was used that \(\mathrm {\Pi }^\dagger \mathrm {\Pi }=\mathrm {\Pi }^2 =1\) and \(S^\dagger (\xi ) S(\xi ') = S(\xi '{-}\xi ) \). Since \(S(\xi '{-}\xi ) \) is unitary, one obtains that \(|\langle \phi | S(\xi '{-}\xi ) \phi \rangle | \le 1\). Finally,

$$\begin{aligned} \quad \Vert \mathrm {\Pi }_{\text {BJ}} \phi \Vert _{L^2}^2 \le \tfrac{1}{16} \iint \textrm{sech}(\xi /2) \textrm{sech}(\xi '/2) \,\textrm{d} \xi \, \textrm{d} \xi ' = \pi ^2/4. \end{aligned}$$

5.2 Spectral Decomposition of the Born–Jordan Parity Operator

We will now adapt results for generalized spectral decompositions, refer to [26, 27, 57, 95]. This will allow us to solve the generalized eigenvalue equation for parity operators and to determine their spectral decompositions.

Recall the following from Sect. 2.1: The distributional pairing for smooth, well-behaved functions \(\psi (x) \in \mathcal {S}(\mathbb R)\) in with respect to tempered distributions \(a\in \mathcal S'(\mathbb R)\) (such as functions of slow growth a(x)) extends to \(L^2\)-scalar products of the form \( \langle a, \psi \rangle = \int _{\mathbb R}a^*(x)\psi (x)\,\textrm{d}x\), which corresponds to a rigged Hilbert space [26, 57] or the Gelfand triple \(\mathcal {S}(\mathbb {R})\subset L^2 (\mathbb {R}) \subset \mathcal S'(\mathbb R)\). This rigged Hilbert space allows us to specify the generalized spectral decomposition of the Born–Jordan parity operator with generalized eigenvectors in \(\mathcal S'(\mathbb R)\) as functions of slow growth.

It was shown in the previous section that the Born–Jordan parity operator \(\mathrm {\Pi }_{\text {BJ}}\) is a composition of a coordinate reflection and a squeezing operator. We now recapitulate results on the spectral decomposition of the squeezing operator from [19, 26, 27], up to minor modifications. Recall that the squeezing operator forms a unitary, strongly continuous one-parameter group \(S(\xi ) = e^{-i \xi {H}}\) with \(\xi \in \mathbb R\) that is generated by the (unbounded) self-adjoint Hamiltonian

$$\begin{aligned} {H}:= - \tfrac{1}{2} (\hat{x} \hat{p} {+} \hat{p} \hat{x}) = - \tfrac{1}{2} [ \hat{a}^2 {-} (\hat{a}^\dagger )^2]. \end{aligned}$$

This Hamiltonian admits a purely continuous spectrum \(E \in \mathbb {R}\) and satisfies generalized eigenvalue equations

$$\begin{aligned} \langle H \psi | \psi ^E_{\pm } \rangle = \langle \psi | H \psi ^E_{\pm } \rangle = E \langle \psi | \psi ^E_{\pm } \rangle \end{aligned}$$

for every \(\psi \in \mathcal {S}(\mathbb {R})\), where the last equation is equivalent to \( H | \psi ^E_{\pm } \rangle = E | \psi ^E_{\pm }\rangle \). The Gelfand–Maurin spectral theorem [26, 57, 95] results in a spectral resolution of

$$\begin{aligned} H = \int _{-\infty }^{\infty } E \, | \psi ^E_{\pm } \rangle \langle \psi ^E_{\pm } | \, \textrm{d} E. \end{aligned}$$

Here, the generalized eigenvectors are specified in terms of their coordinate representations as slowly increasing functions, i.e., \(\psi _\pm ^E(x):= \langle x | \psi ^E_\pm \rangle \in {S}'(\mathbb {R})\) with

$$\begin{aligned} \psi _+^E(x) = \tfrac{1}{2\sqrt{\pi }} |x|^{-(iE + \tfrac{1}{2})} \; \text { and } \; \psi _-^E (x) = \tfrac{1}{2\sqrt{\pi }} \textrm{sgn}(x) |x|^{-(iE + \tfrac{1}{2})}, \end{aligned}$$
(48)

refer to [26, 27, 57, 95] and “Appendix F” for more details. Note that \(\psi _\pm ^E(x)\) are generalized eigenfunctions: They are not square integrable, but the integral \(\int _{\mathbb R} [\psi _\pm ^E(x)]^* \phi (x)\,\textrm{d}x\) exists as a distributional pairing for every \( \phi \in \mathcal S(\mathbb R)\). Also, note that these generalized eigenvectors can be decomposed into the number-state basis with finite expansion coefficients that decrease to zero for large n, refer to “Appendix F.” The spectral decomposition of the squeezing operator is then given by

$$\begin{aligned} S(\xi ) = \int _{-\infty }^{\infty } e^{-i E \xi } [ \, | \psi ^E_{+} \rangle \langle \psi ^E_{+} | +| \psi ^E_{-} \rangle \langle \psi ^E_{-} | \, ] \, \textrm{d} E, \end{aligned}$$

refer to Eq. 6.12 in [26] and Eq. 2.14 in [27]. Note that these eigenvectors are also invariant under the Fourier transform (e.g., as Hermite polynomials).

It immediately follows that the squeezing operator satisfies the generalized eigenvalue equation

$$\begin{aligned} S(\xi ) | \psi ^E_{\pm } \rangle = e^{-i \xi E} \, | \psi ^E_{\pm } \rangle , \end{aligned}$$
(49)

which can be easily verified using the explicit action \(S(\xi ) \psi _\pm ^E(x) = e^{\xi /2} \psi _\pm ^E(e^{\xi } x) = e^{-iE \xi } \psi _\pm ^E(x)\). One can now specify the Born–Jordan parity operator using its spectral decomposition.

Theorem 5

Generalized eigenvectors \(| \psi ^E_{\pm } \rangle \) of the squeezing operator from (48) are also generalized eigenvectors of the Born–Jordan parity operator which satisfy

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \, | \psi ^E_{\pm } \rangle = \pm \tfrac{\pi }{2} \, \textrm{sech}(\pi E) \, | \psi ^E_{\pm } \rangle \end{aligned}$$

for all \(E\in \mathbb R\). The parity operator \(\mathrm {\Pi }_{\text {BJ}}\) therefore admits the spectral decomposition

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} = \tfrac{\pi }{2} \, \int _{-\infty }^{\infty } \textrm{sech}(\pi E) \, [ \, | \psi ^E_{+} \rangle \langle \psi ^E_{+} | - | \psi ^E_{-} \rangle \langle \psi ^E_{-} | \, ] \, \textrm{d}E, \end{aligned}$$

where \( \langle \psi ^E_{\pm } | \mathrm {\Pi }= \pm \langle \psi ^E_{\pm } |\) has been used.

Proof

The generalized eigenvalues can be computed via

$$\begin{aligned} \mathrm {\Pi }_{\text {BJ}} \, | \psi ^E_{\pm } \rangle = \tfrac{1}{4} \int _{-\infty }^{\infty } \, \textrm{sech}(\xi /2) S(\xi ) \,\textrm{d} \xi \; \mathrm {\Pi }\, | \psi ^E_{\pm } \rangle , \end{aligned}$$

where \(\mathrm {\Pi }\, | \psi ^E_{\pm } \rangle = \pm | \psi ^E_{\pm } \rangle \). Using (49), one obtains

$$\begin{aligned} \quad \mathrm {\Pi }_{\text {BJ}} \, | \psi ^E_{\pm } \rangle = \pm \tfrac{ 1}{4} \int _{-\infty }^{\infty } \, \textrm{sech}(\xi /2) e^{-iE \xi } \,\textrm{d} \xi \; | \psi ^E_{\pm } \rangle = \pm \tfrac{\pi }{2} \, \textrm{sech}(\pi E) \, | \psi ^E_{\pm } \rangle . \end{aligned}$$

\(\square \)

Remark 7

Recall that \(\mathrm {\Pi }_{\text {BJ}}\) is a bounded (by Proposition 2) and self-adjoint operator. Consequently, the usual spectral theorem in multiplication operator form [66, Thm. 7.20] yields a \(\sigma \)-finite measure space \((X,\mu )\), a bounded, measurable, real-valued function h on X, and unitary \(U:L^2(\mathbb R)\rightarrow L^2(X,\mu )\) such that

$$\begin{aligned}{}[U\mathrm {\Pi }_{\text {BJ}}U^{-1}(\psi )](\lambda )=h(\lambda )\psi (\lambda ) \end{aligned}$$

for all \(\psi \in L^2(X,\mu )\) and \(\lambda \in X\). While this undoubtedly is a nice representation, the spectral decomposition in Theorem 5 is more readily determined with the help of the Gelfand–Maurin spectral theorem [26, 57, 95]. In particular, the said theorem lets us directly work with the generalized eigenfunctions in Eq. (48), even though they are not square integrable.

5.3 Geometric Interpretation of Parity Operators

While, above, we have comprehensively explored analytic properties of the Born–Jordan and other practically important parity operators, here we relate these mathematical objects to geometric transformations. Even the rather complex Born–Jordan parity operator admits a surprisingly simple decomposition into two elementary geometric transformations. Equation (47) decomposes the Born–Jordan parity operator into an ordinary reflection of the wave function’s coordinate followed by a weighted average of squeezing operations as \(\mathrm {\Pi }_{\text {BJ}} =[\, \tfrac{1}{4} \int _{-\infty }^{\infty } \, \textrm{sech}(\xi /2) S(\xi ) \,\textrm{d} \xi \,] \, \mathrm {\Pi }. \) As such, the action on any wave function \(\psi (x) \in L^2(\mathbb {R})\) can be summarized as the reflected, squeezed function \(S(\xi )\mathrm {\Pi }\, \psi (x) = e^{\xi /2} \psi (- e^{\xi } x)\) averaged over all parameters \(\xi \in (-\infty , \infty )\) with respect to the rapidly decaying weight function \(\textrm{sech}(\xi /2)\).

It is not only the Born–Jordan parity operator that admits a simple geometric interpretation; rather this seems to hint at a universal property, at least in the classes of practically important phase-space representations. In particular, we now state that both \(\mathrm {\Pi }_{\tau }\) and the pivotal parity operator \(\mathrm {\Pi }_s\)—which contains the most popular variants of Wigner, Husimi, and Glauber P phase-space functions as special cases—can be decomposed into elementary geometric transformations.

Remark 8

Applying the substitution \(e^{\xi }:= \tau /(1{-}\tau )\), the parity operator \(\mathrm {\Pi }_\tau \) from (41) can be decomposed for \(0<\tau < 1\) into

$$\begin{aligned} \mathrm {\Pi }_\tau = \cosh (\xi /2) S(\xi )\, \mathrm {\Pi }\end{aligned}$$

which consists of a coordinate reflection and a squeezing.

Consequently, the parity operator \(\mathrm {\Pi }_\tau \) admits a spectral decomposition

$$\begin{aligned} \mathrm {\Pi }_\tau = \cosh (\xi /2) \, \int e^{-i E \xi } [ \, | \psi ^E_{+} \rangle \langle \psi ^E_{+} | -| \psi ^E_{-} \rangle \langle \psi ^E_{-} | \, ] \, \textrm{d} E, \end{aligned}$$

where \( \langle \psi ^E_{\pm } | \mathrm {\Pi }= \pm \langle \psi ^E_{\pm } |\) has been used.

Remark 9

The parity operator with \(\kappa _s:=\ln [(1{+}s)/(1{-}s)]\) and \(-1<s < 1\) is the composition

$$\begin{aligned} \mathrm {\Pi }_s=(1{-}s)^{-1} e^{\kappa _s \hat{a}^\dagger \hat{a}} \, \mathrm {\Pi }\end{aligned}$$

from (40) of the usual coordinate reflection \(\mathrm {\Pi }\) followed by a positive semi-definite operator. In particular, \(\mathrm {\Pi }_{-1}= \tfrac{1}{2} \, | 0 \rangle \langle 0 | \, \mathrm {\Pi }\). Note that the positive semi-definite operator \(e^{\kappa _s \hat{a}^\dagger \hat{a}}\) describes the effective phenomenon of photon loss for \(s < 0\), refer to [89]. Of course domain restrictions might need to be considered for \(s>0\) as discussed earlier.

6 Explicit Matrix Representation of the Born–Jordan Parity Operator

Recall that the s-parametrized parity operators \(\mathrm {\Pi }_s\) are diagonal in the Fock basis and their diagonal entries can be computed using the simple expression in (40). This enables the experimental reconstruction of distribution functions from photon-count statistics [7, 14, 46, 93] in quantum optics.

Remark 10

The Born–Jordan parity operator \(\mathrm {\Pi }_{\text {BJ}}\) is not diagonal in the number-state basis, as its filter function \(K_{\text {BJ}}(\mathrm {\Omega })\) is not invariant under arbitrary phase-space rotations, refer to Property 4. The filter function \(K_{\text {BJ}}(\mathrm {\Omega })\) is, however, invariant under \(\pi /2\) rotations in phase space, and therefore, only every fourth off-diagonal is nonzero.

We now discuss the number-state representation of the parity operator \(\mathrm {\Pi }_{\text {BJ}}\), which provides a convenient way to calculate (or, more precisely, approximate) Born–Jordan distributions.

Theorem 6

The matrix elements \([ \mathrm {\Pi }_{\text {BJ}} ]_{mn}:= \langle m | \mathrm {\Pi }_{\text {BJ}}\, n \rangle \) of the Born–Jordan parity operator in the Fock basis can be calculated in the form of a finite sum

$$\begin{aligned} {[} \mathrm {\Pi }_{\text {BJ}} ]_{mn} = \sum _{k=0}^n \sum _{\begin{array}{c} \ell =0 \\ \ell \text { even} \end{array} }^{m-n} d_{mn}^{k\ell } \; \mathrm {\Phi }_{(m-n-\ell )/2,\ell /2}^k \, , \end{aligned}$$
(50)

for \(m \in \{n, n+4, n+8, \dots \}\) and \([ \mathrm {\Pi }_{\text {BJ}} ]_{mn} = [ \mathrm {\Pi }_{\text {BJ}} ]_{nm}\) with the coefficients

$$\begin{aligned} d_{mn}^{k\ell }&:= (-1)^{(\ell +m-n)/2} \sqrt{\frac{n!}{m!}} 2^{(2 k+m-n)/2} \left( {\begin{array}{c}m\\ n{-}k\end{array}}\right) \left( {\begin{array}{c}m{-}n\\ \ell \end{array}}\right) /k!, \end{aligned}$$
(51)
$$\begin{aligned} \mathrm {\Phi }_{ab}^k&:= [ \partial ^k_{\mu } [\partial ^{a}_{\lambda }\partial ^{b}_{\mu } f(\lambda ,\mu )]|_{\lambda =\mu } ]|_{\mu =1}. \end{aligned}$$
(52)

Here, \(\mathrm {\Phi }_{ab}^k\) denotes the ath and bth partial derivatives of the function \(f(\lambda ,\mu ) = \textrm{arcsinh}[1/\sqrt{\lambda \mu }]\) with respect to its variables \(\lambda \) and \(\mu \), respectively, evaluated at \(\lambda =\mu \), then differentiated again k times and finally its variable is set to \(\mu =1\).

Refer to “Appendix G” for a proof. The derivatives in (52) can be calculated in the form of a finite sum

$$\begin{aligned} \mathrm {\Phi }_{ab}^k = \sqrt{2} (-1)^{a+b} \, 2^{-2a -2b -k} \sum _{j=0}^{a+b+k-1} \xi _j^{a{}b{}k}, \end{aligned}$$
(53)

where \(a+b+k\ge 1\) and \(\xi _j^{a{}b{}k}\) are recursively defined integers, refer to (68) in “Appendix H.” Substituting \(\ell \) for \(2\ell \) in (50), the matrix elements \([ \mathrm {\Pi }_{\text {BJ}}]_{mn}\) then depend only on these integers \(\xi _j^{a{}b{}k}\) via the finite sum

$$\begin{aligned} {[} \mathrm {\Pi }_{\text {BJ}}]_{mn} = \gamma _{mn} \sum _{k=0}^n \sum _{\ell =0}^{(m-n)/2} [\sum _j \xi _j^{a{}\ell {}k} ] \, (-1)^{\ell +m -n} \left( {\begin{array}{c}m\\ n{-}k\end{array}}\right) \left( {\begin{array}{c}m{-}n\\ 2 \ell \end{array}}\right) /k! , \end{aligned}$$

where \(\gamma _{mn}:= 2^{[-(m-n)+1]/2} \sqrt{{n!}/{m!}} \) for \(m \in \{n, n{+}4, n{+}8, \dots \}\) and \(a = (m{-}n{-} 2\ell )/2\).

Fig. 1
figure 1

a The first \(8 \times 8\) matrix elements of the Born–Jordan parity operator from Eq. (50) where \(\textrm{ash}\) denotes \(\textrm{arcsinh}(1)\). b The first \(6 \times 6\) elements of the recursive sequence that defines the Born–Jordan parity operator via Eq. (55). Orange and blue colors represent positive and negative values, respectively, while the color intensity reflects the absolute value of the corresponding numbers

Figure 1a shows the first \(8 \times 8\) entries of \([ \mathrm {\Pi }_{\text {BJ}}]\). One observes the following structure: Only every fourth off-diagonal is nonzero, the matrix is real and symmetric, and the entries along every diagonal and off-diagonal decrease in their absolute value. In particular, the diagonal elements of \(\mathrm {\Pi }_{\text {BJ}}\) admit the following special property.

Proposition 3

For every \(n\in \lbrace 0,1,2,\ldots \rbrace \), the diagonal entries of \(\mathrm {\Pi }_{\text {BJ}}\) in the Fock basis are

$$\begin{aligned} {[} \mathrm {\Pi }_{\text {BJ}}]_{nn}={\text {arcsinh}}(1)-\sqrt{2} \sum _{k=0}^{n-1} \frac{(-1)^k}{k{+}1} \sum _{m=0}^{\lfloor \frac{k}{2}\rfloor } \left( {\begin{array}{c}2m\\ m\end{array}}\right) \Big (\frac{-1}{4}\Big )^{m}. \end{aligned}$$
(54)

In particular, \([ \mathrm {\Pi }_{\text {BJ}}]_{nn}\rightarrow 0\) as \(n\rightarrow \infty \). For a proof, we refer to “Appendix I.”

Note that the sum of these decreasing diagonal entries results in a trace \(\textrm{Tr}[\mathrm {\Pi }_{\text {BJ}}] = 1/2\) (Property 6) in the number-state basis. However, this trace does not necessarily exist in an arbitrary basis, as \(\mathrm {\Pi }_{\text {BJ}}\) is not a trace-class operator.

Remark 11

Let us emphasize that boundedness of \(\mathrm {\Pi }_{\text {BJ}}\) (Proposition 2) guarantees that using a (large enough) finite block of \(\mathrm {\Pi }_{\text {BJ}}\) for computations yields a good approximation of \(F_\rho (\mathrm {\Omega }, \text {BJ})\), refer to “Appendix C” for details.

In the following, we specify a more convenient form for the calculation of these matrix elements, i.e., a direct recursion without summation, which is based on the following conjecture (see “Appendix J”).

Conjecture 1

The nonzero matrix elements

$$\begin{aligned} {[} \mathrm {\Pi }_{\text {BJ}} ]_{k+4\ell ,k} = [ \mathrm {\Pi }_{\text {BJ}} ]_{k, k+4\ell } = \Gamma _{k\ell }\, [ M_{k\ell } + \delta _{\ell 0} \, \textrm{arcsinh}(1)/\sqrt{2} ], \end{aligned}$$
(55)

of the Born–Jordan parity operator are determined by a set of rational numbers \(M_{k\ell }\) where \(\Gamma _{k\ell } = 2^{-2\ell + {1}/{2}} \sqrt{{k!}/{(k{+}4\ell )!}}\) and \(k,\ell \in \{0,1,2, \dots \}\). The indexing is specified relative to the diagonal (where \(\ell =0\)) and \(\delta _{\ell {}m}\) is the Kronecker delta. The rational numbers \( M_{k\ell } \) can be calculated recursively using only 8 numbers as initial conditions, refer to “Appendix J” for details. Unlike (50), this form does not require a summation.

Figure 1b shows the first \(6 \times 6\) elements of the recursive sequence of rational numbers \(M_{k\ell }\). The first column of \(M_{k0}\) corresponds to the diagonal of the matrix \([\mathrm {\Pi }_{\text {BJ}}]_{mn}\) from Fig. 1a. For example, for \(k=5\) one obtains \(M_{5,0} = -{43}/{60}\), which corresponds to \([\mathrm {\Pi }_{\text {BJ}}]_{5,5}=\Gamma _{5,0} [ M_{5,0} + \delta _{0,0} \, \textrm{arcsinh}(1)/\sqrt{2} ]\) and \(\Gamma _{5,0}= \sqrt{2}\), and therefore \([\mathrm {\Pi }_{\text {BJ}}]_{5,5}=-{ 43 \sqrt{2}}/{60} +\textrm{arcsinh}(1) \) as detailed in Fig. 1a.

The direct recursion in Conjecture 1 enables us to conveniently and efficiently calculate the matrix elements \([\mathrm {\Pi }_{\text {BJ}}]_{mn}\), and we have verified the correctness of this approach for up to 6400 matrix elements, i.e., by calculating a matrix representation of size \(80 \times 80\). This facilitates an efficient calculation and plotting of Born–Jordan distributions for harmonic oscillator systems, such as in quantum optics [60, 88, 94]. Note that a recursively calculated \(80 \times 80\) matrix representation, which we have verified with exact calculations, is sufficient for most physical applications, i.e., Figs. 2 and 3 were calculated using \(30 \times 30\) matrix representations. However, a matrix representation of size \(2000 \times 2000\) can be easily calculated on a current notebook computer using the recursive method. Numerical evidence shows that the matrix representation of \(\mathrm {\Pi }_{\text {BJ}}\) can be well-approximated by a low-rank matrix, i.e., diagonalizing the matrix \(\mathrm {\Pi }_{\text {BJ}}\) reveals only very few significant eigenvalues. In particular, the sum of squares of the first 9 eigenvalues corresponds to approximately \(99.97\%\) of the sum of squares of all the eigenvalues of a \(2000 \times 2000\) matrix representation.

Fig. 2
figure 2

a, b Wigner, Born–Jordan, and Husimi Q phase-space plots of the number states \(|0\rangle \) and \(|1\rangle \) (which are eigenstates of the quantum harmonic oscillator). c Corresponding phase-space plots of the Schrödinger cat state \(|\text {cat}\rangle = (|0\rangle {+} |1\rangle )/\sqrt{2}\) (which is a superposition of the previous states)

7 Example Quantum States

Matrix representations of parity operators are used to conveniently calculate phase-space representations for bosonic quantum states via their associated Laguerre polynomial decompositions. The s-parametrized distribution functions of Fock states \(|n\rangle \) are sums

$$\begin{aligned} F_{|n\rangle } (\mathrm {\Omega },s) = \sum _{\mu =0} |[\mathcal {D}(\mathrm {\Omega })]_{n \mu }|^2 \, [\mathrm {\Pi }_s]_{\mu \mu } \end{aligned}$$

of the associated Laguerre polynomials from (10), which are weighted by their parity operator elements. The corresponding phase-space functions are radially symmetric as \(|[\mathcal {D}(\mathrm {\Omega })]_{n \mu }|^2\) depends only on the radial distance \(x^2+p^2\). The Wigner functions in Fig. 2a, b are radially symmetric and show strong oscillations, which are sometimes regarded as a quantum-mechanical feature [88].

Fig. 3
figure 3

a Wigner and Born–Jordan phase-space plots of the number state \(|4\rangle \). b The Born–Jordan distribution is decomposed into functions that correspond to diagonal and off-diagonal entries of the parity operator matrix

In contrast, the Born–Jordan parity operator is not diagonal in the number-state representation and it can be written in terms of projectors as \( \mathrm {\Pi }_{\text {BJ}}= \sum _{\mu =0}[ \mathrm {\Pi }_{\text {BJ}} ]_{\mu \mu } |\mu \rangle \langle \mu |+ \sum _{\mu =0} \sum _{\nu =1} [ \mathrm {\Pi }_{\text {BJ}} ]_{\mu , 4\nu } \, ( |\mu {+} 4 \nu \rangle \langle \mu | + |\mu \rangle \langle \mu {+} 4 \nu | )\). The Born–Jordan distribution of number states \(|n\rangle \) is given by

$$\begin{aligned} F_{|n\rangle } (\mathrm {\Omega },\text {BJ})&= \sum _{\mu =0} |[\mathcal {D}]_{n \mu }|^2 \, [\mathrm {\Pi }_{\text {BJ}} ]_{\mu \mu } \end{aligned}$$
(56)
$$\begin{aligned}&+ \sum _{\mu =0} \sum _{\nu =1} 2 \, \Re ( [\mathcal {D}]_{n \mu } ([\mathcal {D}]_{n, \mu + 4\nu })^* ) [ \mathrm {\Pi }_{\text {BJ}} ]_{\mu ,\mu +4\nu }. \end{aligned}$$
(57)

The Born–Jordan distribution of coherent states, i.e., the displaced vacuum states, closely matches the Wigner functions, see Fig. 2a. The first part in Eq. (56) contains the diagonal elements of the parity operator which correspond to the radially symmetric part of \(F_{|n\rangle } (\mathrm {\Omega },\text {BJ})\), see Fig. 3b (left). The second part in Eq. (57) results in a radially non-symmetric function, see Fig. 3b (right). The radially symmetric parts are quite similar to the Wigner function and have \(n+1\) wave fronts enclosed by the Bohr–Sommerfeld band [48, 88], i.e., the ring with radius \(\sqrt{2n{+}1}\). The radially non-symmetric functions have \(n+1\) local maxima along the outer squares, i.e., along phase-space cuts at the Bohr–Sommerfeld distance \(x,p \propto \sqrt{2n{+}1}\). The sum of these two contributions is the Born–Jordan distribution, and it is not radially symmetric for number states, see Fig. 3a.

8 Conclusion

We have introduced parity operators \(\mathrm {\Pi }_\theta \) which give rise to a rich family of phase-space distribution functions of quantum states. These phase-space functions have been previously defined in terms of convolutions, integral transformations, or Fourier transformations. Our approach using parity operators is both conceptually and computationally advantageous and now allows for a direct calculation of phase-space functions as quantum-mechanical expectation values. This approach therefore averts the necessity for the repeated and expensive computation of Fourier transformations. We motivate the name “parity operator” by the fact that parity operators \(\mathrm {\Pi }_\theta = A_\theta \circ \mathrm {\Pi }\) are composed of the usual parity operator and some specific geometric or physical transformation. We detailed the explicit form of parity operators for various phase spaces and, in particular, for the Born–Jordan distribution. We have also obtained a generalized spectral decomposition of the Born–Jordan parity operator, proved its boundedness, and explicitly calculated its matrix representation in the number-state basis. We conjecture that these matrix elements are determined by a proposed recursive scheme which allows for an efficient computation of Born–Jordan distribution functions. Moreover, large matrix representations of the Born–Jordan parity operator can be well approximated using rank-9 matrices. All this will be useful to connect our results with applications in (e.g.,) quantum optics, where techniques such as squeezing operators and the number-state representation are widely used.