On Mon, 02 May 2022 17:12:47 +0200
Matthias Welwarsky time-nuts@welwarsky.de wrote:

So far so good, but you are missing one major component that is
important for a GPSDO: temperature dependence. All OCXO have
a slight temperature dependence that can be anything from a
few 1e-12/K to tens of 1e-9/K. It's often also not linear
and, generally, not well specified in OCXO data sheets.
There were plenty of discussions about this on time-nuts
and especially Bob Camp has shared a lot of good information.

Random walk is not flicker noise. Random walk has a PSD
that is proportional to 1/f^2 while flicker noise has a PSD
that is proportional to 1/f.

Though, at the time scales that are interesting for a GPSDO
(i.e. (mili-)seconds to minutes) the dominant noise would be flicker
frequency with a (phase) PSD of 1/f^3 or even random walk
frequency with a (phase) PSD of 1/f^4.

The model is ok, though I would definitely add the model for
temperature dependence. Also, for a GPSDO you will need to model
more than just a day. A lot of effects, not just temperature,
have periodes of 12-24h. If you simulate only a single day, you'll
not properly simulate all those long term effects. Not to mention
that temperature looks weird when looking at it on a larger scale.

Attached are two representative plots of what the temperature in
a normal (small) office looks like, when nobody is allowed to enter.
One graph is for Jan-Sept of 2021 and one for May (units in °C).
You will see that there is quite some structure and "weirdness"
to it. Especially note that here is a not-quite daily spike in
temperature when the heating comes on somewhen between 18:00
and 22:00 and shuts off again between 21:00 and 00:00.

The noise in the temperature is mostly due to convection in
the room and not due to the noise of the temperature sensor.

I would also recommend that you create your noise vector with
something like sigma-theta by François Vernotte[1,2].
The bruiteur tool can generate all types of noise directly and
give you a single vector to read and process. This should be
faster than building your own power law noise simulator.

		Attila Kinali

[1] https://theta.obs-besancon.fr/spip.php?article103&lang=fr
[2] https://sourcesup.renater.fr/www/sigmatheta/

--
The driving force behind research is the question: "Why?"
There are things we don't understand and things we always
wonder about. And that's why we do research.
-- Kobayashi Makoto

Hi

By far the best approach is to use actual data. Grab a pair of OCXO’s and
compare them. A single mixer setup is one easy ( = cheap ) way to do it. You
will get the sum of the two devices, but for simulation purposes, it will be much
closer to reality than anything you can brew up with a formula.

Bob

On Mon, May 2, 2022 at 10:01 PM Bob kb8tq kb8tq@n1k.org wrote:

But a programmatic way of generating plausible OCXO noise lets you
generally millions of times more data, over a much larger spectrum of
plausible operating conditions-- so that you can test the stability of
an algorithm over a wider collection of conditions.

It's not a complete replacement for using real data-- you should do
that too.  But it can be a much more comprehensive test.

Hi

….. except that having done this for many decades on hundreds of
designs, , a single data set from a real OCXO is likely to show you
things that millions of simulations from a formula will somehow miss ….

Bob

Matthias,

On 2022-05-02 17:12, Matthias Welwarsky wrote:

What you have shown is 2 out of the standard 5 noise-types. The Leeson
model describes how 4 noise-types can be generated, but in addition to
that the Random Walk Frequency Modulation is included. For standard, see
IEEE Std 1139. For Leeson model, see David Leeson articles from special
issue of Feb 1966, as I believe you can find on IEEE UFFC site. The
Wikipedia Allan Deviation article may also be of some guidance.

Noise-types:

White Phase Modulation - your thermal noise on phase

Flicker Phase Modulation - your 1/f flicker noise on phase

Because we have oscillators that integrate inside the resonators
bandwidth, we also get their integrated equivalents

White Frequency Modulation - your thermal noise on frequency - forms a
Random Walk Phase Modulation

Flicker Frequency Modulation - Your 1/f flicker noise on frequency

Random Walk Frequency Modulation - Observed on disturbed oscillators, a
random walk in frequency, mostly analyzed for finding rare faults.

You have only modelled the first two. This may or may not be relevant
depending on the bandwidth of your control loop and the dominant noise
there. It may not be relevant. The phase-noise graph will illustrate
well where the power of the various noise-types intercept and thus
provide a range over frequency for which one or the other is dominant.
You control loop bandwidth will high-pass filter this for your locked
oscillator, and low-pass filter for your reference oscillator/source.
The end result will be a composit of both. The Q-value / damping factor
will control just how much jitter peaking occurrs at the cut-over
frequency, and the basic recommendation is to keep it well damped at all
times.

The ADEV variant of this is similar.

Now, to simulate flicker you have a basic problem, whatever processing
you do will end up needing the full time of your simulation data to
maintain the slope. You can naturally cheat and do a much reduced setup
that only provides the 1/f slope of PSD at and about the area where it
dominates. Notice that this can occurr both for the FPM and FFM cases.
Also notice that if we respect the model, these should be completely
independent.

Simulation wise, you can turn WPM into WFM by integration, and FPM into
FFM by integration. Similarly the WFM becomes RWFM through a second
integration. Just source each of these independent to respect the model.

The Leeson model for an oscillator does not have these fully
independent, but for most uses, simulate fully independent, and you will
not make more fool of yourself than anyone else usually does, present
company included.

As you see, flicker and random-walk is not the same thing. Often "random
walk" implies Random Walk Frequency Modulation.

Another thing. I think the rand function you use will give you a normal
distribution rather than one being Gaussian or at least pseudo-Gaussian.
A very quick-and-dirty trick to get pseudo-Gaussian noise is to take 12
normal distribution random numbers, subtract them pair-wise and then add
the six pairs. The subtraction removes any bias. The 12 samples will
create a normalized deviation of 1.0, but the peak-to-peak limit is
limited to be within +/- 12, so it may not be relevant for all noise
simultations. Another approach is that of Box-Jenkins that creates much
better shape, but comes at some cost in basic processing.

As I've stated above, I think it misses key noise-types. I really wonder
if the flicker noise model you have is producing real flicker noise.

Using Bob's suggestion of using real data is not bad, noise simulation
is hard and a research field in itself. My quick and dirty methods get
you started at relatively low cost, except for flicker, flicker always hurt.

Cheers,
Magnus

On 5/2/22 6:09 PM, Magnus Danielson via time-nuts wrote:

rand() gives uniform distribution from [0,1). (Matlab's doc says (0,1),
but I've seen zero, but never seen 1.) What you want is randn(), which
gives a zero mean, unity variance Gaussian distribution.

https://www.mathworks.com/help/matlab/ref/randn.html

That would be for uniform distribution. A time-honored approach from the
IBM Scientific Subroutine Package.

Hi Jim,

Thanks for the corrections. Was way to tired to get the uniform and
normal distributions right.

rand() is then by classical UNIX tradition is generated as a unsigned
integer divided by the suitable (32th) power of two, so the maximum
value will not be there, and this is why a small bias is introduced,
since 0 can be reached but not 1.

In practice the bias is small, but care is taken never the less.

Cheers,
Magnus

On 2022-05-03 03:43, Lux, Jim wrote:

On 5/2/22 7:03 PM, Magnus Danielson via time-nuts wrote:

I'll bet it's "pre-Unix" -

System/360 Scientific Subroutine Package Programmer's Manual, Version
III, DOI: 10.3247/SL2Soft08.001 says 1968 version, but I'm pretty sure
the SSP is older (it is Version 3 after all)

I used it on an IBM 1130 as a mere youth in 1969

http://media.ibm1130.org/1130-106-ocr.pdf

    SUBROUTINE RANDU(X,IY,YFL)
    IY = IX*899
    if (IY) 5,6,6
5   IY=IY+32767+1
6   YFL=IY
    YFL=YFL/32767.
    RETURN
    END

GAUSS does normal distribution, with this comment:

Y approaches a true normal distribution asympototically as K approaches
infinity. For this subroutine, K was chosen as 12 to reduce execution time.

It also helps that the variance of a uniform distribution is 1/12, so
summing 12 numbers produces a distribution with a variance of 1.

But it's older than that.. I found a reference to it in some
documentation for 7090 from 1961.  Since Unix wasn't even a name until
1970...

Yes, that's a clever technique.

And the less said about the actual "randomness" of generators from that
era, the better.

Dear all,

thanks for your kind comments, corrections and suggestions. Please forgive if
I don't reply to all of your comments individually. Summary response follows:

Attila - yes, I realize temperature dependence is one key parameter. I model
this meanwhile as a frequency shift over time.

Bob - I agree in principle, real world data is a good reality check for any
model, but there are only so few datasets available and most of the time they
don't contain associated environmental data. You get a mix of effects without
any chance to isolate them.

Magnus, Jim - thanks a lot. Your post encouraged me to look especially into
flicker noise an how to generate it in the time domain. I now use randn() and
a low-pass filter. Also, I think I understood now how to create phase vs
frequency noise.

I've some Timelab screenshots attached, ADEV and frequency plot of a data set
I generated with the following matlab function, plus some temperature response
modeled outside of this function.

function [phase] = synth_osc(samples,da,wpn,wfn,fpn,ffn)
# low-pass butterworth filter for 1/f noise generator
[b,a] = butter(1, 0.1);

# aging
phase = (((1:samples)/86400).^2)*da;
# white phase noise
phase += (randn(1, samples))*wpn;
# white frequency noise
phase += cumsum(randn(1, samples))*wfn;
# 1/f phase noise
phase += filter(b,a,randn(1,samples))*fpn;
# 1/f frequency noise
phase += cumsum(filter(b,a,randn(1,samples))*ffn);

end

osc = synth_osc(400000, -50e-6, 5e-11, 1e-11, 5e-11, 5e-11);

Thanks.

On Montag, 2. Mai 2022 17:12:47 CEST Matthias Welwarsky wrote: