MATLAB's $\tt cpsd$ and $\tt pwelch$ - different results for cross spectral power density I am wondering why the following code does not yield the same results for Sxy1 and Sxy2, where Sxy1 is calculated using the MATLAB's pwelch function, and Sxy2 is calculated using cpsd function.

fs = 1000;
t = 0:1/fs:5-1/fs;
win=hann(500);
overlap=300;
nfft=500;
fs=1000;   

x = cos(2*pi*100*t)+randn(size(t));
y=sin(2*pi*160*t)+cos(2*pi*150*t)+randn(size(t));

[Sxx,f] = pwelch(x,win,overlap,nfft,fs);
[Syy,f] = pwelch(y,win,overlap,nfft,fs);

Sxy1=sqrt(Sxx.*Syy); %******

[Sxy2,f] = cpsd(x,y,win,overlap,nfft,fs);

Sxy2=abs(Sxy2); %******
Whilst I have not worked out how to calculate CPSD from power spectra, I have confirmed the following two methods are equivalent (using amplitude spectra, rather than power spectra or cross-correlation):

Fs = 1000;
t = 0:1/Fs:1-1/Fs;
x = cos(2*pi*100*t) + randn(size(t));
y=sin(2*pi*160*t)+cos(2*pi*150*t)+randn(size(t));
N = length(x);
nfft=500;
lfft=2*nfft;
overlap=0;
win=rectwin(N);

xdft1 = fft(x,2*nfft);  &lt;outputs two-sided complex amplitude spectra
ydft1 = fft(y,2*nfft);
Sxy1a=xdft1.*conj(ydft1);   
Sxy1b=Sxy1a./N^2;
Sxy1c=Sxy1b(1:lfft/2+1);
Sxy1c(2:end-1)=2*Sxy1c(2:end-1);
Sxy1=abs(Sxy1c)';

%CPSD Method 2
Sxy2a=cpsd(x,y,win,overlap,2*nfft,Fs);   &lt; uses the cpsd function with Welch’s overlapped method
Sxy2=abs(Sxy2a);