Phase-Resolved Numerical Lock-In Analysis

• Phase-resolved numerical lock-in analysis is a computational method that decomposes periodic signals into phase-sensitive components using dual-phase demodulation.
• It applies to diverse fields such as thermography, synchronization, and quantum sensing by enabling precise parameter estimation and robust feedback control.
• Its numerical protocols integrate time- and frequency-domain techniques to invert and classify phase data, thereby enhancing measurement accuracy and system stability.

Phase-resolved numerical lock-in analysis denotes a class of computational procedures in which a periodic response is referenced to a known modulation, decomposed into phase-sensitive components, and then interpreted for estimation, classification, or feedback. In the cited literature, the term is used for both dual-phase demodulation of measured signals and numerical analysis of lock-in behavior in nonlinear dynamical systems. Representative implementations include pixelwise Lock-In Thermography for thickness estimation (Tian et al., 2019), phase-based infrared thermography for thermal conductivity extraction (Scott et al., 10 Jun 2026), software lock-in spectral interferometry for ultrafast nonlinear optics (Pandey et al., 7 Mar 2025), FPGA-based quadrature processing in ultrasound pulse-echo experiments (Galeski et al., 2022), and phase-plane computation of lock-in ranges in phase-locked loops (Aleksandrov et al., 2016).

1. Core numerical formalism

In lock-in thermography, an object is periodically heated and the temperature signal at each pixel is recorded over time. The recorded signal $I_n(r,c)$ is multiplied by reference sine and cosine functions at the excitation frequency,

$$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$

and summed to form in-phase and quadrature channels,

$$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$

The lock-in amplitude and phase are then

$$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$

These quantities are the basic phase-resolved observables for subsequent inversion or classification (Tian et al., 2019).

The same dual-phase structure recurs in other implementations. In software lock-in spectral interferometry, each spectral pixel is demodulated to produce

$$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$

with amplitude and phase given by $A=\sqrt{X^2+Y^2}$ and $\phi=\arctan2(Y,X)$; an analogous quadrature demodulation is implemented in FPGA-based ultrasound experiments, where synchronous mixing and boxcar averaging are applied within user-defined echo windows (Pandey et al., 7 Mar 2025, Galeski et al., 2022). In infrared thermography for thermal conductivity measurements, the measured thermal response is represented as

$$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$

so that phase becomes the fitted quantity rather than a purely descriptive output (Scott et al., 10 Jun 2026).

A frequency-domain variant appears in all-phase FFT based frequency locking. There, weighted recombination of $2N-1$ samples followed by an $N$-point FFT yields a center-aligned phase estimate that is described as unbiased with respect to frequency-bin offset, and the resulting phase difference across successive intervals is converted into a frequency-deviation estimate for digital feedback control (Zheng et al., 16 Jun 2026). This establishes that phase-resolved numerical lock-in analysis is not restricted to time-domain mixing; it also includes frequency-domain estimators designed for feedback-grade phase accuracy.

2. Thermographic inversion and parameter extraction

A central thermographic application is the use of Lock-In Thermography to infer object geometry. One study addresses the limitation that LIT “fails to estimate the thickness at a point on the tested object,” and therefore “unable to figure out the 3-dimensional geometry of an object,” by introducing two techniques that relate lock-in observables to local thickness (Tian et al., 2019). In Technique I, an empirical functional relation is built between thickness $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$0, lock-in phase $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$1, and a derived parameter $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$2: $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$3 which is then refined to

$$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$4

where $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$5 is a degree-5 polynomial fit to mean phase versus thickness and $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$6 is an exponential fit capturing the offset between mean phase and intercept. Thickness is estimated by inverting this relation. The same source states that the method “works well in regions with sufficient phase variation (not on extended flat/thick regions)” and that the inversion “cannot always be uniquely or robustly performed” (Tian et al., 2019).

Technique II replaces an explicit inversion with database retrieval. Each entry is $$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$7, constructed from experimental or simulated data, and estimation proceeds by minimizing the Euclidean distance

$$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$8

Similar entries within a preset resolution are averaged, and the database may be reduced by Principal Component Analysis using the normalized covariance matrix

$$K_0(t_n)=\sin(2\pi f t_n), \qquad K_{-90}(t_n)=-\cos(2\pi f t_n),$$9

Evaluation is performed by measuring root-mean-square deviation and calculating successful rate with different tolerances. The same work further proposes Stochastic Gradient Descent for determining “the time when sufficient data have been collected from LIT measurement to generate the estimated geometry accurately” (Tian et al., 2019).

A more recent thermographic formulation targets thermal conductivity rather than thickness. In “phase-based lock-in thermography,” the spatial phase distribution around a modulated laser spot is radially averaged, compared with a multilayer cylindrical heat-diffusion model, and fitted by least squares to extract $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$0 (Scott et al., 10 Jun 2026). The fit minimizes

$$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$1

and uncertainty is evaluated by Monte Carlo simulations that perturb spot size, pixel scaling, layer thicknesses, absorber properties, and interface conductance. Sensitivity is defined by

$$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$2

The study reports measurements on materials spanning “over three orders of magnitude (approximately $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$3 to $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$4),” notes that the method is “insensitive to surface roughness,” and states that a transducer layer is “not strictly required,” although a removable adhesive absorber can improve signal quality (Scott et al., 10 Jun 2026).

These thermographic studies show two distinct numerical roles for phase: as a feature used in direct lookup or inversion, and as a model-constrained observable used in parameter fitting. This suggests that in thermography the lock-in phase is not merely a contrast image but a quantitative state variable.

3. Lock-in ranges, cycle slipping, and nonlinear phase-space analysis

In synchronization theory, phase-resolved numerical lock-in analysis is often formulated in phase space rather than through quadrature channels. For a PLL with sinusoidal phase detector characteristic $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$5 and active PI filter $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$6, the autonomous model is

$$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$7

A rigorous definition of the lock-in range is then given in terms of global asymptotic stability and a lock-in domain containing the relevant equilibria, with cycle slipping defined by

$$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$8

The numerical procedure computes separatrices at $$S_0(r,c)=\sum_{n=1}^{N} I_n(r,c)K_0(t_n), \qquad S_{-90}(r,c)=\sum_{n=1}^{N} I_n(r,c)K_{-90}(t_n).$$9 and obtains the lock-in frequency from

$$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$0

while for the active PI filter the pull-out frequency satisfies $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$1 (Aleksandrov et al., 2016).

The same line of work addresses Gardner’s problem of defining a unique lock-in frequency. A general nonlinear PLL model is written in signal phase space,

$$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$2

and the lock-in range is distinguished from the pull-in range by requiring convergence after abrupt frequency changes without cycle slipping, where cycle slipping occurs if

$$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$3

For active PI filters the pull-in range is described as infinite, while for passive lead-lag filters it is finite; the paper supplies analytical estimates, numerical algorithms, and universal diagrams, and states that the framework extends to two-phase PLLs, Costas loops, and optical Costas loops (Kuznetsov et al., 2017).

A related large-signal formulation appears in DC/AC inverter analysis. There the lock-in domain is defined as the set of initial conditions from which the system converges to the origin without the phase variable $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$4 crossing $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$5. The method combines a quadratic Lyapunov function for a 4-dimensional current controller with a numerically constructed Lyapunov function for a 2-dimensional nonlinear PLL, forming a vector Lyapunov function and a forward-invariant set established through the comparison principle and LaSalle’s invariance principle (Ponomarev et al., 26 May 2025). The emphasis is again on phase-resolved safe convergence, not only asymptotic stability.

The mathematical scope is broader than electronic synchronization. For non-autonomous ODEs on the two-torus,

$$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$6

phase-lock areas are the level sets of the rotation number $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$7 with non-empty interiors. For $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$8, phase-lock areas occur only at rotation numbers that are integer multiples of $$A(r,c)=\frac{\Delta t}{T}\sqrt{S_0(r,c)^2+S_{-90}(r,c)^2}, \qquad \phi(r,c)=\arctan\left(\frac{-S_{-90}(r,c)}{S_0(r,c)}\right).$$9; for other analytic $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$0 with at least two nonzero non-opposite Fourier harmonics, there exists an analytic $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$1 such that phase-lock areas occur for all rational rotation numbers (Glutsyuk et al., 2015). This establishes a rigorous dynamical-systems context in which “lock-in” is a geometrical property of parameter space.

Digital frequency locking provides a high-precision computational endpoint for this tradition. In an FPGA implementation of all-phase FFT with frequency-domain unbiased phase estimation, the standard deviation of frequency fluctuations is reported to decrease from $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$2 rms in the free-running state to $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$3 rms after locking, while the Allan deviation at $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$4 decreases from $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$5 to $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$6 (Zheng et al., 16 Jun 2026). Here, phase-resolved numerical lock-in analysis functions directly as the sensing stage of a feedback loop.

4. Mechanical and acoustic implementations

In vortex-induced vibrations of a circular cylinder near a wall, the lock-in region is defined as the regime in which vibration frequency synchronizes with vortex shedding. Two-dimensional numerical simulations at $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$7, $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$8, and gap ratios $$X=\frac{2}{N}\sum_{k=1}^{N} I_k \cos(\omega_{\text{chop}} t_k+\phi_0),\qquad Y=\frac{2}{N}\sum_{k=1}^{N} I_k \sin(\omega_{\text{chop}} t_k+\phi_0),$$9 to $A=\sqrt{X^2+Y^2}$0 show that “as the gap ratio reduces, the maximum transverse vibration amplitude reduces, and the lock-in region widens” (P et al., 2022). The phase angle $A=\sqrt{X^2+Y^2}$1 between lift and transverse displacement is obtained by Fourier analysis; within lock-in, the response is divided into two branches, and the transition is marked by a “sudden jump in the phase angle $A=\sqrt{X^2+Y^2}$2 between the lift and transverse displacement from 0 to a value slightly greater than $A=\sqrt{X^2+Y^2}$3.” As $A=\sqrt{X^2+Y^2}$4 decreases, Branch II widens, Branch I narrows, and at $A=\sqrt{X^2+Y^2}$5 Branch I vanishes entirely. The study further states that periodic vortex shedding with a single “S” vortex street is observed for all gap ratios and that wall proximity affects the mean lift coefficient not only in lock-in but also in pre and post lock-in regions (P et al., 2022).

In ultrasound pulse-echo experiments, a digital FPGA-based lock-in amplifier performs phase-resolved analysis of multiple echoes in a single acquisition. A numerically controlled oscillator in the FPGA generates the RF excitation; synchronized time windows isolate echoes by time of flight; and digital quadrature demodulation computes

$A=\sqrt{X^2+Y^2}$6

followed by

$A=\sqrt{X^2+Y^2}$7

The phase shift of an echo is related to sound velocity by

$A=\sqrt{X^2+Y^2}$8

The system is reported to resolve relative changes in velocity with sensitivity “better than 1 part in $A=\sqrt{X^2+Y^2}$9,” to allow simultaneous analysis of several echoes, and to support continuous automated operation over $\phi=\arctan2(Y,X)$0 hours without glitches (Galeski et al., 2022).

Taken together, these studies show that the same phase-resolved numerical logic can serve very different purposes: identification of branch transitions in fluid-structure interaction, and precision extraction of minute propagation changes in solid-state ultrasound.

5. Optical and quantum realizations

In femtosecond nonlinear optical spectroscopy, phase-resolved numerical lock-in analysis is embedded in spectral interferometry. A chopper modulates the nonlinear signal, camera frames are time-tagged to laser triggers, and a software lock-in amplifier demodulates each spectral pixel. The measured spectrum is Fourier processed to recover the spectral phase of the weak nonlinear field, while a multidimensional lock-in space separates modulated and unmodulated interference terms. Slow interferometric drifts are corrected numerically through

$\phi=\arctan2(Y,X)$1

where the background-reference phase provides the common-mode drift estimate. The reported consequence is “an order of magnitude or more improvement” in signal-to-noise ratio over simple averaging, together with phase-resolved recovery of weak $\phi=\arctan2(Y,X)$2–$\phi=\arctan2(Y,X)$3 signals without active stabilization (Pandey et al., 7 Mar 2025).

Quantum analogs extend the lock-in principle to state evolution. The quantum double lock-in amplifier uses two quantum mixers driven by orthogonal pulse sequences, PDD and CP, which play the roles of sine and cosine references in a classical double lock-in amplifier. By combining the two outputs, the complete characteristics of a target alternating signal—amplitude, frequency, and phase—can be recovered even when the initial phase is unknown. Numerical calculations are reported to show robustness against finite pulse length and stochastic noise, including signal-to-noise ratios as low as $\phi=\arctan2(Y,X)$4 dB in weak-signal simulations, and the protocol is illustrated with a five-level double-$\phi=\arctan2(Y,X)$5 coherent population trapping system in $\phi=\arctan2(Y,X)$6 atoms (Chen et al., 2023).

A different quantum route combines phase-locked amplification with spin squeezing. In an atomic ensemble with Hamiltonian

$\phi=\arctan2(Y,X)$7

one-axis twisting supplies squeezing, while a train of $\phi=\arctan2(Y,X)$8-pulses performs phase locking to the target signal. The study derives optimal $\phi=\arctan2(Y,X)$9-pulse and $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$0-pulse schemes and states that the resulting phase sensitivity is enhanced and the usable detection window for phase locking is widened (Zhang et al., 3 Jul 2025).

Rydberg atom sensing replaces the external local oscillator with a system-internal quantum reference generated by a closed-loop interferometric scheme. The loop enforces

$$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$1

and detuning one RF field produces an intermediate frequency suitable for lock-in detection of the EIT response. The method is reported to provide full $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$2 phase resolution, to demodulate a four phase-state QPSK signal at $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$3 Hz symbol rate, and to achieve sensitivity within a factor of $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$4–$$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$5 of a traditional LO-based Rydberg mixer (Berweger et al., 2022).

6. Recurring constraints, limitations, and interpretive themes

Several constraints recur across these literatures. In Lock-In Thermography, Technique I “must be empirically calibrated” and “cannot always be uniquely or robustly performed,” while Technique II depends on the chosen database resolution, with coarser resolution trading accuracy for efficiency (Tian et al., 2019). In phase-based thermal conductivity measurements, radial region selection is critical, the fit is typically localized to $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$6 mm, and the sensitivity analysis indicates that $$T(r,\omega)=X(r,\omega)+iY(r,\omega), \qquad \theta(r,\omega)=\tan^{-1}\!\left(\frac{Y(r,\omega)}{X(r,\omega)}\right),$$7 and heat capacity cannot be simultaneously fit because the phase sensitivity to heat capacity is opposite in sign; finite-size and edge effects must also be avoided by selecting appropriate modulation frequency and region of interest (Scott et al., 10 Jun 2026).

The literature also distinguishes sharply between amplitude information and phase information. In vortex-induced vibrations, the “phase-jump” is described as the primary indicator of branch transition and is “not always apparent in amplitude plots for small gap ratios” (P et al., 2022). In PLL and inverter studies, ordinary stability is insufficient when the engineering requirement is synchronization without cycle slipping; separatrix geometry, forward-invariant sets, and explicit phase boundaries are therefore central objects of computation (Aleksandrov et al., 2016, Ponomarev et al., 26 May 2025).

Reference construction is another common concern. Software lock-in spectroscopy corrects slow interferometric drifts by measuring both signal-reference and background-reference phases within the same multidimensional acquisition (Pandey et al., 7 Mar 2025). The quantum double lock-in amplifier resolves unknown initial phase by using orthogonal quantum mixers (Chen et al., 2023). Rydberg closed-loop interferometry creates an internal LO-equivalent reference from the atomic structure itself (Berweger et al., 2022). APFFT-based frequency locking addresses phase bias introduced by ordinary FFT processing when the signal frequency does not coincide with an FFT bin (Zheng et al., 16 Jun 2026).

These recurring features suggest that phase-resolved numerical lock-in analysis is best understood not as a single algorithm but as a family of synchronous estimation and phase-space computation methods. Across imaging, spectroscopy, mechanics, synchronization theory, and quantum sensing, its reliability depends on reference integrity, numerical treatment of ambiguity and drift, and the adequacy of the model that links phase-resolved observables to the underlying physical quantity.