function [f,I11,I22,I12]=crossp(x1,x2,dt)
%
% [f,I11,I22,I12]=crossp(x1,x2,dt)  computes periodograms I11, I22 and
%                                   (complex) cross periodogram I12 for
%                                   time series in vectors x1,x2. 
%                                   Time step is dt and f is a frequency
%                                   vector (units=cycles/units-of-dt).
%
%   I11 = dt/N |fft(x1)|^2 , I22 sim., I12 = dt/N fft(x1).*conj(fft(x2))
%
% For coherence calculations, see PCOH.M.
%
%
% MIM 25/5/92 - based on CROSSP.FOR
%
% Note from BRS, 24/1/96: cross spectrum produced by smoothing crossp's I12 
% by rsmooth, etc. is  C - iQ , for co- and quadrature-spectra, C and Q, 
% as defined by Jenkins & Watts (1968) .

N1=length(x1); N2=length(x2); 
if N1~=N2; error('CROSSP.M requires length(x1)=length(x2)');end;
N=N1;

X1=fft(x1(:)-mean(x1)); X2=fft(x2(:)-mean(x2));
jnyq=fix(N/2)+1;X1=X1(1:jnyq);X2=X2(1:jnyq);
I11=dt/N * real(X1.*conj(X1));
I22=dt/N * real(X2.*conj(X2));
I12=dt/N * X1.*conj(X2);          
I12=conj(I12);%%%makes exactly same as crossp.for for testing purposes
f=[0:jnyq-1]'/(N*dt);