[time-nuts] Modified Total Deviation calculation

[Magnus Danielson]

Dear fellow time-nuts,

In my effort to implement the suite of ADEV and friends, I have been 
implementing various forms of them along with the 1000 frequency sample 
test sequence out of NIST SP1065:
http://tf.nist.gov/timefreq/general/pdf/2220.pdf

However, I seem to be unable to perfectly match the numbers for the MTotDev.

My numbers so far:

         m = 1        m = 10       m = 100
Max     0.9957452943 0.7003371204 0.5489367785
Min     0.0013717599 0.2545924150 0.4533354120
Avg     0.4897744629 0.4897744629 0.4897744629
Sdev    0.2884663647 0.0929635201 0.0320665644
NAdev   0.2922318781 0.0996573606 0.0389780433
OSAdev  0.2922318781 0.0915995342 0.0324134303
OAdev   0.2922318781 0.0915995342 0.0324134303
MAdev   0.2922318781 0.0617237638 0.0217092091
Tdev    0.1687201535 0.3563623166 1.2533817739
Hdev    0.2943883291 0.1052754194 0.0391086056
OHdev   0.2943883291 0.0958108317 0.0323763825
MHdev   0.2942275231 0.0621023549 0.0213087110
TOTdev  0.2922318781 0.0913474326 0.0340653025
MTOTdev 0.2303857898 0.0555288598 0.0195467513

The MTOTdev numbers according to page 118 (of PDF, page number 108 
according to the printed pagenumbers) should be

Modified Total Dev 2.418528e-01 6.499161e-02 2.287774e-02

Do anyone happend to have an implementation (in source) of MTOTdev at 
hand giving the NIST SP1065 numbers?

Please note that W. Riley has about the same document in his Handbook, 
and it reflects the same number. I would also suspect that a STABLE32 
run would give those numbers.

I have been using both SP1065 and the original article as reference:
http://tf.nist.gov/general/pdf/1369.pdf

I've put care into ensuring that I implement it as close to these as 
possible, but with no luck in fixing the numbers.

I use the averaging of the two half-ranges of the 3*m block tecnique 
rather than the minimum square estimator as recommended. The formulation 
given may seem strange, but it is the 3/2*m sample average of the

f(i)=(x[n+3/2*m+i)-x[n+i])/ (3/2*m)

frequency estimation where i varies from 0 to 3/2*m-1.

I have already verified my test-sequence and the numbers produced by the 
other algorithms have been able to match after removing various bugs. An 
independent implementation may be a good clue.

An alternative may be that the published numbers is incorrect for some 
reason, but I don't have sufficient proof for that.

I have however made two different implementation variants that crunch 
out the same numbers.

Cheers,
Magnus