[time-nuts] AVAR <-> S_Y conversion

Sat Mar 14 06:12:37 EDT 2015

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-wallner at 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 at 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 at 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 at febo.com
>> To unsubscribe, go to https://www.febo.com/cgi-bin/mailman/listinfo/time-nuts
>> and follow the instructions there.