11institutetext: Dipartimento di Biotecnologie Mediche e Medicina
Traslazionale, Università degli Studi di Milano, via F.lli Cervi 93, 20090
Segrate, Italy 22institutetext: Division of Physical Chemistry, Department of
Chemistry, Lund University, Lund, Sweden 33institutetext: University of
Vienna, Faculty of Physics, Boltzmanngasse 5, 1090 Vienna, Austria

# Probing roto-translational diffusion of small anisotropic colloidal
particles with a bright-field microscope

Fabio Giavazzi Corresponding author<EMAIL_ADDRESS>Antara Pal 22
Roberto Cerbino 1133

(Received: date / Revised version: date)

###### Abstract

Soft and biological materials are often composed of elementary constituents
exhibiting an incessant roto-translational motion at the microscopic scale.
Tracking this motion with a bright-field microscope becomes increasingly
challenging when the particle size becomes smaller than the microscope
resolution, a case which is frequently encountered. Here we demonstrate
Squared-Gradient Differential Dynamic Microscopy (SG-DDM) as a tool to
successfully use bright-field microscopy to extract the roto-translational
dynamics of small anisotropic colloidal particles, whose rotational motion
cannot be tracked accurately in direct space. We provide analytical
justification and experimental demonstration of the method by successful
application to an aqueous suspension of peanut-shaped particles.

##  1 Introduction

Over a century ago, Jean Perrin conducted a series of experiments to mark the
birth of quantitative microscopy perrin2013brownian ; cerbino2018quantitative
. Among his many intriguing results, the simultaneous measurement of the
translational and rotational diffusion coefficients of micron-sized particles
provided the first convincing experimental demonstration of the energy
equipartition between the translational and rotational degrees of freedom, one
of the most important predictions of the kinetic theory. Since then, the
rotational and the roto-translational motion of colloidal objects has been
repeatedly used as the tool of choice to investigate fundamental physical
processes like mass transport or phase transitions han2006brownian ;
edmond2012decoupling , as a probe of material properties at the micro-scale
PhysRevLett.90.018304 ; colin2014rotational , as a biosensor chen2020trace ,
to monitor the motility pattern of active swimmers and their interactions with
passive particles saragosti2012modeling ; peng2016diffusion , to characterize
complex transport dynamics in cells gao2017seeing , just to mention a few
recent examples.

Nowadays, optical microscopy remains a key tool for probing roto-translational
particle motion, especially in combination with specific (for example
fluorescent) markers. For low particle density, single particle tracking (SPT)
represents a rather straightforward way to extract a quantitative information
on the translational dynamics of well-resolved particles by following their
center of mass; by contrast, the accurate determination of their three-
dimensional rotational dynamics with SPT is usually much more challenging and
time consuming anthony2015tracking . Moreover, most of the SPT algorithms are
developed and optimized for spherical particles.

The study of the roto-translational dynamic of small anisotropic particles
and/or of particles at high densities is thus the realm of ensemble-averaging
techniques, such as depolarized dynamic light scattering
(DDLS)berne2000dynamic ; PhysRevE.52.2707 , polarized fluorescence recovery
after photobleaching (pFRAP) velez1988polarized ; lettinga2004rotational , and
nuclear magnetic resonance (NMR). However, the main limitation or these
approaches is that they cannot be used in complex, highly heterogeneous
environments, which is a common scenario in biology and material science, and
one for which microscopy is particularly utilized.

An intermediate approach is provided by a family of methods known as Digital
Fourier microscopy (DFM) giavazzi2014digital : as in microscopy, temporal
image sequences of the sample are acquired in direct space but the sample
intermediate scattering function $f(q,\tau)$ is extracted via image
correlation in the reciprocal (Fourier) space, as a function of the scattering
wave-vector $q$ and of the time delay $\tau$. One of the most known
implementation of DFM is Differential Dynamic Microscopy (DDM), a technique
that has found vast application with a variety of soft and biological matter
systems cerbino2008differential ; cerbino2017perspective . DDM has been
successfully used for studying the translational dynamics of shape- or
functionally-anisotropic particles reufer2012differential ;
PhysRevLett.121.078001 ; nixon2019differential ; pal2020anisotropic . In
combination with a suitable imaging mode, like polarized
giavazzi2016simultaneous ; edera2017differential  or dark-field
cerbino2017dark  microscopy, DDM enables the accurate measurement of the roto-
translational dynamics of optically and/or shape-anisotropic particles.

While measuring the translational motion in DDM is based on the simple fact
that a translating particle occupies different portions of an image, the key
for DDM to be able to capture the particle rotational motion is to convert the
latter into an intensity fluctuation of the particle image. For example, the
image of optically anisotropic, uniaxial particle observed between cross-
polarizers appears darker or brighter, according to the angle formed between
its optical axis and the axes of the polarizers. In this condition, rotational
diffusion, by promoting the random reorientation of the particle, leads to an
intermittent ”blinking” of the intensity associated to the particle image. If
DDM analysis is now performed on an image sequence collected in these
conditions, a $q$-independent decay is observed in $f(q,\tau)$, which is
particularly evident in the low-$q$ regime. With proper fitting of $f(q,\tau)$
the characteristic correlation time of this decay can be measured and used to
determine the associated rotational diffusion coefficient
giavazzi2016simultaneous . A similar effect, although realized via a different
optical mechanism, is exploited with dark-field microscopy to measure the
characteristic reorientation time of micrometer-scale, shape-anisotropic
particles cerbino2017dark .

Using the same strategy for small anisotropic objects imaged in bright-field
or in non-polarized fluorescence microscopy remains much more challenging: an
in-plane rotation of the object simply leads to a rotated image, while out-of-
plane rotations can introduce subtle shape changes, but with no effect on the
overall associated intensity. In reciprocal space, this corresponds to the
fact that the low-$q$ portion of the image power spectrum is not affected by
rotations. In theory, rotational diffusion contributions become important in
scattering at large wave-vectors $q$, where the object rotation can induce a
fluctuation in the scattered intensity berne2000dynamic ; pecora1968spectral .
However, the relaxation of the intermediate scattering function for such large
$q$ is usually dominated by the translational term whose relaxation rate, in
the case of a Brownian particle, rapidly grows as $D_{T}q^{2}$, where $D_{T}$
is the translational diffusion coefficient. To the best of our knowledge, the
only study reporting a DDM-based measurement of the rotational dynamics of
colloidal particles in bright-field is Ref. wittmeier2015rotational , in which
the particles translational motion was suppressed by partially tethering them
to a solid surface. A general recipe for using bright field microscopy to
simultaneously characterize the translational and rotational diffusive
dynamics of Brownian particles remains thus unavailable.

In this work, we address this issue by combining the standard DDM analysis of
bright-field movies with a simple image preprocessing step, whose function is
to digitally reproduce the intensity change induced by particle rotation that
is obtained by optical means in polarized or dark-field microscopy: for each
image, we generate a local-orientation map whose intensity level encodes the
local particle orientation at those image pixels where a particle is present;
we then perform a standard DDM analysis on these local orientation maps and
obtain the accurate characterization of both translational and rotational
diffusion. A detailed description of the implementation of the method is
reported in Sec. 2.3, while a simple analytic justification is discussed in
Appendix A.

##  2 Materials and methods

###  2.1 Sample preparation and imaging

The sample is a highly diluted aqueous suspension of hematite-silica core-
shell peanut-shaped particles, synthesized according to the protocol described
in detailed Refs. sugimoto1993preparation ; pal2018anomalous . The length and
diameter of the lobes of these particles are $(1723\pm 50)$ nm and $(740\pm
50)$ nm, respectively. Before measurement, the sample is loaded in a glass
capillary of thickness $100$ $\mu m$, which is then sealed with vacuum grease.
After about one hour at room temperature ($(T=22\pm 1)$ oC), the sample
appears to be sedimented close to the bottom of the container as the
sedimentation length $l_{g}\sim 175$ $nm$ of these particles is smaller than
their diameter kamal2020path . The final number density is about $5\cdot
10^{3}$ $mm^{-2}$, as obtained from particle counting (see Fig. 2).

The sample is observed with an inverted Nikon Ti-E bright-field microscope
equipped with a Hamamatsu Orca Flash 4.0 camera (pixel size $6.5$ $\mu m$).
The condenser diaphragm (numerical aperture NA=0.52) is kept completely open
to achieve incoherent illumination of the sample. Two different objectives are
used, a 40X, 0.60 NA objective and a 10X 0.25 NA objective. The size of the
imaged region is $6.9\cdot 10^{3}$ $\mu m^{2}$ and $2.8\cdot 10^{4}$ $\mu
m^{2}$ for the higher and the lower magnification objective, respectively.
Image sequences of $10^{4}$ frames each were recorded at 25 frame/s.

###  2.2 Differential dynamic microscopy

The translational dynamics of the particles is characterized by analyzing the
sample movies with the standard DDM analysis scheme cerbino2008differential ;
giavazzi2009scattering . In short, we calculate the difference $\Delta
I(\mathbf{x},t,\tau)=I(\mathbf{x},t+\tau)-I(\mathbf{x},t)$ between two images
acquired at times $t$ and $t+\tau$. By averaging the spatial Fourier power
spectrum of $\Delta I(\mathbf{x},t,\tau)$ obtained for the same $\tau$ but
different reference times $t$ we obtain the so called image structure function

$D(\mathbf{q},\tau)=\left\langle\left|\Delta\hat{I}(\mathbf{x},t,\tau)\right|^{2}\right\rangle_{t}$
(1)

that captures the sample dynamics as a function of the two-dimensional
scattering wavevector $\mathbf{q}$ and of the lag time $\tau$. The symbol
$\hat{\cdot}$ indicates the two dimensional digital Fourier transform, usually
performed with a Fast Fourier Transform algorithm. In a last step, we take
advantage of the circular symmetry of the problem to also perform an azimuthal
average of $D(\mathbf{q},\tau)$, which leads to the one-dimensional function
$D(q,\tau)$ of the radial wavevector $q=\sqrt{q_{x}^{2}+q_{y}^{2}}$. The
obtained structure function is connected to the intermediate scattering
function (ISF) $f(q,\tau)$ berne2000dynamic  by the relation

$D(q,\tau)=A(q)\left[1-f(q,\tau)\right]+B(q)$ (2)

where the term $B(q)$ accounts for the camera noise and the amplitude $A(q)$
depends on the scattering properties of the sample and the transfer function
of the imaging system giavazzi2014digital .

###  2.3 Squared-gradient differential dynamic microscopy

Figure 1:  ”Squared-gradient” maps for a single anisotropic particle whose
long axis lays on the image plane. First row: simulated images of the same
particle with bivariate Gaussian intensity profile and different orientations.
Second row: ”gradient” maps. Each map correspond to the partial derivative
along the vertical direction of the image above it. Third row: ”squared-
gradient” maps. Each map corresponds to the square of the map above it. In
each row, maps are shown with the same grayscale.

In order to probe also the rotational dynamics of the particles we used a
digital pre-processing step of the bright-field images to reveal information
about the particles orientation. A simple implementation of this concept,
which we call ”squared-gradient” (SG), involves the calculation of the partial
derivative of the image with respect to a given component $\mu$, which is then
squared to give

$c_{\mu}(\mathbf{x},t)=\left[\partial_{\mu}I(\mathbf{x},t)\right]^{2}.$ (3)

While other implementations are also possible DIP , in this work, the partial
derivatives in Eq. 3 are calculated simply as the difference between the
original image and its copy translated along one of the main axes by a single
pixel. As it can be appreciated in Fig. 1, where a single elongated particle
is considered, the obtained map has the required property, as the overall
intensity changes according to the particle’s orientation. Indeed, when the
long axis of the particle is oriented along the $\mu$ axis (the vertical
direction in Fig. 1), the contrast in the gradient map is lower compared to
the case where the particle is oriented in the perpendicular direction. In the
square-gradient map $c_{\mu}(\mathbf{x},t)$ the particle is replaced by two
smaller lobes. The most evident effect of the particle rotation on
$c_{\mu}(\mathbf{x},t)$ is a modulation of the intensity associated with the
two lobes. As shown in Appendix A, a similar intensity modulation in
$c_{\mu}(\mathbf{x},t)$ is also produced by changes in the angle formed by the
long axis of the particle with optical axis.

Once a sequence of SG-maps are obtained from the original images, the standard
DDM analysis is performed to calculate the corresponding structure function

$D_{SG}(q,\tau)=\sum_{\mu=1,2}\langle|\hat{c}_{\mu}(\mathbf{q},t+\tau)-\hat{c}_{\mu}(\mathbf{q},t)|^{2}\rangle_{t,|\mathbf{q}|=q}.$
(4)

where we also average over the two orthogonal directions. As it is shown
explicitly in Appendix A for the special case $q=0$, $D_{SG}(q,\tau)$ can be
expressed in terms of the normalized ISF via Eq. 2 and we employ the general
expression

$f(q,\tau)=\alpha(q)f_{RT}(q,\tau)+[1-\alpha(q)]f_{T}(q,\tau),$ (5)

where $f_{T}$ is the translational ISF, $f_{RT}$ is the roto-translational one
and $\alpha(q)$ is the $q$-dependent relative amplitude of the roto-
translational term. In Eq.5, we neglect a term describing the coupling between
orientational and translational diffusion, which is expected to be of little
significance han2006brownian . For monodisperse Brownian particles, the ISFs
take the well know expressions $f_{T}(q,\tau)=e^{-D_{T}q^{2}\tau}$ and
$f_{RT}(q,\tau)=e^{-6D_{R}\tau-D_{T}q^{2}\tau}$, where $D_{T}$ and $D_{R}$ are
the average translational and rotational diffusion coefficients, respectively
berne2000dynamic . It is worth stressing that the taking the square in Eq. 3
represents an essential step of the procedure, as it enables ”translating” the
orientation-dependent spatial modulation introduced by the gradient operation
into a global intensity variation, which affects in a $q$-independent fashion
all the Fourier components of the SG map.

###  2.4 Single particle tracking

To validate SG-DDM, we analyze the same images also with single particle
tracking (SPT), by using the MATLAB code developed by the Kilfoil group at the
University of Massachusetts PhysRevLett.102.188303  and available at
https://github.com/dsseara/microrheology. We reconstruct individual
trajectories of the center of mass of the particles and we calculate their
mean-squared displacement $\langle\Delta r^{2}(\tau)\rangle$ as a function of
the delay time $\tau$.

To track over time the orientation of the particles, we use a custom MATLAB
code implementing a procedure similar to the one described in Refs.
PhysRevE.80.011403 ; zheng2010self , which is described in detail in Appendix
B. We applied the described procedure, which provides an independent estimate
of the rotational diffusion coefficient, only to the image sequence recorded
with the higher magnification (40X), as we were not able to successfully
perform the same orientational tracking-based analysis on the image sequences
acquired with a 10X objective.

##  3 Results and discussion

We show in Fig. 3 representative ISFs obtained from both standard DDM (panel
a) and SG-DDM (panel b) analysis of the same sequence of bright-field images
(magnification 40X) of the sample. The ISFs obtained from standard DDM
analysis exhibit a single exponential decay $e^{-\Gamma(q)\tau}$, whose
relaxation rate $\Gamma(q)$ displays a clean $\sim q^{2}$ scaling (Fig. 4,
gray triangles).

Figure 2:  a) Representative bright-field image of anisotropic peanut-shaped
particles diffusing close to the bottom of the glass capillary (see text for
details) taken with a 40X, 0.60 NA objective. The size of the image
corresponds to 39 $\mu$m in real space. b) Square-gradient map obtained from
the image shown in panel (a). Each particle appears to be replaced by two
bright spots whose intensity depends on the particle orientation.

A fit of $\Gamma(q)$ with a quadratic model $\Gamma(q)=D_{T}q^{2}$ provides an
estimate for the translational diffusion coefficient $D_{T}=0.24\pm 0.02$ $\mu
m^{2}/s$ of the particles, which is in excellent agreement with the value
$D_{T,PT}=0.23\pm 0.01$ $\mu m^{2}/s$ obtained with SPT (Fig. 5a) via a linear
fit $\langle\Delta r^{2}(\tau)\rangle=4D_{T}\tau$ of the MSD of the particles.

Figure 3:  a) Representative ISFs obtained from standard DDM analysis of an
image sequence of Brownian anisotropic particles for different $q$-values,
logarithmically spaced in the interval $[0.15,3.4]$ $\mu m^{-1}$. Continous
black curves are best fitting curves to the data with a simple exponential
model. b) Representative ISFs obtained from SG-DDM analysis of the same image
sequence and for the same $q$-values considered in panel (a). Continous black
curves are best fitting curves to the data with a double exponential model.

Inspection of the ISFs obtained from SG-DDM (Fig. 3b) reveals instead the
existence of a double decay, which is particularly evident for the lowest
values of $q$. Fitting a double exponential model to the data allows to
reliably estimate, at least for $q\lesssim 2.5$ $\mu m^{-1}$, two well
distinct decay rates: a fast decay with an almost $q$-indpendent rate
$\Gamma_{1}$, and a slow diffusive decay with rate $\Gamma_{2}\sim q^{2}$ that
is fully compatible with the relaxation rate obtained from standard DDM (Fig.
4). For $q>2.5$ $\mu m^{-1}$, the two relaxation rates are too close to be
reliably separated and the fitting procedure provides a single relaxation
rate.

According the discussion in Sec.2.3, we attribute the fast decay to the roto-
traslational diffusion of the particles, $\Gamma_{1}(q)=6D_{R}+D_{T}q^{2}$,
and we obtain the estimate $D_{R}=0.30\pm 0.05$ $s^{-1}$. This result agrees
with those obtained with SPT according to the procedures described in Section
2.4, in particular with the estimate $D_{R,PT_{1}}=0.28\pm 0.05$ $s^{-1}$
obtained from a linear fit of the angular MSD
$\langle\Delta\phi^{2}(\tau)\rangle=2D_{R,PT_{1}}\tau$ (Fig. 5b), and with the
estimate $D_{R,PT_{2}}=0.29\pm 0.05$ $s^{-1}$ obtained from an exponential fit
of the effective aspect ratio autocorrelation function
$C_{\epsilon}(\tau)=e^{-6D_{R,PT_{2}}\tau}$ (Fig. 5c).

We perform the same analyses on an image sequence of the same sample, acquired
with a 10X objective (NA=0.25). In this imaging condition, the particle size
is very close to the diffraction limit ($d\sim\overline{\lambda}/(2NA)\simeq
1.2$ $\mu m$, assuming $\overline{\lambda}\simeq 0.6$ $\mu m$) and the effect
of the convolution with the point-spread-function is to markedly reduce the
apparent anisotropy in the particle image, that is now barely appreciable
(inset of Fig. 4). As a consequence, the uncertainty in the frame-by-frame
determination of single particle orientation is so large that the angular
tracking procedure described in Section 2.4 no longer provides reliable
results. On the contrary, the results of both DDM and SG-DDM are in excellent
agreement with the one obtained with the 40X objective (NA=0.60), as shown in
Fig. 4.

As a further validity check for the proposed approach, we compare the obtained
diffusion coefficients with the ones predicted by the theory for a rigid rod,
which represents a fair model for our peanut-shape particles. If we consider a
cylinder of length $1723$ nm and diameter $740$ nm, in a Newtonian fluid of
viscosity $\eta=1.0$ mPa at temperature $T=22$ oC, by using the theoretical
expressions reported in Ref. PhysRevE.50.1232  (Eqs. 1-4) we obtain the
following values . $D_{R,th}=0.44$ $s^{-1}$, $D_{T,th}=0.25\cdot$ $\mu
m^{2}/s$ for the rotational and the translational diffusion coefficient,
respectively. We note that, while $D_{T,th}$ is in very good agreement with
the corresponding values obtained with both SG-DDM and SPT, the calculated
diffusion coefficient $D_{R,th}$ is off by about 30% with respect the
experimental ones. This moderate discrepancy could be explained as a
systematic effect of the adopted cylindrical approximation. However, it can be
also a genuine effect due the hydrodynamic interactions between particles and
the bottom plate of the container, which can slow down the particles’
diffusivity compared to the bulk Yang .

Figure 4:  Roto-translational dynamics of anisotropic Brownian particles
obtained from DDM and SG-DDM. Panels (a-b) and (c-d) show results obtained
from the analysis of image sequences of the same sample, acquired with a 40X,
0.60 NA objective and a 10X, 0.25 NA objective, respectively. a) Gray
triangles: relaxation rate $\Gamma(q)$ obtained from the single exponential
fit of the structure functions obtained from standard DDM analysis. Orange
circles and blue squares: relaxation rates $\Gamma_{1}(q)$ and $\Gamma_{2}(q)$
obtained from a double exponential fit of the structure functions obtained
from SG-DDM analysis on the same image sequence. The continous line
$\Gamma=D_{T}q^{2}$ is a best fitting curve to the blue squares. The analogous
curve obtained from the gray triangle is indistinguishable. The dashed curve
is a best fitting curve $\Gamma_{1}=6D_{R}+D_{T}q^{2}$ to the orange circles.
The vertical dashed line marks the $q$-value above which the two relaxation
rates $\Gamma_{1}$ and $\Gamma_{2}$ cannot be resolved any more. The small
inset shows a representative ROI of size 10.4 $\mu m$. b) Orange circles:
relative amplitude $\alpha(q)$ of the roto-translational contribution to the
ISF (see Eq.5), obtained from a double exponential fit of the SG-DDM structure
functions. The description of panels c) and d) is identical to the one of
panels a) and b), respectively

Figure 5:  Roto-translational dynamics of anisotropic Brownian particles
obtained from single particle tracking. a) Symbols: MSD of the center of mass
of the particles. The continous line a best fitting curve to the data with a
linear model $\langle\Delta r^{2}(\tau)\rangle=4D_{T}\tau$, with
$D_{T}=(0.23\pm 0.01)$ $\mu m^{2}s^{-1}$ b) Symbols: polar angular MSD of the
particles. The continuous line a best fitting curve to the data with a linear
model $\langle\Delta\phi^{2}(\tau)\rangle=2D_{R}\tau$, with $D_{R}=(0.28\pm
0.05)$ $s^{-1}$ . c) Symbols: normalized temporal correlation function of the
apparent aspect ratio $\epsilon(t)=b^{2}(t)/a^{2}(t)$ of each particle, as
obtained from a bivariate Gaussian fit of the associated intensity
distribution (see main text for details). The continuous line a best fitting
curve to the data with an exponential decay
$C_{\epsilon}(\tau)=\exp(-\gamma\tau)$, with $\gamma=(1.75\pm 0.2)$ $s^{-1}$

##  4 Conclusions

In this work, we introduced SG-DDM, a simple procedure enabling the tracking-
free determination of the rotational and translational dynamics of anisotropic
particles in microscopy images. SG-DDM was here demonstrated for bright-field
microscopy experiments, but we trust that it could be applied without
modifications to a broad range of imaging modes, including wide-field or
confocal fluorescence microscopy.

Compared to SPT-based approaches, we highlight a number of significant
advantages. First of all, at a variance with most SPT methods, the proposed
approach does not require user-dependent parameters and does not rely on a
specific model for the particle shape. Moreover, by exploiting the intrinsic
ensemble-averaging capability of the method, we proved that SG-DDM is able to
provide a quantitative information on the rotational dynamics even when the
size of the particles is very close to the diffraction limit and a low signal-
to-noise ratio does not allow a reliable frame-by-frame determination of
single-particle orientation in real space. We thus expected the method to be
compatible with a wide range of shape anisotropic particles, such as for
example metallic nanowires farain2018universal , as well as with spherical,
but optically anisotropic, objects, like Janus particles
wittmeier2015rotational . Computationally, SG-DDM is as effective as DDM,
whose use becomes advantageous in particular when a large number of particles
are considered. In fact, the computation time does not depend on the number of
imaged particles, so that thousands of particles can be analyzed in parallel
in a single movie.

Being based on a non-linear pre-processing step, a condition under which the
method is expected to provide the correct translational dynamics of the sample
is the absence of a significant overlap between images of the different
particles giavazzi2009scattering ; cerbino2017dark . While this condition is
almost automatically fulfilled for rigid objects in two dimensions, it could
represent a more serious limitation on the particle concentration in the three
dimensional case.

In addition, if the particles deviate significantly from a simple uniaxial
symmetry, such as the one investigated here, higher other contributions
(corresponding to high order terms in the spherical harmonics expansion of the
particles’ shape) become relevant in the ISF. In those cases, a more refined
system-dependent model for the ISF, including multiple exponential decays,
should be used to properly interpret the SC-DDM signal.

Finally, also the presence of a strongly inhomogeneous optical background can
affect the quality of the results, as the non-linear processing can introduce
a coupling between moving objects of interest and static features. A
mitigation of this problem could be provided by a background subtraction
before calculating the SG-map. For example, if the observation time window is
long enough to enable a complete decorrelation in the particles positions, a
good estimate of the background intensity could be obtained as the time
average of the whole sequence.

In all the above mentioned cases, a quite compelling internal quality test is
provided by the comparison between the results obtained from standard DDM and
SG-DDM on the same sequence, as the exact same translational dynamics is
expected to be measured in DDM and SG-DDM. Any lack of correspondence between
the translational ISFs obtained with the two methods is a proxy of a potential
artifact introduced by the SG mapping, offering a simple and sensitive check
of the validity of the mapping procedure.

Finally, beyond offering a novel way to determine the rotational dynamics of
small anisotropic entities, SG-DDM represents also the first example of non-
linear pre-processing of microscope images prior to DDM analysis, which we
believe will expand the range of applications of DDM and other DFM techniques
beyond the current frontier.

## Appendix A. Analytical calculation of $D_{SG}(q=0,\tau)$

In this Appendix we provide analytical justification of SG-DDM as a probe of
the rotational diffusion of anisotropic particles. In particular, we show
that, in the $q\rightarrow 0$ limit, the structure function $D_{SG}(q,\tau)$
is dominated by a $q$-independent relaxation mirroring the rotational
dynamics.

We assume a simple model where each ”particle” described by a 3D Gaussian
profile with variances $(a^{2},a^{2},b^{2})$ along its three main axes
($a<b$). Let $(x_{0},y_{0},z_{0})$ be the coordinate of its center of mass,
while the orientation of the long axis is determined by the polar angle
$\theta$ (i.e. the angle formed with the $z$-axis, corresponding to the
optical axis) and the azimuthal angle $\phi$.

In what follows, we consider for simplicity the case of a single particle, the
generalization to a collection of identical, non overlapping, particles being
straightforward. We assume that the intensity distribution $I(x,y)$ on the
image plane is given by the projection on the $x-y$ plane of the object
density distribution (good for example in the case of incoherent bright-field
or wide-field fluorescence microscopy of a particle close to the object plane
giavazzi2009scattering ). Under this hypothesis, $I(x,y)$ is a 2D Gaussian
distribution $I(x,y)=1/(\sigma_{x}\sigma_{y})\exp\left(-{x^{\prime
2}}/{2\sigma_{x}^{2}}-{y^{\prime 2}}/{2\sigma_{y}^{2}}\right)$, where
$x^{\prime}=(x-x_{0})\cos\phi+(y-y_{0})\sin\phi$,
$y^{\prime}=(y-y_{0})\cos\phi-(x-x_{0})\sin\phi$, $\sigma_{x}=a$ and
$\sigma_{y}^{2}=a^{2}\cos^{2}\theta+b^{2}\sin^{2}\theta$.

The total intensity $i_{0}=\int
d^{2}\mathbf{x}c_{\mu}(\mathbf{x},t)=\hat{c}_{\mu}(\mathbf{q}=0,t)$ associated
with the gradient map $c_{\mu}(\mathbf{x},t)$ can be explicitly calculated. We
find
$i_{0}(\theta,\phi)=c\left(\frac{\sigma_{x}^{2}+\sigma_{y}^{2}}{\sigma_{x}^{3}\sigma_{x}^{3}}\right)\left[1+\left(\frac{\sigma_{y}^{2}-\sigma_{x}^{2}}{\sigma_{x}^{2}+\sigma_{y}^{2}}\right)\cos
2\phi\right]$, where $c$ is a dimensionless constant.

We initially consider the case of a moderately anisotropic particle (i.e. one
whose eccentricity $e\equiv\sqrt{1-a^{2}/b^{2}}$ is small), for which we
obtain the following expression, valid to the first order in $e$:

$i_{0}(\theta,\phi)=c^{\prime}\left[1+e^{2}\sin^{2}\theta\left(1+(1/2)\cos
2\phi\right)\right],$ (6)

where $c^{\prime}=(2c/a^{4})$.

It is worth expressing $i_{0}$ in terms of the spherical harmonics
$Y_{l}^{m}(\theta,\phi)$

$i_{0}(\theta,\phi)=a_{00}Y_{0}^{0}+a_{20}Y_{2}^{0}+a_{22}[Y_{2}^{-2}+Y_{2}^{2}],$
(7)

where $a_{00}=2\sqrt{2}(1+2e^{2}/3)$, $a_{20}=(4/3)e^{2}\sqrt{\pi/5}$ and
$a_{22}=e^{2}\sqrt{2\pi/15}$.

Combining Eq. 7 with the identity berne2000dynamic

$\langle
Y_{l}^{m}Y_{l^{\prime}}^{m^{\prime}}\rangle=\frac{1}{4\pi}e^{-l(l+1)D_{R}\tau}\delta_{l,l^{\prime}}\delta_{m,m^{\prime}},$
(8)

we obtain the following expression for $D_{SG}(q=0,\tau)=2\left(\langle
i^{2}_{0}(t)\rangle-\langle i_{0}(t)i_{0}(t+\tau)\rangle\right)$

$D(q=0,\tau)=A_{0}(1-e^{-6D_{r}\tau}).$ (9)

This simple expression is valid under the hypothesis of small shape
anisotropy. Within the same approximation, the effective aspect ratio
$\epsilon=\sigma_{y}/\sigma_{x}$ is given by
$\epsilon=1+\frac{1}{2}e^{2}\sin^{2}\theta$ and the corresponding time
autocorrelation function
$C_{\epsilon}(\tau)=\left[\langle\epsilon(t+\tau)\epsilon(t)\rangle-\langle\epsilon(t)\rangle^{2}\right]/\left[\langle\epsilon(t)^{2}\rangle-\langle\epsilon(t)\rangle^{2}\right]$
can be easily calculated as

$C_{\epsilon}(\tau)=e^{-6D_{r}\tau}.$ (10)

If the eccentricity is large, higher order spherical harmonics are expected to
appear in Eq.7, leading to multiple exponential decay of the correlation
functions berne2000dynamic .

## Appendix B. Rotational tracking algorithm

Figure 6:  Schematic representation of the rotational tracking algorithm. The
coordinates $(x_{m},y_{m})$ of the centroid of each particle are determined
$via$ a standard particle tracking algorithm PhysRevLett.102.188303 . The
shape and orientation of the particle’s image are determined via a bivariate
Gaussian fit to the local intensity profile. This enables estimating the
length $\sigma_{x}$ and $\sigma_{y}$ of the short and of the long axis,
respectively, and the azimuthal angle $\phi$ formed by the particle’s long
axis with the horizontal direction. The effective aspect ratio
$\epsilon=\sigma_{y}/\sigma_{x}$ depends of the polar angle $\theta$ between
the particle’s axis and the vertical direction, as described in the text.

In this Appendix, we present a detailed description of the algorithm used to
track the rotational motion of the particles in real space. As described in
Par. 2.4, we first apply a standard particle tracking algorithm
PhysRevLett.102.188303  to determine, frame by frame, the position of the
centroid of each particle. For each tracked particle $i$ and for each time
$t$, a square ROI of size 16x16 pixel (corresponding to area of about 27 $\mu
m^{2}$ when the 40X objective is used) is obtained from the original image,
centered on the particle’s centroid; a bivariate Gaussian distribution with
arbitrary orientation of the main axes $f(x,y)=A\exp{-[x^{\prime
2}/2\sigma_{x}^{2}-y^{\prime 2}/2\sigma_{y}^{2}]}$, with
$x^{\prime}=[(x-x_{m})\cos\phi+(y-y_{m})\sin\phi]$ and
$y^{\prime}=[(x-x_{m})\sin\phi-(y-y_{m})\cos\phi]$, is fitted to the ROI
intensity distribution, which provides an estimate for the variances
$\sigma_{x}^{2}(t)$ and $\sigma_{y}^{2}(t)$ along the main axes, and for the
azimuthal angle $\phi_{i}(t)$ identifying the direction of projection of the
long axis on the $x,y$ plane. The values obtained for $\sigma_{x}(t)$ and
$\sigma_{y}(t)$ are quite consistent along each trajectory and display time-
correlated fluctuation around their mean values, which we attribute to random
variations in the polar angle $\theta$ formed by the long axis of the particle
with the vertical direction (see Fig. 7).

Figure 7:  a) Blue (orange) curve: time-dependent length of the major (minor)
axis $\sigma_{y}(t)$ ($\sigma_{x}(t)$) of a single particle, as obtained from
a bivariate Gaussian fit of its intensity distribution. The black continuous
line corresponds to the mode $\overline{\sigma_{y}}$ of the frequency
distribution of $\sigma_{y}(t)$

, while the black dashed line corresponds to threshold value
$0.95\overline{\sigma_{y}}$ (see main text for details). b) Frequency
distribution of $\sigma_{y}(t)$. The vertical line corresponds to the mode
$\overline{\sigma_{y}}=(0.90\pm 0.02)$ $\mu m$ of the distribution.

We then calculate the autocorrelation function
$C_{\epsilon}(\tau)=\left[\langle\epsilon(t+\tau)\epsilon(t)\rangle-\langle\epsilon(t)\rangle^{2}\right]/\left[\langle\epsilon(t)^{2}\rangle-\langle\epsilon(t)\rangle^{2}\right]$
of the effective aspect ratio $\epsilon(t)=\sigma_{y}(t)/\sigma_{x}(t)$. For
an elongated particle performing rotational Brownian motion,
$C_{\epsilon}(\tau)$ is expected to display a simple exponential relaxation
$C_{\epsilon}(\tau)=e^{-\gamma\tau}$, with $\gamma=6D_{R}$, as shown in
Appendix A. This provides a robust means to measure the rotational dynamics of
the particles, namely by fitting an exponential function to the obtained
$C_{\epsilon}(\tau)$.

In order to extract a similar information from the azimuthal degree of
freedom, which is encoded in $\phi$, the simultaneous knowledge of the polar
angle $\theta$ is also needed. This is due to the fact that the same change in
the azimuthal angle $\delta\phi$ can correspond to actual displacements of
very different amplitudes according to the value of $\theta$
anthony2015tracking . To account for this, a different variable $\psi$ is
often considered, which is defined by the relation $d\psi=d\phi\sin\theta$
colin2014rotational . For a particle performing rotational diffusion, the mean
squared value of $\psi$ follows the simple relation
$\langle\Delta\psi^{2}(\tau)\rangle=2D_{R}\tau$.

In our experiments, although we could in principle calculate $\theta(t)$ from
$\sigma_{x}(t)$, $\sigma_{y}(t)$, we found that this direct inversion was too
noisy to provide reliable results. As a consequence, we adopt a simplified yet
more robust approach. For each particle, we consider the histogram of the
values assumed by the long axis length $\sigma_{y}(t)$ and we identify the
value $\overline{\sigma_{y}}$ corresponding to the right peak of the
distribution (see Fig.6). According to the model described in the Appendix,
$\overline{\sigma_{y}}$ corresponds to the length of the long axis when it
lies on the $x-y$ plane, i.e. when $\theta=\pi/2$. For each particle, we
identify the portions of the trajectory such that

$(\overline{\sigma_{y}}-\sigma_{y}(t))/\overline{\sigma_{y}}<0.05.$ (11)

This corresponds to $|\theta_{i}(t)-\pi/2|\lesssim\pi/8$. Within each of these
portions, $\Delta\phi\simeq\Delta\psi$ and thus
$\langle\Delta\phi^{2}(\tau)\rangle\simeq 2D_{R}\tau$, where the average is
performed only over those portions for which Eq. 11 is satisfied.

## Acknowledgments

We acknowledge funding from from the Associazione Italiana per la Ricerca sul
Cancro (AIRC) - Project MFAG # 22083, and from the Italian Ministry of
University and Scientific Research (MIUR) - Project RBFR125H0M.

## Author contribution statement

FG and RC designed research. FG performed experiments and analyzed the data.
AP synthesized the particles and prepared the samples. All authors discussed
the experimental results. FG and RC wrote the paper.

## References

  * (1) J. Perrin, _Brownian movement and molecular reality_ (Courier Corporation, 2013)
  * (2) R. Cerbino, Current Opinion in Colloid & Interface Science 34, 47 (2018)
  * (3) Y. Han, A.M. Alsayed, M. Nobili, J. Zhang, T.C. Lubensky, A.G. Yodh, Science 314, 626 (2006)
  * (4) K.V. Edmond, M.T. Elsesser, G.L. Hunter, D.J. Pine, E.R. Weeks, Proceedings of the National Academy of Sciences 109, 17891 (2012)
  * (5) Z. Cheng, T.G. Mason, Phys. Rev. Lett. 90, 018304 (2003)
  * (6) R. Colin, L. Chevry, J.F. Berret, B. Abou, Soft Matter 10, 1167 (2014)
  * (7) W.L. Chen, H.S. Chuang, Analytical Chemistry 92, 12996 (2020)
  * (8) J. Saragosti, P. Silberzan, A. Buguin, PloS one 7, e35412 (2012)
  * (9) Y. Peng, L. Lai, Y.S. Tai, K. Zhang, X. Xu, X. Cheng, Physical review letters 116, 068303 (2016)
  * (10) Y. Gao, Y. Yu, L. Sanchez, Y. Yu, Micron 101, 123 (2017)
  * (11) S.M. Anthony, Y. Yu, Analytical Methods 7, 7020 (2015)
  * (12) B.J. Berne, R. Pecora, _Dynamic light scattering: with applications to chemistry, biology, and physics_ (Courier Corporation, 2000)
  * (13) V. Degiorgio, R. Piazza, R.B. Jones, Phys. Rev. E 52, 2707 (1995)
  * (14) M. Velez, D. Axelrod, Biophysical journal 53, 575 (1988)
  * (15) M. Lettinga, G. Koenderink, B. Kuipers, E. Bessels, A. Philipse, The Journal of chemical physics 120, 4517 (2004)
  * (16) F. Giavazzi, R. Cerbino, Journal of Optics 16, 083001 (2014)
  * (17) R. Cerbino, V. Trappe, Physical review letters 100, 188102 (2008)
  * (18) R. Cerbino, P. Cicuta, The Journal of chemical physics 147, 110901 (2017)
  * (19) M. Reufer, V.A. Martinez, P. Schurtenberger, W.C. Poon, Langmuir 28, 4618 (2012)
  * (20) C. Kurzthaler, C. Devailly, J. Arlt, T. Franosch, W.C.K. Poon, V.A. Martinez, A.T. Brown, Phys. Rev. Lett. 121, 078001 (2018)
  * (21) R. Nixon-Luke, G. Bryant, Journal of Physics: Condensed Matter 32, 115102 (2019)
  * (22) A. Pal, V.A. Martinez, T.H. Ito, J. Arlt, J.J. Crassous, W.C. Poon, P. Schurtenberger, Science Advances 6, eaaw9733 (2020)
  * (23) F. Giavazzi, C. Haro-Pérez, R. Cerbino, Journal of Physics: Condensed Matter 28, 195201 (2016)
  * (24) P. Edera, D. Bergamini, V. Trappe, F. Giavazzi, R. Cerbino, Physical Review Materials 1, 073804 (2017)
  * (25) R. Cerbino, D. Piotti, M. Buscaglia, F. Giavazzi, Journal of Physics: Condensed Matter 30, 025901 (2017)
  * (26) R. Pecora, The Journal of Chemical Physics 48, 4126 (1968)
  * (27) A. Wittmeier, A. Leeth Holterhoff, J. Johnson, J.G. Gibbs, Langmuir 31, 10402 (2015)
  * (28) T. Sugimoto, M.M. Khan, A. Muramatsu, Colloids and Surfaces A: Physicochemical and Engineering Aspects 70, 167 (1993)
  * (29) A. Pal, T. Zinn, M.A. Kamal, T. Narayanan, P. Schurtenberger, Small 14, 1802233 (2018)
  * (30) M.A. Kamal, A.V. Petukhov, A. Pal, The Journal of Physical Chemistry B 124, 5754 (2020)
  * (31) F. Giavazzi, D. Brogioli, V. Trappe, T. Bellini, R. Cerbino, Physical Review E 80, 031403 (2009)
  * (32) R.C. Gonzales, R.E. Woods, _Digital Image Processing_ (Addison Wesley, 1992)
  * (33) V. Pelletier, N. Gal, P. Fournier, M.L. Kilfoil, Phys. Rev. Lett. 102, 188303 (2009)
  * (34) Y. Han, A. Alsayed, M. Nobili, A.G. Yodh, Phys. Rev. E 80, 011403 (2009)
  * (35) Z. Zheng, Y. Han, The Journal of chemical physics 133, 124509 (2010)
  * (36) H. Löwen, Phys. Rev. E 50, 1232 (1994)
  * (37) Y. Yang, M.A. Bevan, The Journal of Chemical Physics 147, 054902 (2017), https://doi.org/10.1063/1.4995949
  * (38) K. Farain, A. Esfandiar, A.Z. Moshfegh, Applied Physics Letters 113, 063101 (2018)


# Modulation-free Laser Stabilization Technique Using Integrated Cavity-
Coupled Mach-Zehnder Interferometer

Mohamad Hossein Idjadi<EMAIL_ADDRESS>Kwangwoong Kim
Nokia Bell Labs, 600 Mountain Ave, Murray Hill, NJ 07974, USA.

###### Abstract

Stable narrow-linewidth light sources play a significant role in many
precision optical systems. Electro-optic laser frequency stabilization
systems, such as the well-known Pound-Drever-Hall (PDH) technique, have been
key components of stable laser systems for decades. These control loops
utilize an optical frequency noise discriminator (OFND) to measure frequency
noise and convert it to an electronic servo signal. Despite their excellent
performance, there has been a trade-off between complexity, scalability, power
consumption, and noise measurement sensitivity. Here, we propose and
experimentally demonstrate a modulation-free laser stabilization technique
using an integrated cavity-coupled Mach-Zehnder interferometer (MZI) as an
OFND. The proposed architecture maintains the sensitivity and performance of
the PDH architecture without the need for any modulation. This significantly
improves overall power consumption, simplifies the architecture, and makes it
easier to miniaturize into an integrated photonic platform. An on-chip
microring resonator with a loaded quality factor of 2.5 million is used as the
frequency reference. The implemented chip suppresses the frequency noise of a
semiconductor laser by 4 orders of magnitude. The integral linewidth of the
free-running laser is suppressed from 6.1 MHz to 695 KHz. The passive
implemented photonic integrated circuit occupies an area of 0.456 mm2 and is
integrated on AIM Photonics 180 nm silicon-on-insulator process.

##  I Introduction

Precise laser frequency control is a crucial requirement in various
applications, including optical communication Al-Taiy et al. (2014);
Blumenthal et al. (2020), optical atomic clocks Ludlow et al. (2015),
microwave photonic Li et al. (2013), and sensing Marra et al. (2018) which
makes stable and narrow linewidth lasers indispensable part of precision
optical experiments. Researchers have explored various methods to suppress
unwanted laser frequency noise using optical feedbackLiang et al. (2015);
Dahmani et al. (1987), electro-optic feedback Drever et al. (1983); Idjadi and
Aflatouni (2020), and electro-optic feed-forward Aflatouni and Hashemi (2012)
techniques. The Optical feedback method relies on the Rayleigh scattering
Gorodetsky et al. (2000) within the cavity and can be an effective and
promising way to suppress the laser frequency noise. This mandates the co-
design of laser system and optical signals with an ultra-high quality factor
(Q-factor) cavity. Electro-optic techniques, on the other hand, leverage
state-of-the-art and mature electronic devices and systems to control the
laser frequency precisely. This will push some of the challenges from optical
domain into the electrical one where control and manipulation of the signals
and systems are comparatively more manageable and cost-effective.
A key building block in electro-optic laser stabilization techniques is an
optical frequency noise discriminator. An OFND measures frequency fluctuations
by comparing to a frequency reference and translates it into an electronic
signal that can be processed in the electrical domain. Different OFND
configurations have been explored such as ”squash” locking technique Diorico
et al. (2022); Chabbra et al. (2021), the PDH laser stabilization method
Drever et al. (1983), and the unbalanced MZI Sorin et al. (1992); Idjadi and
Aflatouni (2020). The PDH loop stands out as the most well-known precision
laser instrumentation technique among the extensively utilized OFNDs Liu et
al. (2022); Spencer et al. (2016); Wang et al. (2020); Kelleher et al. (2023).
Using the PDH technique, a sharp asymmetric error signal can be generated,
which can then be utilized as a servo signal to stabilize the laser frequency.
Despite excellent performance, the PDH scheme requires electro-optic phase
modulator and relatively fast and complex electronics for modulation and
demodulation which increases the power consumption and area for an integrated
PDH chip Idjadi and Aflatouni (2019, 2017). Alternatively, a passive-only
unbalanced MZI can serve as an OFND where the two arms of the MZI are phase
locked at the quadrature point Bagheri et al. (2009). Although the unbalanced
MZI has a simple architecture to implement on a chip, achieving high-frequency
detection sensitivity comparable to the PDH method requires either a large
optical delay line or a substantial electronic gain that comes at the cost of
chip area of overall system power consumption.
Here, we propose and experimentally demonstrate a modulation-free laser
stabilization technique using a cavity-coupled MZI on a silicon photonic chip
as an OFND. The proposed frequency noise discrimination technique utilizes a
high Q-factor cavity coupled to a MZI which breaks the trade-off between
sensitivity, complexity, chip area, and power consumption.

With a careful design of the on-chip high Q-factor cavity coupled to an MZI, 4
orders of magnitude suppression in frequency noise of a commercially available
distributed feedback (DFB) laser is achieved. On-chip thermal tuners are
implemented for the potential trimming of fabrication-induced errors and also
for facilitating the broadband operation of the OFND. The proposed
architecture combines the advantages of a passive-only structure, which offers
simplicity and less electronic power consumption for error signal processing,
with the high sensitivity of the widely utilized PDH technique. The proof-of-
concept integrated photonic chip is fabricated in AIM Photonic commercially
available 180 nm silicon-on-insulator (SOI) process. The photonic chip
occupies an area of 0.456 mm2 and consumes only about 50 $\mu$W power for
reverse biasing the balanced photo detector. The proposed architecture offers
a promising solution for achieving sensitive, simple, and low power laser
frequency stabilization systems and sets the stage for the development of low-
cost, scalable, and stable integrated lasers.

##  II Results

###  II.1 The principle of operation

Fig. 1: The cavity-coupled MZI frequency discriminator. (a) The block diagram
of an electro-optic laser frequency stabilization using an optical frequency
noise discriminator (OFND). The OFND response is asymmetric around the
frequency reference point,$f_{ref}$, and hence small frequency fluctuations
translate into electrical signal by the OFND gain. (b) The conceptual diagram
of the proposed cavity-coupled MZI. $T(\omega)$ is the transfer function of a
frequency reference (e.g. optical cavity) coupled into the top arm of the MZI.
(c) Numerical analysis comparing the normalized error signal for the proposed
cavity-coupled MZI, a conventional unbalanced MZI and the PDH architecture.
$\Delta f$ and $\Delta\nu$ are the laser offset frequency compared to
$f_{ref}$ and the free-spectral range of the cavity.

Figure 1(a) shows the block diagram of an electro-optic laser frequency noise
reduction loop. As mentioned earlier, the key part of this loop is an OFND
that senses the frequency fluctuations of the incoming laser signal, compares
it with an optical frequency reference ($f_{ref}$ ), and generates an
electronic signal whose amplitude is proportional to the intensity of
frequency fluctuations. This error signal is amplified in the electrical
domain and fed back into the laser to stabilize its frequency. The details of
the close-loop operation and the linearized block diagram of the control loop
is discussed comprehensively in Supplementary Note 1. The OFND can
mathematically be represented by an asymmetric transfer function where the
small frequency perturbation around $f_{ref}$ is amplified by a gain
($K_{FD}$) and converted into an electronic error signal. The slope of the
transfer function also indicates the sensitivity of the OFND in measuring
frequency noise.
Figure 1(b) shows the proposed cavity-coupled MZI as an OFND that maintains a
simple passive structure without any need for fast optical phase modulation.
In this method, the incoming laser intensity is split equally into two MZI
branches using a broadband Y-junction. The top branch of the MZI is coupled
into an optical frequency reference (e.g. an optical resonator or a cavity).
The amplitude and phase of the electric field at the output of the cavity is
affected by that of the frequency reference, $T(\omega)$. The bottom arm of
the MZI is used to interfere with the frequency reference output electric
field using a directional coupler. a balanced photodetector is used to
photodetect the output of the MZI and by subtracting the currents, the error
signal ($i_{error}$) is generated. In other words, the proposed architecture
is a coherent detector that uses the input laser signal (bottom arm of the
MZI) to down-convert the signal at the output of the cavity. In this way the
sharp asymmetric phase transition in the transfer function of the cavity at
frequency of $f_{ref}$ will translate in a sharp electrical error signal that
can be used to lock the laser in a feedback loop. As discussed in
Supplementary Note 2, the error signal can be written as

$i_{error}(\omega)=RP_{0}|T(\omega)|\times Sin(\psi(\omega)-\phi),$ (1)

where $|T(\omega)|,\psi(\omega),\phi,P_{0}$, and $R$ are the amplitude and
phase of the optical reference at the frequency of $\omega$, phase difference
between arms of the MZI controlled by a thermal phase shifter, intensity of
the electric field at the input of the MZI, and the responsivity of the
photodetectors, respectively. Figure 1(c) is the numerical analysis of the
normalized error signal of the PDH, a conventional unbalanced MZI, and the
proposed cavity-coupled MZI structure using Eq. (1). To ensure a fair
comparison, assuming a fixed area to implement the OFND, a 1 mm circumference
ring resonator is used as the frequency reference in both the PDH and the
cavity-coupled MZI. The length mismatch between the arms of the unbalanced MZI
is also set to 1 mm . The waveguide loss is assumed 0.2 dB/cm. The details of
the numerical comparison is presented in Supplementary Note 3. As shown in
Fig. 1(c), the error signal of the cavity-coupled MZI is significantly sharper
than the conventional unbalanced MZI for a same given setting, and indeed is
comparable to that of the PDH offering the same level of sensitivity but much
simpler architecture.

Fig. 2: The modulation-free laser stabilization scheme. (a) The proposed
cavity-coupled MZI is used as an OFND in a feedback loop for laser frequency
noise suppression. A high Q-factor Euler microring resonator is utilized as a
frequency reference. Half of the on-chip laser intensity is injected to the
microring resonator and coherently interfered with the bottom branch and
converted to electrical signal using balanced Germanium photodetectors. A
small portion of the ring resonator output is used to monitor its response.
The generated error signal, $i_{error}$, is then amplified, and fed into the
PID controller. (b) The micro-photograph of the integrated photonic chip in
AIM Photonics 180 nm silicon-on-insulator (SOI) process. The size of the
photonic integrated circuit measures 0.95 mm $\times$ 0.48 mm. PID:
proportional-integral-derivative, TIA: trans-impedance amplifier.

###  II.2 Modulation-free laser frequency noise suppression system

Figure 2(a) shows the block diagram of the implemented modulation-free laser
stabilization scheme using the implemented cavity-coupled MZI as an OFND. As
illustrated in Fig. 2(a), the laser output is coupled into the chip via a
grating coupler and divided in half by a Y-junction. In the top branch, a high
Q-factor microring resonator is used as an optical frequency reference which
filters the amplitude and phase of the incoming light. The phase response of
the ring resonator exhibits rapid and asymmetric change around its resonance.
When combined with the amplitude response, this characteristic provides
sufficient information to determine the offset between the laser frequency and
$f_{ref}$, as well as whether the laser frequency is higher or lower than
$f_{ref}$. The output of the MZI uses a 2x2 adiabatic broadband couplerChen et
al. (2017) terminated to balanced photodetectors. As described in detail in
Supplementary Note 2, the output error signal is asymmetric around $f_{ref}$
and can be used as a servo signal. The error signal is amplified and converted
to voltage using an off-chip trans-impedance amplifier (TIA). The voltage
signal is then fed into a proportional-integral-derivative (PID) controller,
which modulates the laser current and corrects any frequency error. To
compensate for potential fabrication-induced errors and adjust the microring
resonance frequency and the optical phase of the bottom MZI branch, thermal
phase shifters that are thermally isolated by a deep trench are utilized.
Figure 2(b) shows the micro-photograph of the photonic integrated circuit
implemented in AIM Photonics 180 nm SOI process. The silicon photonic chip
area is 0.456 mm2.

###  II.3 High Q-factor silicon microring resonator

Utilization of a microring resonator as an optical frequency reference plays a
crucial role in OFND performance. A stable high Q-factor microring resonator ,
when coupled to the MZI, effectively enhances the sensitivity of the OFND.
This, in turn, directly impacts the closed-loop operation and the ultimate
laser frequency noise suppression. It is important to highlight that the
proposed architecture can be implemented not only on various material
platforms, such as low-loss silicon nitride, but also using bench-top ultra-
low expansion and stable etalons. Choosing silicon as a platform to implement
the optical frequency reference is a trade-off between different design
considerations such as potential co-integration with CMOS electronics, chip
area, scalability, cost, and ultimate frequency stability.
In order to achieve a high Q-factor microring cavity, it is necessary to
minimize intra-cavity losses including silicon nanophotonic waveguides. The
propagation loss is influenced by several factors, with interface scattering
and bend radiation loss being the most significant in a state-of-the-art
silicon photonic foundry process Bauters et al. (2011); Vlasov and McNab
(2004); Fahrenkopf et al. (2019). The top-bottom surface roughness of a
waveguide is well-controlled by the foundry and it is not a design parameter.
However, careful design of waveguide width can significantly improve the
propagation loss. Figure 3 illustrates the implemented microring resonator in
silicon. Utilizing multi-mode waveguides can reduce TE-mode overlap with side-
wall roughness and greatly enhances the waveguide transmission Bauters et al.
(2011). The implemented multi-mode waveguide is 2.2 $\mu$m wide and the
theoretical fit model applied to the measured microring response suggests a
propagation loss of approximately 0.2 dB/cm. Although wide multi-mode
waveguides can greatly reduce interfacial scattering loss, bend radiation loss
increases significantly if a tight bend is used Jiang et al. (2018),
especially for a highly multi-mode waveguide. An Euler bend is employed to
achieve a compact multi-mode bend with minimal excitation of higher-order
modes and mode cross-talk. The ring resonator is designed to achieve critical
coupling to maximize OFND sensitivity, however, the fabricated ring are
slightly under-coupled due to the potential fabrication-induced errors.
Supplementary Note 4 discusses OFND gain sensitivity to waveguide loss and
ring coupling ratio.

Fig. 3: High Q-factor silicon microring resonator. In order to minimize the
interface scattering loss due to waveguide edge roughness, a wide multi-mode
waveguide is utilized. Moreover, to avoid excitation of higher order modes and
Euler bend is used. As shown by the Finite-Difference Time-Domain simulation
of the waveguide cross-sections (y-z and x-y planes), the fundamental TE-mode
is well preserved inside the cavity. The higher order mode excitation is less
than -23 dB.

As illustrated in Fig. 3, the fundamental TE-mode remains preserved within the
bent multi-mode waveguide. The implemented microring resonator has a
circumference of about 950 $\mu$m which corresponds to a free spectral range
(FSR) of about 83.5 GHz with a loaded Q-factor of about 2.5 million at the
resonance wavelength of 1550.73 nm.

###  II.4 The open-loop operation: device characterization and error signal

Fig. 4: The open loop response. (a) Measurement setup for the open loop and
microring resonator characterization. (b) The microring resonator response is
measured via the top photodetector sniffer. The ring has a measured extinction
ratio (ER) of 4.8 dB and a loaded Q-factor of 2.5 million. (c) The measured
asymmetric error signal. Heaters are off during the test. CW: continuous wave,
OSA: optical spectrum analyzer, GC: grating coupler, BPD: balanced
photodetector.

Figure 4(a) shows the schematic of the measurement setup to characterize the
open-loop performance and the error signal. A tunable continues wave laser
(TOPTICA CTL 1550) with a wavelength of 1550.7 nm is coupled into the chip via
the on-chip grating coupler. To ensure linear operation and thermal stability
of the high Q-factor microring resonator, laser power coupled to the grating
coupler is set to 0.7 mW using a variable optical attenuator (VOA). The
silicon chip temperature is stabilized at 27oC. The laser frequency is
continuously scanned within the range of 30 GHz. A calibrated fiber-based MZI
with a FSR of 20 MHz is used for time-frequency translation. As shown in Fig.
4(a), a sniffer photodetector with 5% coupling ratio is utilized to monitor
the resonance response. The resonance response of the ring and the error
signal are measured simultaneously using an oscilloscope . Figure 4(b) shows
the response of the ring where the Q-factor is about 2.5 million. The measured
extinction ratio is about 4.8 dB that suggests the fabricated microring is
slightly under coupled, likely due to fabrication error in the coupling gap.
Figure 4(c) shows the measured normalized asymmetric error signal that
suggests OFND sensitivity of 1.3$\times$10-8 Hz-1 which agrees well with
analytical models. Since the ring and MZI are optimized for 1550 nm, there was
no need to utilize Heaters 1 and 2 during the measurement.

###  II.5 The closed-loop operation and laser frequency noise suppression

Fig. 5: The closed-loop operation. (a) Measurement setup for closed-loop
operation and laser locking experiment. The small portion of the laser output
power is used in the delayed self-heterodyne (DSH) for frequency noise
measurement. Heaters are off during this measurement.(b) The power spectral
density of frequency noise of the AeroDiode DFB laser under free-running and
closed-loop operation. The highlighted regions indicate frequency band which
contribute to the laser linewidth. The frequency noise of the DFB laser is
suppressed by 40 dB at low Fourier frequencies and is limited to the silicon
microring thermorefractive noise (TRN). VOA: variable optical attenuator, TIA:
trans-impedance amplifier, CW: continues wave, GC: grating coupler.

As mentioned previously, the asymmetric error signal can serve as the servo
signal to precisely adjust and lock the laser frequency to the resonance
frequency of the cavity. Figure 5 (a) shows the block diagram of the closed-
loop measurement setup. As shown in Fig. 5(a), a DFB laser (AeroDiode
1550LD-2-0-0-1) with a free-running integral linewidth of 6.1 MHz is used. The
laser output is directed into the chip by utilizing a 90%/10% coupler, which
is then followed by a VOA. This arrangement allows for the adjustment of the
coupled power into the chip to about 0.5 mW. The generated error signal is
amplified using a low noise TIA with a gain of 5 K$\Omega$ and total input-
referred current noise of 3.4 pA/$\sqrt{Hz}$ which is three orders of
magnitude lower than the thermorefractive noise (TRN) of the microring. Part
of the laser power is used in a fiber-based delayed self-heterodyne (DSH)
interferometer and the output signal is digitized and sent to a computer for
processing Yuan et al. (2022). Figure 5(b) shows the frequency noise
measurement results for the free-running and locked laser. As shown in Fig.
5(b), the cavity-coupled MZI OFND can effectively suppress the frequency noise
of the free-running DFB laser by 40 dB at low Fourier frequencies. The
frequency noise reduction bandwidth is about 300 KHz which is mainly limited
by the frequency modulation response of the laser Vankwikelberge et al. (1989)
and the modulation bandwidth of the current driver ($<$1 MHz). The
$\beta$-separation line indicates highlighted regions that contribute to the
laser linewidth Di Domenico et al. (2010). As shown by the highlighted regions
in Fig. 5(b), the integral linewidth of the free-running laser is reduced from
6.1 MHz to 695 KHz under closed-loop condition. It is worth mentioning that
the suppressed frequency noise is limited to the TRN of the silicon microring
resonator Huang et al. (2019); Panuski et al. (2020).

##  III Discussion

Ideally, a cavity that supports a larger mode volume and is made from a
temperature insensitive material can enable a lower TRN limit, resulting in a
more stable laser. Ultra-low loss silicon nitride and micro-fabricated mirrors
Jin et al. (2022) present themselves as promising technology platforms for
fulfilling these requirements. However, this improvement comes at the cost of
chip area, packaging complexity, and availability of the technology for large-
scale and robust integrated photonic chip manufacturing. On the other hand,
silicon photonic platforms are evolving rapidly and while silicon may not be
the most suitable choice for achieving highly stable micro-cavities, it
certainly serves as an excellent alternative in numerous applications where
strict level of laser stability may not be required. The combination of robust
large-scale silicon photonic manufacturing and the ability to monolithically
integrate with mature CMOS electronics provides integrated silicon photonics
with a significant advantage over other alternatives. It is important to note
that the choice of technology depends on the specific application. Our
proposed architecture can be implemented on a suitable technology platform to
effectively meet the requirements of an application.

In conclusion, we have experimentally demonstrated the first modulation-free
laser frequency noise stabilization using integrated cavity-coupled MZI. The
integrated cavity-coupled MZI consists of a high-Q Euler microring resonator
with a loaded quality factor of about 2.5 million and extinction ratio of
about 4.8 dB which suggests normalized OFND gain of 1.3$\times$10-8 Hz-1. The
cavity-coupled MZI is used to suppress the frequency noise of a commercially
available DFB laser by 40 dB at 1 KHz offset frequency. This corresponds to
the integral linewidth reduction of the free-running laser from 6.1 MHz to 695
KHz. The implemented chip occupies 0.456 mm2 and is integrated on AIM
Photonics 180 nm process.

##  IV Methods

###  IV.1 Photonic chip implementation

All photonic devices are monolithically integrated on AIM Photonics
commercially available 180 nm SOI process. The laser is coupled into the chip
vertically using a grating coupler. A Y-junction with excess loss of less than
0.5 dB is utilized at the input of the MZI to split the laser intensity
equally. A high Q-factor microring resonator is coupled to the top arm of the
MZI. The ring resonator is made with wide multi-mode waveguides and Euler
bends to both reduce the interfacial scattering loss and avoid excitation of
higher-order modes. Thermal phase shifters are integrated into the photonic
chip for potential adjustment of phase or compensation for fabrication-induced
errors. The two arms of the MZI are combined using an adiabatic broadband
directional coupler. A balanced photodetector is used to generate the error
signal. The responsivity and dark current of the on-chip Germanium
photodetector at 2.5 V reverse bias voltage are 1.16 A/W (at 1550 nm) and 40
nA, respectively. The photonic chip occupies area of 0.456 mm2.

###  IV.2 The closed-loop operation in presence of noise sources

Different noise sources contribute to the ultimate frequency noise limit. The
two main sources of noise are the TRN of the reference cavity and the total
electronic noise that includes input referred noise of electronics and shot
noise of balanced photodetector. The closed-loop operation and minimum
achievable frequency noise in the presence of these noise sources are
discussed in Supplementary Note 1.

###  IV.3 Delayed self-heterodyne frequency noise measurement

A fiber-based delayed self-heterodyne setup is built for laser frequency noise
measurement. The first null of the MZI is at 8.7 MHz, corresponding to a
single-mode fiber length of about 23 m. An acousto-optic frequency shifter is
used in the second arm of the MZI to up-convert the laser phase noise to 27
MHz which helps to mitigate electronic flicker noise for better measurement
sensitivity. The two arms of the MZI are cross-coupled and fed into two
independent low noise TIAs. The outputs of the TIAs are digitized at 250 MS/s
and recorded for a duration of 5 msec and sent to a computer for processing.

###  IV.4 The OFND gain sensitivity analysis

The design of the high Q-factor microring resonator in the OFND directly
impacts the frequency noise detection gain and laser frequency noise
suppression. Supplementary Note 4 discusses the sensitivity of OFND gain to
both the coupling ratio of the ring and the waveguide propagation loss.

## Data availability

The data that support the plots within this paper and other findings of this
study are available from the corresponding author upon reasonable request.

## Acknowledgements

We thank Nicolas K. Fontaine and Andrea Blanco-Redondo for helpful discussions
and support.

## Author contributions

M.H.I. conceived the project idea. M.H.I. designed, simulated and laid out the
integrated photonic circuits and conducted measurements. M.H.I designed the
printed circuit board and the on-board electronics. K.K. packaged the
integrated photonic chip. M.H.I. wrote the manuscript.

## Competing interests

The authors declare no competing interests.

## References

  * Al-Taiy et al. (2014) H. Al-Taiy, N. Wenzel, S. Preußler, J. Klinger, and T. Schneider, Optics Letters 39, 5826 (2014).
  * Blumenthal et al. (2020) D. J. Blumenthal, H. Ballani, R. O. Behunin, J. E. Bowers, P. Costa, D. Lenoski, P. A. Morton, S. B. Papp, and P. T. Rakich, Journal of Lightwave Technology 38, 3376 (2020).
  * Ludlow et al. (2015) A. D. Ludlow, M. M. Boyd, J. Ye, E. Peik, and P. O. Schmidt, Reviews of Modern Physics 87, 637 (2015).
  * Li et al. (2013) J. Li, H. Lee, and K. J. Vahala, Nature communications 4, 2097 (2013).
  * Marra et al. (2018) G. Marra, C. Clivati, R. Luckett, A. Tampellini, J. Kronjäger, L. Wright, A. Mura, F. Levi, S. Robinson, A. Xuereb, et al., Science 361, 486 (2018).
  * Liang et al. (2015) W. Liang, V. Ilchenko, D. Eliyahu, A. Savchenkov, A. Matsko, D. Seidel, and L. Maleki, Nature communications 6, 7371 (2015).
  * Dahmani et al. (1987) B. Dahmani, L. Hollberg, and R. Drullinger, Optics letters 12, 876 (1987).
  * Drever et al. (1983) R. W. Drever, J. L. Hall, F. V. Kowalski, J. Hough, G. Ford, A. Munley, and H. Ward, Applied Physics B 31, 97 (1983).
  * Idjadi and Aflatouni (2020) M. H. Idjadi and F. Aflatouni, Nature Photonics 14, 234 (2020).
  * Aflatouni and Hashemi (2012) F. Aflatouni and H. Hashemi, Optics letters 37, 196 (2012).
  * Gorodetsky et al. (2000) M. L. Gorodetsky, A. D. Pryamikov, and V. S. Ilchenko, JOSA B 17, 1051 (2000).
  * Diorico et al. (2022) F. Diorico, A. Zhutov, and O. Hosten, arXiv preprint arXiv:2203.04550 (2022).
  * Chabbra et al. (2021) N. Chabbra, A. R. Wade, E. R. Rees, A. J. Sutton, A. Stochino, R. L. Ward, D. A. Shaddock, and K. McKenzie, Optics Letters 46, 3199 (2021).
  * Sorin et al. (1992) W. V. Sorin, K. W. Chang, G. A. Conrad, and P. R. Hernday, Journal of Lightwave Technology 10, 787 (1992).
  * Liu et al. (2022) K. Liu, N. Chauhan, J. Wang, A. Isichenko, G. M. Brodnik, P. A. Morton, R. O. Behunin, S. B. Papp, and D. J. Blumenthal, Optica 9, 770 (2022).
  * Spencer et al. (2016) D. T. Spencer, M. L. Davenport, T. Komljenovic, S. Srinivasan, and J. E. Bowers, Optics express 24, 13511 (2016).
  * Wang et al. (2020) Z. Wang, H. Wei, Y. Li, R. Kan, and W. Ren, Optics Letters 45, 1148 (2020).
  * Kelleher et al. (2023) M. L. Kelleher, C. A. McLemore, D. Lee, J. Davila-Rodriguez, S. A. Diddams, and F. Quinlan, Optics Express 31, 11954 (2023).
  * Idjadi and Aflatouni (2019) M. H. Idjadi and F. Aflatouni, in _2019 IEEE Radio Frequency Integrated Circuits Symposium (RFIC)_ (IEEE, 2019), pp. 319–322.
  * Idjadi and Aflatouni (2017) M. H. Idjadi and F. Aflatouni, Nature communications 8, 1209 (2017).
  * Bagheri et al. (2009) M. Bagheri, F. Aflatouni, A. Imani, A. Goel, and H. Hashemi, Optics letters 34, 2979 (2009).
  * Chen et al. (2017) G. F. Chen, J. R. Ong, T. Y. Ang, S. T. Lim, C. E. Png, and D. T. Tan, Scientific reports 7, 7246 (2017).
  * Bauters et al. (2011) J. F. Bauters, M. J. Heck, D. John, D. Dai, M.-C. Tien, J. S. Barton, A. Leinse, R. G. Heideman, D. J. Blumenthal, and J. E. Bowers, Optics express 19, 3163 (2011).
  * Vlasov and McNab (2004) Y. A. Vlasov and S. J. McNab, Optics express 12, 1622 (2004).
  * Fahrenkopf et al. (2019) N. M. Fahrenkopf, C. McDonough, G. L. Leake, Z. Su, E. Timurdogan, and D. D. Coolbaugh, IEEE Journal of Selected Topics in Quantum Electronics 25, 1 (2019).
  * Jiang et al. (2018) X. Jiang, H. Wu, and D. Dai, Optics express 26, 17680 (2018).
  * Yuan et al. (2022) Z. Yuan, H. Wang, P. Liu, B. Li, B. Shen, M. Gao, L. Chang, W. Jin, A. Feshali, M. Paniccia, et al., Optics Express 30, 25147 (2022).
  * Vankwikelberge et al. (1989) P. Vankwikelberge, F. Buytaert, A. Franchois, R. Baets, P. Kuindersma, and C. Fredriksz, IEEE journal of quantum electronics 25, 2239 (1989).
  * Di Domenico et al. (2010) G. Di Domenico, S. Schilt, and P. Thomann, Applied optics 49, 4801 (2010).
  * Huang et al. (2019) G. Huang, E. Lucas, J. Liu, A. S. Raja, G. Lihachev, M. L. Gorodetsky, N. J. Engelsen, and T. J. Kippenberg, Physical Review A 99, 061801 (2019).
  * Panuski et al. (2020) C. Panuski, D. Englund, and R. Hamerly, Physical Review X 10, 041046 (2020).
  * Jin et al. (2022) N. Jin, C. A. McLemore, D. Mason, J. P. Hendrie, Y. Luo, M. L. Kelleher, P. Kharel, F. Quinlan, S. A. Diddams, and P. T. Rakich, Optica 9, 965 (2022).