Here is (5) from the paper roughly translated into Mathematica.
2 DeltaT/Nw Sum[Abs[1/Sqrt[Nd]Sum[w[[n]]x[[n]]^l E^((-I 2 Pi k(n-1))/Nd)]^2, {l, 1, Nw}], {n, 1, Nd}]
To use this you need to assign appropriate constant values to DeltaT, Nw, Nd, choose and initialize a windowing function in list w, choose and initialize some plausible time domain signal list x.
Then you should probably try using you data in both Matlab and Mathematica so you can compare the results and we start the process of trying to figure out why the results are completely different, what mistakes have been made, how to fix those, etc, etc, etc.
I would start with really really simple data. A single simple sine wave that exactly fills a single segment with complete cycles is always a nice start. Then the sum of two different sine waves. Then... you get the idea. Gradually ratchet it up with multiple windows and overlapping windows and more complicated time domain data, exterminating every small error before making the next step.
If Matlab anywhere in their documentation describes exactly what they do then we might avoid some of the confusion by trying to port the Matlab method to Mathematica, rather than me just grabbing the first random page that Google showed me and asking you if this was it.
