# Unified learning-based lossy and lossless JPEG recompression††thanks:
*Corresponding author.

###### Abstract

JPEG is still the most widely used image compression algorithm. Most image
compression algorithms only consider uncompressed original image, while
ignoring a large number of already existing JPEG images. Recently, JPEG
recompression approaches have been proposed to further reduce the size of JPEG
files. However, those methods only consider JPEG lossless recompression, which
is just a special case of the rate-distortion theorem. In this paper, we
propose a unified lossly and lossless JPEG recompression framework, which
consists of learned quantization table and Markovian hierarchical variational
autoencoders. Experiments show that our method can achieve arbitrarily low
distortion when the bitrate is close to the upper bound, namely the bitrate of
the lossless compression model. To the best of our knowledge, this is the
first learned method that bridges the gap between lossy and lossless
recompression of JPEG images.

Index Terms—  JPEG recompression, quantization table

##  1 Introduction

With the rapid development of the internet and smart devices, multimedia data
is growing exponentially. Due to the expansive cost of transmission and
storage, data compression is indispensable for data centers, cloud storage,
and network filesystems. JPEG [1] is a commonly used method of lossy
compression for digital images, which was first published in 1992. After
decades of technological development, a variety of new image compression
algorithms are emerging, like JPEG2000 [2], PNG [3], WEBP [4], BPG [5], the
intra-coding of VVC/H.266 [6] and learning based methods [7, 8, 9, 10, 11,
12]. However, these subsequent methods only focus on uncompressed original
image data, ignoring the large number of JPEG images that already exist on the
internet. According to the survey [13] of Dropdox, JPEG files take up about
35% of the total storage space.

Fig. 1: Objective quality comparison on Kodak dataset when recompressing JPEG
images of QP 75 with different methods: JPEG-XL [14], MLCC [15] and ELIC [12].

Due to the demand to further reduce the size of JPEG files, several methods
for JPEG recompression are proposed, e.g. Lepton [13], JPEG-XL [14] and MLCC
[15]. However, these methods only consider the lossless recompression for JPEG
files, which is a special case of the rate-distortion theorem. In other words,
those methods cannot provide variable bitrates for different distortion
requirements, which means lossy recompression of JPEG files.

Intuitively, those new lossy compression methods [4, 5, 14, 7, 8, 9, 10, 11]
can also be used for lossy JPEG recompression. Specifically, the JPEG file is
firstly decompressed to RGB data, and we can use those new lossy compression
methods to compress the RGB data. However, this straightforward method will
induce the Reverse-Compression problem when the desired distortion is below a
certain value. As shown in Fig. 1, whether we use hand-designed compression
methods or learnable compression methods to recompress JPEG files in this way,
there is still non-negligible compression distortion when the bitrate exceeds
the bitrate of lossless recompression, which should be the upper bound for any
lossy recompression. According to the rate-distortion theory [16], it means
there is a gap between current lossless recompression methods and lossy
recompression methods for JPEG images.

In this paper, we propose a unified JPEG lossy and lossless recompression
framework, which consists of one learned quantization table, one learned
inverse quantization table and a DCT-domain lossless compression model
implemented by Markovian hierarchical variational autoencoder. The
contributions of this paper are summarized as follows:

  1. 1.

We proposed an end-to-end learned framework for unified lossy and lossless
JPEG recompression. To the best of our knowledge, this is the first approach
that bridges the gap between lossy recompression and lossless recompression
for JPEG images.

  2. 2.

Experiments show that our method for lossy recompression can achieve
arbitrarily low distortion while keeping the bitrate always below the upper
bound, $i.e.$ the bitrate of lossless recompression.

##  2 Related Work

###  2.1 Overview of JFIF standard

As one of the most common image formats, JFIF [1] defines supplementary
specifications for the JPEG algorithm. During the compression phase, JFIF
first converts RGB signals into YCbCr color space. Then, since the human
visual system is more sensitive to the luma component compared to the chroma
component, commonly used YCbCr 4:2:0 subsampling is applied to reduce data
redundancy. After that, each component is divided into $8\times 8$ blocks. For
each $8\times 8$ block, discrete cosine transformation (DCT) is applied to get
frequency data, so-called dct coefficients. Subsequently, these three
components are quantized with two quantization tables. Y component is
quantized with one quantization table, while CbCr components share another
quantization table. Finally, quantized dct coefficients are compressed with
lossless Huffman coding. The decompression process is the reverse of the
compression process. For the sake of clarity, in the following, we denote the
QP 75 as the meaning of using the example table defined by Annex K of the JPEG
standard [1] with scale factor QP=75 according to the IJG scaling algorithm.

###  2.2 JPEG recompression methods

Lepton [13] was proposed by Horn et al., which mainly focuses on the
optimization of the entropy model and symbol representations. Lepton achieves
more than $20\%$ storage saving for JPEG lossless recompression.
JPEG-XL [14] is a new image coding standard, which supports both lossless and
lossy compression. For JPEG images, lossless recompression is also supported
by JPEG-XL. JPEG-XL achieves better compression ratio by a series of
optimization, like variable-size DCT block, local adaptive quantization, more
DC coefficient prediction modes, and using Asymmetric Numeral Systems in place
of Huffman coding.
MLCC [15] proposed by Guo et al. is the first learned lossless JPEG
recompression method. They proposed a Multi-Level Cross-Channel entropy model
for lossless recompression of JPEG images, which achieves state-of-the-art
performance. For JPEG images with QP 75, MLCC obtains about $30\%$ compression
savings.

###  2.3 Learned image compression

Lossy Compression In recent years, since Balle et al. [7] first proposed an
end-to-end image compression model based on autoencoder architecture, a series
of subsequent methods [8, 8, 9, 10, 11, 12] are proposed to improve the
compression ratio. They mainly focus on the optimization of rate estimation
with non-differential quantization and entropy estimation by designing a more
accurate entropy model. The newest method [12] has outperformed the intra-
coding of VVC/H.266 [6], which is the best hand-crafted method.
Lossless Compression Lossless compression is highly relevant to the
probabilistic model. Theoretically, any probabilistic model combined with an
entropy encoder can be used for lossless compression. With the development of
probabilistic models, lossless compression methods [17, 18, 19, 20] have also
made great progress. The recent method [19] has outperformed the popular
engineered lossless compression codecs, like PNG [3], WebP [4] and JPEG2000
[2].

However, all of these methods only consider the compression of RGB image data,
and cannot directly be applied to JPEG image data [15].

Fig. 2: Rate-distortion of lossless compression models and lossy compression
models for JPEG QP 75 recompression.

##  3 Method

###  3.1 Motivation

For the rate-distortion problem of JPEG recompression, without losing
generality, we consider bpp as the rate metric and mean squared error as the
distortion metric. Assuming the range of each pixel value is $[0,1]$. Let’s
take the latest lossy compression method ELIC [12] and lossless compression
method for JPEG recompression MLCC [15] as an example. For the sake of
convenience, we denote the rate of MLCC as $H$. As shown in Fig. 2, there is a
remarkable gap between the rate-distortion curve of ELIC and the rate of MLCC.
According to rate-distortion theory [16], there must be one theoretical convex
rate-distortion curve (marked as Theoretical RD in Fig. 2) for the MLCC-based
lossy compression method, which can achieve arbitrarily low distortion when
the rate is close to $H$. In other words, theoretically, there is a better
method for JPEG recompression compared to ELIC when the rate is close to $H$.

To sum up, as shown in Fig. 2, neither hand-designed lossy compression methods
like JPEG-XL nor learnable lossy compression methods like ELIC performs well
for JPEG lossy recompression. Once the rate exceeds $H$, it is better to use
MLCC rather than those lossy compression methods for JPEG recompression.

Fig. 3: Overall Architecture of our method. ED denotes entropy decoder followd
by inverse quantization operation, EC denotes entropy encoder. $Q_{t}$ and
$Q_{t}^{{}^{\prime}}$ represent the learned quantization table and learned
inverse quantization table respectively. The lossless compression model for
quantized DCT coefficients is adapted from [8], where Q denotes quantization,
AE and AD denote arithmetic coding and arithmetic decoding respectively. The
architecture of $G_{a}$, $G_{s}$, $H_{a}$, $H_{s}$ is in Tab. 1.

###  3.2 Framework

The framework of our method is shown in Fig. 3. For a given JPEG file with QP
75, it is first decoded to DCT coefficients with the entropy decoder followed
by inverse quantization operation. Then, the learned quantization table is
applied to quantize the DCT coefficients. Quantized DCT coefficients are
losslessly compressed into the bitstream by a lossless compression model with
Markovian hierarchical variational autoencoder architecture, which is adapted
from [8]. During the decompression phase, the bitstream is first decompressed
to quantized DCT coefficients. And then, it is compressed into a JPEG file
together with the learned inverse quantization table by the entropy encoder.

In our framework, the lossless compression model is used to optimize the rate,
the learned quantization table and inverse quantization table are used to
optimize the rate and distortion. We adopt iterative optimization way to
optimize each of these two items separately. Given a JPEG image as input, we
denote its DCT coefficients as $I$, the learned quantization table as $Q_{t}$,
the learned inverse quantization table as $Q_{t}^{{}^{\prime}}$ and the
parameters of the lossless model as $\theta$, respectively. Then the input of
the lossless compression model as shown in Fig. 3 is
$\boldsymbol{x}=\left\lfloor\frac{I}{Q_{t}}\right\rceil$, $\lfloor\cdot\rceil$
denotes straight-through estimator. And there are two hidden variables
$\boldsymbol{y},\boldsymbol{z}$ in the framework, where
$\boldsymbol{y}=G_{a}(\boldsymbol{x})$ and
$\boldsymbol{z}=H_{a}(\boldsymbol{y})$. The loss function for optimizing the
lossless compression model is formulated as follows:

$\displaystyle\mathcal{L}_{m}=R_{\theta}(\boldsymbol{x})$ (1)

where $Q_{t}$ is fixed during optimizing $\theta$. The loss function for
optimizing the learned quantization table and inverse quantization table is
formulated as follows:

$\displaystyle\mathcal{L}_{q}=\lambda_{r}R_{\theta}(\boldsymbol{x})+\lambda_{d}D(I,\boldsymbol{x}\cdot
Q_{t}^{{}^{\prime}})$ (2)

where $\lambda_{r}$ and $\lambda_{d}$ are hyper-parameters used to control the
rate-distortion trade-off. The parameters $\theta$ is fixed during optimizing
$Q_{t}$ and $Q_{t}^{{}^{\prime}}$. Following the setting in [8], the expected
negative entropy is supervised as the $R_{\theta}$ term and the MSE measure on
RGB domain is used as the $D$ term. Specifically, $R_{\theta}$ term and $D$
term are formulated as follows, respectively:

$\displaystyle R_{\theta}(\boldsymbol{x})=$
$\displaystyle-\mathbf{E}\left[\log{p(\boldsymbol{x},\hat{\boldsymbol{y}},\hat{\boldsymbol{z}})}\right]$
(3) $\displaystyle=$
$\displaystyle-\mathbf{E}\left[\log{p(\boldsymbol{x}|\hat{\boldsymbol{y}})}+\log{p(\hat{\boldsymbol{y}}|\hat{\boldsymbol{z}})}+\log{p(\hat{\boldsymbol{z}})}\right]$
$\displaystyle D(I,\boldsymbol{x}\cdot Q_{t}^{{}^{\prime}})$
$\displaystyle=\mathbf{E}\left[\left(f_{idct}(I)-f_{idct}(\boldsymbol{x}\cdot
Q_{t}^{{}^{\prime}})\right)^{2}\right]$

where
$\hat{\boldsymbol{y}}=\lfloor\boldsymbol{y}\rceil,\hat{\boldsymbol{z}}=\lfloor\boldsymbol{z}\rceil$,
$f_{idct}$ denotes inverse discrete cosine transformation,
$p(\boldsymbol{x}|\hat{\boldsymbol{y}})$,
$p(\hat{\boldsymbol{y}}|\hat{\boldsymbol{z}})$ are respectively modeled as
$i.i.d.$ gaussian distribution with mean 0, variance
$\sigma_{\boldsymbol{x}}=G_{s}(\hat{\boldsymbol{y}})$,
$\sigma_{\hat{\boldsymbol{y}}}=H_{s}(\hat{\boldsymbol{z}})$ and
$p(\hat{\boldsymbol{z}})$ is modeled as a fully factorized density model [8].

$G_{a}$ | $G_{s}$
---|---
Conv $5\\!\\!\times\\!\\!5$, s2 $\mid$ ReLU | Conv $3\\!\\!\times\\!\\!3$, s1 $\mid$ ReLU
Conv $5\\!\\!\times\\!\\!5$, s2 $\mid$ ReLU6 | Conv $3\\!\\!\times\\!\\!3$, s1 $\mid$ ReLU
Conv $5\\!\\!\times\\!\\!5$, s2 $\mid$ ReLU6 | Conv $3\\!\\!\times\\!\\!3$, s1
Conv $5\\!\\!\times\\!\\!5$, s2 | ConvTranspose $4\\!\\!\times\\!\\!4$, s2 $\mid$ ReLU
| ConvTranspose $4\\!\\!\times\\!\\!4$, s2 $\mid$ ReLU
| ConvTranspose $4\\!\\!\times\\!\\!4$, s4

$H_{a}$ | $H_{s}$
---|---
Conv $3\\!\\!\times\\!\\!3$, s1 $\mid$ ReLU | ConvTranspose $4\\!\\!\times\\!\\!4$, s2 $\mid$ ReLU
Conv $5\\!\\!\times\\!\\!5$, s2 $\mid$ ReLU | ConvTranspose $4\\!\\!\times\\!\\!4$, s2 $\mid$ ReLU
Conv $5\\!\\!\times\\!\\!5$, s2 | ConvTranspose $3\\!\\!\times\\!\\!3$, s1

Table 1: Details of $G_{a}$, $G_{s}$, $H_{a}$, $H_{s}$ modules in our lossless
compression model.

##  4 Experiments

###  4.1 Setup

Datasets Our training dataset consists of $8000$ images chosen from ImageNet
[21] validation set, we test all methods on Kodak [22] dataset. During
training, the RGB images are first randomly cropped into $256\times 256$ size
and then compressed into JPEG files with QP 75. During testing, the RGB images
are directly compressed into JPEG files with QP 75.
Implementations details We implement our experiments on PyTorch. We adopt a
two-step training strategy. We first train our lossless compression model to
converge on quantized DCT coefficients quantized by standard QP 75
quantization table. Then, for specific $\lambda_{r},\lambda_{d}$ pair in
equation 2, we iteratively finetune the lossless compression model and learned
quantization/inverse quantization tables which are initialized with standard
QP 75 quantization table. We use Adam optimizer with $\beta_{1}=0.9$,
$\beta_{2}=0.999$ and set batch-size to $16$ for training. In the first stage,
we set the initial learning rate to $10^{-4}$, training epochs to $2000$, and
decay the learning rate by $\gamma=0.9$ after every $500$ epochs. In the
second stage, we finetune the lossless compression model by setting the
initial learning rate to $10^{-3}$, training epochs to $50$, and decaying the
learning rate by $\gamma=0.94$ after every epoch. For learned
quantization/inverse quantization tables, we finetune it by setting the
initial learning rate to $5e-3$, training epochs to $30$, and decaying the
learning rate by $\gamma=0.8$ after every epoch. We iteratively finetune the
lossless compression model and the learned quantization/inverse quantization
tables for total of $8$ times. By setting
$(\lambda_{r},\lambda_{d})=\\{(0.1,100),(0.2,100),(0.4,100),(0.6,100),(0.8,100),(1,\\\
100),(2,100),(4,100),(8,100),(16,100),(32,100),(48,50),\\\
(64,20),(64,10),(64,5)\\}$, we obtain a series of models under different rate-
distortion trade-offs.

###  4.2 Performance

We compare the proposed model against the state-of-art learned based image
lossy compression method ELIC [12] and handcrafted image lossy compression
method JPEG-XL [14]. As shown in Fig. 4, our model can achieve arbitrarily
high PSNR and MS-SSIM while the rate never exceeds the rate of the lossless
compression model in the proposed framework. It means that the proposed model
unifies the lossless compression and lossy compression for JPEG images. $D=0$,
namely PSNR equal to $+\infty$ or MS-SSIM equal to $1$, denotes the lossless
compression while $D>0$ denotes the lossy compression. For the proposed model,
there is no longer a gap between lossless compression and lossy compression
for JPEG images. In the high rate region, namely near the lossless compression
rate, our model significantly outperforms ELIC [12] and JPEG-XL [14].

###  4.3 Ablation study

We test the proposed model with a fixed quantization table and inverse
quantization table, which are set to the JPEG quantization table with QP from
$10$ to $74$. As shown in Fig. 4, the proposed model with learned $Q_{t}$,
$Q_{t}^{{}^{\prime}}$ significantly outperforms the one with fixed $Q_{t}$,
$Q_{t}^{{}^{\prime}}$, which verifies the effectiveness of learned
quantization table and inverse quantization table for downstream lossless
compression model.

Fig. 4: Objective quality comparison on Kodak dataset when recompressing JPEG
images of QP 75.

##  5 Conclusion

In this paper, we reveal that there is a gap between lossless compression and
lossy compression for JPEG recompression. To solve this problem, we propose an
efficient JPEG recompression framework, which firstly unifies the lossless
compression and lossy compression for JPEG images. Our framework consists of
one learned quantization table, one learned inverse quantization table, and a
learned lossless compression model in the DCT domain. The proposed model can
achieve arbitrarily low distortion while the rate never exceeds the upper
bound, namely the rate of the lossless compression model in the framework.

A major limitation of this study is that the rate-distortion performance falls
too fast as the rate goes down. In the future, we will dive deeper into this
issue and extend our method to other specialized data compression and
transcoding problem.

##  6 ACKNOWLEDGMENTS

This work is supported by Baidu Inc. through Apollo-AIR Joint Research Center.

## References

  * [1] William B Pennebaker and Joan L Mitchell,  JPEG: Still image data compression standard,  Springer Science & Business Media, 1992.
  * [2] Majid Rabbani,  “Book review: Jpeg2000: Image compression fundamentals, standards and practice,” 2002.
  * [3] Thomas Boutell,  “Png (portable network graphics) specification version 1.0,”  Tech. Rep., 1997.
  * [4] Li Lian and Wei Shilei,  “Webp: A new image compression format based on vp8 encoding,”  Microcontrollers & Embedded Systems, vol. 3, 2012.
  * [5] Fabrice Bellard,  “Bpg image format,” https://bellard.org/bpg, 2015,  Accessed: 2022-11-16.
  * [6] Jens-Rainer Ohm and Gary J Sullivan,  “Versatile video coding–towards the next generation of video compression,”  in Picture Coding Symposium, 2018, vol. 2018.
  * [7] Johannes Ballé, Valero Laparra, and Eero P. Simoncelli,  “End-to-end optimized image compression,”  in International Conference on Learning Representations, 2017.
  * [8] Johannes Ballé, David Minnen, Saurabh Singh, Sung Jin Hwang, and Nick Johnston,  “Variational image compression with a scale hyperprior,”  in International Conference on Learning Representations, 2018.
  * [9] David Minnen, Johannes Ballé, and George D Toderici,  “Joint autoregressive and hierarchical priors for learned image compression,”  Advances in neural information processing systems, vol. 31, 2018\.
  * [10] Zhengxue Cheng, Heming Sun, Masaru Takeuchi, and Jiro Katto,  “Learned image compression with discretized gaussian mixture likelihoods and attention modules,”  in Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2020, pp. 7939–7948.
  * [11] Ge Gao, Pei You, Rong Pan, Shunyuan Han, Yuanyuan Zhang, Yuchao Dai, and Hojae Lee,  “Neural image compression via attentional multi-scale back projection and frequency decomposition,”  in Proceedings of the IEEE/CVF International Conference on Computer Vision, 2021, pp. 14677–14686.
  * [12] Dailan He, Ziming Yang, Weikun Peng, Rui Ma, Hongwei Qin, and Yan Wang,  “Elic: Efficient learned image compression with unevenly grouped space-channel contextual adaptive coding,”  in Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2022, pp. 5718–5727.
  * [13] Daniel Reiter Horn, Ken Elkabany, Chris Lesniewski-Lass, and Keith Winstein,  “The design, implementation, and deployment of a system to transparently compress hundreds of petabytes of image files for a file-storage service,”  in 14th USENIX Symposium on Networked Systems Design and Implementation (NSDI 17), 2017, pp. 1–15.
  * [14] Jyrki Alakuijala, Ruud Van Asseldonk, Sami Boukortt, Martin Bruse, Iulia-Maria Comșa, Moritz Firsching, Thomas Fischbacher, Evgenii Kliuchnikov, Sebastian Gomez, Robert Obryk, et al.,  “Jpeg xl next-generation image compression architecture and coding tools,”  in Applications of Digital Image Processing XLII. SPIE, 2019, vol. 11137, pp. 112–124.
  * [15] Lina Guo, Xinjie Shi, Dailan He, Yuanyuan Wang, Rui Ma, Hongwei Qin, and Yan Wang,  “Practical learned lossless jpeg recompression with multi-level cross-channel entropy model in the dct domain,”  in Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2022, pp. 5862–5871.
  * [16] Raymond W Yeung,  Information theory and network coding,  Springer Science & Business Media, 2008.
  * [17] Scott Reed, Aäron Oord, Nal Kalchbrenner, Sergio Gómez Colmenarejo, Ziyu Wang, Yutian Chen, Dan Belov, and Nando Freitas,  “Parallel multiscale autoregressive density estimation,”  in International Conference on Machine Learning. PMLR, 2017, pp. 2912–2921.
  * [18] Friso Kingma, Pieter Abbeel, and Jonathan Ho,  “Bit-swap: Recursive bits-back coding for lossless compression with hierarchical latent variables,”  in International Conference on Machine Learning. PMLR, 2019, pp. 3408–3417.
  * [19] Fabian Mentzer, Eirikur Agustsson, Michael Tschannen, Radu Timofte, and Luc Van Gool,  “Practical full resolution learned lossless image compression,”  in Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, 2019, pp. 10629–10638.
  * [20] Rianne van den Berg, Alexey A. Gritsenko, Mostafa Dehghani, Casper Kaae Sønderby, and Tim Salimans,  “Idf++: Analyzing and improving integer discrete flows for lossless compression,”  in International Conference on Learning Representations, 2021.
  * [21] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei,  “Imagenet: A large-scale hierarchical image database,”  in 2009 IEEE conference on computer vision and pattern recognition. Ieee, 2009, pp. 248–255.
  * [22] Eastman Kodak,  “Kodak lossless true color image suite (photocd pcd0992),” , 1993\.


# Multi-mode Heterodyne Laser Interferometry Realized via Software Defined
Radio

X. Lin1  M. T. Hartman1  S. Zhang1,2  S. Seidelin3  B. Fang1 and Y. Le Coq  1
1LNE-SYRTE, Observatoire de Paris, Université PSL, CNRS, Sorbonne Université,
Paris, France
2currently with Menlo Systems GmbH, Bunsenstraße 5, D-82152 Martinsried,
Germany
3Univ. Grenoble Alpes, CNRS, Grenoble INP and Institut Néel, 38000 Grenoble,
France<EMAIL_ADDRESS>
††journal: opticajournal††articletype: Research Article

The agile generation and control of multiple optical frequency modes combined
with the realtime processing of multi-mode data provides access to
experimentation in domains such as optomechanical systems, optical information
processing, and multi-mode spectroscopy. The latter, specifically spectroscopy
of spectral-hole burning (SHB), has motivated our development of a multi-mode
heterodyne laser interferometric scheme centered around a software-defined
radio platform for signal generation and processing, with development in an
entirely open-source environment. A challenge to SHB is the high level of shot
noise due to the laser power constraint imposed by the spectroscopic sample.
Here, we have demonstrated the production, detection, and separation of
multiple optical frequency modes to the benefit of optical environment sensing
for realtime phase noise subtraction as well as shot noise reduction through
multi-mode averaging. This has allowed us to achieve improved noise
performance in low-optical-power interferometry. Although our target
application is laser stabilization via SHB in cryogenic temperature rare-earth
doped crystals, these techniques may be employed in a variety of different
contexts.

##  1 Introduction

Optical heterodyne detection is a well developed technique, having seen
innumerable applications in science and industry, where the measurement of a
relatively weak optical signal is achieved by interference with a relatively
strong local oscillator to produce a detectable beat note. Implementations of
the technique vary widely through the domains of industrial application,
including fiber communication [1], LiDAR and atmospheric sciences [2],
astronomy such as radiometry [3] and astrophysical spectrometry [4],
instrumentation and metrology such as laser frequency stabilization [5, 6, 7],
and fundamental physics, for example gravitational-wave detection [8] and the
search for axion-like particles [9, 10], just to name a few.

In many use cases of spectroscopy the probe beam is constrained to low optical
powers, either by the sources including astronomical applications or to avoid
damaging the sample, for instance in many biological applications [11, 12].
Spectral hole burning (SHB) is a spectroscopic technique of imprinting
spectral lines in an inhomogeneously-broadened absorption spectrum [13]. This
is accomplished by using a pump-laser to selectively excite an absorptive dye
or dopant atoms/ions in the host to a higher state, and, upon returning to the
ground state, they occupy a different hyperfine level which is off-resonant
with pump light. The result is a narrow transmissive spectral line at the
pump-laser frequency that has been ‘burned’ into the absorption spectrum. When
probing these imprinted spectral features, the laser power used is constrained
due to the risk of ‘overburning’, where the use of relatively high power
causes, after some time, a broadening of the spectral hole leading to a
degradation in the spectroscopic sample [14].

Rare-earth ion doped crystals present an attractive medium for spectral-hole
spectroscopy by providing optical transitions with large coherence times
(order ms) and hyperfine splitting at the ground state which provide long-
lived population states at cryogenic temperatures. Spectral-hole spectroscopy
of rare-earth ions can be, and has been, used in many applications ranging
from classical [15, 16] and quantum [17, 18, 19, 20] information processing to
optomechanics [21, 22] and high-precision laser measurements and stabilization
[23, 24, 25, 26]. The domain where our experimental group concerns itself is a
metrology-purposed SHB heterodyne spectroscopy, namely its application to the
stabilization of laser frequency to spectral holes burned into a rare-earth
ion doped crystal, $\mathrm{Eu}^{3+}:\mathrm{Y}_{2}\mathrm{Si}\mathrm{O}_{5}$.
Each of these applications requires precise and agile manipulation of laser
frequency for the burning and probing of spectral holes. For our purpose, we
developed a laser control system based on a Software Defined Radio (SDR)
platform which provides a flexible, easily customizable, and greatly
extensible field-programmable gate array (FPGA)-based tool for acquisition and
control. We have implemented this platform successfully in our experiment,
demonstrating its use in laser stabilization to, as well as physical
characterization of, spectral holes [6, 7, 27, 28, 29, 30] .

In the context of laser frequency stabilization via SHB, the laser power
constraint (of order $10\,\mathrm{nW}$ for continuous probing of a single
spectral hole) to avoid hole overburning leads to a substantial level of shot
noise in the heterodyne probe signal resulting in laser frequency error,
random in time, in the detection process. In this paper, we address this
challenge by extending our SDR-based control system to produce, detect, and
demodulate multiple heterodyne beat signals, and we demonstrate the subsequent
reduced detection noise obtained by averaging over these modes. This novel
technique allows us to simultaneously probe several spectral holes and
synthesize an error signal from any combination of resulting heterodyne
spectroscopy signals. In addition to a simple average of the narrow spectral
hole signals, we demonstrate the production of additional monitor modes,
outside the band around the narrow spectral holes, which can measure and
subtract unwanted differential optical phase noise between a reference
heterodyne path and the spectroscopic sample path in real time. We believe
this technique and platform will benefit not only our ultra-stable laser
experiment, but any spectroscopy or heterodyne interferometry scheme where
either agile frequency control is required, or the power per frequency-mode is
constrained, or an appropriate use of additional probing modes may help reduce
the effect of technical noise sources, including other laser stabilization
experiments (for example, iodine-based laser stabilization [31]). This article
is organized as follows: a first section describes the general optical
apparatus that is used in conjunction with the SDR-based platform. In a second
part, we describe the main building blocks of the platform and the way it can
be used to perform various experimental tasks pertaining to spectral hole
burning in rare-earth-doped crystals at cryogenic temperatures. In a third
part, we present the results of sensing noise measurements exploiting the
fixed, low-power multi-mode heterodyne detection technique using our SDR-based
system. A last section discusses the advantages and the limitations of this
platform and explores ways to go beyond what is currently achieved.

##  2 Experimental Apparatus

Past iterations of our apparatus have been described succinctly at its various
development stages in our publications [6, 7, 27]. Our experimental apparatus
is primarily built for exploring the use of spectral hole burning in
$\mathrm{Eu}^{3+}:\mathrm{Y}_{2}\mathrm{Si}\mathrm{O}_{5}$ in the realization
of an ultra-stable laser of unprecedented short-term (near $1\,\mathrm{s}$)
frequency stability [7]. The ultimate goal of this experiment is to produce a
continuous wave (CW) laser whose performance allows optical lattice clocks
[32] to reach their fundamental frequency stability limit, namely the so-
called quantum projection noise limit [33]. While this limit can be breached
by quantum-engineered protocols (spin-squeezing, NOON states[34]…), those
implementations likewise necessitate a probe laser with ultra-high short term
frequency stability. A secondary target of our experiment is to perform high-
precision spectroscopy of micro-fabricated mechanical resonators realized via
rare-earth-ion-doped crystals, in order to probe and control their
oscillations down to the quantum level [22]. Both scientific targets utilize
the same apparatus and SDR-based platform, only the rare-earth-ion-doped
crystal is changed between the two applications: a bulk centimeter-scale
crystal for the ultra-stable laser application in contrast to various micro-
fabricated samples for the micro-resonator application. Here we give an
overview of the current setup and summarize parts that are particularly
important to understand the use of our SDR-based platform in the next
sections.

Figure 1:  Illustration of our current experimental setup. ppLN = periodically
poled lithium niobate optical frequency doubler; PDH = Pound–Drever–Hall; DM =
dichroic mirror; AOM = acousto-optic modulator; PLL = phase-lock-loop; PII =
proportional+double integrator; $\mathrm{\lambda/2}$ = half-wave plate; PM/SM
fiber = polarization-maintaining single-mode fiber; PD = photodetector.

An illustration of the experimental apparatus can be seen in Fig. 1. The
Master and Slave lasers are extended cavity diode lasers (Toptica DL Pro)
operating nominally at $1160\,\mathrm{nm}$. The master laser, which is
frequency locked via the Pound-Drever-Hall technique [35, 36] to a commercial
metrological-grade Fabry-Perot cavity (Stable Laser Systems), serves only as a
frequency reference and to provide an initial frequency stability to the slave
laser, which is offset phase-locked to the master via a phase-lock loop (PLL).
Adjusting the PLL offset allows us to coarse-tune the slaver laser frequency
to the interesting spectroscopic sites which occur between the resonances of
the reference cavity (Free Spectral Range $\approx 1.5\,\mathrm{GHz}$). The
slave laser, which is the main science laser for spectroscopy, is frequency
doubled in fiber-coupled waveguide ppLN crystals (NTT Electronics) to the
spectroscopic sites (nominally $580\,\mathrm{nm}$) of our typical sample,
$\mathrm{Eu}^{3+}:\mathrm{Y}_{2}\mathrm{Si}\mathrm{O}_{5}$.

The yellow light of the slave laser is frequency shifted in a double-pass
configuration [37] through a high-speed acousto-optic modulator (AA Opto-
electronics, MT350-A0.12-VIS) to the exact desired spectroscopic frequencies.
The high-speed and large bandwidth of the acousto-optic modulator (AOM) are
required to provide a high-fidelity reproduction of the wide range of
frequency components in its radio-frequency (RF) drive signal. It is at this
point where we apply the fine frequency control in the generation of the
multiple frequency modes in the slave laser beam, namely the spectroscopic
probe modes, an optical path length monitor mode, and a local oscillator mode
to beat the probe/monitor modes down to detectable RF heterodyne frequency.
The production of the AOM RF drive signal is realized through the output of
the Ettus USRP X300 FPGA-Based Software Defined Radio (SDR). The electrical RF
signal containing the desired frequency modes is programmed via the GNU Radio
environment in the controlling PC. This process is explained in great detail
in Section 3. The signal is Digitally Up-Converted (DUC) from base-band to the
SDR output center frequency, $f_{\mathrm{DUC}}=80\,\mathrm{MHz}$, which is
then up-mixed with a function generator (a $270\,\mathrm{MHz}$ sinusoid) to
reach the center drive frequency of the AOM,
$f_{\mathrm{AOM}}=350\,\mathrm{MHz}$, which is then amplified before driving
the AOM.

The frequency shifted light is typically composed of several optical spectral
modes with near identical spatial modes, thanks to the double-path AOM
configuration. It is passed through a polarization-maintaining single-mode
(PM/SM) fiber before being spatially separated into a reference path (‘Ref’,
which provides the reference heterodyne beat note after free-space) and sample
path (‘Xtal’, which, in spectral-hole experiments, contains the sample crystal
housed in a cryostat (MyCryoFirm Optidry) within a vacuum system). In laser
stabilization experiments, the desired signal is an optical-frequency-
dependant phase-shift, $\phi_{\mathrm{s}}$, of the optical mode along the Xtal
path which, due to the Kramers-Kronig relation, is proportional to the Hilbert
transform of the absorption in the crystal. When the mode is close to the
center of a spectral hole, this phase shift is simply proportional to the
frequency error, $\nu_{\mathrm{err}}\,[\mathrm{Hz}]$, between the laser probe-
mode frequency and the spectral hole center frequency with the proportionality
coefficient defining the frequency discriminant of the spectral hole,
$D_{\phi/\nu}\,[\mathrm{rad}\,\mathrm{Hz}^{-1}]$,

$\phi_{\mathrm{s}}=\nu_{\mathrm{err}}D_{\phi/\nu}.$ (1)

The reference and signal beat-notes are detected at the
$\mathrm{PD}^{\mathrm{Ref}}$ and $\mathrm{PD}^{\mathrm{Xtal}}$ avalanche
photodetectors (Thorlabs, APD410A/M), respectively. The beat signals are
digitally down converted in the SDR and the signal phase is recovered via
digitally extracting the phases of the various spectral components and
comparing these quantities between the Xtal path and the Ref path within the
computer control program. The clocks for the analog-to-digital converters
(ADCs) (ie. the clock of the SDR) and the clock in the external upconverting
function generators are all synchronized to a $10\,\mathrm{MHz}$ lab-wide
distribution of a reference clock signal. In a laser stabilization experiment,
this real-time complex data treatment provides the error signal which acts on
the slave laser PLL frequency offset actuating the slave laser towards the
spectral hole center frequency. In this paper we demonstrate a multi-mode
heterodyne interferometric readout where the simultaneous acquisition and
real-time analysis of several optical frequency channels in the SDR platform
is detailed in Section 3. The detection sensitivity measurement results are
discussed in Section 4.

##  3 The Software Defined Radio Based Platform

Agile RF signal generation and treatment are necessary to allow the creation
and use of multi-spectral-hole patterns. An SDR platform is particularly
suited for this application as it provides such capability both for emission
(TX) and reception (RX) channels, at a moderate cost and with relative ease of
programming and reconfiguration, where only knowledge of Python and, in some
rare cases, C/C++ programming languages is required. It is important to note
that some optimizations, such as the correct synchronization between the RX
and TX channels and latency minimization, are needed for the application of an
SDR platform to our research; significant instrumentation notes on the subject
are added in Appendix A.

For our SDR platform, we use commercially available open hardware (Ettus
Research X310) with open source software development tools (GNU Radio). The
heart of this hardware lies in a pre-programmed Field Programmable Gate Array
(FPGA), interfaced, on one side, by an ethernet port to a computer, and, to
the other side, to analog-to-digital and digital-to-analog converters (ADC and
DAC), operating at 200 megasamples per second (MSPS). These ADC and DAC are
themselves connected to analog frontends, called “daugtherboards”, which can
be interchanged, typically for use with signals of various carrier
frequencies. We use in our applications exclusively the BasicRX and BasicTX
daughterboards, which simply provide near direct access to the ADC and DAC
through a balun transformer. We note that many other, more complex,
daughterboards exist, allowing in particular upconversion to the GHz range by
analog in-phase/quadrature (I/Q) mixing with an on-chip programmable sinewave
synthesizer, however, in the context of our experiments, compactness is less
relevant than direct access to the signal throughout the entire detection
chain.

An ensemble of GNU Radio digital processing blocks, connected together through
a data flowgraph and controlled by a graphical user interface (GUI), allows
the realization of all the signal generation (through the SDR transmission
channels, TX) and recording (through the SDR reception channels, RX) as well
as the necessary signal processing for our experimental protocols.

###  3.1 Generation of Multiple Optical Frequency Modes:
Multi-channel TX in SDR

In TX, the SDR receives a stream of data generated by a GNU Radio-made data
flow algorithm from the control computer. This stream rate determines the
bandwidth of the desired spectral pattern to be generated (in our use case,
typically 20 MSPS). This data stream is sent to the FPGA which makes a Digital
Up Conversion (DUC), to a programmable carrier frequency, $f_{\mathrm{DUC}}$.
In essence, the FPGA is receiving a stream of samples from the computer,
processes it through an interpolating filter to reach 200 MSPS, and mixes the
result at the modulation frequency, generated by a numerically controlled
oscillator (NCO). Note the modulation/demodulation frequency can be set to
zero, should no actual up/down conversion be required. In principle, in this
implementation the SDR allows us to generate an arbitrary waveform pattern
centered around a programmable frequency in the 0-100MHz range, with a
bandwidth equal to the Nyquist frequency of the GNU Radio program sample rate.

####  3.1.1 TX data processing blocks

We use two distinct methods, depending on the intended experimental goal, for
generating the control RF signals of the experiment.

The first method is, in essence, time-based, and produces a single tone signal
whose frequency can be swept in a controlled way. This function is not
available in the standard libraries of GNU Radio; as such, we wrote a custom
data processing block “sweep-frequency” in C++. Ample documentation exists
about writing GNU Radio Blocks, a process which is almost entirely automatized
thanks to the gr-modtool utility which is part of the GNU Radio standard
implementation. The sweep-frequency block is, at its core, a numerically
controlled oscillator which generates a constant amplitude complex stream of
samples, for which the phase of each sample is equal to that of the previous
one plus an increment, and this increment is itself incremented at each sample
by a value proportional to the sweep rate of the frequency expressed in Hz/s.
Using such a block, the platform can easily generate frequency ramps of
controlled amplitude, sweep rate, and initial frequency and phase.

The second method is frequency based. It produces an ensemble of single tone
modes of controlled amplitudes and phases. It is implemented by using the GNU
Radio built-in inverse Fast Fourier Transform (iFFT) block, which, when fed
with the spectral pattern to realize, generates the corresponding complex
time-domain signal data stream. Our typical use of such an approach is for
realizing the illumination of the rare-earth-doped crystal with multiple
optical frequency tones in parallel. In such applications, the waveform to
generate is most easily described, in frequency domain, by a discrete Fourier
series composed of zero everywhere except at the places where a mode is to be
generated.

It should be noticed that these two approaches are complementary, as one is
based in the time domain and the other in the frequency domain. We combine the
two signals by numerically multiplying the time samples produced by these two
functions resulting in the convolution of the spectra in the frequency domain,
whereas adding them simply provides the sum of the two spectra. With the
desired spectrum programmed, the resulting RF signal generated by the TX
channel of the SDR is sent to drive the AOM in our setup (Fig. 1) to produce
the desired optical spectrum; the production of the optical frequency modes
contains some caveats, described in Appendix B.

Figure 2:  Measurements of optical spectra produced in our experiment. Here,
in addition to $1\times$LO-Mode, the RF signal driving the AOM produces
$1\times$Probe-mode (a), $3\times$Probe-modes (b), and $5\times$Probe-mode
(c). The spectra observed here are the result of beatnotes between the optical
LO-mode and the various optical probe modes.

Using these techniques, we demonstrate multi-mode optical spectra generated by
our experiment in Fig. 2. In this measurement, we simultaneously generate a
single fixed-frequency local-oscillator (LO) mode located
$f_{\mathrm{LO}}=2.1\,\mathrm{MHz}$ from the center of the IFFT generated
multi-mode spectrum to beat the multi-mode spectrum down to a signal
detectable by our APDs. Within each spectrum, the RF signals produced by the
individual desired beat notes vary by less than $0.15\,\mathrm{dBm}$ from
minimum to maximum translating to a variation in beat-note optical power of
less than $2\,\%$. The origin of the additional sidebands is unknown and
effects such as AOM response linearity and RF-chain/detection saturation have
been investigated. Noting the log scale, the sidebands are negligibly small,
and, without regard to their magnitude, do not contribute to the error signal,
which is generated by demodulation at the specific desired frequencies
followed by aggressive filtering.

###  3.2 Reception of Multiple Optical Frequency Modes:
Multi-channel RX in SDR

In the SDR, the input signal is digitized and then digitally down-converted
(DDC) internally by a digitally synthesized demodulation signal of
programmable frequency in the 0-100MHz range, $f_{\mathrm{DDC}}$, referenced
to an external 10 MHz clock. We set the digital down-conversion frequency to
frequency separation between the center of the IFFT multi-mode spectrum and
the local oscillator mode so that $f_{\mathrm{DDC}}=f_{\mathrm{LO}}$. The
resulting baseband signal is streamed to the connected control computer, where
a driver (called Universal Hardware Driver “UHD”) receives it and sends it to
a GUI-designed/python-programmed data flow treatment algorithm realized in the
GNU Radio framework. We typically use a 2 MSPS data rate from the FPGA to the
computer, which is sufficient for our applications, but have also successfully
worked at higher sampling rates.

####  3.2.1 RX data processing blocks

After DDC, the decomposition of the signals generated by the two photodiodes
into their separate frequency channels is realized on the control computer via
poly-phase channelizer blocks, a standard component of the GNURadio library
[38]. Typically, our current implementation uses the input stream of data at 2
MSPS from a photodiode to generate 10 streams of data at 200 kSPS per channel,
equidistant in frequency-space, separated typically by
$f_{\mathrm{Sp}}=200\,\mathrm{kHz}$. These streams of data correspond to
individual probing channels, and will allow simultaneous probing of several
narrow or large spectral holes previously photo-imprinted on the crystal. The
poly-phase channelizer structure allows efficient algorithmic implementation
with a large control over the filtering of channels so as to minimize cross-
talk and distortion. As the base of our polyphase channelizer, we typically
use a finite impulse response equiripple low-pass filter [39] with pass band
of 10 kHz (less than 1 dB ripple in-band), stopband at 40 kHz and 120 dB
attenuation out-of-band, which offers sufficient performance for our
applications and is implemented in less than 250 taps. Other filters with
different performances may of course be used depending on the application.

For the comparison of the sample heterodyne phase to the reference heterodyne
phase, the complex data stream of each individual frequency channel (10 in our
current setup) coming from the channelized DDC Xtal signal is divided by the
corresponding frequency channel from the similarly channelized DDC Ref signal
and the result is separated into amplitude and phase components. In our
typical application to multi-spectral-hole spectroscopy, for a frequency
channel where light was effectively applied (i.e., from the appropriate TX
channel signal), these quantities correspond to the absorption and phase
retardation corresponding to the Kramers-Kronig relations resulting from
transmission through a narrow spectral hole. Frequency channels where no light
was actually applied exhibit only noise-related content and should be
discarded.

##  4 Multi-channel Low-power Heterodyne Detection

In the previous sections we discussed the production (using an SDR-driven AOM)
and detection (SDR-digitized/channelized signals from two APDs) of multiple
optical frequency modes. Here we apply this to an interferometric detection
scenario where we demonstrate the noise reduction capabilities of this multi-
mode interferometry.

###  4.1 Avalanche-enhanced Shot Noise

Figure 3:  Measured noise power at different incident laser powers for
avalanche multiplication gain $G_{\mathrm{A}}=10$ (‘$\Diamond$’ Markers) and
for $G_{\mathrm{A}}=100$ (‘$+$’ Markers). The calculated noise contributions
for $G_{\mathrm{A}}=10$ (Solid Lines) and $G_{\mathrm{A}}=100$ (Dotted Lines)
are also plotted. The photodetector (NEP) noise contribution is plotted as the
Black Line. The levels of photon shot noise are plotted in Blue Lines and as
avalanche-enhanced shot noise in Green Lines. The calculated sum total noises
are plotted in Red Lines.

In the experiment, we have selected avalanche photodetectors for fast
detection of low optical powers. For this demonstration, we used the Thorlabs
APD410A amplified photoreceiver, which has an adjustable avalanche current
gain (from $G_{\mathrm{A}}=10\rightarrow 100$) and uses a Hamamatsu S12023-10
avalanche photodiode. This photoreceiver has an electronic noise at the output
that varies greatly with measurement frequency. Around our measurement beat
frequency (typically 1 MHz), the manufacturer specifications indicate an
output noise spectral density of $2\times 10^{-14}\,\mathrm{V^{2}Hz^{-1}}$ or
$-95\,\mathrm{dBm}$ into $50\,\mathrm{\Omega}$ with a $1\,\mathrm{kHz}$
measurement bandwidth.

For low-power precision interferometry, shot noise is a primary concern.
Additionally, avalanche photodetectors manifest a higher level of measured
shot noise due to the additional statistical process of the avalanche effect
[40]. This so-called excess noise factor,

$F_{\mathrm{x}}={\kappa}{G_{\mathrm{A}}}+(1-{\kappa})\left(2-\frac{1}{G_{\mathrm{A}}}\right),$
(2)

depends on the ratio of the hole to electron ionization rates, $\kappa$. For
the Hamamatsu S12023-10 in the avalanche diode, this value is $\kappa=0.4$ at
our measurement wavelength. For an average photon-incidence rate $\bar{n}$ on
a detector with quantum efficiency $\eta$, the single sided power spectral
density of avalanche-enhanced current shot noise (ASN) is given by

$I_{ASN}=2q_{\mathrm{e}}^{2}F_{\mathrm{x}}G_{\mathrm{A}}^{2}\eta\bar{n},$ (3)

where $q_{\mathrm{e}}$ is the electron charge. These modeled noise
contributions, alongside their sums, are plotted in curves as a function of DC
incident power in Fig. 3. These are in close agreement with the measured noise
power plotted with ‘$\Diamond$’ markers (avalanche multiplication set to
$G_{\mathrm{A}}=10$) and with ‘$+$’ markers ($G_{\mathrm{A}}=100$). As seen in
the figure, when going from a $10\times$ to $100\times$ avalanche gain, the
$20\,\mathrm{dB}$ increase in signal is accompanied by a $\geq
23\,\mathrm{dB}$ increase in noise across the power measurement range implying
the best SNR case is achieved using the detector adjusted to
$G_{\mathrm{A}}=10$ avalanche gain.

###  4.2 Illustrative Case: Single Channel Low-power Heterodyne Detection

In application to the probing of delicate spectral structures (e.g. burned
spectral holes) these features necessitate the use of low optical-power
spectroscopy. In a simple illustrative case, we generate a single probe mode
($\mathrm{p}_{1}$), at optical frequency $\nu_{\mathrm{p}_{1}}$, such that an
optical power $P_{\mathrm{p}_{1}}=10\,\mathrm{nW}$ is incident on each of the
two photodetectors, $\mathrm{PD}^{\mathrm{Xtal}}$ and
$\mathrm{PD}^{\mathrm{Ref}}$. This low-power mode would probe a narrow
spectral hole and accumulate a phase proportional to its frequency difference
with the spectral-hole center frequency. Alongside the probe mode, we generate
the local-oscillator mode (LO) with power $P_{\mathrm{LO}}=1\,\mathrm{\mu{W}}$
at optical frequency $\nu_{\mathrm{LO}}$, spaced a distance
$f_{\mathrm{LO}}=\nu_{\mathrm{LO}}-\nu_{\mathrm{p}_{1}}$ from the probe mode,
to produce a beat-note at $f_{\mathrm{LO}}$, within the bandwidth of our
photodetectors. The LO passes through a broad-band transmissive hole making
its phase sensitivity to frequency shifts orders of magnitude lower than the
probe narrow-hole signal. An illustration of this scenario can be seen in Fig.
4(a). In this configuration (as well as all other configurations described
later on), the DDC chains of the detection path are set to demodulate against
the $f_{\mathrm{LO}}$ frequency, therefore putting the beatnotes between the
LO mode and the probe mode at baseband in the GNU Radio data treatment. In our
experimental setup, we derive our interferometric signal from the difference
between the reference path ${\mathrm{LO}}/{\mathrm{p}_{1}}$ heterodyne beat
note phase ($\phi^{\mathrm{Ref}}_{\mathrm{p}_{1}}$) and the crystal path
${\mathrm{LO}}/{\mathrm{p}_{1}}$ beat note phase
($\phi^{\mathrm{Xtal}}_{\mathrm{p}_{1}}$),

$\phi_{\mathrm{p}_{1}}(t)=\phi^{\mathrm{Xtal}}_{\mathrm{p}_{1}}(t)-\phi^{\mathrm{Ref}}_{\mathrm{p}_{1}}(t).$
(4)

In this simple configuration, we have seen that technical noise sources begin
to exceed the shot noise sensing level at long time scales
($>10\,\mathrm{s}$). A more complicated probing scheme, which is described in
section 4.4, was demonstrated to effectively subtract these noise sources.

We test the conversion of amplitude to phase noise by intentionally applying,
in software, a fixed amplitude modulation of index $m=0.1$ for both the LO and
the probe mode in the TX chain, and measuring the resulting phase modulation
in the RX detection chain. In this way, we measured the coefficient of
amplitude to phase fluctuations to be $7.1\times 10^{-4}$ rad per 100%
amplitude modulation ($m=1$) in our operation regime. This very high level of
rejection of amplitude noise, owing largely to the digital approach to our
differential measurements and from the naturally small signals on the
avalanche photodiodes (thus operating far from saturation) makes this
detection scheme largely immune to laser power fluctuations.

###  4.3 Novel Case: Multi-probe-mode Averaging for Reduction of Incoherent
Noise Sources

(a)

(b)

Figure 4:  Example spectral hole configurations, illustrated in Blue,
alongside the measurement laser modes, drawn in Yellow, at the corresponding
optical frequencies ($\nu$). Figure 4(a) consists of one spectral probe mode
($\mathrm{p}_{1}$) alongside a local-oscillator mode ($\mathrm{LO}$) to form
the heterodyne detection frequency, $f_{\mathrm{LO}}$, as described in Section
4.2. Figure 4(b) consists of two spectral probe modes ($\mathrm{p}_{1}$,
$\mathrm{p}_{2}$) plus one phase noise monitor mode ($\mathrm{m}_{1}$)
alongside the local-oscillator ($\mathrm{LO}$) as described in Section 4.4.
Note the monitor and local oscillator modes transmit through broad
transmission bands to avoid spectroscopic interaction with the sample.

In the cases where optical beatnote power is restricted by the system being
probed (e.g. low power necessary in the probing of delicate spectral
structures, such as burned spectral holes), we can improve the phase detection
sensitivity by implementing a multi-mode heterodyne detection scheme. In this
case, along side one relatively high power local oscillator (which is removed
in frequency space from the spectral features we wish to measure), we create
multiple low-power probe modes, separated in frequency space by
$f_{\mathrm{Sp}}$, chosen to occur at the center frequencies of several
frequency-separated spectral features.

Over $N$ probe modes, for each individual probe mode, $i$, the phase
difference between the reference and sample paths, $\phi_{\mathrm{p}_{i}}$, is
calculated individually and then averaged evenly to find the mean differential
phase between the reference and sample paths,

$\widebar{\phi_{\mathrm{p}}}(t)=\frac{1}{N}\sum_{i=1}^{N}\phi_{\mathrm{p}_{i}}(t).$
(5)

The purpose of this exercise is to decrease the detection noise without
increasing the power of the individual probe modes. For some uncorrelated
noise power, $S_{\phi_{\mathrm{p}_{i}}}$, in individual differential phase
signals, $\phi_{\mathrm{p}_{i}}$, the total noise power in the mean
differential phase is the sum

$S_{\widebar{\phi_{\mathrm{p}}}}=\sum_{i=1}^{N}\left(\frac{\partial{\widebar{{\phi_{\mathrm{p}}}}}}{\partial{\phi_{\mathrm{p}_{i}}}}\right)^{2}S_{\phi_{\mathrm{p}_{i}}}=\frac{1}{N^{2}}\sum_{i=1}^{N}S_{\phi_{\mathrm{p}_{i}}}.$
(6)

For the typical case where
$S_{\phi_{\mathrm{p}_{1}}}\\!\\!\\!=\\!S_{\phi_{\mathrm{p}_{2}}}\\!\\!\\!=\\!\ldots\\!=\\!S_{\phi_{\mathrm{p}_{N}}}\\!\\!\\!\equiv\\!S_{\phi_{\mathrm{p}}}$
(e.g. shot noise in $N$ equal optical-power probe modes) the expression
becomes:

$S_{\widebar{{\phi_{\mathrm{p}}}}}=\frac{S_{\phi_{\mathrm{p}}}}{N}.$ (7)

We demonstrate this noise reduction by calculating the power spectral
densities of the recorded time series of detection noise, that is,
differential phase noise between the Ref and Xtal paths. This phase
sensitivity measurement was made in the condition where no spectroscopic
sample is in place (i.e. Ref path is through freespace and the Xtal path is
through the cryostat, but with no crystal installed). For the case of a single
probe mode ($N=1$), we measured a mean power spectral density of detection
noise dominated by avalanche-enhanced shot noise at
$S_{\phi_{\mathrm{p}}}=-93.9\,\mathrm{dB({rad^{2}}{Hz^{-1}})}$ in the Fourier
frequency range from $10\,\mathrm{Hz}$ to $1\,\mathrm{kHz}$. The cases of
$N=3$ and $N=5$ probe modes are expected to result in a
$10\log_{10}(\frac{1}{3})\approx-4.8\,\mathrm{dB}$ and
$10\log_{10}(\frac{1}{5})\approx-7.0\,\mathrm{dB}$ reduction in detection
noise compared to the single probe-mode case, respectively. These are
confirmed by our measurements, plotted in Fig. 5.

Figure 5:  The PSD of the detection phase noise. Averaging the digitally
demodulated signals produces the predicted $-4.8\,\mathrm{dB}$ (for $N=3$,
red) and $-7.0\,\mathrm{dB}$ (for $N=5$, green) reduction in shot-noise
compared to the single mode case (black).

###  4.4 Complete Case: Reducing Differential Path Noise Using a Monitor Mode

The real power of our platform is revealed when utilizing its multiple
detection channels in a more complex way. Of particular interest is the use of
dedicated channels for the suppression of common-mode differential path noise
(CMDN), that is, phase noise which is differential between the two separate
optical paths but common amongst the modes within each path. As sources of
CMDN are present and correlated in each of the probe modes, they sum
coherently as we average over several probe modes, and subsequently do not
reduce when averaging over probe modes. For this we take a different approach,
adding additional “Monitor Modes” ($\mathrm{m}_{i}$) monitoring differential
path noise. We produce the monitor modes using the same method as the probe
modes, and, like the probe modes, they are digitally demodulated in the SDR
after detection. By adding monitor modes at frequencies distanced from the
spectral features of interest, we are able to measure undesired differential
phase fluctuations between the reference and sample paths independent from the
desired signal, the spectroscopic phase contribution. An illustration of the
generated laser spectrum (yellow) is laid against an illustration of a
spectral hole pattern (blue) in Fig. 4(b). In this example, two probe modes
are generated at the location of two narrow line spectral holes, while a
single monitor mode is produced between them at the location of a broad
transmissive hole. Narrow holes are pre-burned by selectively pumping the ions
at a particular resonance frequency using a relatively high-power narrow-
linewidth laser, frequency stabilized to a cavity, over a short time (of order
$1\,\mathrm{{\mu}{W}}$ over $1\,\mathrm{{s}}$); a broad spectral hole is
produced by sweeping the hole-burning laser over the desired frequency range
to burn. A local oscillator, which brings these from optical frequencies to a
measurable RF frequency, transmits through a broad spectral hole at a distance
$f_{\mathrm{LO}}=\nu_{\mathrm{LO}}-\nu_{\mathrm{m}_{1}}$ from the monitor
mode.

Sources of CMDN include differential optical path length noise and phase noise
resulting from unwanted on-axis interference from undesired reflections at
surfaces at normal incidence. In our experiment, the sample path contains
several opportunities for these unwanted Fabry-Perot effects as the sample is
housed within several thermal shields inside a vacuum chamber, each with a
glass window to transmit the sample beam. These parasitic etalon effects are
of particular concern as the resulting phase noise has a frequency dependence,
causing an imperfect cancellation when demodulating against the local-
oscillator mode alone, which is limited to a cancellation factor set by the
ratio of the beat frequency to optical frequency. This prompts us to use a
monitor mode centered in frequency space between probe modes. When using an
even number of probe modes, this symmetric configuration in frequency space of
the probe/monitor modes allows for a more complete first-order cancellation of
the phase noise resulting from this parasitic interference, a detail discussed
further in Appendix C. Furthermore, a common monitor mode measures and allows
us to subtract additional differential phase noise between the channels, such
as timing differences between the two ADC channels and electrical signal delay
between the photodetectors.

Figure 6:  A simplified diagram of the RX signal detection in the experiment
for the case of 4 modes: $2\times$probe-mode, $1\times$monitor-mode,
$1\times$LO-mode case illustrated in Fig. 4(b). The details are described in
the text.

An illustration of the signal flow for this case (Fig. 6) can demonstrate the
noise subtraction using a monitor mode. Following this signal flow, the laser
field in the reference path accumulates phase according to its optical path
before being detected at $\mathrm{PD}^{\mathrm{Ref}}$, while the laser field
in the sample path additionally accumulates the spectroscopic signal phase and
is detected at $\mathrm{PD}^{\mathrm{Xtal}}$. These signals are digitized
separately and digitally demodulated at $f_{\mathrm{DDC}}=f_{\mathrm{LO}}$,
bringing the signals into baseband around the center frequency signal, i.e.
the monitor mode. The I-Q signals are passed to the control program on the
computer where ‘Xtal’ path signals and ‘Ref’ path signals are filtered by the
Polyphase Channelizer into their frequency components. The phase difference,
$\phi_{x}=\phi^{\mathrm{Xtal}}_{x}-\phi^{\mathrm{Ref}}_{x}$, is calculated for
each corresponding component
($x=\mathrm{p}_{1},\mathrm{m}_{1},\mathrm{p}_{2}$). Using the channelized
monitor phase signals, the mean over $M$ accumulated differential monitor-mode
phases ($\phi_{\mathrm{m}_{i}}$) can be subtracted from the simple
differential probe-mode phase average (previous case given in Eqn. (5)) in
realtime, leaving only the average of the desired signal phases:

$\widebar{\phi_{\mathrm{s}}}(t)=\frac{1}{N}\sum_{i=1}^{N}\phi_{\mathrm{p}_{i}}(t)-\frac{1}{M}\sum_{i=1}^{M}\phi_{\mathrm{m}_{i}}(t).$
(8)

Figure 7:  Overlapping Allan deviation measurement of phase fluctuations
projected as frequency error w.r.t. average center frequency of narrow
spectral holes with typical discriminant
$D_{\phi/\nu}=300\,\mathrm{{\mu}{rad}}\,\mathrm{Hz}^{-1}$. Compared to the
single mode case (‘$\times$’ Markers), the inclusion of a monitor mode to
subtract differential path length noise reduces sensing noise by an order of
magnitude at $10^{3}\,\mathrm{s}$ of integration time (‘$\Diamond$’ Markers).
For visual reference, a line $\propto\tau^{-1/2}$ has been included on the
plot.

We see the result of CMDN best in long phase noise measurements. To
demonstrate the baseline measurement noise, a three hour long Allan deviation
of the Ref-Xtal differential phase, which included the vacuum chamber and
thermal shields (but no spectroscopy samples) in the sample path, is plotted
in Fig. 7. Here we projected the measured phase noise as fractional frequency
difference between a single probe mode and the center frequency of a
hypothetical spectral hole with a typical frequency discriminant
$D_{\phi/\nu}=300\,\mathrm{{\mu}{rad}}\,\mathrm{Hz}^{-1}$. A measurement of
phase noise using a single probe mode beat down with a local-oscillator mode
is plotted with ‘$\times$’ markers (the case described in Section 4.2). Here
we see detection white noise at high frequencies dominates the short
averaging-time ($\tau$) measurement points, rolling off as $\tau^{-1/2}$ until
about $10\,\mathrm{seconds}$ where sources of CMDN take over at low
frequencies/long integration times. A measurement that utilized a monitor mode
is plotted with ‘$\Diamond$’ markers, this result was obtained using one local
oscillator mode, one monitor mode, and two probe modes (the example
illustrated in Fig 4(b)). Using the monitor mode to subtract CMDN, we obtain
an order of magnitude reduction in sensing noise at long timescales.

####  4.4.1 Comments on detection noise using a monitor mode

Note in Eqn. (8) that, for successful subtraction of CMDN we subtract the
monitor mode outside of the averaging of the probe modes, thus it becomes an
additional contribution to shot-noise. The total noise power can be described
by the sum of uncorrelated noise powers of the individual modes,

$\displaystyle S_{\widebar{{\phi_{\mathrm{s}}}}}$
$\displaystyle=\sum_{i=1}^{N}\left(\frac{\partial{\widebar{{\phi_{\mathrm{s}}}}}}{\partial{\phi_{\mathrm{p}_{i}}}}\right)^{2}S_{\phi_{\mathrm{p}_{i}}}+\sum_{i=1}^{M}\left(\frac{\partial{\widebar{{\phi_{\mathrm{s}}}}}}{\partial{\phi_{\mathrm{m}_{i}}}}\right)^{2}S_{\phi_{\mathrm{m}_{i}}}$
(9)
$\displaystyle=\frac{1}{N^{2}}\sum_{i=1}^{N}S_{\phi_{\mathrm{p}_{i}}}+\frac{1}{M^{2}}\sum_{i=1}^{M}S_{\phi_{\mathrm{m}_{i}}}.$
(10)

For the case where probe modes have equal noise power
($S_{\phi_{\mathrm{p}_{1}}}\\!\\!\\!=\\!S_{\phi_{\mathrm{p}_{2}}}\\!\\!\\!=\\!\ldots\\!=\\!S_{\phi_{\mathrm{p}_{N}}}\\!\\!\\!\equiv\\!S_{\phi_{\mathrm{p}}}$)
and monitor modes have equal noise power
($S_{\phi_{\mathrm{m}_{1}}}\\!\\!\\!=\\!S_{\phi_{\mathrm{m}_{2}}}\\!\\!\\!=\\!\ldots\\!=\\!S_{\phi_{\mathrm{m}_{M}}}\\!\\!\\!\equiv\\!S_{\phi_{\mathrm{m}}}$)
the expression becomes:

$S_{\widebar{{\phi_{\mathrm{s}}}}}=\frac{S_{\phi_{\mathrm{p}}}}{N}+\frac{S_{\phi_{\mathrm{m}}}}{M}.$
(11)

For the case of shot-noise limited detection, the power spectral density of
phase noise in each individual probe or monitor mode,

$S_{\phi}=F_{\mathrm{x}}\frac{{h}{\nu}}{4{\eta}{\bar{P}}},$ (12)

is inversely proportional to the mean power, $\bar{P}$, in that mode. If we
have individual probe modes of equal power
($\widebar{P_{\mathrm{p}_{1}}}\\!\\!\\!=\\!\widebar{P_{\mathrm{p}_{2}}}\\!\\!\\!=\\!\ldots\\!=\\!\widebar{P_{\mathrm{p}_{N}}}\\!\\!\\!\equiv\\!\widebar{P_{\mathrm{p}}}$)
and monitor modes of equal power
($\widebar{P_{\mathrm{m}_{1}}}\\!\\!\\!=\\!\widebar{P_{\mathrm{m}_{2}}}\\!\\!\\!=\\!\ldots\\!=\\!\widebar{P_{\mathrm{m}_{M}}}\\!\\!\\!\equiv\\!\widebar{P_{\mathrm{m}}}$)
the total noise spectral density is then given by

$S_{\widebar{{\phi_{\mathrm{s}}}}}=F_{\mathrm{x}}\frac{{h}{\nu}}{4{\eta}}\left(\frac{1}{{N}{\widebar{P_{\mathrm{p}}}}}+\frac{1}{{M}{\widebar{P_{\mathrm{m}}}}}\right).$
(13)

The shot noise contribution from monitor mode can be made arbitrarily small by
fixing ${M}{\widebar{P_{\mathrm{m}}}}\gg{N}{\widebar{P_{\mathrm{p}}}}$.
However, if we have a total working power limit, $P_{\mathrm{t}}$, set by
either the available laser power for mode production or by the saturation
limit of the detectors, to be split amongst all the probe/monitor modes such
that
$P_{\mathrm{t}}={N}{\widebar{P_{\mathrm{p}}}}+{M}{\widebar{P_{\mathrm{m}}}}$,
we see that Eqn. (13) is minimized when
${N}{\widebar{P_{\mathrm{p}}}}={M}{\widebar{P_{\mathrm{m}}}}$. Additionally,
in the application to SHB, we wish to limit the power in the probe and monitor
modes to avoid ‘overburning’, or the spectral distortion of the spectral
holes.

Returning to the example of probing burned spectral holes, it is impractical
to make several broad transmissive spectral holes to allow multiple monitor
modes, and, in fact, it is mostly the probe-mode power that is limited by the
delicate narrow spectral features, much more so than a monitor mode which
interacts with a broad spectral feature. Thus it is possible to fix the number
of monitor modes to $M=1$ and minimize Eqn. (13) with
${\widebar{P_{\mathrm{m}}}}=N{\widebar{P_{\mathrm{p}}}}$.

Figure 8:  The PSD of the detection phase noise using a monitor mode for CMDN
subtraction. The Black curve represents Case I, the Red curve represents Case
II and the Green curve represents Case III, described in the text.

We demonstrate these results in measurements of the noise power spectral
density for three cases:
Case I: One single probe mode (N=1) of power
$\widebar{P_{\mathrm{p}}}=10\,\mathrm{nW}$ and no monitor mode (black)
Case II: Two probe modes (N=2) and one monitor mode
($\widebar{P_{\mathrm{m}}}=2\widebar{P_{\mathrm{p}}}=20\,\mathrm{nW}$) (red)
Case III: Four probe modes (N=4) and one monitor mode
($\widebar{P_{\mathrm{m}}}=4\widebar{P_{\mathrm{p}}}=40\,\mathrm{nW}$) (green)
It is clear from Eqn. (11) and Eqn. (13) that moving from Case I to Case II
the shot-noise contribution from the probe modes are halved, however we add an
equal contribution from a monitor mode resulting in the same shot-noise level
for both cases. Moving to Case III we would then expect a
$10\log_{10}(\frac{2}{4})\approx-3\,\mathrm{dB}$ reduction in shot-noise.
These are confirmed by our measurements plotted in Fig. 8. For the long-term
Allan deviation measurement (Fig. 7) the cases correspond roughly to Case I (1
probe + 0 monitor) and Case II (2 probe + 1 monitor), however, with the
difference being for that measurement we did not increase the power in the
monitor mode (i.e.
$\widebar{P_{\mathrm{m}}}=\widebar{P_{\mathrm{p}}}=10\,\mathrm{nW}$) resulting
in a slightly higher level of shot noise ($\propto\tau^{-1/2}$) in the 2 probe
+ 1 monitor mode measurement in Fig. 7.

##  5 Conclusion

In the context of our SHB laser stabilization experiment, we have demonstrated
a scheme for flexible multi-mode heterodyne interferometry. The scheme’s
agility and ease of implementation derive from the use of open-source hardware
(Ettus USRP SDR) and software (GNU Radio) platforms. We have shown an
adjustable interrogating setup and detection noise reduction that allows us to
work around the power-limit restrictions of spectral hole probing. We have
analysed the detection noise performances in the system and the results
suggest that the simultaneous probing ability and the noise performance of
this novel technique show the potential to meet the requirement of the SHB-
based ultra-stable lasers. Additionally, the strategic utilization of
additional monitor modes enables the subtraction of technical noise sources,
such as parasitic interferometric effects. We are confident that these
techniques and platform will bring benefits not only to our ultra-stable laser
experiment but also to any spectroscopy or heterodyne interferometry schemes
that require agile frequency control or encounter power-per-frequency mode
limitations (e.g., laser locking to iodine[31]). The relative ease of use of
the SDR-based platform has allowed a full noise characterization of the system
and tests of spectral hole locking with a unity gain frequency $\approx
100\,\mathrm{Hz}$. Future work will involve investigating ways to decrease
signal delay in the system to improve locking bandwidth.

##  Appendix A Appendix: Synchronisation Between RX and TX, SDR-based
Servoloops, and Other Optimizations

For direct analysis of the transmission heterodyne signal the SDR and
GNURadio-based platform is very well suited as, in addition to saving data for
post processing, it allows direct access to the multiple channels of
information. However, it should be mentioned, that with no special care taken,
the time offset and thus relative phase between the two RX channels is random
(although constant over a given run of the program). In order to provide a
minimal and reproducible phase offset between the two inputs of the SDR
platform, it is necessary to synchronize their DDC chains. This is done simply
in software at the beginning of the program, by using timed commands to set
the two input channels’ DDC chain parameters simultaneously. The corresponding
Python code is:

⬇

now = self.USRP_RX.get_time_now()

cmd_time = now + uhd.time_spec(1)

self.USRP_RX.set_command_time(cmd_time)

self.USRP_RX.set_center_freq(freq, 0)

self.USRP_RX.set_center_freq(freq, 1)

self.USRP_RX.clear_command_time()

where the `USRP_RX` object is instantiated in the flowgraph by a “UHD: USRP
Source” block and `freq` is the demodulation frequency of the DDC chain of the
Ettus X310.

Similarly, in closed-loop implementations where the error signal that is
produced by the multichannel heterodyning is used to generate a correction
signal to be emitted by one of the TX channels, it is desirable to minimize
the time of propagation of the signal between the RX and TX channels.
Furthermore, making this propagation time repeatable is also desirable as
otherwise the in-loop servo parameters would need to be fine tuned at each run
of the program for optimal performance (defined in particular as maximal gain
and bandwidth while maintaining servo-loop stability). The GNURadio framework
is not particularly well suited for this application as it is largely designed
to optimize throughput, not latency. However, it is possible to implement a
closed-loop application with acceptable performance, at least for testing and
prototyping, provided special care is taken to force the scheduler to
accommodate this non-standard requirement.

By default, the GNURadio implementation we use (version 3.8.4.0 compiled from
the GitHub sources) sets a 100 ms delay time between the start of the RX
digitization and that of the TX generation. As a consequence, this imposes a
maximal locking bandwidth of the order of 2.5 Hz, no matter how optimized the
flowgraph may be. It is therefore necessary to modify this default behavior,
which can be achieved easily with the following piece of code run at the time
of initialization:

⬇

now = self.USRP_RX.get_time_now()

# prepare RXTX sync 1s from now

t0 = now + uhd.time_spec(1)

dt = uhd.time_spec(delay_time)

self.USRP_RX.set_start_time(t0)

self.USRP_TX.set_start_time(t0+dt)

where, as before, the `USRP_RX` and `USRP_TX` objects are instantiated in the
flowgraph by the “UHD: USRP Source” and “UHD: USRP Sink” blocks, respectively.
The `delay_time` parameter should be set to the smallest value possible before
triggering failures of the scheduler to meet the timing requirement (such
failure is readily observed by repeated occurrence of late packet events).
Note that it is possible to set the `delay_time` parameter to zero, by
providing an extra delay block of fixed delay value in the flowgraph. This
somewhat surprising behavior is related to specifics of the GNURadio framework
and delay block programming, and may not necessarily be expected to be true
for all versions.

Once this is done, it then becomes useful to optimize the working of the
flowgraph to allow decreasing the latency as much as available computing power
allows. We use a computer equipped with an Intel Xeon W-2275 multi-core
processor. We found it beneficial to de-activate hyperthreading and force the
operating clock frequency of all cores to the maximum availability (4.0 GHz in
our case, using the Linux cpu-freq-set package). Additionally, using the isol-
cpu Linux package, we isolate 12 out of the 14 cores, making them available
only to the GNURadio blocks involved in the closed loop part of the flowgraph.
We have also found it beneficial to assign each GNURadio block within the
servo-loop signal path to a specific processor core (whereas in normal
operation, the Linux scheduler would be free to assign and re-assign
dynamically the blocks’ threads to the various cores at run-time).

Utilizing this approach, with careful optimization of core use by each
GNURadio block, we are able to achieve a 100 Hz feedback bandwidth with our
complete multichannel flowgraph (10 optical frequency channels), consistent
with the minimal value of 2 ms of the `delay_time` without significant package
loss or late arrival. Future work in our experiment involves the
implementation of a custom purpose-built FPGA system, in place of the GNU
Radio/SDR system presented here, in order to reduce signal latency. We will
soon implement this change in our setup where we expect the latency reduction
to translate into a higher locking bandwidth in our spectral-hole-based ultra-
stable laser.

##  Appendix B Appendix: Managing Multiple Frequency Components in a Double
Pass Acousto-optic Modulator

We imprint the frequency modes on the laser beam using an AOM in a double-pass
configuration to produce an optical beam whose geometrical properties do not
depend on the RF spectral pattern which drives the AOM. However, an RF driving
signal with multiple frequency components results in an intermodulation
between the modes upon the second pass through the AOM. Fortunately, as
described in [41], a simple algorithmic solution exists to counteract these
effects. This can be seen mathematically if one writes, for a given
monochromatic input laser field ($A_{\rm in}\exp{[{\mathrm{i}}2\pi\nu_{0}t]}$)
of amplitude $A_{\rm in}$ and frequency $\nu_{0}$, entering the double-pass
AOM driven by an RF field ($A_{\rm RF}(t)\exp{[{\mathrm{i}}(2\pi
f_{0}t+\phi_{\rm RF}(t))]}$) at a carrier frequency $f_{0}$ and time-varying
amplitude $A_{\rm RF}(t)$ and phase $\phi_{\rm RF}(t)$, the resulting field,
$A_{\mathrm{out}}$, at the output of the first-order diffracted double pass
AOM is proportional to the square of the driving field, $A_{\mathrm{RF}}$:

$A_{\mathrm{out}}\propto A_{\rm
in}{\mathrm{e}}^{{\mathrm{i}}2\pi\nu_{0}t}\cdot A_{\rm
RF}^{2}(t){\mathrm{e}}^{{\mathrm{i}}[4\pi f_{0}t+2\phi_{\rm RF}(t)]}.$ (14)

Note that we have omitted the position dependence of the travelling wave here.

The production of unwanted cross-terms can therefore be avoided by simply
taking the complex-square root of each time sample before sending it to the
AOM RF driving field, a simple-enough process in an SDR platform. Note,
however, that this simple mathematical expression requires the AOM to operate
in its linear regime (i.e. driven at an RF power significantly lower than that
giving maximum diffraction efficiency in CW operation). This is an acceptable
constraint in our application where we require only a very low optical power
(of order $\mathrm{\mu W}$) for probing the spectral holes, three orders of
magnitude less than our available source laser. Note, also, that we need to
use a complex-number implementation of the square root function. For this, we
have chosen to decompose the complex stream into its amplitude and phase; in
parallel we divide the phase by 2 and take the real square root of the
amplitude, finally recombining the resulting values to obtain the complex
square root of the input sample.

With more complex spectral patterns being desired for experimentation, an
extra subtlety appears. The phase data in a digital system is typically
wrapped in the range of $-\pi$ to $\pi$. However, in our data treatment, after
passing through the complex square root function, the phase value will become
half of the original one, which will be in the range $-\pi/2$ to $\pi/2$. In
the complex plane, two successive samples
$A\exp{[{\mathrm{i}}(\pi-\epsilon)]}$ and
$A\exp{[{\mathrm{i}}(\pi+\epsilon)]}=A\exp{[{\mathrm{i}}(-\pi+\epsilon)]}$,
which are continuous in phase for very small $\epsilon$, will therefore
produce, after complex-square-rooting, two successive samples
$A\exp{[{\mathrm{i}}(\pi/2-\epsilon/2)]}$ and
$A\exp{[{\mathrm{i}}(-\pi/2+\epsilon/2)]}$. Although the squaring effect of
the double pass AOM is mathematically canceled out, the physical signal
produced will exhibit a sudden jump in phase at this point (except if the
amplitude happens to be exactly zero at that moment, which is for example the
case for simple spectra which are symmetric with reference to the carrier
frequency). Such a sudden jump is not perfectly realizable by the AOM which
effectively acts as a low pass filter, leading to distortion of the signal. To
overcome this problem, we wrote a “phase re-wrap” function to ensure that the
phase data wrapped in the $-\pi$ to $\pi$ interval is re-wrapped in the range
of $-2\pi$ to $2\pi$ before the calculation of the complex square root.

##  Appendix C Appendix: Cancellation at First Order of Parasitic Fabry-Perot
Cavity Effects for Appropriate Multi-mode Probing

In a spectral hole configuration described in Fig. 4(b) and section 4.4, the
particular use of narrow and broad spectral holes allows generating an error
signal which is largely immune to parasitic Fabry-Perot etalon effects:
effectively, the optical beatnote between the LO mode and the monitor mode
provides a measurement of the phase retardation due to the parasitic Fabry-
Perot effect, while the two symmetric probe modes, subject to a similar
parasitic Fabry-Perot etalon effect, are also sensitive to the phase change
due to the narrow spectral holes. An appropriate combination of the phases of
these different beatnotes provides a signal suitable to frequency lock a laser
onto the central frequency of the narrow holes with removal of the parasitic
Fabry-Perot etalon effect. The beatnote signal between the LO mode and any of
the other three modes exhibits a component related to the dispersion of the
spectral hole pattern in addition to a component related to parasitic Fabry-
Perot effect.

Considering, for simplicity but without loss of generality, a single Fabry-
Perot parasitic cavity of length $L$ in vacuum, it comes, from the well known
complex transmission coefficient $T/(1-Re^{i.4{\pi}{\nu}{L}/c})$ (where $\nu$
is the optical frequency, $c$ the speed of light and $T$ and $R$ are the
combined transmissivity and reflectivity of the mirrors - assumed to be
independent of $\nu$):

$\phi_{i}-\phi_{\rm LO}=\arctan\left[\frac{R\sin(4\pi
L\nu_{i}/c)}{1+R\cos(4\pi L\nu_{i}/c)}\right]-\arctan\left[\frac{R\sin(4\pi
L\nu_{\rm LO}/c)}{1+R\cos(4\pi L\nu_{\rm LO}/c)}\right],$ (15)

where $i=\\{{\rm p_{1}},{\rm m_{1}},{\rm p_{2}}\\}$ for the case of two probe
modes and one monitor mode (ie. the case used in Fig. 7, the measurement
parasitic Fabry-Perot induced noise cancellation). For a parasitic effect, we
can consider $R\ll 1$, leading to

$\phi_{i}\simeq R\sin(4\pi L\nu_{i}/c).$ (16)

Therefore, combining the measured phase of the beatnotes between the LO-mode
and the $i$-modes which are spectrally separated by $\Delta\nu$, with relative
weighing coefficients (1/2, -1, 1/2) which provides an error signal suitable
to lock the probing laser on the central frequencies of the narrow spectral
holes, we get:

$\displaystyle\epsilon_{\rm FP}$ $\displaystyle=$
$\displaystyle\frac{1}{2}(\phi_{\rm p_{1}}-\phi_{\rm LO})-(\phi_{\rm
m_{1}}-\phi_{\rm LO})+\frac{1}{2}(\phi_{\rm p_{2}}-\phi_{\rm LO})$ (17)
$\displaystyle\simeq$ $\displaystyle R\sin(4\pi L\nu_{\rm m}/c)[\cos(4\pi
L\Delta\nu/c)-1].$ (18)

Because $\Delta\nu\ll c/L$ for practical experimental implementations (in a
table top experiment the largest $L$ may be of order $1\,\mathrm{m}$ whereas
in our SDR-based system, $\Delta\nu$ is typically of the order of a few 100
kHz), this strongly reduces the effect of the parasitic Fabry-Perot cavity
compared to single channel probing (in fact, it cancels out at first order in
$\Delta\nu$). Note that using a more complex (non-linear) combination of the
phases of probe and monitor modes (with an increased number of monitor modes),
it is also possible to cancel it at higher orders as well (effectively using
the phase of the monitor modes to fit a polynomial estimation of the parasitic
Fabry Perot effect in the vicinity of the probe modes and subtracting this
effect from the measured probe-mode phases), should the necessity arise in a
particular highly-demanding application.

Funding The project has received financial support from Ville de Paris
Emergence Program, the Région Ile de France DIM C’nano and SIRTEQ, the LABEX
Cluster of Excellence FIRST-TF (ANR-10-LABX-48-01) within the Program
“Investissement d’Avenir” operated by the French National Research Agency
(ANR), and the 15SIB03 OC18 and 20FUN08 NEXTLASERS projects from the EMPIR
program cofinanced by the Participating States and from the European Union’s
Horizon 2020 research and innovation program, and the UltraStabLaserViaSHB
(GAP-101068547) from Marie Skłodowska-Curie Actions (HORIZON-TMA-MSCA-PF-EF)
from the European Commission Horizon Europe Framework Programme (HORIZON).

Acknowledgments We thank D. Nicolodi for the use of his custom GNURadio
blocks, on which we based our sweep-frequency block.

Disclosures The authors declare no conflicts of interest.

Data availability Data underlying the results presented in this paper are
available in [42]

## References

  * [1] T. Okoshi, “Recent advances in coherent optical fiber communication systems,” Journal of Lightwave Technology 5, 44–52 (1987).
  * [2] G. G. Goyer and R. Watson, “The laser and its application to meteorology,” Bulletin of the American Meteorological Society 44, 564 – 570 (1963).
  * [3] J. H. McElroy, “Infrared heterodyne solar radiometry,” Appl. Opt. 11, 1619–1622 (1972).
  * [4] M. M. Abbas, M. J. Mumma, T. Kostiuk, and D. Buhl, “Sensitivity limits of an infrared heterodyne spectrometer for astrophysical applications,” Appl. Opt. 15, 427–436 (1976).
  * [5] J. Eichholz, D. B. Tanner, and G. Mueller, “Heterodyne laser frequency stabilization for long baseline optical interferometry in space-based gravitational wave detectors,” Phys. Rev. D 92, 022004 (2015).
  * [6] O. Gobron, K. Jung, N. Galland, K. Predehl, R. L. Targat, A. Ferrier, P. Goldner, S. Seidelin, and Y. Le Coq, “Dispersive heterodyne probing method for laser frequency stabilization based on spectral hole burning in rare-earth doped crystals,” Opt. Express 25, 15539–15548 (2017).
  * [7] N. Galland, N. Lučić, S. Zhang, H. Alvarez-Martinez, R. L. Targat, A. Ferrier, P. Goldner, B. Fang, S. Seidelin, and Y. Le Coq, “Double-heterodyne probing for an ultra-stable laser based on spectral hole burning in a rare-earth-doped crystal,” Opt. Lett. 45, 1930–1933 (2020).
  * [8] B. P. Abbott and The LIGO Scientific Collaboration, “LIGO: the laser interferometer gravitational-wave observatory,” Reports on Progress in Physics 72, 076901 (2009).
  * [9] M. Diaz Ortiz, J. Gleason, H. Grote, A. Hallal, M. Hartman, H. Hollis, K.-S. Isleif, A. James, K. Karan, T. Kozlowski, A. Lindner, G. Messineo, G. Mueller, J. Põld, R. Smith, A. Spector, D. Tanner, L.-W. Wei, and B. Willke, “Design of the ALPS II optical system,” Physics of the Dark Universe 35, 100968 (2022).
  * [10] A. Hallal, G. Messineo, M. D. Ortiz, J. Gleason, H. Hollis, D. Tanner, G. Mueller, and A. Spector, “The heterodyne sensing system for the ALPS II search for sub-ev weakly interacting particles,” Physics of the Dark Universe 35, 100914 (2022).
  * [11] N. Kuhar, S. Sil, T. Verma, and S. Umapathy, “Challenges in application of raman spectroscopy to biology and materials,” RSC Adv. 8, 25888–25908 (2018).
  * [12] E. I. Nagaev, I. V. Baimler, A. S. Baryshev, M. E. Astashev, and S. V. Gudkov, “Effect of laser-induced optical breakdown on the structure of bsa molecules in aqueous solutions: An optical study,” Molecules 27 (2022).
  * [13] M. Hercher, “An analysis of saturable absorbers,” Appl. Opt. 6, 947–954 (1967).
  * [14] E. S. Maniloff, F. R. Graf, H. Gygax, S. B. Altner, S. Bernet, A. Renn, and U. P. Wild, “Power broadening of the spectral hole width in an optically thick sample,” Chemical Physics 193, 173–180 (1995).
  * [15] P. Berger, Y. Attal, M. Schwarz, S. Molin, A. Louchet-Chauvet, T. Chanelière, J.-L. Le Gouët, D. Dolfi, and L. Morvan, “RF spectrum analyzer for pulsed signals: Ultra-wide instantaneous bandwidth, high sensitivity, and high time-resolution,” J. Lightwave Technol. 34, 4658–4663 (2016).
  * [16] C. Venet, M. Bocoum, J.-B. Laudereau, T. Chaneliere, F. Ramaz, and A. Louchet-Chauvet, “Ultrasound-modulated optical tomography in scattering media: flux filtering based on persistent spectral hole burning in the optical diagnosis window,” Opt. Lett. 43, 3993–3996 (2018).
  * [17] M. Nilsson and S. Kröll, “Solid state quantum memory using complete absorption and re-emission of photons by tailored and externally controlled inhomogeneous absorption profiles,” Optics Communications 247, 393–403 (2005).
  * [18] F. Bussières, C. Clausen, A. Tiranov, B. Korzh, V. B. Verma, S. W. Nam, F. Marsili, A. Ferrier, P. Goldner, H. Herrmann, C. Silberhorn, W. Sohler, M. Afzelius, and N. Gisin, “Quantum teleportation from a telecom-wavelength photon to a solid-state quantum memory,” Nature Photonics 8, 775–778 (2014).
  * [19] A. Walther, L. Rippe, Y. Yan, J. Karlsson, D. Serrano, A. N. Nilsson, S. Bengtsson, and S. Kröll, “High-fidelity readout scheme for rare-earth solid-state quantum computing,” Phys. Rev. A 92, 022319 (2015).
  * [20] N. Maring, P. Farrera, K. Kutluer, M. Mazzera, G. Heinze, and H. de Riedmatten, “Photonic quantum state transfer between a cold atomic gas and a crystal,” Nature 551, 485–488 (2017).
  * [21] K. Mølmer, Y. Le Coq, and S. Seidelin, “Dispersive coupling between light and a rare-earth-ion–doped mechanical resonator,” Phys. Rev. A 94, 053804 (2016).
  * [22] S. Seidelin, Y. Le Coq, and K. Mølmer, “Rapid cooling of a strain-coupled oscillator by an optical phase-shift measurement,” Phys. Rev. A 100, 013828 (2019).
  * [23] B. Julsgaard, A. Walther, S. Kröll, and L. Rippe, “Understanding laser stabilization using spectral hole burning,” Opt. Express 15, 11444–11465 (2007).
  * [24] Q.-F. Chen, A. Troshyn, I. Ernsting, S. Kayser, S. Vasilyev, A. Nevsky, and S. Schiller, “Spectrally narrow, long-term stable optical frequency reference based on a ${\mathrm{Eu}}^{3+}:{\mathrm{Y}}_{2}{\mathrm{SiO}}_{5}$ crystal at cryogenic temperature,” Phys. Rev. Lett. 107, 223202 (2011).
  * [25] M. J. Thorpe, L. Rippe, T. M. Fortier, M. S. Kirchner, and T. Rosenband, “Frequency stabilization to $6\times 10^{-16}$ via spectral-hole burning,” Nature Photonics 5, 688–693 (2011).
  * [26] M. J. Thorpe, D. R. Leibrandt, and T. Rosenband, “Shifts of optical frequency references based on spectral-hole burning in ${\mathrm{Eu}}^{3+}:{\mathrm{Y}}_{2}{\mathrm{SiO}}_{5}$,” New Journal of Physics 15, 033006 (2013).
  * [27] N. Galland, N. Lučić, B. Fang, S. Zhang, R. Le Targat, A. Ferrier, P. Goldner, S. Seidelin, and Y. Le Coq, “Mechanical tunability of an ultranarrow spectral feature of a rare-earth-doped crystal via uniaxial stress,” Phys. Rev. Applied 13, 044022 (2020).
  * [28] S. Zhang, N. Galland, N. Lučić, R. Le Targat, A. Ferrier, P. Goldner, B. Fang, Y. Le Coq, and S. Seidelin, “Inhomogeneous response of an ion ensemble from mechanical stress,” Phys. Rev. Research 2, 013306 (2020).
  * [29] S. Zhang, N. Lučić, N. Galland, R. Le Targat, P. Goldner, B. Fang, S. Seidelin, and Y. Le Coq, “Precision measurements of electric-field-induced frequency displacements of an ultranarrow optical transition in ions in a solid,” Applied Physics Letters 117, 221102 (2020).
  * [30] S. Zhang, S. Seidelin, R. Le Targat, P. Goldner, B. Fang, and Y. Le Coq, “First-order thermal insensitivity of the frequency of a narrow spectral hole in a crystal,” Phys. Rev. A 107, 013518 (2023).
  * [31] C. Philippe, R. Le Targat, D. Holleville, M. Lours, T. Minh-Pham, J. Hrabina, F. Du Burck, P. Wolf, and O. Acef, “Frequency tripled $1.5\,\mathrm{{\mu}m}$ telecom laser diode stabilized to iodine hyperfine line in the $10^{-15}$ range,” in _2016 European Frequency and Time Forum (EFTF),_ (2016), pp. 1–3.
  * [32] T. L. Nicholson, M. J. Martin, J. R. Williams, B. J. Bloom, M. Bishof, M. D. Swallows, S. L. Campbell, and J. Ye, “Comparison of two independent sr optical clocks with $1\mathbf{\times{}}{10}^{-17}$ stability at ${10}^{3}\text{ }\text{ }\mathbf{s}$,” Phys. Rev. Lett. 109, 230801 (2012).
  * [33] W. M. Itano, J. C. Bergquist, J. J. Bollinger, J. M. Gilligan, D. J. Heinzen, F. L. Moore, M. G. Raizen, and D. J. Wineland, “Quantum projection noise: Population fluctuations in two-level systems,” Phys. Rev. A 47, 3554–3570 (1993).
  * [34] I. Afek, O. Ambar, and Y. Silberberg, “High-NOON states by mixing quantum and classical light,” Science 328, 879–881 (2010).
  * [35] R. W. P. Drever, J. L. Hall, F. V. Kowalski, J. Hough, G. M. Ford, A. J. Munley, and H. Ward, “Laser phase and frequency stabilization using an optical resonator,” Applied Physics B 31, 97–105 (1983).
  * [36] E. D. Black, “An introduction to Pound–Drever–Hall laser frequency stabilization,” American Journal of Physics 69, 79–87 (2001).
  * [37] N. A. Riza, “Scanning heterodyne optical interferometers,” Review of Scientific Instruments 67, 2466–2476 (1996).
  * [38] “Polyphase Channelizer,” https://wiki.gnuradio.org/index.php/Polyphase_Channelizer. [Online; published (at time of access) 31 December 2021].
  * [39] “GNU Radio Manual and C++ API Reference: Filter Signal Processing Blocks,” https://www.gnuradio.org/doc/doxygen/page_filter.html. [Online reference manual version 3.9.4.0, as accessed on 07 June 2023].
  * [40] R. McIntyre, “Multiplication noise in uniform avalanche diodes,” IEEE Transactions on Electron Devices ED-13, 164–168 (1966).
  * [41] P. Jobez, N. Timoney, C. Laplane, J. Etesse, A. Ferrier, P. Goldner, N. Gisin, and M. Afzelius, “Towards highly multimode optical quantum memory for quantum repeaters,” Phys. Rev. A 93, 032327 (2016).
  * [42] X. Lin, M. Hartman, S. Zhang, S. Seidelin, B. Fang, and Y. Le Coq, “Dataset for: Multi-mode heterodyne laser interferometry realized via software defined radio,” http://dx.doi.org/10.5281/zenodo.8390942 (2023).