how to use a MATLAB function called MC estimator to determine frequency?

the following is the matlab function for the MC estimator

 function w = mc(x);
 % w = mc(x) is used to estimate frequency
 % the MC method from the vector x
 % x is supposed to be a noisy single-tone sequence
 % w is the estimated frequency in radian
 N=max(size(x));
 t1=0;
 t2=0;
 for n=3:N
 t1=t1+x(n-1)*(x(n)+x(n-2));
 t2=t2+2*(x(n-1))^2;
 end
 r = t1/t2;
 if (r>1)
 r=1;
 end
 if (r<-1)
 r=-1;
 end
 w=acos(r);

attempted to determine the frequency estimate of the data, which can be found at the following link,using the MC function.

 http://solarscience.msfc.nasa.gov/greenwch/spot_num.txt 
(please note that the first row and rows after July 2012 are removed)

this is my attempt where the data of thrid column are used only (observing the data, sampling interval is 1 month, that is Fs=1)

 load spot_num.txt
 ssn = spot_num(:,3);
 ssn = ssn-mean(ssn);
 x=ssn;
 N=max(size(x));
 t1=0;
 t2=0;
 for n=3:N
 t1=t1+x(n-1)*(x(n)+x(n-2));
 t2=t2+2*(x(n-1))^2;
 end
 r = t1/t2;
 if (r>1)
 r=1;
 end
 if (r<-1)
 r=-1;
 end
 w=acos(r);

w should be close to 0.0076 according to model answer, but somehow i just couldnt get the right answer? could anyone help, thx!

Answer 1

[Wayne King]

It may be that this method is just not very robust for noise-corrupted data.

For example:

 t = 0:159;
 x = cos(pi/4*t);
 w = mc(x);

I get 0.7854, which is pi/4 and the method does a perfect job, but if I add some noise

 xnew = x+0.5*randn(size(t));
 plot(xnew)
 w = mc(xnew);

then the estimate is way off.

Why not use the periodogram, which will get the frequency correct for your sunspot data?

   [Pxx,F] = periodogram(ssn,rectwin(length(ssn)),length(ssn),1);
   [~,idx] = max(Pxx);
   F(idx)

[modified covariance]

thank you Wayne King, the periodogram u suggested certainly gets the answer which is 0.0076.

unfortunately, i must have to use the following code

 load spot_num.txt
 ssn = spot_num(:,3);
 ssn = ssn-mean(ssn);

together with MC estimator function to get the frequency "As Required by Assignment Guidelines"...

but anyway, thx a lot!