function c=correlation(x,y)
% c=correlation(x,y)

  [nx,nt]=size(x);
  xm=mean(x,2);
  ym=mean(y,2);

  x=x-xm(:,ones(1,nt));
  y=y-ym(:,ones(1,nt));

  c=zeros(nx,1);
  xstd=sqrt(mean(x.^2,2));
  ystd=sqrt(mean(y.^2,2));

  ok=(xstd>0)&(ystd>0);
  c(ok)=mean(x(ok,:).*y(ok,:),2)./xstd(ok)./ystd(ok);

  c(~ok)=NaN;
