Analytically exploiting noise correlations inside the feedback loop to improve locked-oscillator performance
AAnalytically exploiting noise correlations inside the feedback loop to improvelocked-oscillator performance
J. Sastrawan, C. Jones, I. Akhalwaya, H. Uys, and M.J. Biercuk ∗ ARC Centre for Engineered Quantum Systems, School of Physics, The University of Sydney, NSW 2006 Australia and
National Measurement Institute, West Lindfield, NSW 2070 Australia National Laser Centre, Council for Scientific and Industrial Research, Pretoria, South Africa (Dated: October 15, 2018)We introduce concepts from optimal estimation to the stabilization of precision frequency stan-dards limited by noisy local oscillators. We develop a theoretical framework casting various measuresfor frequency standard variance in terms of frequency-domain transfer functions, capturing the ef-fects of feedback stabilization via a time-series of Ramsey measurements. Using this frameworkwe introduce a novel optimized hybrid predictive feedforward measurement protocol which employsresults from multiple past measurements and transfer-function-based calculations of measurementcovariance to improve the accuracy of corrections within the feedback loop. In the presence of com-mon non-Markovian noise processes these measurements will be correlated in a calculable manner,providing a means to capture the stochastic evolution of the LO frequency during the measurementcycle. We present analytic calculations and numerical simulations of oscillator performance undercompeting feedback schemes and demonstrate benefits in both correction accuracy and long-termoscillator stability using hybrid feedforward. Simulations verify that in the presence of uncompen-sated dead time and noise with significant spectral weight near the inverse cycle time predictivefeedforward outperforms traditional feedback, providing a path towards developing a new class ofstabilization “software” routines for frequency standards limited by noisy local oscillators.
High-performance passive frequency standards play amajor role in technological applications such as networksynchronization and GPS [1] as well as many fields ofphysical inquiry, including radioastronomy (very-long-baseline interferometry) [2], tests of general relativity [3],and particle physics [4]. Atomic clocks exploiting thestability of Cs [5–8] or other atomic references [9–13]to stabilize an oscillator are known as the most precisetimekeeping devices available, but constant performancegains are sought for technical and scientific applications.In many settings, such as miniaturized deployable fre-quency standards or in GPS-denied environments, a ma-jor performance limitation aries from the quality of thelocal oscillator (LO) that probes and is locked to theatomic transition. The LO frequency may evolve ran-domly in time due to intrinsic noise processes in the un-derlying hardware [10, 11], leading to time-varying devi-ations of the LO frequency from that of the stable atomicreference. These instabilities are partially compensatedthrough use of a feedback protocol designed to transferthe stability of the reference to the LO, but their effectscannot be mitigated completely.Early work characterizing the so-called Dick effect [14]demonstrated that no matter how good the reference be-comes, LO noise will still produce residual instabilities inthe locked LO (LLO) through the feedback protocol it-self. The dominant mechanism for this is evolution of theLO’s frequency on timescales rapid compared with theshortest measurement and feedback cycle. Major con-tributors to this phenomenology relate to the presence ∗ [email protected] of uncompensated LO evolution during initialization andreadout stages of the measurement cycle (dead time), aswell as aliasing of LO noise at harmonics of the feedback-loop period – the Dick effect [14–16]. Accordingly, signif-icant research focus in the frequency standards commu-nity has been placed on improving LO performance, usinge.g. ultra-low-phase-noise cryogenic sapphire oscillatorsor similar [17, 18], with concomitant increases in hard-ware infrastructure requirements and complexity. Otherapproaches to mitigating the impact of LO instabilitiesinvolve significant modification of the relevant referencehardware, for instance employing multiple atomic refer-ences [10, 19].In this Manuscript we devise and analyze a methodby which both the accuracy of the LLO relative to theatomic reference, and the stability of the composite pas-sive frequency standard, can be improved without theneed for hardware modification. We develop new an-alytic tools casting time-domain statistical measures offrequency-standard performance in terms of analyticallycalculable transfer functions [20], exploiting recent re-lated work in quantum information [21–24]. This ap-proach reveals opportunities to exploit non-Markovianityin the dynamics of LO frequency fluctuations in order toimprove feedback stabilization by bringing optimal esti-mation inside the feedback loop of the LO.Our method expresses the properties of the LLO interms of the statistics of the unlocked LO at differenttimes as well as correlations between those measure-ments. We present the relevant transfer functions fortime-series measurements of arbitrary-duration Ramseymeasurements, and introduce the pair-covariance trans-fer function explicitly capturing correlations betweenmeasurement outcomes at different times. Thus, given a r X i v : . [ qu a n t - ph ] A p r statistical knowledge of the LO noise characteristics, wecraft a new form of hybrid feedforward stabilization in-corporating the results of an arbitrary number of pastmeasurements with variable duration to calculate an im-proved correction to the LO. This approach shares con-cepts with techniques of optimal estimation [25] com-monly used in engineering to predict the evolution of adynamical system – here the noisy LO.In cases where dead time is significant and there issubstantial uncompensated LO evolution, we use numer-ical simulations to show that this approach allows cor-rections of improved accuracy to be applied to the LO.Simulations demonstrate that long-term stability of theLLO is improved through a moving-average correctionscheme, where corrections are made based on weight-ing values determined analytically in the same hybridfeedforward approach. The method described here isa technology-independent software-oriented approach toimproving the performance of frequency standards de-rived from locked local oscillators. It may be freely usedin conjunction with hardware modifications targeted atreducing the same limitations identified, such as inter-leaving the cycles of two clocks to reduce dead time[10, 19].The remainder of this manuscript is organized as fol-lows. In Section I we provide an analytic descriptionof the deleterious effects of LO noise on frequency stan-dards, introducing the relevant metrics for performanceof interest. This includes presentation of novel ana-lytic expressions explicitly capturing the effects of feed-back stabilization on the aggregate system performancethrough a recursive formulation. Section II demonstrateshow to convert these time-domain statistical measures offrequency-standard performance to the Fourier domain,introducing both transfer functions for individual mea-surements and the pair-covariance transfer function cap-turing the correlations between arbitrary-duration Ram-sey measurements conducted at arbitrary times. We thenexploit these tools in Sec. III in order to devise a newhybrid-feedforward correction scheme similar in spirit toconcepts from optimal estimation in order to maximizethe accuracy of corrections applied to the LO. We demon-strate improvements in correction accuracy and LLO sta-bility via this approach using numerical simulations withrealistic LO noise power spectra. Finally, we concludewith a summary and discussion in Sec. IV. I. THE EFFECT OF LOCAL OSCILLATORNOISE ON FREQUENCY STANDARDSTABILITY
Our primary objective is to suppress the impact ofLO noise on the ultimate performance of the locked
LO,which is stabilized to an (in general atomic) reference.Accordingly, throughout this analysis we do not considersystematic shifts or uncertainties in the reference and ex-plicitly assume that the reference is perfect. -50050 y ( t ) ( a . u . ) Correction Cycle g ( t ) g ( t ) t (a.u.) C ( n ) k C ( n ) k +1 T c T R T D t s t s t s t s t s t e t e t e t e t e C ( n ) k +2 LO LLO y ( LO ) k y ( LLO ) k y ( LO ) ( t ) y ( LLO ) ( t ) FeedbackHybrid Feedforward (b)(c) M ea s u r e m en t D ead T i m e FIG. 1. Effect of LO noise on the performance of a locked os-cillator. Simulated evolution for a noisy LO, unlocked (black)and locked with traditional feedback (red). The dotted hor-izontal bars indicate the measurement outcomes ( samples )over each cycle, ¯ y k , which are applied as correction at the endof the cycle, indicated by the bent arrow in the first cycle.Measurement period of duration T R (white background) isfollowed by dead time with duration T D (grey background).Total cycle time T c = T D + T R , and here we represent a 50%duty factor, d . Undetected evolution of the LO during thedead time leads corrections to incompletely cancel frequencyoffsets at the time of correction. The arrows on the far rightschematically indicate how locking reduces the variance of y ( t ) though it does not eliminate it. We represent the fractional frequency offset of theLO relative to an ontologically perfect reference y ( t ) ≡ ( ν ( t ) − ν ) /ν , where ν is the reference frequency and ν ( t ) is the LO frequency. This limit provides a reason-able approximation to the performance of many deploy-able frequency standards where LO stability is far worsethan that of the associated atomic reference. A. Time-domain description of Ramseymeasurements and feedback stabilization
In such a setting, Ramsey spectroscopy provides ameans to determine the fractional frequency offset of theLO relative to the reference over a period T R . Point-likerealisations of the stochastic process y ( t ) cannot be ob-tained experimentally; instead, the LO frequency errorproduces integrated samples , denoted ¯ y k and indexed intime by k : ¯ y k ≡ T ( k ) R (cid:90) t ek t sk y ( t ) g ( t − t sk ) dt (1)where T ( k ) R ≡ t ek − t sk , [ t sk , t ek ] is the time interval overwhich the k th sample is taken, and g ( t ) is a sensitivityfunction capturing the extent to which LO fluctuationsat some instant t contribute to the measured outcomefor that sample [26]. The range of g ( t ) is [0 ,
1] and itsdomain is t ∈ [0 , T ( k ) R ]. The ideal case is the rectangularwindow case, where g ( t ) = (cid:40) t ∈ [0 , T ( k ) R ]0 otherwise (2)in which case ¯ y k reduces to the time-average of y ( t ) overthe interval [ t sk , t ek ].In traditional feedback stabilization, the samples, ¯ y k ,are used to determine corrections to be applied to theLO in order to reduce frequency differences from the ref-erence (Fig. 1). Consider the trajectory of the same fre-quency noise realisation y ( t ) in the cases of no correction, y LO ( t ) and correction, y LLO (t). The relation betweenthese two cases of y ( t ) is y LLO ( t ) = y LO ( t ) + n (cid:88) k =1 C k (3)where C k refers to the value of the k th frequency correc-tion applied to the LO, n of which have occurred beforetime t .Under traditional feedback stabilization, each correc-tion is directly proportional to the immediately preced-ing measurement outcome: C k = w k ¯ y LLOk , where w k iscorrection gain. Since ¯ y LLOk is calculated by convolv-ing y LLO ( t ) with a sensitivity function pertaining to themeasurement parameters, (3) is a recursive equation ingeneral. It is possible to cancel all but one of the re-cursive terms by setting the correction gain equal to theinverse of the average sensitivity ¯ g k ≡ (cid:82) T ( k ) R g ( t ) /T ( k ) R dt of the preceding measurement, i.e. w k = − ¯ g − k , wherethe minus sign indicates negative feedback. With thisconstraint we can write¯ y LLOk = ¯ y LOk − ¯ g k ¯ g k − ¯ y LOk − (4)and for a Ramsey interrogation and measurement withnegligibly short pulses, ¯ g k = 1. Applying feedback cor-rections sequentially after measurements is able to effec-tively reduce y ( t ) over many cycles, improving long-termstability.In the limit of a static offset, a single (perfect) cor-rection will set the frequency offset error of the LLO tozero; however, such perfect correction is in general notachieved. The primary reason for this in the limit of per-fect measurements and corrections is dynamic evolution of the LO on timescales rapid compared to the measure-ments which cannot be fully compensated by the feed-back loop.In Fig. 1 we demonstrate how evolution of the LO fre-quency during T R leads the feedback protocol to incom-pletely correct the offset y ( t ). From the formalism pre-sented above we see that incomplete feedback arises be-cause the corrections are based only on the average valueof the frequency offset as measured over the k th period, ¯ y k (horizontal solid lines in Fig. 1), rather than the in-stantaneous value of the LO frequency offset at the timeof correction (here the end of a cycle) which cannot beknown . The difference between these two values leadsto incomplete compensation of time-varying frequencyoffsets, and hence residual fractional instability in thequantity y ( LLO ) ( t ).The impact of these effects on the ultimate stability ofthe LLO is exacerbated in circumstances where there isnonzero dead time , T D , during which the LO may evolve,but this evolution is not captured by a measurement.Dead time arises due to e.g. the need to reinitializethe reference between measurements, or perform classicalprocessing of the measurement outcome before a correc-tion can be applied.The net impact of this uncompensated evolution is areduction in the long-term stability of the locked local os-cillator. We now move on to describe the relevant quan-titative metrics for LLO variance in both free-runningand feedback-locked settings. B. Measures of frequency standard stability forunlocked and locked LOs
The performance of the frequency standard is statis-tically characterized by various time-domain measurescapturing the evolution of LO frequency as a functionof time.The variance of ¯ y k , denoted σ y ( k ) and often called truevariance [26] is, σ y ( k ) = Var[¯ y k ] = (cid:0) (cid:104) ¯ y k (cid:105) − (cid:104) ¯ y k (cid:105) (cid:1) → E[¯ y k ] (5)= E (cid:20)(cid:18) T ( k ) R (cid:90) t ek t sk y ( t ) g ( t − t sk ) dt (cid:19) (cid:21) (6)where in the first line we assume that the true varianceis simply equal to the expected value of ¯ y k , since y ( t ) isassumed to be a zero-mean process. The true variancecaptures the spread of measurement outcomes due to dif-ferent noise realizations in a single timestep. However,in a measurement context one does not have immediateaccess to an infinite ensemble of noise realizations, butrather a single series of measurement outcomes conductedsequentially over a single noise realization. As a resultwe rely on a measure more conducive to this setting, the sample variance σ y [ N ] = 1 N − N (cid:88) k =1 (¯ y k − N N (cid:88) l =1 ¯ y l ) (7)for N sequential finite-duration measurements { ¯ y k } [26].In this work we will rely on such measures of frequencystability, rather than the more commonly employed Al-lan variance, in line with recent experiments [27]. TheAllan variance is calculated by finding the variance ofthe difference between consecutive pairs of measurementoutcomes: A σ y ( y ) = 12 (cid:104) (¯ y k +1 − ¯ y k ) (cid:105) (8)where ¯ y k is the k th measurement outcome and (cid:104)· · · (cid:105) mayindicate a time average or an ensemble average, depend-ing on whether y ( t ) is assumed to be ergodic. Our deci-sion to avoid the Allan variance is deliberate, as its form– effectively a moving average – specifically masks the ef-fect of LO noise components with long correlation times.In fact the Allan variance is employed by the communityin part because it does not diverge at long integrationtimes τ due to LO drifts, as would the sample or true variance [26, 28–30]. In the limit where the stability ofa frequency reference is dominated by LO noise (and thereference can be treated as perfect) this approach givesphysically meaningful results.The standard measures for oscillator performance con-sider either a free-running LO or provide a means onlyto statistically characterize measurement outcomes underblack-box conditions. We may derive explicit analyticforms for different measurements of variance in the pres-ence of feedback locking in order to provide insights intoopportunities to improve net LLO performance throughmodification of the stabilization protocol.We write time-domain expressions for variance usingthe relevant definitions provided above and the link be-tween corrections in feedback and the history of theLLO’s evolution. For the true variance we substituteEq. 4 to find σ yLLO ( k ) = Var[¯ y LLOk ] (9)= σ yLO ( k ) + (cid:18) ¯ g k ¯ g k − (cid:19) σ yLO ( k − − g k ¯ g k − σ (¯ y LOk − , ¯ y LOk ) (10)and calculate the expected value of the LLO sample variance in a similar manner using Eq. 3E[ σ yLLO [ N ]] = 1 N − N (cid:88) k =1 (cid:26)(cid:18) σ yLO ( k ) + ¯ g k k − (cid:88) r =1 k − (cid:88) s =1 σ ( C r , C s ) − g k k − (cid:88) u =1 σ (¯ y LOk , C u ) (cid:19) + 1 N N (cid:88) p =1 N (cid:88) q =1 (cid:18) σ (¯ y LOp , ¯ y LOq ) + ¯ g p ¯ g q p − (cid:88) w =1 q − (cid:88) x =1 σ ( C x , C y ) (cid:19) − N N (cid:88) l =1 (cid:18) σ (¯ y LOk , ¯ y LOl ) + ¯ g k ¯ g l k − (cid:88) y =1 l − (cid:88) z =1 σ ( C y , C z ) (cid:19)(cid:27) (11)We see that the characteristics of the locked LO canbe expressed in terms of the unlocked LO and the co-variance covariance between two quantities, σ ( x, y ), cap-turing correlations between them. This may include thecovariance of different measurement outcomes on the LO,or different corrections applied to the LO. It is this ob-servation – that we may express relevant statistical quan-tities surrounding the performance of locked local oscil-lators in terms of measurement covariances – that willprovide a path towards the development of new stabi-lization routines exploiting temporal correlations in theLO noise (and hence measurement outcomes). II. PERFORMANCE MEASURES FORFREQUENCY STANDARDS IN THE FOURIERDOMAIN
We require an efficient theoretical framework in whichto capture these effects, and hence transition to the fre-quency domain, making use of the power spectral density of the LO, S y ( ω ), in order to characterize average perfor-mance over a hypothetical statistical ensemble. In thisdescription residual LLO instability persists because thefeedback is insensitive to LO noise at high frequenciesrelative to the inverse measurement time. Additional in-stability due to the Dick effect comes from aliasing ofnoise at harmonics of the loop bandwidth.We may analytically calculate the effects of measure-ment, dead time, and the feedback protocol itself on fre-quency standard performance in the frequency domain asfollows. Defining a normalised, time-reversed sensitivityfunction ¯ g ( t mk − t ) = g ( t − t sk ) /T ( k ) R , where g ( t ) is assumedto be time-reversal symmetric about t mk , the midpoint of[ t sk , t ek ], we can express, for instance, the true varianceas a convolution σ y ( k ) = E (cid:20)(cid:18) (cid:82) ∞−∞ y ( t )¯ g ( t mk − t ) dt (cid:19) (cid:21) .Expanding this expression gives σ y ( k ) = (cid:90) ∞−∞ (cid:90) ∞−∞ E[ y ( t ) y ( t (cid:48) )]¯ g ( t mk − t )¯ g ( t mk − t (cid:48) ) dt (cid:48) dt (12)= (cid:90) ∞−∞ (cid:90) ∞−∞ R T Syy (∆ t )¯ g ( t mk − t )¯ g ( t mk − t (cid:48) ) dt (cid:48) dt (13)where R T Syy (∆ t ) is the two-sided autocorrelation functionand ∆ t ≡ t (cid:48) − t . Using the Wiener-Khinchin theorem wewrite R T Syy (∆ t ) = F − { S T Syy ( ω ) } , relating the autocor-relation function to the Fourier transform of the powerspectral density of the LO noise. Defining the Fouriertransform of ¯ g ( t mk − t ): G k ( ω ) ≡ (cid:90) ∞−∞ ¯ g ( t mk − t ) e iωt dt (14)We may then express the true variance σ y ( k ) = 12 π (cid:90) ∞ S y ( ω ) | G k ( ω ) | dω (15)where the substitution of the one-sided PSD S y ( ω ) ispossible because | G k ( ω ) | is even. This result is similarto the convolution theorem, which states that F{ f (cid:63)g } = F{ f } · F{ g } , where (cid:63) denotes a convolution and f and g are Fourier-invertible functions.Here | G k ( ω ) | is called the transfer function for the k th sample, describing the spectral properties of themeasurement protocol itself. For measurements per-formed using Ramsey interrogation with π/ | G k ( ω ) | =(sin ( ωT ( k ) R / / ( ωT ( k ) R / . This framework has recentlyseen broad adoption in the quantum information com-munity where time-varying dephasing noise is a majorconcern for the stability of quantum bits [21–24, 31–34].Recalling that statistical measures of LLO variancerely not only on expressions for the true variance overnoise ensembles, but also of covariances between mea-surements or corrections, we must equivalently expressthe covariance in terms of transfer functions. Using theidentity σ ( A ± B ) = σ ( A ) + σ ( B ) ± σ ( A, B ), we de-fine a sum and a difference sensitivity function: g + k,l ( t )and g − k,l ( t ), with respect to two measurements indexed k and l . These expressions are general functions of timewith two regions of high sensitivity corresponding to theindividual measurement periods. g ± k,l ( t ) ≡ g ( t − t sk ), for t ∈ [ t sk , t ek ] ± g ( t − t sl ), for t ∈ [ t sl , t el ]0, otherwise (16) These time-domain sum and difference sensitivity func-tions have their corresponding frequency-domain transferfunctions, defined as their Fourier transforms normalisedby T ( k,l ) R : G ± k,l ( ω ) ≡ (cid:90) ∞−∞ (cid:18) g ( t mk − t ) T ( k ) R ± g ( t ml − t ) T ( l ) R (cid:19) e iωt dt (17)Substituting this and the form of the true variance (15)into the variance identity above and rearranging termsgives the covariance of the two measurement outcomes σ (¯ y k , ¯ y l ) = 12 π (cid:90) ∞ S y ( ω )4 (cid:18) (cid:12)(cid:12)(cid:12) G + k,l ( ω ) (cid:12)(cid:12)(cid:12) − (cid:12)(cid:12)(cid:12) G − k,l ( ω ) (cid:12)(cid:12)(cid:12) (cid:19) dω (18) ≡ π (cid:90) ∞ S y ( ω ) G k,l ( ω ) dω (19)whereby G k,l ( ω ) is defined to be the pair covariancetransfer function . For the case of flat-top Ramsey mea-surements over the intervals [ t sk,l , t ek,l ] this term takes theform G k,l ( ω ) = ( ω T ( k ) R T ( l ) R ) − (cid:104) cos ( ω ( t sl − t sk )) + cos ( ω ( t el − t ek )) − cos ( ω ( t el − t sk )) − cos ( ω ( t sl − t ek )) (cid:105) . (20)This is a generalization of the transfer function previouslyderived for the special case of periodic, equal-durationRamsey interrogations [26, 29], and allows effective es-timation of y ( t ) for any t and for any set of measuredsamples ¯y k .We thus see that this approach allows expression oftime-domain LO variances as overlap integrals between S y ( ω ) and the transfer functions capturing the effects ofthe measurement and feedback protocol, including cor-relations between measurements or corrections in time.Through this formalism we may incorporate arbitrarymeasurement protocols (e.g. arbitrary and dynamicRamsey periods and dead times): the underlying physicsof e.g. changing linewidth of the measurement is ex-plicitly captured through the form and implicit time-dependence of the transfer function used to characterizethe measurement protocol. III. EXPLOITING NOISE CORRELATIONS TOIMPROVE FEEDBACK STABILIZATION
Recasting variance metrics for the stability of LOs interms of transfer functions is particularly powerful be-cause it provides a path to craft new measurement feed-back protocols designed to reduce residual variance mea-sures for the LLO by modifying the protocol’s spectralresponse. Our key insight is that the non-Markovianityof dominant noise processes in typical LOs – capturedthrough the low-frequency bias in S y ( ω ) [26, 29] – impliesthe presence of temporal correlations in y ( t ) that may beexploited to improve feedback stabilization. These corre-lations are captured in the set of measurement outcomes ¯y k ; accordingly future evolution of y ( t ) may be predictedbased on a past set of measurements within ¯y k , so long asthe past measurements and point of prediction fall withinthe characteristic correlation time for the LO noise givenby S y ( ω ). This approach provides a direct means to ac-count for LO evolution that is normally not compensatedduring dead time in the measurement process. A. Optimal estimator for corrections
The formal basis of our analytic approach, in summary,is to calculate a covariance matrix in the frequency do-main via transfer functions to capture the relative corre-lations between sequential measurement outcomes of anLLO, and use this matrix to derive a linear predictor ofthe LLO frequency offset at the moment of correction.Under appropriate conditions this predictor provides acorrection with higher accuracy than that derived from asingle measurement, allowing us to improve the ultimateperformance of the LLO. Since the predictor is found us-ing information from previous measurements (feedback)and a priori statistical knowledge of the LO noise to pre-dict the evolution of the LO (feedforward), we call thescheme hybrid feedforward .This approach shares common objectives with applica-tion of optimal control techniques such as Kalman filter-ing in the production of composite frequency standardsfrom an ensemble of physical clocks [35], or in compen-sating for deterministic frequency shifts due to e.g. agingor changes in the ambient temperature of a clock [36, 37].The primary advance of this work is the insight that stochastic evolution of the LO can be predicted and com-pensated using optimal control protocols inside the feed-back loop .In hybrid feedforward, results from a set of n pastmeasurements are linearly combined with weighting co-efficients c k optimized such that the k th correction, C k ,provides maximum correlation to y ( t ck ) at the instant ofcorrection t ck (Fig. 1c). Assuming that the LO noise isGaussian, the optimal least minimum mean squares es-timator (MMSE) is linear, and the optimal value of thecorrection is given by C k = c k · ¯y k : the dot product ofa set of correlation coefficients c k derived from knowl-edge of S y ( ω ) and a set of n past measured samples, ¯y k = { ¯ y k, , · · · , ¯ y k,n } . We define an ( n + 1) × ( n + 1) co-variance matrix where the ( n + 1)th term represents anideal zero-duration sample at t ck and in the second line we write the covariance matrix in block form:Σ k ≡ σ (¯ y k, , ¯ y k, ) · · · σ (¯ y k, , y ( t ck )) σ (¯ y k, , ¯ y k, ) · · · σ (¯ y k, , y ( t ck )) · · · · · · · · · σ ( y ( t ck ) , ¯ y k, ) · · · σ ( y ( t ck ) , y ( t ck )) (21) ≡ (cid:20) M k F k F T k σ ( y ( t ck ) , y ( t ck )) (cid:21) . (22)In this form the matrix M k describes correlations be-tween measurement outcomes while the vector F k de-scribes correlations between each measurement and theLLO at the time of correction. The MMSE optimalitycondition is then fulfilled for c k = F k (cid:113) F T k M k F k w k π (cid:90) ∞ S y ( ω ) dω (23)where w k is an overall correction gain. The covariancematrix elements are calculated as defined above in termsof the LO noise power spectrum.In the practical setting of a frequency standard exper-iment, we wish to improve both the accuracy of each cor-rection, by maximising the correlation between C k and y ( t ck ), and the long-term stability of the LLO output,captured by the metrics of frequency variance, samplevariance, and Allan variance.Although the LLO frequency variance under hybridfeedforward for more than a single cycle cannot be ex-pressed in a closed non-recursive form, a consideration ofa single cycle can provide a value for (cid:104) y LLO ( t ck ) (cid:105) in termsof covariance matrix elements. This in turn provides g ( t ) t (a.u.) C ( n ) k C ( n ) k +1 T R T D t s t s t s t s t s t e t e t e t e t e C ( n ) k +2 Hybrid Feedforward: Overlapping corrections
FIG. 2. Schematic diagram of hybrid feedforward withan example protocol using n = 3. Start and end times ofmeasurements are defined arbitrarily permitting non-uniform-duration measurements, although measurements are illus-trated as uniform for clarity. Corrections C ( n =3) k are appliedin either non-overlapping blocks of three measurements or asa moving average (depicted here). In the latter case, the co-variance matrix must be recalculated to correctly account forany variations in measurement duration. Dashed red arrowsindicate the first corrections performed without full calcula-tion of the covariance matrix. This effect vanishes for k > n . a metric for the correction accuracy for hybrid feedfor-ward, defined as the extent to which a correction brings y LLO ( t ) → t = t ck A k ≡ (cid:104) y LO ( t ck ) (cid:105)(cid:104) y LLO ( t ck ) (cid:105) (24)= (cid:18) w k − w k | F k | (cid:113) F Tk M k F k (cid:19) − (25)We can gain insights into the performance of the cor-rection protocol by considering limiting cases. For in-stance, in the limit of white noise with negligible cor-relations, M k → I ,the identity matrix. In this limitthe rightmost term in Eq. 24 reduces to w k | F k | , whichis small (there are negligible correlations between mea-surement outcomes and y ( t ck )). In this limit, accuracy A k → / (1 + w k ), and is maximized by setting w k = 0(not performing feedback at all) as corrections are un-correlated with y ( t ck ). By contrast with perfect corre-lations all elements of the covariance matrix take valueunity. Standard feedback works perfectly by selectingunity gain and selecting the number of measurements tobe combined, n = 1, to correct based on a single mea-surement.In intermediate regimes the ensemble-averaged accu-racy of the hybrid feedforward correction is determinedin by a balance of covariance between elements of ¯y k andcovariance between ¯y k and y ( t ck ), the LO noise at thetime of correction. Achieving correction which improvesLLO variance requires setting the term in parentheses inEq. 24 to less than unity. This in turn places a conditionon the correlations in the system (cid:113) F Tk M k F k < | F k | w k (26)We can interpret the effect of M k as an effective rotationmatrix, reducing the magnitude of the left-hand side ofthe expression above by effectively maximizing the “an-gle” between M k F k and F k . While it is unphysical toreduce this to zero based on the limiting cases discussedabove, it is possible to appropriately select k , based oncharacteristics of S y ( ω ) in order to improve correctionaccuracy.In all slaved frequency standards we rely on repeatedmeasurements and corrections to provide long-term sta-bility , a measure of how the output frequency of the LLOdeviates from its mean value over time. We study thisby calculating the sample variance of a time-sequenceof measurement outcomes averaged over an ensemble ofnoise realizations, (cid:104) σ y [ N ] (cid:105) . A “moving average” style ofhybrid feedforward provides improved long-term stabil-ity, as the correction C k will depend on the set of mea-surement outcomes ¯y k = { ¯ y k − n +1 , · · · , ¯ y k } , among whichprevious corrections have been interleaved, as illustratedin Fig. 2. In this case the covariance matrix must be updated to reflect the action of each correction. See Ap-pendix for a detailed form of the Sample Variance in thecase of this form of stabilization. B. Numerical Simulations
In order to test the general performance of hybrid feed-forward in different regimes we perform numerical simu-lations of noisy LOs with user-defined statistical proper-ties, characterized by S y ( ω ). We produce a fixed numberof LO realizations in the time domain and then use theseto calculate measures such as the sample variance over asequence of “measurement” outcomes with user-definedRamsey measurement times, dead times, and the like. Inthese calculations we may assume that the LO is free run-ning, experiencing standard feedback, or employing hy-brid feedforward, and then take an ensemble average overLO noise realizations. Our calculations include variousnoise power spectra, with tunable high-frequency cutoffs,including common ‘flicker frequency’ ( S y ( ω ) ∝ /ω ), and‘random walk frequency’ ( S y ( ω ) ∝ /ω ) noise, as appro-priate for experiments incorporating realistic LOs,Tunability in the hybrid feedforward protocol comesfrom the selection of n , in determining { C k } as well as theselected Ramsey periods, permitting an operator to sam-ple different parts of S y ( ω ). As an example, we fix ourpredictor to consider n = 2 sequential measurements andpermit the Ramsey durations to be varied as optimizationparameters. A Nelder-Mead simplex optimization overthe measurement durations finds that a hybrid feedfor-ward protocol consisting of a long measurement periodfollowed by a short period maximizes correction accu-racy (Fig. 3). This structure ensures that low-frequencycomponents of S y ( ω ) are sampled but the measurementsampling the highest frequency noise contributions aremaximally correlated with y ( t ck ). With S y ( ω ) ∝ /ω and S y ( ω ) ∝ /ω we observe increased accuracy under hy-brid feedforward while the rapid fluctuations in y ( t ) aris-ing from a white power spectrum mitigate the benefits ofhybrid feedforward, as expected. In the parameter rangeswe have studied numerically we find that correction ac-curacy is maximized for n = 2 to 3, with diminishingperformance for larger n . Again, this is determined bythe relevant correlation time of the LO noise.In Fig. 4b we demonstrate the resulting normalizedimprovement in (cid:104) σ y [ N ] (cid:105) up to N = 100 measurements,calculated using feedback and hybrid feedforward with n = 2, and assuming uniform T R . We observe clearimprovement (reduction) in (cid:104) σ y [ N ] (cid:105) through the hybridfeedforward approach, with benefits of order 5 − (cid:104) σ y [ N ] (cid:105) relative performance improvement over stan-dard measurement feedback. We present data for differ-ent functional forms of S y ( ω ), including low-frequencydominated flicker noise ( ∝ /ω ), and power spectra( ∝ /ω / ) with more significant noise near T − c . Thebenefits of our approach are most significant in the longterm when high-frequency noise reduces the efficacy of g ( t ) N o r m a li z ed A cc u r a cy Ratio of Ramsey Measurements
White noise 1/ ω noise 1/ ω noise FIG. 3. Calculated correction accuracy of the first correctionfor hybrid feedforward normalized to feedback (accuracy =1), under different forms of S y ( ω ) as a function of the ratioof Ramsey periods between the two measurements employedin constructing C (2) k . Correction accuracy for feedback is cal-culated assuming the minimum Ramsey time; thus for theratio of Ramsey measurements taking value unity on the x -axis, the hybrid feedforward scheme takes twice as long asfeedback. Inset: depiction of the form of C (2) k used in hy-brid feedforward, depicting the “slower” measurement beingperformed first. standard feedback. Notably, because of well known re-lationships between LO phase noise and LO frequencynoise [28], significant high-frequency weight in S y ( ω ) iscommonly encountered.In Fig. 4c we calculate the expectation value of thesample variance at a fixed value of N = 20 for a LLO sta-bilized using either traditional feedback or hybrid feed-forward. The sample variances are normalized by that forthe free-running LO, meaning that values of this metricless than unity demonstrate improvement due to stabi-lization, and smaller values indicate better stabilization.On the horizontal axis we vary the duty factor d , definedas the ratio of the interrogation time to total cycle time: d ≡ T R /T c from 1% to unity (no dead time), and wecompare S y ( ω ) ∝ /ω and S y ( ω ) ∝ /ω / . These powerspectra are conservative but inspired by typical LO phasenoise specifications weighted to enhanced high-frequencycontent due to the conversion between phase and fre-quency instability [28].This improvement provided by hybrid feedforward ismost marked for low duty factor d . As d → S y ( ω ) withfrequency weight near T c , feedback efficacy diminishesdue to uncompensated evolution of the LO during thedead time.In this regime knowledge of correlations in the noiseallows hybrid feedforward to provide metrologically sig-nificant gains in stability relative to traditional feedback.In Fig. 4d we further demonstrate that in the presence ω S a m p l e V a r i an c e ( a . u . ) Number of Samples
Duty Cycle (T R /T C ) ω ω a) b) c) d) LO FB HFF ω ω FB, 1/ ω HFF, 1/ ω FB, 1/ ω + spurs HFF, 1/ ω + spurs 〈 σ y [ ] 〉 ( H FF / FB ) / 〈 σ y [ ] 〉 ( L O ) 〈 σ y [ N ] 〉 ( H FF ) / 〈 σ y [ N ] 〉 ( FB ) FIG. 4. (a, b) Calculated sample variance for an unlockedLO, feedback, and hybrid feedforward, as a function of mea-surement number N , for different power spectra (indicatedon graphs). Calculations assume S y ( ω ) ∝ /ω , with a high-frequency cutoff ω c / π = 100 /T c and S y ( ω ) ∝ /ω / witha cutoff frequency ω c / π = 1 /T c , demonstrating the impor-tance of high-f noise near ω/ π = T − c . PSDs with different ω -dependences are normalised to have the same value at ω low =1 / T c . (c) Normalized sample variance data from panels (a)and (b) presented as the ratio of (cid:104) σ y [ N ] (cid:105) ( HFF ) / (cid:104) σ y [ N ] (cid:105) ( FB ) in order to demonstrate improvement due to hybrid feedfor-ward (numbers less than unity indicate smaller sample vari-ance under hybrid feedforward). (d) Calculated (cid:104) σ y [ N ] (cid:105) for N = 20 as a function of duty factor, normalized to the samplevariance for the free-running LO. Data above red dashed lineindicate that the standard feedback approach produces insta-bility larger than that for the free-running oscillator. Bothdata sets assume S y ( ω ) ∝ /ω , with ω c / π = 100 /T c . Crossesrepresent data with ten noise spurs superimposed on S y ( ω ),starting at ω/ π = 1 . T − c , and increasing linearly with stepsize 0 . T − c . of a typical 1 /ω power spectrum, the inclusion of noisespurs near ω/ π = T − c results in certain regimes wherestandard feedback makes long-term stability worse thanapplying no feedback at all, while feedforward providesuseful stabilization. This significant difference arises be-cause even though the noise processes are random, knowl-edge of the statistical properties of the noise provides ameans to effectively model the average dynamical evolu-tion of the system, and accurately predict how the systemwill evolve in the future. Exact performance depends sen-sitively on the form and magnitude of S y ( ω ), but resultsdemonstrate that systems with high-frequency noise con-tent around ω/ π ≈ T − c benefit significantly from hy-brid feedforward. IV. CONCLUSION
In summary, we have presented a set of analytical toolsdescribing LLO performance in the frequency domain forarbitrary measurement times, durations, and duty cycles.We have employed these generalized transfer functions todevelop a new software approach to LO feedback stabi-lization in slaved passive frequency standards, bringingoptimal estimation techniques inside the feedback loop.This technique leverages a series of past measurementsand statistical knowledge of the noise to improve the ac-curacy of feedback corrections and ultimately improvethe stability of the slaved LO. We have validated thesetheoretical insights using numerical simulations of noisylocal oscillators and calculations of relevant stability met-rics.The results we have presented have not by any meansexhausted the space of modifications to clock protocolsavailable using this framework. For instance we havenumerically demonstrated improved correction accuracyusing nonuniform-duration T R over a cycle, as well aslong-term stability improvement using only the simplestcase of uniform T R . These approaches may be combinedto produce LLOs with improved accuracy relative to thereference at the time of correction and improved long-term stability. In cases where the penalty associated withincreasing T R is modest (lower high-frequency cutoff),such composite schemes can provide substantial bene-fits as well, improving both accuracy of correction to theLLO and overall frequency standard stability. Other ex-pansions may leverage the basic analytic formalism wehave introduced; we have introduced the transfer func-tions, | G ( ω ) | and G k,l ( ω ), but have assumed only thesimplest form for the time-domain sensitivity functionand fixed overall gain. However, it is possible to crafta measurement protocol to yield | G ( ω ) | that suppressesthe dominant spectral features of the LO noise. We haveobserved that through such an approach one may reduce the impact of aliasing on clock stabilization, indicatinga path for future work on reducing of the so-called Dicklimit in precision frequency references.In the parameter regimes we have studied the relativeperformance benefits of the hybrid feedforward approachare of metrological significance - especially consideringthey may be gained using only “software” modificationwithout the need for wholesale changes to the clock hard-ware. We believe the approach may find special signif-icance in tight-SWAP (size, weight, and power) appli-cations such as space-based clocks where significantlyaugmenting LO quality is generally impossible due tosystem-level limitations. Overall, we believe that thiswork indicates clear potential to improve passive fre-quency standards by incorporation of optimal estima-tion techniques in the feedback loop itself. Note:
Whilepreparing this manuscript we became aware of relatedwork seeking to employ covariance techniques to improvemeasurements of quantum clocks [38].
Acknowledgements : The authors thank H. Ball, D.Hayes, and J. Bergquist for useful discussions. Thiswork partially supported by Australian Research Coun-cil Discovery Project DP130103823, US Army ResearchOffice under Contract Number W911NF-11-1-0068, andthe Lockheed Martin Corporation.
APPENDIXVariances for locked local oscillators with hybridfeedforward
The standard measures for oscillator performance con-sider either a free-running LO or provide a means onlyto statistically characterize measurement outcomes underblack-box conditions. Here we present explicit analyticforms for different measurements of variance in the pres-ence of feedback locking.The expected value of the LLO sample variance can befound by substituting (3) into the definition of the samplevariance, producing a generic expression for traditionalfeedback (one measurement per correction cycle) and hy-brid feedforward (multiple measurements per cycle):0E[ σ yLLO [ N ]] = 1 N − N (cid:88) k (cid:48) =1 (cid:26) σ yLLO ( k (cid:48) ) + 1 N N (cid:88) p (cid:48) =1 N (cid:88) q (cid:48) =1 σ (¯ y LLOp (cid:48) , ¯ y LLOq (cid:48) ) − N N (cid:88) l (cid:48) =1 σ (¯ y LLOk (cid:48) , ¯ y LLOl (cid:48) ) (cid:27) (27)= 1 N − N (cid:88) k (cid:48) =1 (cid:26)(cid:18) σ yLO ( k (cid:48) ) + ¯ g k (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) r =1 (cid:98) k (cid:48) /n (cid:99) (cid:88) s =1 σ ( C r , C s ) − g k (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) u =1 σ (¯ y LOk (cid:48) , C u ) (cid:19) + 1 N N (cid:88) p (cid:48) =1 N (cid:88) q (cid:48) =1 σ (¯ y LOp (cid:48) + ¯ g p (cid:48) (cid:98) p (cid:48) /n (cid:99) (cid:88) p =1 C p , ¯ y LOn + ¯ g q (cid:48) (cid:98) q (cid:48) /n (cid:99) (cid:88) q =1 C q ) − N N (cid:88) l (cid:48) =1 σ (¯ y LOk (cid:48) + ¯ g k (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) u =1 C u , ¯ y LOl (cid:48) + ¯ g l (cid:48) (cid:98) l (cid:48) /n (cid:99) (cid:88) v =1 C v ) (cid:27) (28)= 1 N − N (cid:88) k (cid:48) =1 (cid:26)(cid:18) σ yLO ( k (cid:48) ) + ¯ g k (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) r =1 (cid:98) k (cid:48) /n (cid:99) (cid:88) s =1 σ ( C r , C s ) − g k (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) u =1 σ (¯ y LOk (cid:48) , C u ) (cid:19) + 1 N N (cid:88) p (cid:48) =1 N (cid:88) q (cid:48) =1 (cid:18) σ (¯ y LOp (cid:48) , ¯ y LOq (cid:48) ) + ¯ g p (cid:48) ¯ g q (cid:48) (cid:98) p (cid:48) /n (cid:99) (cid:88) p =1 (cid:98) q (cid:48) /n (cid:99) (cid:88) q =1 σ ( C p , C q ) (cid:19) − N N (cid:88) l (cid:48) =1 (cid:18) σ (¯ y LOk (cid:48) , ¯ y LOl (cid:48) ) + ¯ g k (cid:48) ¯ g l (cid:48) (cid:98) k (cid:48) /n (cid:99) (cid:88) k =1 (cid:98) l (cid:48) /n (cid:99) (cid:88) l =1 σ ( C k , C l ) (cid:19)(cid:27) (29)where in the case of hybrid feedback, N is defined to betotal number of measurements and n is the number ofmeasurements per cycle. The summation signs with un-primed indices are sums over whole cycles (of which thereare (cid:98) N/n (cid:99) ) and the primed indices are sums over all N measurements. In general, E[ σ yLLO [ N ]] contains recur-sive terms that cannot be concisely expressed in termsof the LO PSD S y ( ω ) and covariance transfer function G ( ω ).The Allan variance, the conventional measure of fre-quency standard instability, can be expressed analo-gously A σ y ( y ) = 12 π (cid:90) ∞ S y ( ω ) (cid:12)(cid:12) A G ( ω ) (cid:12)(cid:12) dω (30)where the transfer function, for ideal Ramsey interroga-tion, is (cid:12)(cid:12) A G ( ω ) (cid:12)(cid:12) = 2 sin ( ωT R / ωT R / (31) where T R lacks an index because the definition of the Al-lan variance assumes equal-duration interrogation bins[26]. The Allan variance calculated via this frequency-domain approach can be compared to its value via thetime-domain approach, which consists of finding the vari-ance of the difference between consecutive pairs of mea-surement outcomes: A σ y ( y ) = 12 (cid:104) (¯ y k +1 − ¯ y k ) (cid:105) (32)where ¯ y k is the k th measurement outcome and (cid:104)· · · (cid:105) mayindicate a time average or an ensemble average, depend-ing on whether y ( t ) is assumed to be ergodic.The LLO Allan variance can be found by substituting(4) into the definition of the Allan variance (32):1 A σ yLLO ( k ) = 12 E[(¯ y LLOk +1 − ¯ y LLOk ) ] (33)= 12 E (cid:20)(cid:18) ¯ y LOk +1 − ¯ g k +1 ¯ g k ¯ y LOk − ¯ y LOk + ¯ g k ¯ g k − ¯ y LOk − (cid:19) (cid:21) (34)= 12 (cid:18) σ yLO ( k + 1) + (cid:18) g k +1 ¯ g k (cid:19) σ yLO ( k ) + (cid:18) ¯ g k ¯ g k − (cid:19) σ yLO ( k − g k ¯ g k − σ (¯ y LOk +1 , ¯ y LOk − ) − (cid:18) g k +1 ¯ g k (cid:19) σ (¯ y LOk , ¯ y LOk +1 ) − g k + ¯ g k +1 )¯ g k − σ (¯ y LOk , ¯ y LOk − ) (cid:19) (35) [1] H. Katori, Nat. Photon. , 203 (2011).[2] N. Nand, J. Hartnett, E. N. Ivanov, and G. Santarelli,IEEE Trans. Microwave Thy and Tech. (2011).[3] C. W. Chou, D. B. Hume, T. Rosenband, and D. J.Wineland, Science , 1630 (2010).[4] S. Blatt, A. D. Ludlow, G. K. Campbell, J. W. Thomsen,T. Zelevinsky, M. M. Boyd, J. Ye, X. Baillard, M. Fouch´e,R. Le Targat, A. Brusch, P. Lemonde, M. Takamoto, F.-L. Hong, H. Katori, and V. V. Flambaum, Phys. Rev.Lett. , 140801 (2008).[5] J. Gu´ena, M. Abgrall, D. Rovera, P. Laurent, B. Chupin,M. Lours, G. Santarelli, P. Rosenbusch, M. Tobar, R. Li,K. Gibble, A. Clairon, and S. Bize, IEEE Trans. Ultra-son., Ferroelec., and Freq. Control (2012).[6] O. I. de la Convention du M`etre, “The international sys-tem of units (si),” (2006).[7] D. B. Sullivan, J. C. Bergquist, J. J. Bollinger, R. E.Drullinger, W. M. Itano, S. R. Jefferts, W. D. Lee,D. Meekhof, T. E. Parker, F. L. Walls, and D. J.Wineland, J. Res. Natl. Inst. Stand. Technol. , 47(2001).[8] C. Audoin and B. Guinot, The Measurement of Time ,1st ed. (Cambridge University Press, 2001).[9] P. T. H. Fisk, Rep. Prog. Phys. , 761 (1997).[10] N. Hinkley, J. Sherman, N. Phillips, M. Schioppo,N. Lemke, K. Beloy, M. Pizzocaro, C. Oates, and A. Lud-low, Science (2013).[11] N. Huntemann, M. Okhapkin, B. Lipphardt, S. Weyers,C. Tamm, and E. Peik, phys. rev. lett. (2012).[12] C. W. Chou, D. B. Hume, J. C. J. Koelemeij, D. J.Wineland, and T. Rosenband, Phys. Rev. Lett. ,070802 (2010).[13] T. Rosenband, D. B. Hume, P. O. Schmidt, C. W. Chou,A. Brusch, L. Lorini, W. H. Oskay, R. E. Drullinger,T. M. Fortier, J. E. Stalnaker, S. A. Diddams, W. C.Swann, N. R. Newbury, W. M. Itano, D. J. Wineland,and J. C. Bergquist, Science (2008).[14] G. Dick, Proc. 19th Annual Precise Time and Time In-terval (PTTI) Appls. and Planning Meeting (1987).[15] Y. Y. Jiang, A. D. Ludlow, N. D. Lemke, R. W. Fox, J. A.Sherman, L.-S. Ma, and C. W. Oates, Nat. Photon. ,158 (2011).[16] M. Takamoto, T. Takano, and H. Katori, Nature Pho-tonics , 288 (2011).[17] J. Hartnett and N. Nand, IEEE Trans. Microwave Thy and Tech. , 3580 (2010).[18] C. Hagemann, C. Grebing, T. Kessler, S. Falke,N. Lemke, C. Lisdat, H. Schnatz, F. Riehle, and U. Sterr,IEEE Trans. Instr. and Meas. , 1556 (2013).[19] G. Biedermann, K. Takase, X. Wu, L. Deslauriers, S. Roy,and M. Kasevich, Phys. Rev. Lett. (2013).[20] B. Girod, R. Rabenstein, and A. Stenger, Signals andSystems (Wiley, New York, 2001).[21] A. Kofman and G. Kurizki, Phys. Rev. Lett. , 270405(2001).[22] T. J. Green, H. Uys, and M. J. Biercuk, Phys. Rev. Lett. , 020501 (2012).[23] T. J. Green, J. Sastrawan, H. Uys, and M. J. Biercuk,New J. Phys. , 095004 (2013).[24] A. Soare, H. Ball, D. Hayes, J. Sastrawan, M. C. Jar-ratt, J. J. McLoughlin, X. Zhen, T. J. Green, and M. J.Biercuk, Nat. Phys. (2014).[25] R. F. Stengel, Optimal Control and Estimation (Dover,Mineola, New York, 1994).[26] J. Rutman, Proceedings of the IEEE , 1048 (1978).[27] W. H. Oskay, S. A. Diddams, E. A. Donley, T. M. Fortier,T. P. Heavner, L. Hollberg, W. M. Itano, S. R. Jefferts,M. J. Delaney, K. Kim, F. Levi, T. E. Parker, and J. C.Bergquist, Phys. Rev. Lett. , 020801 (2006).[28] Characterization of Clocks and Oscillators , Tech. Rep.1337 (National Institute of Standards and Technology,1990).[29] J. Barnes, A. Chi, L. Cutler, D. Healey, D. Leeson, T. Mc-Gunigal, J. Mullen, W. Smith, R. Sydnor, R. C. Vessot,and G. Winkler, IEEE Trans. Instr. and Meas.
IM-20 ,105 (1971).[30] C. Greenhall, D. Howe, and D. Percival, IEEE Trans.Ultrason., Ferroelec., and Freq. Control , 1183 (1999).[31] G. Uhrig, Phys. Rev. Lett. , 100504 (2007).[32] L. Cywinski, R. M. Lutchyn, C. P. Nave, and S. D.Sarma, Phys. Rev. B , 174509 (2008).[33] M. J. Biercuk, H. Uys, A. P. VanDevender, N. Shiga,W. M. Itano, and J. J. Bollinger, Nature , 996 (2009).[34] M. J. Biercuk, A. C. Doherty, and H. Uys, J. Phys. B , 154002 (2011).[35] C. A. Greenhall, Metrologia , S335 (2003).[36] B. Penrod, in Frequency Control Symposium, 1996. 50th.,Proceedings of the 1996 IEEE International. (1996) pp.980–987.[37] C. Nichols, P. Wu, G. Carleton, and S. Beaudin, (2006). [38] M. Mullan and E. Knill, arXiv1404.3810