Wolfgang,

Remember to scale the double-integral right:

d(t) = D
y(t) = integrate(d(t),t) = y_0 + Dt
x(t) = integrate(y(t),t) = x_0 + y_0*t + D/2*t^2

Did you miss the 1/2 factor somewhere?
That would make sense for the Random Walk phase noise.

Cheers,
Magnus

On 03/09/2015 04:57 PM, Wolfgang Wallner wrote:


On 03/06/2015 10:29 PM, Magnus Danielson wrote:
I have checked several sources, and they match up with the IEEE 1139 in
this regard.

I have also evaluated the equation for Allan variance for the random
walk noise, and it matches up with the references and what I put here:
https://en.wikipedia.org/wiki/Allan_variance#Power-law_noise

Thanks a lot for your effort!

So, the A formula you have matches up.

You will need to find another source of the mismatch.

I will take a step back and describe the overall picture of what I'm doing.
Maybe someone can help me spot where I do something wrong.
(As stated later, the part where I'm quite unsure what I'm doing is the PSD 
estimation part.)

My main goal is to simulate powerlaw noise. I then analyze the generated noise 
to check if my simulation is reasonable.

So the basic workflow would be the following:

1) Generate noise
2) Analyze the noise in the time and frequency domain
3) See that everything agrees and be happy :)

Step 1: Noise generation
-----------------------------------------------------------

I generate powerlaw noise with the method described by Kasdin and Walter in [1].
So basically I generate white noise and apply a filter as described in [1] to 
get a PSD shape corresponding to the different values of alpha.
The part of the PSD that will have the correct shape depends on the filter 
length and the simulated sampling frequency.
Basically: the length of a simulation I would like to carry out specifies a 
lower bound on the filter length to get correct results.

For WPM, WFM and RW noise I can use a shortcut: for these types of noise the 
filter coefficients are basically a discrete derivative, an identity filter and 
a cumulative sum.
This is expected, as it agrees with [2], which states that integration of 
powerlaw noise decreases alpha by 2 (chapter 3.4 in [2]).
Thus for even values of alpha I can even skip the expensive convolution to 
apply the filter and implement the filters directly.

As input white noise I use a Gaussian distribution, mainly because that is what 
is used in the original paper.
(I have also found another implementation [3] that optionally provides a 
uniform distribution).

I'm quite confident that the noise generation part works as expected.
However, even if I do something wrong here, it should not influence the 
analyzing part.

Step 2: Analyzing noise
-----------------------------------------------------------

2.1 Time domain

To analyze powerlaw noise in the time domain, I use a Matlab script called 
'allan' [4], which calculates the Allan Deviation.
I also found another Matlab tool called 'Stability Analyzer' [5], which can 
also calculate ADEV values.
These two tools are developed by different authors and expect different input 
formats, but their results agree for any noise example I have tried so far.
Thus I would say both of them can be trusted to work as expected.

2.2 Frequency domain

IEEE 1139[6] defines S_y as: "frequency spectrum Sy(f): One-sided spectral density 
of the normalized frequency fluctuations, as defined in normalized frequency fluctuations 
y(t)."

However, I'm not sure how to calculate this measure for a given noise sample.
Anything I describe below is just based on 'I think this might work'.

If anyone knows a better way of calculating S_y, or tools that can be used for 
this task, I would be glad to hear about it :)

As already stated in the earlier mail I use the method described in [7] to 
estimate the one-sided PSD of my noise data in FFD format.
These plots are quite noisy, and to improve the graphical presentation I use 
the averaging method described in [8].
I split the noise vector in parts of equal length, calculate the individual 
PSDs and average over them.
Using this averaging method, the PSD plots converge to lines on a log-log plot 
with the expected slopes.
I have an example figure attached to the mail that shows the effect of the 
averaging (PSD_Average.png).

Step 3: Comparing time and frequency domain results
-----------------------------------------------------------

At this point I have plots for both the Allan Deviation and the FFD-PSD, and 
would like to compare them.
As first step I estimate h_alpha from the Allan Deviation plot (I'm aware that I need 
to take care for the Allan Deviation <-> Allan Variance conversion).
Then I try to estimate the expected PSD values and compare them with my actual 
plot using the formulas from IEEE 1139.

However, at this point a see that RW noise behaves unexpected :(

Numerical Example:
-----------------------------------------------------------

Suppose the figure attached as 'Numeric_example.png':

At Tau = 0.1s the ADEV plot has a value of 0.005849, so de AVAR would be 
3.4211e-05 at this point.
The constant A is 2 * pi^2/3 = 6.5797.

Thus the value of h_-2 could be roughly estimated as AVAR / (Tau * A) = 
~5.2e-05.
This would lead to an expected S_y value at a frequency f = 10Hz of

h_-2 * f = 5.2000e-07, or -62.84dB

The actual plot value is at -59.83, so its ~3dB too larger than expected.

regards, Wolfgang


[1] Kasdin and Walter, Discrete Simulation of Power Law noise, 1992
[2] Riley, NIST SP 1065: Handbook of Frequency Stability Analysis, 2008
[3] http://people.sc.fsu.edu/~jburkardt%20/c_src/cnoise/cnoise.html
[4] http://de.mathworks.com/matlabcentral/fileexchange/13246-allan
[5] 
http://de.mathworks.com/matlabcentral/fileexchange/31319-stability-analyzer-53230a
[6] IEEE 1139
[7] http://de.mathworks.com/help/signal/ug/psd-estimate-using-fft.html
[8] http://www.dspguide.com/ch9/1.htm


Cheers,
Magnus

On 03/06/2015 11:04 AM, Wolfgang Wallner wrote:


On 03/05/2015 07:23 PM, Attila Kinali wrote:
Servus!

Servus :)

On Thu, 05 Mar 2015 14:35:51 +0100
Wolfgang Wallner <wolfgang-wall...@gmx.at> wrote:

For the random walk noise the expected line is off by a factor of
exactly 2 from the calculated plot, and I don't know how to explain
this
behavior.

I'm probably the wrong one to answer, as I have never done any noise
simulation or even read up the relevant papers, but...
A factor of 2 sounds like the difference you would get between one sided
and two sided noise PSD's.


I calculate the one-sided PSD of the FFD data as described in [1] (first
paragraph), so the code looks like this:

    xdft = fft(x);
    xdft = xdft(1:N/2+1);
    psdx = (1/(Fs*N)) * abs(xdft).^2;
    psdx(2:end-1) = 2*psdx(2:end-1);

Remark: Before calculating the PSD, I split the data into parts of equal
size, calculate the PSD for each one, and average over the set of PSDs.
This improves the graphical visualization a lot.

As the result matches my expectation exactly for 4 different kinds of
noise, I would have assumed that this PSD calculation approach is quite
reasonable.

As I see the unexpected behavior only with random walk noise, and the
main difference in the calculation is the term A, I would suspect that
it has something to do with it.

However, I'm a novice in this field, so any hint is very appreciated.

regards, Wolfgang


[1] http://de.mathworks.com/help/signal/ug/psd-estimate-using-fft.html
_______________________________________________
time-nuts mailing list -- time-nuts@febo.com
To unsubscribe, go to
https://www.febo.com/cgi-bin/mailman/listinfo/time-nuts
and follow the instructions there.

_______________________________________________
time-nuts mailing list -- time-nuts@febo.com
To unsubscribe, go to
https://www.febo.com/cgi-bin/mailman/listinfo/time-nuts
and follow the instructions there.


_______________________________________________
time-nuts mailing list -- time-nuts@febo.com
To unsubscribe, go to https://www.febo.com/cgi-bin/mailman/listinfo/time-nuts
and follow the instructions there.
_______________________________________________
time-nuts mailing list -- time-nuts@febo.com
To unsubscribe, go to https://www.febo.com/cgi-bin/mailman/listinfo/time-nuts
and follow the instructions there.

Reply via email to