# Stimulated parametric down-conversion for spatiotemporal metrology

Filippus S. Roux<EMAIL_ADDRESS>National Metrology Institute of South
Africa, Meiring Naudé Road, Brummeria 0040, Pretoria, South Africa

###### Abstract

A detailed analysis of the stimulated parametric down-conversion (PDC) process
is performed to investigate the effects of the spatiotemporal degrees of
freedom. The analysis provides information that would be useful for PDC-based
metrology applications. Using a Wigner functional approach, we obtain the
parametric down-converted state as the Bogoliubov transformed input state, in
terms of Bogoliubov kernel functions. The result is used to consider the case
for a coherent state seeding stimulated PDC. We also compute the background
which is obtained from spontaneous PDC.

##  I Introduction

Parametric down-conversion (PDC) Hong and Mandel (1985) has become one of the
most widely studied optical phenomena, thanks to its relevance for photonic
quantum information systems. Various quantum information applications are
based on the PDC process, including the preparation of entangled quantum
states Mair et al. (2001); Torres et al. (2003); Miatto et al. (2011); Walborn
et al. (2010); Roux (2020a), preparation of squeezed states Lvovsky (2015);
Kolobov (1999); Wasilewski et al. (2006); Roux (2021), and SU(1,1)
interferometry Yurke et al. (1986); Seyfarth et al. (2020); Ferreri et al.
(2021), to name a few.

Here, we focus on metrology applications. Quantum metrology Giovannetti et al.
(2011); Polino et al. (2020) is usually associated with estimating quantities
better than the standard quantum limit Giovannetti et al. (2004). However, the
association between metrology and quantum physics goes beyond that. While
quantum information science can be applied to enhance metrology, standard
metrology can also enhance quantum information technology. The latter may play
a more significant role in the near future than the former. The reason is that
many of the subsystems in quantum information systems require
characterization, calibration and standardization to ensure the successful
operation of the complete system. Various techniques in quantum optics such as
PDC can be used to provide such metrology applications.

One such application is a scheme to measure the absolute intensity of a
source. Such a source would be used as the seed optical field in a stimulated
PDC experiment, from which the output intensity would be measured. Together
with the measured output intensity from the spontaneous PDC in the same
experiment — the output intensity without the seed field — one should be able
to determine the absolute intensity of the seed field. What makes this
proposal challenging is the fact that the outputs that are obtained in these
two scenarios give very different intensity distributions that depend on the
spatiotemporal degrees of freedom in the experiment. Here, we investigate the
output intensity distributions that are obtained from stimulated PDC
experiments to determine how they depend on the various degrees of freedom.

PDC has been analyzed with diverse formalisms. While many investigations focus
on individual photons Mair et al. (2001); Arnaut and Barbosa (2000); Brambilla
et al. (2004, 2010); Gatti et al. (2012); Romero et al. (2012); Miatto et al.
(2012); Perina Jr. (2019), recent advances include bright squeezed states
Chekhova et al. (2015); Barral et al. (2017); Sharapova et al. (2020) and non-
Gaussian states Kolobov (1999); Biswas and Agarwal (2007); Lvovsky et al.
(2020); Walschaers et al. (2020). The use of PDC for metrology applications
have been considered in multiphoton cases, with several studies focussing on
the temporal frequency aspects of the output state Brecht and Silberhorn
(2013); Triginer et al. (2020). In contrast, few such studies provide analyses
of the full spatiotemporal degrees of freedom of the output state. To address
the more challenging analyses of such multiphoton states, a continuous
variable formalism is often used Braunstein and Van Loock (2005); Weedbrook et
al. (2012); Adesso et al. (2014). It leads to a discretization and truncation
of the spatiotemporal degrees of freedom to facilitate the use of Bloch-
Messiah reductions Bloch and Messiah (1962); Horoshko et al. (2019).

Here, we use a Wigner functional approach Roux (2020a, 2021), which is based
on the incorporation of the spatiotemporal degrees of freedom with the
particle-number degrees of freedom Roux (2018, 2020b); Mrowczynski and Mueller
(1994) to produce a functional phase space for a Moyal representation
Groenewold (1946); Moyal (1949); Curtright and Zachos (2012) of all quantum
optical states. This approach was recently used to derive an evolution
equation for the PDC process in a second order nonlinear crystal and to
determine the output state obtained from a spontaneous PDC process Roux
(2020c). In the current investigation, we consider the case where a seed field
enters the nonlinear crystal together with the pump, thus producing a
stimulated PDC process, often called parametric amplification or difference
frequency generation, respectively depending on whether the focus is the
signal field or the idler field. We’ll consider the complete process
incorporating both the signal and idler fields. A solution is obtained for an
arbitrary quantum state entering as seed field by assuming that the PDC
process can be represented as a Bogoliubov transformation of the seed field.
We thus find expressions for the Bogoliubov kernels that serve to produce a
solution for the output state valid under general conditions.

It is often helpful to use a thin-crystal approximation, which is applicable
under typical experimental conditions. We consider two cases: Under some
conditions, we can remove all propagation distance dependencies from the
kernels, which we call the thin-crystal limit. In other cases, we retain up to
subleading order terms in the propagation distance in the argument of the
exponent and the prefactors of the kernels. It is called the thin-crystal
approximation. In the thin crystal limit, we can compute the contributions in
the expansions of the kernels to all orders in the squeezing parameter. It
shows how successive orders contribute to the output intensity distribution,
but some information about the spatial distribution is lost in this limit. The
thin-crystal approximation provides a more accurate description of the shape
of the output distribution, but, although one can calculate higher order
contributions, they become progressively more challenging.

Using the solution for the stimulated PDC process that we obtain as the
Bogoliubov transformation of the seed field, we consider the example of a
coherent state with a Gaussian parameter function acting as the seed field.
This example is studied in the thin crystal limit to see how the different
orders contribute to the output intensity distribution. We then use a thin-
crystal approximation to determine the conditions for the most efficient
difference frequency generation. We also compute the spontaneous PDC
background that is produced concurrently with the stimulated PDC output for
comparison. With the aid of all these contributions to the output intensity
distribution represented in terms of analytic expressions, we then propose a
procedure to measure the absolute intensity of the seed field.

##  II Solving the evolution equation

Assuming that the pump is a coherent state that remains unentangled with the
down-converted light, we can apply the semi-classical approximation. Since the
contributions beyond the semi-classical term require un-enhanced vertices,
they are severely suppressed relative to the semi-classical term Roux (2020c).

Under the semi-classical approximation, the evolution equation for PDC Roux
(2020c) can be written as

$\displaystyle\partial_{z}W_{\hat{\rho}}=$
$\displaystyle\frac{1}{2}\alpha^{*}\diamond H^{*}\diamond\frac{\delta
W_{\hat{\rho}}}{\delta\alpha}+\frac{1}{2}\frac{\delta
W_{\hat{\rho}}}{\delta\alpha^{*}}\diamond H\diamond\alpha,$ (1)

where $W_{\hat{\rho}}$ is the down-converted state’s Wigner functional, $H$ is
the bilinear semi-classical kernel function for the PDC process, $\alpha$ is
the functional phase space field variable, and the $\diamond$-contraction is
defined by

$f\diamond H\diamond g\equiv\int
f(\mathbf{k})H(\mathbf{k},\mathbf{k}^{\prime})g(\mathbf{k}^{\prime})\
d\bar{}\hskip 1.00006ptk\ d\bar{}\hskip 1.00006ptk^{\prime},$ (2)

with

$d\bar{}\hskip
1.00006ptk\equiv\frac{\text{d}^{2}k\text{d}\omega}{(2\pi)^{3}}.$ (3)

The semi-classical approximation implies that the phase space field variable
for the pump field is replaced by its parameter function. As a consequence,
the vertex function for the second-order parametric process is always
contracted with this parameter function, so that it becomes the bilinear
kernel $H$ in the down-converted degrees of freedom only. It is assumed that
the pump parameter function is given by a Gaussian function in the Fourier
domain:

$\zeta(\mathbf{k})=\sqrt{2\pi}\zeta_{0}w_{\text{p}}h(\omega-\omega_{\text{p}},\delta_{\text{p}})\exp\left(-\tfrac{1}{4}w_{\text{p}}^{2}|\mathbf{K}|^{2}\right),$
(4)

where $\zeta_{0}=|\zeta_{0}|\exp(i\varphi)$ is a complex amplitude,
$w_{\text{p}}$ is the beam waist radius, $\mathbf{K}$ is the two-dimensional
transverse part of the wave vector $\mathbf{k}$, and $h(\omega)$ is a
normalized real-valued spectral function, with a bandwidth
$\delta_{\text{p}}$, and a center frequency $\omega_{\text{p}}$. The magnitude
of the pump profile function is $\|\zeta(\mathbf{k})\|^{2}=|\zeta_{0}|^{2}$.
Under monochromatic conditions,
$h^{2}(\omega-\omega_{\text{p}},\delta_{\text{p}})\rightarrow
2\pi\delta(\omega-\omega_{\text{p}})$, for $\delta_{\text{p}}\rightarrow 0$.

The expression for the bilinear semi-classical kernel function, obtained by
contracting the pump parameter function with the vertex for the PDC, reads
Roux (2020c)

$\displaystyle H(\mathbf{k}_{1},\mathbf{k}_{2},z)=$
$\displaystyle\frac{-{i}4}{\hbar}\int\zeta^{*}(\mathbf{k})T(\mathbf{k}_{1},\mathbf{k}_{2},\mathbf{k},z)\
d\bar{}\hskip 1.00006ptk$ $\displaystyle=$
$\displaystyle-{i}\Omega_{0}\sqrt{\omega_{1}\omega_{2}}h(\omega_{1}+\omega_{2}-\omega_{\text{p}},\delta_{\text{p}})$
$\displaystyle\times\exp\left(-\tfrac{1}{4}w_{\text{p}}^{2}|\mathbf{K}_{1}+\mathbf{K}_{2}|^{2}+{i}\Delta
k_{z}z\right),$ (5)

where $T(\mathbf{k}_{1},\mathbf{k}_{2},\mathbf{k},z)$ is the vertex function
for the second-order parametric process, and

$\Omega_{0}=\frac{4\sqrt{2\pi\omega_{\text{p}}}\zeta_{0}^{*}\sigma_{\text{ooe}}w_{\text{p}}}{c^{2}},$
(6)

with $\sigma_{\text{ooe}}$ being the nonlinear coefficient of the nonlinear
medium for type I phase matching, represented as a scattering cross-section
(with the units of an area), and $c$ being the speed of light. The wave vector
mismatch $\Delta k_{z}$ for noncollinear phase matching, is (see Appen. A)

$\displaystyle\Delta k_{z}=$
$\displaystyle\frac{1}{2}\frac{k_{z}(\omega_{1})k_{z}(\omega_{2})}{k_{z}(\omega_{1})+k_{z}(\omega_{2})}\left|\frac{\mathbf{K}_{1}}{k_{z}(\omega_{1})}-\frac{\mathbf{K}_{2}}{k_{z}(\omega_{2})}\right|^{2}$
$\displaystyle-\frac{k^{2}(\omega_{1})}{2k_{z}(\omega_{1})}-\frac{k^{2}(\omega_{2})}{2k_{z}(\omega_{2})}+\frac{1}{2}k_{z}(\omega_{1})+\frac{1}{2}k_{z}(\omega_{2}),$
(7)

where

$k_{z}(\omega)=\frac{\omega\cos[\theta(\omega)]}{v(\omega)},~{}~{}~{}k(\omega)=\frac{\omega}{v(\omega)},$
(8)

with $\theta(\omega)$ being the PDC angle as a function of the frequency and
$v(\omega)$ being the weakly dispersive phase velocity in the crystal. The PDC
angle is determined by the phase matching condition and is not to be confused
with the incident angle of the seed beam. Note that $k_{z}$ is the
$z$-component of the beam axis and not that of the wave vector of any plane
wave. Therefore, it only depends on the frequency and does not depend on the
transverse wave vectors.

###  II.1 Bogoliubov kernels

If $W_{\hat{\rho}}[\alpha^{*},\alpha](0)$ represents the Wigner functional of
the seed field just before it enters the crystal at $z=0$, then the state at
$z>0$ is assumed to be represented by a Bogoliubov transformation of the
initial seed field. The Bogoliubov transformation changes the arguments of the
Wigner functional by

$\displaystyle\begin{split}\alpha&\rightarrow\bar{\alpha}=U\diamond\alpha+V\diamond\alpha^{*}\\\
\alpha^{*}&\rightarrow\bar{\alpha}^{*}=\alpha^{*}\diamond
U^{{\dagger}}+\alpha\diamond V^{{\dagger}},\end{split}$ (9)

where $U$ and $V$ are Bogoliubov kernels such that

$U^{{\dagger}}\diamond U-V^{{\dagger}}\diamond V=\mathbf{1}.$ (10)

Applied to the initial seed field, the Bogoliubov transformation produces

$W_{\hat{\rho}}[\alpha](z)=W_{\hat{\rho}}[U(z)\diamond\alpha+V(z)\diamond\alpha^{*}](0).$
(11)

The Bogoliubov transformed seed field in Eq. (11) substituted into the
evolution equation in Eq. (1) gives

$\displaystyle\left[\alpha^{*}\diamond\partial_{z}U^{{\dagger}}(z)+\alpha\diamond\partial_{z}V^{{\dagger}}(z)\right]\diamond\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}^{*}}$
$\displaystyle+\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}}\diamond\left[\partial_{z}U(z)\diamond\alpha+\partial_{z}V(z)\diamond\alpha^{*}\right]$
$\displaystyle=$ $\displaystyle\frac{1}{2}\left[\alpha^{*}\diamond
H^{*}(z)\diamond V^{{\dagger}}(z)\diamond\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}^{*}}\right.$
$\displaystyle\left.+\alpha^{*}\diamond H^{*}(z)\diamond
U^{T}(z)\diamond\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}}\right.$
$\displaystyle\left.+\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}^{*}}\diamond
U^{*}(z)\diamond H(z)\diamond\alpha\right.$ $\displaystyle\left.+\frac{\delta
W_{\hat{\rho}}[\bar{\alpha}^{*},\bar{\alpha}]}{\delta\bar{\alpha}}\diamond
V(z)\diamond H(z)\diamond\alpha\right].$ (12)

We can separate the result into four equations, which can then be reduces to
two equations:

$\displaystyle\begin{split}\partial_{z}U(z)=&\tfrac{1}{2}V(z)\diamond H(z),\\\
\partial_{z}V(z)=&\tfrac{1}{2}U(z)\diamond H^{*}(z).\end{split}$ (13)

###  II.2 Consistency with spontaneous process

If the Bogoliubov transformed seed state is a solution for the stimulated PDC
process, then a Bogoliubov transformed vacuum state should be the solution for
the spontaneous PDC process. As a result, the Bogoliubov kernels $U$ and $V$
must combine to produce the squeezed vacuum state kernels $A$ and $B$. It
implies that we should be able to use the differential equations for $U$ and
$V$, given in Eq. (13) to derive differential equations for $A$ and $B$. These
latter equations can then be compared with those in Roux (2020c) to see if
they are consistent. It follows that the squeezed vacuum state kernels are
represented by the following contractions of the Bogoliubov kernels

$\displaystyle\begin{split}A=&U^{{\dagger}}\diamond U+V^{T}\diamond V^{*},\\\
B=&U^{{\dagger}}\diamond V+V^{T}\diamond U^{*},\\\ B^{*}=&U^{T}\diamond
V^{*}+V^{{\dagger}}\diamond U.\end{split}$ (14)

We apply a derivative with respect to $z$, substitute Eq. (13) into the
derivatives of the Bogoliubov kernels, and replace the contractions of
Bogoliubov kernels in terms of $A$ and $B$, using Eq. (14). We also use the
fact that $H$ is symmetric. The resulting equations

$\displaystyle\begin{split}\partial_{z}A(z)=&\tfrac{1}{2}H^{{\dagger}}(z)\diamond
B^{*}(z)+\tfrac{1}{2}B(z)\diamond H(z),\\\
\partial_{z}B(z)=&\tfrac{1}{2}H^{{\dagger}}(z)\diamond
A^{T}(z)+\tfrac{1}{2}A(z)\diamond H^{*}(z),\\\
\partial_{z}B^{*}(z)=&\tfrac{1}{2}H^{T}(z)\diamond
A(z)+\tfrac{1}{2}A^{T}(z)\diamond H(z).\end{split}$ (15)

are the same as those obtained in Roux (2020c), apart from interchanging the
definitions of $B$ and $B^{*}$. It shows that the Bogoliubov solution is
consistent with the previous solution in Roux (2020c).

###  II.3 Solutions for Bogoliubov kernels

To solve the equations in Eq. (13), we can proceed in the same way that was
followed in Roux (2020c): We integrate the two equations with respect to $z$
and then perform progressive back substitutions to obtain an expansion in
terms of integrals of contracted $H$-kernels. The initial conditions for the
expansion are assumed to be $U(0)=\mathbf{1}$ and $V(0)=0$, which give the
Wigner functional of the initial seed field. The resulting expressions for the
Bogoliubov kernels are

$\displaystyle\begin{split}U(z)=&\mathbf{1}+\frac{1}{4}\int_{0}^{z}\int_{0}^{z_{1}}H^{*}(z_{2})\diamond
H(z_{1})\ \text{d}z_{2}\ \text{d}z_{1}\\\
&+\frac{1}{16}\int_{0}^{z}\int_{0}^{z_{1}}\int_{0}^{z_{2}}\int_{0}^{z_{3}}H^{*}(z_{4})\diamond
H(z_{3})\\\ &\diamond H^{*}(z_{2})\diamond H(z_{1})\ \text{d}z_{4}\
\text{d}z_{3}\ \text{d}z_{2}\ \text{d}z_{1}...,\\\
V(z)=&\frac{1}{2}\int_{0}^{z}H^{*}(z_{1})\
\text{d}z_{1}+\frac{1}{8}\int_{0}^{z}\int_{0}^{z_{1}}\int_{0}^{z_{2}}H^{*}(z_{3})\\\
&\diamond H(z_{2})\diamond H^{*}(z_{1})\ \text{d}z_{3}\ \text{d}z_{2}\
\text{d}z_{1}....\end{split}$ (16)

Hence, the representation of the stimulated PDC process as a Bogoliubov
transformation, produces a successful solution for the evolution equation.

###  II.4 Thin-crystal

Often in experimental implementations of the PDC process, the pump beam has a
Rayleigh range that is much longer than the length of the nonlinear crystal.
Under such conditions, we can use a thin-crystal approximation to simplify the
expressions.

In the extreme thin-crystal limit, the $z$-dependence is completely removed
from the bilinear kernel $H$. The integrals over the $z$’s in the expressions
of the Bogoliubov kernels in Eq. (16) can all be evaluated, leading to

$\displaystyle\begin{split}U(z)=&\mathbf{1}+\frac{z^{2}}{2^{2}2!}H^{*}\diamond
H+\frac{z^{4}}{2^{4}4!}H^{*}\diamond H\diamond H^{*}\diamond H...\\\
=&\mathbf{1}+\sum_{m=1}^{\infty}\frac{z^{2m}}{2^{2m}(2m)!}|H|^{2m\diamond}\\\
=&\cosh_{\diamond}\left(\tfrac{1}{2}z|H|\right),\\\
V(z)=&\tfrac{1}{2}H^{*}+\frac{z^{3}}{2^{3}3!}H^{*}\diamond H\diamond
H^{*}...\\\
=&\exp(-{i}\varphi)\sum_{m=1}^{\infty}\frac{z^{2m-1}}{2^{2m-1}(2m-1)!}|H|^{(2m-1)\diamond}\\\
=&\exp(-{i}\varphi)\sinh_{\diamond}\left(\tfrac{1}{2}z|H|\right),\end{split}$
(17)

where $H=|H|\exp({i}\varphi)$, with $\varphi$ being the phase of the complex
amplitude of the pump parameter function, defined beneath Eq. (4). Moreover,
the contractions of the sequences of $H$’s can also be evaluated. The results
for odd and even numbers of contracted $H$’s in the thin-crystal limit are

$\displaystyle H_{m}^{(\text{o})}(\mathbf{k}_{1},\mathbf{k}_{2})=$
$\displaystyle\frac{-{i}M_{0}M_{1}^{m}}{m^{5/4}m!}\exp\left(-\tfrac{1}{4m}w_{\text{p}}^{2}|\mathbf{K}_{1}+\mathbf{K}_{2}|^{2}\right)\omega_{1}^{\frac{m}{2}}$
$\displaystyle\times\omega_{2}^{\frac{m}{2}}h(\omega_{\text{p}}-\omega_{1}-\omega_{2},\sqrt{m}\delta_{\text{p}}),$
(18)

where $m$ is an odd integer, and

$\displaystyle H_{m}^{(\text{e})}(\mathbf{k}_{1},\mathbf{k}_{2})=$
$\displaystyle\frac{M_{0}M_{1}^{m}}{m^{5/4}m!}\exp\left(-\tfrac{1}{4m}w_{\text{p}}^{2}|\mathbf{K}_{1}-\mathbf{K}_{2}|^{2}\right)\omega_{1}^{\frac{m}{2}}$
$\displaystyle\times(\omega_{\text{p}}-\omega_{1})^{\frac{m}{2}}h(\omega_{1}-\omega_{2},\sqrt{m}\delta_{\text{p}}),$
(19)

where $m$ is an even integer. The two quantities in the prefactors are given
by

$\displaystyle\begin{split}M_{0}=&\frac{\pi^{5/4}w_{\text{p}}^{2}}{\sqrt{\delta_{\text{p}}}},\\\
M_{1}=&\frac{4\sqrt{2}L|\zeta_{0}|\sigma_{\text{ooe}}\sqrt{\omega_{\text{p}}\delta_{\text{p}}}}{\pi^{3/4}c^{2}w_{\text{p}}},\end{split}$
(20)

where $L$ is the length of the nonlinear crystal.

Using the expressions for the odd and even orders in Eqs. (18) and (19) with
the additional factors of $\tfrac{1}{2}$, we can express the Bogoliubov
kernels by

$\displaystyle\begin{split}U(\mathbf{k}_{1},\mathbf{k}_{2})=&\mathbf{1}+\sum_{n=1}^{\infty}\frac{1}{4^{n}}H_{2n}^{(\text{e})}(\mathbf{k}_{1},\mathbf{k}_{2}),\\\
V(\mathbf{k}_{1},\mathbf{k}_{2})=&\exp(-{i}\varphi)\sum_{n=1}^{\infty}\frac{2}{4^{n}}H_{2n-1}^{(\text{o})}(\mathbf{k}_{1},\mathbf{k}_{2}).\end{split}$
(21)

##  III Coherent state seed

As an illustrative example, we consider the case where the seed field is a
coherent state, which is sent into the nonlinear crystal, together with the
coherent state pump. We apply a Bogoliubov transformation to the arguments of
the coherent state of the seed field to obtain the Wigner functional for the
stimulated down-converted state exiting the crystal. The result is

$\displaystyle\mathcal{N}_{0}\exp\left(-2\|\alpha-\xi\|^{2}\right)$
$\displaystyle\rightarrow$
$\displaystyle\mathcal{N}_{0}\exp\left(-2\|U\diamond\beta+V\diamond\beta^{*}-U\diamond\zeta-V\diamond\zeta^{*}\|^{2}\right)$
$\displaystyle=$
$\displaystyle\mathcal{N}_{0}\exp\left[-2(\beta^{*}-\zeta^{*})\diamond
A\diamond(\beta-\zeta)\right.$ $\displaystyle-(\beta^{*}-\zeta^{*})\diamond
B\diamond(\beta^{*}-\zeta^{*})$ $\displaystyle\left.-(\beta-\zeta)\diamond
B^{*}\diamond(\beta-\zeta)\right],$ (22)

where $\xi$ is the parameter function of the seed field prior to entering the
nonlinear crystal, and $\zeta$ is a parameter function whose Bogoliubov
transformation produces the original parameter function $\xi$. The
combinations of $U$ and $V$ are expressed in terms of the kernels $A$ and $B$
according to Eq. (14), assuming that $U$ is Hermitian and $V$ is symmetric,
and also that $U\diamond V$ is symmetric. These assumed properties of the
Bogoliubov kernels allow the definition of an inverse Bogoliubov
transformation, which can be used to define $\zeta$ in terms of the original
parameter function $\xi$:

$\displaystyle\begin{split}\zeta&\rightarrow U\diamond\xi-V\diamond\xi^{*}\\\
\zeta^{*}&\rightarrow\xi^{*}\diamond U^{{\dagger}}-\xi\diamond
V^{{\dagger}}.\end{split}$ (23)

Hence, the effect of a stimulated PDC process on a coherent state is to apply
a Bogoliubov transformation on the phase space field variables and an inverse
Bogoliubov transformation on the parameter function.

The expression of the Bogoliubov transformed coherent state in Eq. (22) is
that of a displaced squeezed vacuum state. It can now be used to investigate
the output intensity distribution that is obtained from stimulated PDC, seeded
by a coherent state. For this purpose, we need to consider the measurement
process.

###  III.1 Output intensity distribution

To obtain the output intensity distribution that is produced by stimulated PDC
with a coherent state as seed field, we perform an intensity measurement on
the stimulated PDC state with the aid of a localized number operator that
represents the detection process. The Wigner functional for the localized
number operator is

$W_{\hat{n}}[\beta]=\beta^{*}\diamond
D\diamond\beta-\tfrac{1}{2}\text{tr}\\{D\\},$ (24)

where $D$ is the detector kernel. Placed in the exponent of Eq. (22), and
multiplied by an auxiliary variable $J$, it gives a generating function for
the measurement:

$\displaystyle{\cal W}(J)=$
$\displaystyle\mathcal{N}_{0}\int\exp\left[-2\alpha^{*}\diamond
A\diamond\alpha-\alpha^{*}\diamond B\diamond\alpha^{*}\right.$
$\displaystyle-\alpha\diamond
B^{*}\diamond\alpha+J(\alpha^{*}+\zeta^{*})\diamond D\diamond(\alpha+\zeta)$
$\displaystyle\left.-\tfrac{1}{2}J\text{tr}\\{D\\}\right]\ {\cal
D}^{\circ}[\alpha],$ (25)

where we shifted $\beta\rightarrow\alpha+\zeta$, and ${\cal
D}^{\circ}[\alpha]$ is the functional integration measure over the field
variable $\alpha$. The generic result of a functional integration of this form
is given by

$\displaystyle\mathcal{N}_{0}\int\exp\left(-2\alpha^{*}\diamond
A\diamond\alpha-\alpha^{*}\diamond B\diamond\alpha^{*}-\alpha\diamond
B^{*}\diamond\alpha\right.$ $\displaystyle\left.-\alpha^{*}\diamond
F-F^{*}\diamond\alpha\right)\ {\cal D}^{\circ}[\alpha]$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{\det\\{A\\}\det\\{A-B\diamond(A^{*})^{-1}\diamond
B^{*}\\}}}$ $\displaystyle\times\exp\left\\{\tfrac{1}{4}F^{*}\diamond
A^{-1}\diamond F+\tfrac{1}{4}\left[F^{*}-F\diamond(A^{*})^{-1}\diamond
B^{*}\right]\right.$
$\displaystyle\diamond\left[A-B\diamond(A^{*})^{-1}\diamond B^{*}\right]^{-1}$
$\displaystyle\left.\diamond\left[F-B\diamond(A^{*})^{-1}\diamond
F^{*}\right]\right\\},$ (26)

where $\mathcal{N}_{0}$ is a normalization constant. For the expression in Eq.
(25), we replace $A\rightarrow A-\tfrac{1}{2}JD$,
$F^{*}\rightarrow-J\zeta^{*}\diamond D$ and $F\rightarrow-JD\diamond\zeta$ in
Eq. (26), and multiply it by $\exp(J\zeta^{*}\diamond
D\diamond\zeta-\tfrac{1}{2}J\text{tr}\\{D\\})$. Applying a derivative with
respect to $J$ and setting $J=0$, we obtain

$\displaystyle\langle\hat{n}\rangle=$ $\displaystyle\left.\partial_{J}{\cal
W}(J)\right|_{J=0}$ $\displaystyle=$ $\displaystyle\zeta^{*}\diamond
D\diamond\zeta+\tfrac{1}{2}\text{tr}\\{(A-\mathbf{1})\diamond D\\}$
$\displaystyle+\tfrac{1}{4}\text{tr}\\{A^{-1}\diamond B\diamond A^{*}\diamond
B^{*}\diamond A^{-1}\diamond D\\}$
$\displaystyle-\tfrac{1}{4}\text{tr}\\{B\diamond(A^{*})^{-1}\diamond
B^{*}\diamond D\\},$ (27)

where we used that fact that the argument of a trace can always be transposed
and, for a pure squeezed state,

$A-B\diamond(A^{*})^{-1}\diamond B^{*}=A^{-1}.$ (28)

If we can assume that $A$ is real-valued and commutes with $B$, so that

$A^{-1}\diamond B\diamond A^{*}\diamond B^{*}\diamond A^{-1}=A-A^{-1},$ (29)

then the last two terms in Eq. (27) would cancel. If that is not the case,
then at least the leading contributions cancel, leaving terms that are of
fourth order in the squeezing parameter Note (1). Therefore, we consider the
expression for the average number of photons in the output to be

$\langle n\rangle=\zeta^{*}\diamond
D\diamond\zeta+\tfrac{1}{2}\text{tr}\\{(A-\mathbf{1})\diamond D\\}.$ (30)

The first term represents the intensity distribution of the stimulated PDC
field, where $\zeta$ is given in terms of the parameter function of the seed
field by Eq. (23). The second term is the spontaneous PDC field, which acts as
a background noise term.

###  III.2 Detector kernel

To compute the terms in Eq. (30), we need to specify the detector kernel,
which depends on the details of the experimental setup. For the purpose of a
metrology application, we assume that the quantity of interest is obtained
from the far-field intensity distribution of the PDC field. It can be obtained
with the aid of a $2f$ system, as shown in Fig. 1.

The image of the output intensity distribution is measured with the aid of a
CCD array. The detector kernel represents one small detector in the CCD array
located at a position $\mathbf{X}_{0}$ that can vary over the output plane.
We’ll assume that the detector is smaller than the resolution in the output
plane. The $2f$ system performs a Fourier transform of the field distribution
in the plane of the nonlinear crystal. The latter field distribution is the
inverse Fourier transform of the down-converted state in the Fourier domain.
The combination of these two Fourier transformations produces the output
intensity distribution of the down-converted state in the Fourier domain via a
replacement

$\mathbf{K}\rightarrow\frac{k_{\text{d}}}{f}\mathbf{X}_{0},$ (31)

where $k_{\text{d}}$ is the wavenumber of the down-converted light in the
degenerate case, and $f$ is the focal length of the lens in the $2f$ system.
The detector includes a narrow spectral filter that fixes the angular
frequency to be the degenerate PDC frequency

$\omega\rightarrow\omega_{\text{d}}=\tfrac{1}{2}\omega_{\text{p}}.$ (32)

Figure 1: Diagram of the optical setup for the far-field intensity
measurement.

The replacements in Eqs. (31) and (32) can be implemented via Dirac delta
functions. However, the resulting kernel may lose idempotency. To retain
idempotency, we use a limit process to replace each Dirac delta function with
a Gaussian function:

$\delta(k)=\lim_{w_{\text{D}}\rightarrow\infty}\frac{w_{\text{D}}}{\sqrt{\pi}}\exp(-k^{2}w_{\text{D}}^{2}).$
(33)

The result can be seen as a single-mode detector kernel
$D(\mathbf{k}_{1},\mathbf{k}_{2})=M(\mathbf{k}_{1})M^{*}(\mathbf{k}_{2})$,
where $M(\mathbf{k})$ is a normalized function that is represented by

$\displaystyle M(\mathbf{k})=$ $\displaystyle
2\sqrt{2\pi}w_{\text{D}}h(\omega-\omega_{\text{d}},\delta_{\text{D}})$
$\displaystyle\times\exp\left(-w_{\text{D}}^{2}\left|\mathbf{K}-\frac{k_{\text{d}}}{f}\mathbf{X}_{0}\right|^{2}\right).$
(34)

where $w_{\text{D}}$ is the size of the aperture in the crystal plane that
leads to the resolution size in the output plane, and
$h(\omega-\omega_{\text{p}},\delta_{\text{D}})$ is the spectral filter
function, with a bandwidth $\delta_{\text{D}}$. For a single-mode detector
kernel, it follows that $\text{tr}\\{D\\}=1$. The detector area $\mathcal{A}$
is assumed to be equivalent to the area of the output resolution

$\mathcal{A}=\frac{2\pi f^{2}}{k_{\text{d}}^{2}w_{\text{D}}^{2}}.$ (35)

For a large enough $w_{\text{D}}$ and a small enough $\delta_{\text{D}}$, the
detector kernel will effectively perform the replacements in Eqs. (31) and
(32), so that

$\zeta^{*}\diamond D\diamond\zeta=|M^{*}\diamond\zeta|^{2}\rightarrow
K_{\text{D}}\left|\zeta\left(\frac{k_{\text{d}}}{f}\mathbf{X}_{0},\omega_{\text{d}}\right)\right|^{2},$
(36)

where

$K_{\text{D}}=\frac{(2\pi)^{3}\sqrt{\pi}k_{\text{d}}^{2}\mathcal{A}\delta_{\text{D}}}{f^{2}}.$
(37)

The result in Eq. (36) represents the modulus square of the sum over multiple
orders.

###  III.3 Transformed parameter function

Here, we compute the transformed parameter function given in Eq. (23). The
Bogoliubov transformed coherent state is specified in terms of the transformed
parameter function and the Bogoliubov kernels.

The inverse Bogoliubov transformed parameter function consists of two terms,

$\zeta_{1}=U\diamond\xi,~{}~{}~{}\zeta_{2}=V\diamond\xi^{*},$ (38)

which we’ll compute separately. The complete transformed parameter function is
given by $\zeta=\zeta_{1}-\zeta_{2}$. The original parameter function is
assumed to be a Gaussian spectral function, shifted in the Fourier domain to
represent the angle that the seed beam makes with respect to the pump beam as
they enter the crystal. It is defined as

$\displaystyle\xi(\mathbf{k})=$
$\displaystyle\sqrt{2\pi}\xi_{0}w_{\xi}h(\omega-\omega_{\xi},\delta_{\xi})$
$\displaystyle\times\exp\left(-\tfrac{1}{4}w_{\xi}^{2}|\mathbf{K}-\mathbf{K}_{\xi}|^{2}\right),$
(39)

where $\xi_{0}$ is the complex amplitude, $w_{\xi}$ is the beam width,
$\omega_{\xi}$ is the center frequency, $\delta_{\xi}$ is the bandwidth, and
$\mathbf{K}_{\xi}$ is the shift in the spatial Fourier domain. We’ll assume
that the center frequency is equal to the degenerate frequency
$\omega_{\xi}=\omega_{\text{d}}=\tfrac{1}{2}\omega_{\text{p}}$, The associated
shift in the output plane is then given by

$\mathbf{X}_{\xi}=\frac{f}{k_{\text{d}}}\mathbf{K}_{\xi}.$ (40)

Next, we calculate the contractions given in Eq. (38). To perform explicit
calculations, we’ll assume that the experiment satisfies the conditions for
the thin-crystal approximation. The calculations are done in two different
ways to address different aspects. First, we use the expressions for the
Bogoliubov kernels in the thin-crystal limit, as given in Eq. (21) in terms of
Eqs. (18) and (19), to calculated the contractions. It allows us to compute
the transformed parameter function to all orders, showing how the different
orders contribute to the total intensity distribution and revealing how they
change for increasing orders. However, some detailed information about the
spatial distribution is lost in the thin-crystal limit. For instance, the
background term lacks any spatial degrees of freedom in this limit. For more
accurate calculations, we assume that the leading order in the expansions of
the kernels in terms of the squeezing parameter dominates. Then we apply a
thin crystal approximation that retains some $z$-dependence while allowing a
tractable solution. Thus, we obtain a more detailed description of the spatial
properties of the output intensity distribution. However, in this case the
calculations of higher orders become more challenging. Therefore, we only
consider the leading order terms.

###  III.4 Thin-crystal limit

For $\zeta_{1}$ and $\zeta_{2}$, respectively, we use the definitions of
$U(\mathbf{k}_{1},\mathbf{k}_{2})$ and $V(\mathbf{k}_{1},\mathbf{k}_{2})$ in
Eq. (21) in terms of Eqs. (18) and (19), together with Eq. (39). For the
spectra $h(\cdot)$, we assume normalized Gaussian functions. The integrals are
all readily evaluated. The results are

$\displaystyle\begin{split}\zeta_{1}(\mathbf{k}_{1})=&\int
U(\mathbf{k}_{1},\mathbf{k}_{2})\xi(\mathbf{k}_{2})\ d\bar{}\hskip
1.00006ptk_{2}\\\
=&\xi(\mathbf{k}_{1})+\sum_{m=1}^{\infty}\frac{\sqrt{2\pi}\xi_{0}w_{e}^{2}\Xi^{2m}h(\omega-\omega_{\text{d}},\delta_{e})}{(2m)!w_{\xi}(1+2m\eta)^{1/4}}\\\
&\times\exp\left(-\tfrac{1}{4}w_{e}^{2}|\mathbf{K}_{1}-\mathbf{K}_{\xi}|^{2}\right),\\\
\zeta_{2}(\mathbf{k}_{1})=&\int
V(\mathbf{k}_{1},\mathbf{k}_{2})\xi^{*}(\mathbf{k}_{2})\ d\bar{}\hskip
1.00006ptk_{2}\\\
=&{i}\sum_{m=1}^{\infty}\frac{\sqrt{2\pi}\xi_{0}^{*}w_{o}^{2}\Xi^{2m-1}h(\omega-\omega_{\text{d}},\delta_{o})}{(2m-1)!w_{\xi}(1+2m\eta-\eta)^{1/4}}\\\
&\times\exp(-{i}\varphi)\exp\left(-\tfrac{1}{4}w_{o}^{2}|\mathbf{K}_{1}+\mathbf{K}_{\xi}|^{2}\right),\end{split}$
(41)

where

$\displaystyle\begin{split}w_{e}=&\frac{w_{\xi}w_{\text{p}}}{\sqrt{w_{\text{p}}^{2}+2mw_{\xi}^{2}}},\\\
w_{o}=&\frac{w_{\xi}w_{\text{p}}}{\sqrt{w_{\text{p}}^{2}+(2m-1)w_{\xi}^{2}}},\\\
\delta_{e}=&\delta_{\xi}\sqrt{1+2m\eta},\\\
\delta_{o}=&\delta_{\xi}\sqrt{1+2m\eta-\eta},\\\
\eta=&\frac{\delta_{\text{p}}^{2}}{\delta_{\xi}^{2}},\\\
\Xi=&\frac{L|\zeta_{0}|\sigma_{\text{ooe}}\omega_{\text{p}}^{3/2}\sqrt{\delta_{\text{p}}}}{\sqrt{2}\pi^{3/4}c^{2}w_{\text{p}}}.\end{split}$
(42)

To obtain the output intensity, we perform the replacements in Eqs. (31) and
(32), compute the modulus square of $\zeta_{1}$ and $\zeta_{2}$ separately,
assuming they don’t produce overlapping intensity distributions, and multiply
the result by $K_{\text{D}}$, as shown in Eq. (36).

Figure 2: Amplitude distributions for the different orders from $m=0$ to $m=5$
are plotted alons $x$. The intensity of the sum of all the orders is shown by
the solid line. The experimental parameters are
$\lambda_{\text{d}}=0.8~{}\mu$m, $f=100$ mm, $w_{\text{p}}=1$ mm,
$w_{\xi}=0.6$ mm, $|\mathbf{X}_{\xi}|=0.3$ mm, $\eta=1$, and $\Xi=1.7$.

In Fig. 2, we plot the one-dimensional curves for the amplitude distributions
of the different orders in $\zeta$. We also show the intensity distribution
obtained as the modulus square of the sum of all the amplitude distributions.
The even orders, which contribute to $\zeta_{1}$, appear on the right-hand
side and the odd orders, which contribute to $\zeta_{2}$, are on the left-hand
side. It shows how the higher order distributions become broader due to the
multiple convolutions caused by the multiple contractions of $H$.

Some conclusions can be made from the expressions for the output intensity
distribution of the stimulated PDC field and the curves in Fig. 2. The signal
field $\zeta_{1}$ that represents the parametric amplified seed field contains
the original seed field as a leading contribution. The idler field that
represents the difference frequency generated field does not contain any
photons from the original seed field. Since it is produced with the complex
conjugate of the original parameter function, it represents a natural phase
conjugation process. Its leading contribution is a broadened shifted and phase
conjugated version of the seed field. The signal and idler fields are shifted
in opposite directions due to the opposite signs found in the superposition of
the transverse wave vectors in the expressions of the Bogoliubov kernels that
are responsible for the signal and idler fields, respectively. The higher
order contributions grow progressively broader. The odd numbered higher order
terms all contribute to the idler field and the even numbered higher order
terms all contribute to the signal field. The expressions for all the orders
of the idler field contain global phase factors that are not present in any of
the orders of the signal field. Higher orders also come with higher powers of
the squeezing parameter $\Xi$. Therefore, in the case of high gain, the
contributions of the higher orders will become more prominent, leading to
broader intensity distributions.

###  III.5 Thin-crystal approximation

For the more detailed analysis, we consider only the leading order terms in
the signal and idler beams, respectively. The leading order term for the
signal beam $\zeta_{1}$ is the original seed field. The output intensity
distribution for this field is directly obtained by applying the replacements
to Eq. (39), as in Eq. (36).

The leading order term for the idler beam, generated as the difference
frequency beam, is

$\displaystyle\zeta_{2}(\mathbf{k}_{1})=$
$\displaystyle\frac{1}{2}\int_{0}^{L}\int
H^{*}(\mathbf{k}_{1},\mathbf{k}_{2},z_{1})\xi^{*}(\mathbf{k}_{2})\
d\bar{}\hskip 1.00006ptk_{2}\ \text{d}z_{1},$ (43)

where $H(\mathbf{k}_{1},\mathbf{k}_{2},z_{1})$ is given in Eq. (5). The
integrations over $\mathbf{k}_{2}$ can be evaluated readily, but the
integration over $z_{1}$ is challenging. Therefore, we use the thin-crystal
approximation to expand the prefactor and the argument of the exponent,
respectively, to subleading order in $z_{1}$. The resulting expression is

$\displaystyle\zeta_{2}\approx$
$\displaystyle\frac{{i}\Omega_{1}h\left(\omega-\omega_{\text{d}},\sqrt{\delta_{\text{p}}^{2}+\delta_{\xi}^{2}}\right)}{w_{0}^{2}Lk_{z}(\omega_{\text{d}})}\int_{0}^{L}\exp\left(-\frac{w_{\text{p}}^{2}w_{\xi}^{2}|\mathbf{K}_{a}|^{2}}{4w_{0}^{2}}\right.$
$\displaystyle\left.-\frac{{i}z_{1}|\mathbf{K}_{b}|^{2}}{4k_{\text{d}}w_{0}^{4}}+{i}z_{1}\chi_{\text{d}}\right)\left(k_{\text{d}}w_{0}^{2}-{i}z_{1}\right)\
\text{d}z_{1},$ (44)

where

$\Omega_{1}=\frac{2\Xi\xi_{0}^{*}w_{\xi}w_{\text{p}}^{2}\sqrt{2\pi\delta_{\xi}}}{\left(\delta_{\text{p}}^{2}+\delta_{\xi}^{2}\right)^{1/4}w_{0}^{2}},$
(45)

and

$\displaystyle\begin{split}w_{0}=&\sqrt{w_{\text{p}}^{2}+w_{\xi}^{2}},\\\
\mathbf{K}_{a}=&\mathbf{K}_{1}+\mathbf{K}_{\xi},\\\
\mathbf{K}_{b}=&\left(2w_{\text{p}}^{2}+w_{\xi}^{2}\right)\mathbf{K}_{1}-w_{\xi}^{2}\mathbf{K}_{\xi},\\\
\chi_{\text{d}}=&\frac{k^{2}(\omega_{\text{d}})}{k_{z}(\omega_{\text{d}})}-k_{z}(\omega_{\text{d}})=\frac{k(\omega_{\text{d}})\sin^{2}[\theta(\omega_{\text{d}})]}{\cos[\theta(\omega_{\text{d}})]}.\end{split}$
(46)

Integrated over $z_{1}$, the expression reads

$\displaystyle\zeta_{2}\approx$
$\displaystyle\Omega_{1}h\left(\omega-\omega_{\text{d}},\sqrt{\delta_{\text{p}}^{2}+\delta_{\xi}^{2}}\right)\exp\left(-\frac{w_{\text{p}}^{2}w_{\xi}^{2}|\mathbf{K}_{a}|^{2}}{4w_{0}^{2}}\right)$
$\displaystyle\times\left\\{\frac{(\kappa-\beta)}{\kappa^{2}}[1-\exp(-{i}\kappa)]+\frac{{i}\beta}{\kappa}\exp(-{i}\kappa)\right\\},$
(47)

where

$\displaystyle\begin{split}\kappa=&\frac{L|\mathbf{K}_{b}|^{2}}{4k_{z}(\omega_{\text{d}})w_{0}^{4}}-L\chi_{\text{d}},\\\
\beta=&\frac{L}{k_{z}(\omega_{\text{d}})w_{0}^{2}}.\end{split}$ (48)

We recognize the first exponential function in Eq. (47) as being associated
with the original seed field, but shifted in the opposite direction and
broadened due to the convolution with the kernel. Since the original field is
a Gaussian, the complex conjugation does not have any effect on its shape.
This field is modulated by a function of $\kappa$. The effect of this
modulation is to vary the peak amplitude of the idler field — i.e., the
efficiency of the difference frequency generation. To see how the efficiency
varies, we consider the amplitude at the peak, which is obtained by
substituting

$\mathbf{K}_{1}\rightarrow-\mathbf{K}_{\xi},$ (49)

so that $\mathbf{K}_{a}\rightarrow 0$,

$\mathbf{K}_{b}\rightarrow\mathbf{K}_{b}^{\prime}=-2w_{0}^{2}\mathbf{K}_{\xi},$
(50)

and

$\kappa\rightarrow
a=\frac{L|\mathbf{K}_{\xi}|^{2}}{k_{z}(\omega_{\text{d}})}-L\chi_{\text{d}}.$
(51)

The dimensionless variable $a$ represents the mismatch in the incident angle,
as determined by $\mathbf{K}_{\xi}$, and the PDC angle on which
$\chi_{\text{d}}$ depends and which is determined by the phase matching
condition. The intensity of the peak is governed by the modulus square of the
part of the expression that contains $a$, which is given by the efficiency
function

$\displaystyle f(a)=$
$\displaystyle\frac{\beta^{2}}{a^{2}}+\frac{2(a-\beta)\beta\sin a}{a^{3}}$
$\displaystyle+\frac{2(a-\beta)^{2}(1-\cos a)}{a^{4}}.$ (52)

It shows how the efficiency of the difference frequency generation depends on
the mismatch between the incident angle and the PDC angle.

The curve for $f(a)$, shown in Fig. 3, resembles the shape of a sinc-squared
function. The location of the peak efficiency, where the most efficient
difference frequency generation is obtained, is approximately given by

$a_{\text{peak}}=-\frac{6\beta}{6+\beta^{2}}.$ (53)

We thus obtain a relationship between the phase matching condition in terms of
$\chi(\omega_{\text{d}})$ and the incident angle of the seed beam in terms of
$\mathbf{K}_{\xi}$, given by

$\displaystyle\left|\mathbf{K}_{\xi}\right|^{2}=$ $\displaystyle
k_{z}(\omega_{\text{d}})\chi_{\text{d}}-\frac{k_{z}(\omega_{\text{d}})}{L}\frac{6\beta}{6+\beta^{2}}$
$\displaystyle\approx$ $\displaystyle
k^{2}(\omega_{\text{d}})\sin^{2}[\theta(\omega_{\text{d}})]-w_{0}^{-2}.$ (54)

Figure 3: Plot of the efficiency function in Eq. (52) for $\beta=0.4$ as a
function of the dimensionless variable $a$.

In Fig. 3, a relatively large value for $\beta$ is selected to emphasis the
shift in the location of the peak efficiency. In most experiments, the value
of $\beta$ would be much smaller, moving the peak efficiency closer to the
origin.

###  III.6 Background term

The second term in Eq. (30) is a background term produced by spontaneous PDC.
Usually, the background term would be negligible compared to the term produced
by stimulated PDC, due to the enhancement by the number of photons in the seed
field. However, we can envisage scenarios where the measurement is made for a
small average number of photons in the seed field. Moreover, for absolute
intensity measurements, knowledge of this spontaneous PDC term is required.
Therefore, we also compute the leading order contribution for the background
term.

If an identity operator is used for the detector kernel and the thin-crystal
limit is applied to the kernel function $A$, the trace would produce a
divergent result. With the localized detector kernel discussed in Section
III.2, the result is finite, but all the spatial degrees of freedom are lost
in the thin-crystal limit, giving a constant background. For a more precise
calculation that provides some detail of the spatial behavior of the
background, we use the full expression and then apply a thin-crystal
approximation to alleviate the calculation.

The leading order term in the expansion of $A-\mathbf{1}$ is second order in
the bilinear kernel:

$A-\mathbf{1}\approx\frac{1}{2}\int_{0}^{L}\int_{0}^{L}H^{*}(z_{1})\diamond
H(z_{2})\ \text{d}z_{2}\ \text{d}z_{1}.$ (55)

The contraction of the two bilinear kernels produces

$\displaystyle\int
H^{*}(\mathbf{k}_{1},\mathbf{k}_{2},z_{1})H(\mathbf{k}_{2},\mathbf{k}_{3},z_{2})\
d\bar{}\hskip 1.00006ptk_{2}$ $\displaystyle=$
$\displaystyle\Omega_{2}h(\omega_{1}-\omega_{3},\sqrt{2}\delta_{\text{p}})\exp\left[-R_{1}|\mathbf{K}_{1}-\mathbf{K}_{3}|^{2}\right.$
$\displaystyle\left.-R_{2}|\mathbf{K}_{1}+\mathbf{K}_{3}|^{2}-R_{3}\left(|\mathbf{K}_{1}|^{2}-|\mathbf{K}_{3}|^{2}\right)+R_{4}\right],$
(56)

where

$\Omega_{2}=\frac{16\pi^{5/4}2^{3/4}\Xi^{2}w_{\text{p}}^{2}\omega_{1}(\omega_{\text{p}}-\omega_{1})}{L^{2}\omega_{\text{p}}^{2}\sqrt{\delta_{\text{p}}}\tau(\omega_{1},z_{1},z_{2})},$
(57)

with $\Xi$ given in Eq. (42), and

$\displaystyle\begin{split}R_{1}=&\frac{w_{\text{p}}^{2}}{8}+\frac{{i}(z_{1}-z_{2})k_{z}(\omega_{\text{p}}-\omega_{1})}{8k(\omega_{\text{p}})k_{z}(\omega_{1})}\\\
&+\frac{(z_{1}+z_{2})^{2}}{8w_{\text{p}}^{2}k^{2}(\omega_{\text{p}})\tau(\omega_{1},z_{1},z_{2})},\\\
R_{2}=&\frac{{i}(z_{1}-z_{2})k(\omega_{\text{p}})}{8k_{z}(\omega_{1})k_{z}(\omega_{\text{p}}-\omega_{1})\tau(\omega_{1},z_{1},z_{2})},\\\
R_{3}=&\frac{{i}(z_{1}+z_{2})}{4k_{z}(\omega_{1})\tau(\omega_{1},z_{1},z_{2})},\\\
R_{4}=&{i}\tfrac{1}{2}(z_{1}-z_{2})[\chi(\omega_{\text{p}}-\omega_{1})+\chi(\omega_{1})],\end{split}$
(58)

with

$\displaystyle\begin{split}\tau(\omega_{1},z_{1},z_{2})=&1+\frac{{i}(z_{1}-z_{2})k_{z}(\omega_{1})}{w_{\text{p}}^{2}k(\omega_{\text{p}})k_{z}(\omega_{\text{p}}-\omega_{1})},\\\
\chi(\omega)=&\frac{k^{2}(\omega)}{k_{z}(\omega)}-k_{z}(\omega).\end{split}$
(59)

As discussed in Sec. III.2, the effect of the detector kernel is to perform
the replacements in Eqs. (31) and (32), and multiply the result by
$K_{\text{D}}$. Since both $\mathbf{K}_{1}$ and $\mathbf{K}_{3}$ are replaced
in terms of $\mathbf{X}_{0}$, the terms with $R_{1}$ and $R_{3}$ drop away,
leaving only the $R_{2}$ and $R_{4}$ terms.

The remaining integrations over $z_{1}$ and $z_{2}$ are alleviated by using
the thin-crystal approximation as we did in Sec. III.5. Therefore, we use the
approximation

$\frac{1}{\tau(\omega_{1},z_{1},z_{2})}\approx\tau^{*}(\omega_{1},z_{1},z_{2}).$
(60)

The result is

$\displaystyle\begin{split}\Omega_{2}\rightarrow\Omega_{2}^{\prime}=&\frac{4\pi^{5/4}2^{3/4}\Xi^{2}w_{\text{p}}^{2}}{L^{2}\sqrt{\delta_{\text{p}}}}\tau^{*}(\omega_{1},z_{1},z_{2}),\\\
R_{2}\rightarrow
R_{2}^{\prime}=&\frac{{i}(z_{1}-z_{2})}{4k_{z}(\omega_{\text{d}})},\\\
R_{4}\rightarrow R_{4}^{\prime}=&{i}(z_{1}-z_{2})\chi_{\text{d}},\end{split}$
(61)

where we also discarded a $(z_{1}-z_{2})^{2}$-term in $R_{2}^{\prime}$, and
used the fact that $k(\omega_{\text{p}})=2k_{z}(\omega_{\text{d}})$.

Figure 4: Intensity distribution of the spontaneous PDC background field with
$L=3$ mm, $\lambda_{\text{d}}=0.8~{}\mu$m, $f=100$ mm, and $w_{\text{p}}=0.2$
mm.

After evaluating the $z$-integrations, we obtain

$\displaystyle\tfrac{1}{2}\text{tr}\\{(A-\mathbf{1})\diamond D\\}$
$\displaystyle=$
$\displaystyle\Omega_{3}\left[\frac{\beta_{0}\mathcal{S}}{\left(r^{2}-r_{0}^{2}\right)^{2}}+\frac{2\left(r^{2}-r_{0}^{2}-R^{2}\right)(1-\mathcal{C})}{\left(r^{2}-r_{0}^{2}\right)^{3}}\right],$
(62)

where

$\displaystyle\begin{split}\Omega_{3}=&\frac{2K_{\text{D}}\pi^{3/2}\Xi^{2}w_{\text{p}}^{2}R^{4}}{\beta_{0}^{2}\delta_{\text{p}}},\\\
\mathcal{S}=&\sin\left[\frac{\left(r^{2}-r_{0}^{2}\right)\beta_{0}}{R^{2}}\right],\\\
\mathcal{C}=&\cos\left[\frac{\left(r^{2}-r_{0}^{2}\right)\beta_{0}}{R^{2}}\right],\\\
r=&|\mathbf{X}_{0}|,\\\ r_{0}=&f\sin[\theta(\omega_{\text{d}})],\\\
R=&\frac{f}{k_{\text{d}}w_{\text{p}}},\\\
\beta_{0}=&\frac{L}{w_{\text{p}}^{2}k_{\text{d}}\cos[\theta(\omega_{\text{d}})]}.\end{split}$
(63)

The output intensity distribution in Eq. (62) is a non-negative real-valued
rotationally symmetric function of $r$. The shape of the output intensity
distribution is governed by three parameters composed of the experimental
parameters.

The intensity distribution of the background field is plotted in Fig. 4 along
the $x$-axis for three different PDC angles, representing different phase
matching conditions. In the experiment, the phase matching condition is
changed by rotating the nonlinear crystal to change the angle between the
crystal axis and the pump beam axis. For $\theta(\omega_{\text{d}})=0$
degrees, the intensity distribution is a broadened sinc-function. With
increasing angle, the peak develops a dip in the center, which eventually
produces a ring-shaped intensity distribution.

###  III.7 Combined intensity distribution

The combined intensity distribution is obtained by adding the intensities of
the signal field, idler field and background field. The result is shown in
Fig. 5 for a small number of photons in the seed field and for three different
angles between the seed and pump beams to demonstrate the effect of the
efficiency function. The peaks for the signal field are found on the right-
hand side and those for the idler field are on the left-hand side. All the
peaks are added on top of the spontaneous PDC background.

The small chosen values for the number of photons in the seed field
$|\xi_{0}|^{2}$ and the squeezing parameter $\Xi$ allow us to see the
spontaneous PDC background together with the stimulated PDC fields. For larger
values of these parameters, the stimulated PDC fields would be much larger
than the spontaneous PDC background.

Using the quantities defined in Eq. (63), we can express the location of the
idler peak intensity for the most efficient difference frequency generation in
the output plane, based on Eq. (54), as

$\left|\mathbf{X}_{\text{peak}}\right|^{2}=r_{0}^{2}-\tfrac{1}{2}R^{2}.$ (64)

A parameter $G$ is used to represent the location of the idler peak intensity
in the output plane, according to

$\mathbf{X}_{\xi}=G\mathbf{X}_{\text{peak}}.$ (65)

It can be varied experimentally by changing the incident angle of the seed
beam for fixed phase matching conditions. In Fig. 5, we show the curves for
$G=0.8,1,1.2$. The highest idler peak intensity is obtained for $G=1$, and its
location coincides with the location of the peak intensity of the spontaneous
PDC background, as expected.

Figure 5: Intensity distribution of the combined output field for three
different displacements of the seed field. The experimental parameters: $L=3$
mm, $\lambda_{\text{d}}=0.8~{}\mu$m, $f=100$ mm, $w_{\text{p}}=w_{\xi}=0.2$
mm, $|\xi_{0}|^{2}=4$, and $\Xi=1$.

The effective squeezing parameter $\Xi$, given in Eq. (42), serves as the
expansion parameter for the different orders of the kernels. In strongly
pumped scenarios, the effective squeezing parameter may become quite large.
However, the factorials in the denominators in Eq. (41) ensure that the higher
orders will eventually converge. In Fig. 2, we used a relatively large value
for the effective squeezing parameter so that multiple orders would be
visible. For the case in Fig. 5 on the other hand, we assumed that the
effective squeezing parameter is relatively small so that the leading order
would completely dominate.

##  IV Toward absolute intensity

The fact that the expected output intensity distribution can be reproduced
with analytic expressions for the stimulated and spontaneous PDC fields means
that the measured output distribution can be used to determine the parameters
that govern its shape. By capturing an image of such an output intensity
distribution with a CCD array, we can use the analytic expressions to fit the
intensity distribution with the unknown parameters as fitting parameters.

Consider for example the case where the intensity of the seed field, which is
proportional to the average number of photons $|\xi_{0}|^{2}$, is unknown. The
three fields shown in Fig. 5 — signal, idler and background — respectively
receive enhancement factors $|\xi_{0}|^{2}$, $\Xi|\xi_{0}|^{2}$, and $\Xi$,
with some minor modifications due to the various bandwidths and beam widths.
By fitting the three fields to their respective analytic expressions, the
absolute intensity of the seed field can be determined in terms of the average
number of photons, even if $\Xi$ is unknown. It may require separate
measurements with and without the seed field if the dynamic range of the CCD
array cannot see the background together with the stimulated fields.

##  V Conclusions

The analytic expressions for the output intensity distributions, produced by a
stimulated PDC process are obtained with a Wigner functional approach. It
produces these results without knowledge of the eigenstates of the process,
i.e., without discretizing and truncating the kernel functions of the process.
Moreover, it avoids the need for numerical simulations.

This analysis allows the calculation of higher order contributions. However,
they do become more challenging as the order increases, unless the thin
crystal limit can be applied where all the $z$-dependencies are removed from
the expressions of the bilinear kernel. In the thin crystal limit, the
information about the phase matching condition is lost. If such information is
important, the thin crystal approximation should be used instead, which
retains the subleading order $z$-dependencies in the exponents and the
prefactors of the kernel functions.

Equipped with these analytic expressions, we may be able to use the measured
output intensity distributions from such PDC experiments to determine the
values of unknown parameters. It can therefore serve as a tool for metrology
applications in radiometry and photometry. One such application is the
measurement of the absolute intensity of a light source, which is used as the
seed field. It should be possible to use this method for arbitrary large
intensities, which would require separate measurements of the stimulated
fields and the background field, as well as the calculation of higher order
contributions to the output intensity distribution.

## Acknowledgement

This work was supported in part by funding from the National Research
Foundation of South Africa (Grant Numbers: 118532).

##  Appendix A Noncollinear phase matching

Here, we discuss the derivation of the expression for $\Delta k_{z}$ under
paraxial conditions when the PDC angle $\theta(\omega)$ is too large to be
considered paraxial. We impose the paraxial condition on each of the beams
separately by assuming that the beam divergence angle of each of the three
beams is small, even for arbitrary large $\theta(\omega)$.

The pump beam propagates along the $z$-axis, while the signal and idler beams
propagate at angles of $\theta(\omega)$ with respect to the $z$-axis. The
$z$-components of the wave vectors are given in terms of the angular frequency
and the transverse part of the wave vectors by

$k_{z}=\sqrt{\frac{\omega^{2}}{v^{2}(\omega)}-|\mathbf{K}|^{2}}.$ (66)

To impose the paraxial approximation for the pump beam, it suffices to assume
that $|\mathbf{K}|\ll\omega/v(\omega)$. However, since the other two beams
propagate at different angles that may be large compared to their beam
divergence angles, their situation is more complicated.

For these beams, we replace the transverse part by

$|\mathbf{K}|^{2}\rightarrow\frac{\omega^{2}\sin^{2}[\theta(\omega)]}{v^{2}(\omega)}+\left\\{|\mathbf{K}|^{2}-\frac{\omega^{2}\sin^{2}[\theta(\omega)]}{v^{2}(\omega)}\right\\}\epsilon,$
(67)

where we tagged the part in the expression that is small under paraxial
conditions by an auxiliary parameter $\epsilon$. We can then substitute the
appropriate expressions for the three beams’ $z$-components into

$\Delta k_{z}=k_{z}^{(p)}-k_{z}^{(s)}-k_{z}^{(i)},$ (68)

and expand the result to subleading order in $\epsilon$, after which we set
$\epsilon=1$. For critical phase matching, the leading order term cancels,
leaving the subleading order which is represented by the expression in Eq.
(7).

## References

  * Hong and Mandel (1985) C. K. Hong and L. Mandel, “Theory of parametric frequency down conversion of light,” Phys. Rev. A 31, 2409 (1985).
  * Mair et al. (2001) A. Mair, A. Vaziri, G. Weihs, and A. Zeilinger, “Entanglement of the orbital angular momentum states of photons,” Nature 412, 313 (2001).
  * Torres et al. (2003) J. P. Torres, A. Alexandrescu, and L. Torner, “Quantum spiral bandwidth of entangled two-photon states,” Phys. Rev. A 68, 050301 (2003).
  * Miatto et al. (2011) F. M. Miatto, A. M. Yao, and S. M. Barnett, “Full characterization of the quantum spiral bandwidth of entangled biphotons,” Phys. Rev. A 83, 033816 (2011).
  * Walborn et al. (2010) S. P. Walborn, C. H. Monken, S. Pádua, and P. H. Souto Ribeiro, “Spatial correlations in parametric down-conversion,” Phys. Rep. 495, 87 (2010).
  * Roux (2020a) F. S. Roux, “Quantifying entanglement of parametric down-converted states in all degrees of freedom,” Phys. Rev. Research 2, 023137 (2020a).
  * Lvovsky (2015) A. I. Lvovsky, “Squeezed light,” Photonics: Scientific Foundations, Technology and Applications 1, 121 (2015).
  * Kolobov (1999) M. I. Kolobov, “The spatial behavior of nonclassical light,” Rev. Mod. Phys. 71, 1539 (1999).
  * Wasilewski et al. (2006) W. Wasilewski, A. I. Lvovsky, K. Banaszek, and C. Radzewicz, “Pulsed squeezed light: Simultaneous squeezing of multiple modes,” Phys. Rev. A 73, 063819 (2006).
  * Roux (2021) F. S. Roux, “Spatiotemporal effects on squeezing measurements,” Phys. Rev. A 103, 013701 (2021).
  * Yurke et al. (1986) B. Yurke, S. L. McCall, and J. R. Klauder, “SU(2) and SU(1,1) interferometers,” Phys. Rev. A 33, 4033 (1986).
  * Seyfarth et al. (2020) U. Seyfarth, A. B. Klimov, H. de Guise, G. Leuchs, and L. L. Sánchez-Soto, “Wigner function for SU(1,1),” Quantum 4, 317 (2020).
  * Ferreri et al. (2021) A. Ferreri, M. Santandrea, M. Stefszky, K. H. Luo, H. Herrmann, C. Silberhorn, and P. R. Sharapova, “Spectrally multimode integrated SU(1,1) interferometer,” Quantum 5, 461 (2021).
  * Giovannetti et al. (2011) V. Giovannetti, S. Lloyd, and L. Maccone, “Advances in quantum metrology,” Nature Photon. 5, 222 (2011).
  * Polino et al. (2020) E. Polino, M. Valeri, N. Spagnolo, and F. Sciarrino, “Photonic quantum metrology,” AVS Quantum Sci. 2, 024703 (2020)).
  * Giovannetti et al. (2004) V. Giovannetti, S. Lloyd, and L. Maccone, “Quantum-enhanced measurements: beating the standard quantum limit,” Science 306, 1330 (2004).
  * Arnaut and Barbosa (2000) H. H. Arnaut and G. A. Barbosa, “Orbital and intrinsic angular momentum of single photons and entangled pairs of photons generated by parametric down-conversion,” Phys. Rev. Lett. 85, 286 (2000).
  * Brambilla et al. (2004) E. Brambilla, A. Gatti, M. Bache, and L. A. Lugiato, “Simultaneous near-field and far-field spatial quantum correlations in the high-gain regime of parametric down-conversion,” Phys. Rev. A 69, 023802 (2004).
  * Brambilla et al. (2010) E. Brambilla, L. Caspani, L. A. Lugiato, and A. Gatti, “Spatiotemporal structure of biphoton entanglement in type-II parametric down-conversion,” Phys. Rev. A 82, 013835 (2010).
  * Gatti et al. (2012) A. Gatti, T. Corti, E. Brambilla, and D. B. Horoshko, “Dimensionality of the spatiotemporal entanglement of parametric down-conversion photon pairs,” Phys. Rev. A 86, 053803 (2012).
  * Romero et al. (2012) J. Romero, D. Giovannini, S. Franke-Arnold, S. M. Barnett, and M. J. Padgett, “Increasing the dimension in high-dimensional two-photon orbital angular momentum entanglement,” Phys. Rev. A 86, 012334 (2012).
  * Miatto et al. (2012) F. M. Miatto, H. D. L. Pires, S. M. Barnett, and M. P. van Exter, “Spatial Schmidt modes generated in parametric down-conversion,” Eur. Phys. J. D 66, 263 (2012).
  * Perina Jr. (2019) J. Perina Jr., “Waves in spatio-spectral and temporal coherence of evolving ultra-intense twin beams,” Sci. Rep. 9, 4256 (2019).
  * Chekhova et al. (2015) M. V. Chekhova, G. Leuchs, and M. Żukowski, “Bright squeezed vacuum: Entanglement of macroscopic light beams,” Opt. Commun. 337, 27 (2015).
  * Barral et al. (2017) D. Barral, N. Belabas, L. M. Procopio, V. d’Auria, S. Tanzilli, K. Bencheikh, and J. A. Levenson, “Continuous-variable entanglement of two bright coherent states that never interacted,” Phys. Rev. A 96, 053822 (2017).
  * Sharapova et al. (2020) P. R. Sharapova, G. Frascella, M. Riabinin, A. M. Pérez, O. V. Tikhonova, S. Lemieux, R. W. Boyd, G. Leuchs, and M. V. Chekhova, “Properties of bright squeezed vacuum at increasing brightness,” Phys. Rev. Research 2, 013371 (2020).
  * Biswas and Agarwal (2007) A. Biswas and G. S. Agarwal, “Nonclassicality and decoherence of photon-subtracted squeezed states,” Phys. Rev. A 75, 032104 (2007).
  * Lvovsky et al. (2020) A. I. Lvovsky, P. Grangier, A. Ourjoumtsev, V. Parigi, M. Sasaki, and R. Tualle-Brouri, “Production and applications of non-gaussian quantum states of light,” arXiv preprint arXiv:2006.16985 (2020).
  * Walschaers et al. (2020) M. Walschaers, V. Parigi, and N. Treps, “Practical framework for conditional non-Gaussian quantum state preparation,” PRX Quantum 1, 020305 (2020).
  * Brecht and Silberhorn (2013) B. Brecht and C. Silberhorn, “Characterizing entanglement in pulsed parametric down-conversion using chronocyclic wigner functions,” Phys. Rev. A 87, 053810 (2013)).
  * Triginer et al. (2020) G. Triginer, M. D. Vidrighin, N. Quesada, A. Eckstein, M. Moore, W. S. Kolthammer, J. E. Sipe, and I. A. Walmsley, “Understanding high-gain twin-beam sources using cascaded stimulated emission,” Phys. Rev. X 10, 031063 (2020)).
  * Braunstein and Van Loock (2005) S. L. Braunstein and P. Van Loock, “Quantum information with continuous variables,” Rev. Mod. Phys. 77, 513 (2005).
  * Weedbrook et al. (2012) C. Weedbrook, S. Pirandola, R. García-Patrón, N. J. Cerf, T. C. Ralph, J. H. Shapiro, and S. Lloyd, “Gaussian quantum information,” Rev. Mod. Phys. 84, 621 (2012).
  * Adesso et al. (2014) G. Adesso, S. Ragy, and A. R. Lee, “Continuous variable quantum information: Gaussian states and beyond,” Open Syst. Inf. Dyn. 21, 1440001 (2014).
  * Bloch and Messiah (1962) C. Bloch and A. Messiah, “The canonical form of an antisymmetric tensor and its application to the theory of superconductivity,” Nucl. Phys. 39, 95 (1962).
  * Horoshko et al. (2019) D. B. Horoshko, L. La Volpe, F. Arzani, N. Treps, C. Fabre, and M. I. Kolobov, “Bloch-Messiah reduction for twin beams of light,” Phys. Rev. A 100, 013837 (2019).
  * Roux (2018) F. S. Roux, “Combining spatiotemporal and particle-number degrees of freedom,” Phys. Rev. A 98, 043841 (2018).
  * Roux (2020b) F. S. Roux, “Erratum: Combining spatiotemporal and particle-number degrees of freedom [Phys. Rev. A 98, 043841 (2018)],” Phys. Rev. A 101, 019903(E) (2020b).
  * Mrowczynski and Mueller (1994) S. Mrowczynski and B. Mueller, “Wigner functional approach to quantum field dynamics,” Phys. Rev. D 50, 7542 (1994).
  * Groenewold (1946) H. J. Groenewold, “On the principles of elementary quantum mechanics,” Physica 12, 405 (1946).
  * Moyal (1949) J. E. Moyal, “Quantum mechanics as a statistical theory,” Math. Proc. Camb. Philos. Soc. 45, 99 (1949).
  * Curtright and Zachos (2012) T. L. Curtright and C. K. Zachos, “Quantum mechanics in phase space,” Asia Pacific Physics Newsletter 1, 37 (2012).
  * Roux (2020c) F. S. Roux, “Parametric down-conversion beyond the semiclassical approximation,” Phys. Rev. Research 2, 033398 (2020c).
  * Note (1) Both $A$ and $A^{-1}$ have the identity $\mathbf{1}$ as their leading terms. When we replace $A$ and $A^{-1}$ by $\mathbf{1}$ in the last two terms in Eq. (27), they cancel.