Plotting an Accurate PSD Estimate in MATLAB

% This script plots the power spectral density estimate
% of a signal --Alper Ucar, ucar {at} ieee {dot} org
% data: output stream in time-domain
% fftlen: FFT length
% Fs: sampling rate
% pxf: mean-square power of the signal
function[pxf] = pspec(data,fftlen,Fs)

data = squeeze(data);
N = 2^floor(log2(fftlen));
data = data(1:N);
% hann-windowing
W = hann(N);
x = data.*W;
% Spectral content from 0 to Fs
xf = abs(fft(x))/N;
% One-sided PSD
xf = xf(1:N/2);
% coherent gain factor
% See "On the Use of Windows for Harmonic Analysis
% with the DFT" by Fred J Harris
CG = 8/3;
psdxf = CG*2*xf.^2;
% DC component should not be doubled
psdxf(1) = psdxf(1)/2;
pxf = 10*log10(sum(psdxf));
psdxfdB = 10*log10(psdxf);
v = 0:1/(N-1):0.5;
F = v*Fs;
plot(F,psdxfdB)
xlabel('Frequency [Hz]')
ylabel('PSD [dBV]')
title('Output PSD')
grid on
text(50,-5,sprintf('Signal Power: %4.2f',pxf),...
'BackgroundColor',[1 1 1]);
grid on