function corr = corr(N,Mmax)

%computes g(m)=<xn(i)*xn(i+m)>
%at m=0, g(m)=1/3, m~=0 g(m)=1/4

%initialize
rnd=rand(1,N);
sumold=0;

for m=0:Mmax
    for i=1:N
        if i<=N-m
            sumnew=sumold+rnd(i)*rnd(i+m);
        else
            sumnew=sumold+rnd(i)*rnd(i+m-N); %periodic
        end
    
    sumold=sumnew;
    end
    g(m+1)=sumold/N;
    sumold=0;
end

avg=mean(rnd)
stdev=std(rnd)

for i=1:length(g)
    m(i)=i-1;
end

%length(g)

%plot results
figure(1)
plot(m,g)
figure(2)

%plot histogram to show uniform distribution of rnd
hist(rnd)

%examine value at g(m=0)
g(1)