最近设计了一个8位的DAC，准备仿真它的SFDR，我已经在DAC前面用verilogA写了一个8位理想的ADC，时钟频率Fclk=100MHz，然后信号输入频率Fin=Fclk*127/256，得到了输出曲线，然后通过cadence取出了256个采样点，打算利用matlab分析，但是在网上找的matlab代码跑了之后出现了一些问题，如图所示，代码如下所示，想问问论坛里的各位大神这是什么原因呢？
clear
load DACtest.txt -ascii %cadence导出的数据文件，文件必须以英文开头
A= DACtest; %将测量数据赋给A，此时A为N×2的数组
x=A(:,1); %将A中的第一列赋值给x，形成时间序列
y=A(:,2); %将A中的第二列赋值给y，形成被测量序列
%显示原始数据曲线图（时域）
subplot(2,1,1);
plot(x,y) ; %显示原始数据曲线图
%axis([min(x) max(x) 1.1*floor(min(y)) 1.1*ceil(max(y))]) %优化坐标，可有可无
xlabel('时间 (s)');
ylabel('被测变量y');
title('原始信号(时域)');
grid on;
format long;
%傅立叶变换
y=y-mean(y); %消去直流分量，使频谱更能体现有效信息
fclk=(length(x)-1)/(max(x)-min(x)); %仪器的采样频率
numpt=length(y); %data.txt中的被测量个数，即采样个数
%If no window function is used, the input tone must be chosen to be unique and with
%regard to the sampling frequency. To achieve this prime numbers are introduced and the
%input tone is determined by fIN = fSAMPLE * (Prime Number / Data Record Size).
%To relax this requirement, window functions such as HANNING and HAMING (see below) can
%be introduced, however the fundamental in the resulting FFT spectrum appears 'sharper'
%without the use of window functions.
Doutw=y;
%Doutw=y'.*hanning(numpt);
%Doutw=y'.*hamming(numpt);
%Performing the Fast Fourier Transform
Dout_spect=fft(Doutw);
%Recalculate to dB
Dout_dB=20*log10(abs(Dout_spect));
%Display the results in the frequency domain with an FFT plot
figure; %建立图形
maxdB=max(Dout_dB(1:numpt/2));
%For TTIMD, use the following short routine, normalized to —6.5dB full-scale.
%plot([0:numpt/2-1].*fclk/numpt,Dout_dB(1:numpt/2)-maxdB-6.5);
%plot([0:30/2-1].*fclk/numpt,Dout_dB(1:30/2)-maxdB);
h=stem([0:256/2-1].*fclk/numpt,Dout_dB(1:256/2)-maxdB,'fill',':');
set(h,'MarkerEdgeColor','red','Marker','*')
grid on;
title('FFT PLOT');
xlabel('ANALOG INPUT FREQUENCY (Hz)');
ylabel('AMPLITUDE (dB)');
a1=axis; axis([a1(1) a1(2) -120 a1(4)]);
%Calculate SNR, SINAD, THD and SFDR values
%Find the signal bin number, DC = bin 1
fin=find(Dout_dB(1:numpt/2)==maxdB);
%Span of the input frequency on each side
%span=max(round(numpt/200),5);%span=max(round(numpt/200),5);
%Approximate search span for harmonics on each side
spanh=2;%spanh=2;
%Determine power spectrum
spectP=(abs(Dout_spect)).*(abs(Dout_spect));
%Find DC offset power
%Pdc=sum(spectP(1:span));
%Extract overall signal power
%Ps=sum(spectP(fin-span:fin+span));
%Vector/matrix to store both frequency and power of signal and harmonics
Fh=[];
%The 1st element in the vector/matrix represents the signal, the next element represents
%the 2nd harmonic, etc.
Ph=[];
%Find harmonic frequencies and power components in the FFT spectrum
for har_num=1:10 %har_num谐波总数
%Input tones greater than fSAMPLE are aliased back into the spectrum
%tone=rem((har_num*(fin-1)+1)/numpt,1); %rem（x，y）x除以y的余数 numpt(Number of Points)
%if tone>0.5
%Input tones greater than 0.5*fSAMPLE (after aliasing) are reflected
%tone=1-tone;
%end
%Fh=[Fh tone];
%For this procedure to work, ensure the folded back high order harmonics do not overlap
%with DC or signal or lower order harmonics
%har_peak=max(spectP(round(tone*numpt)-spanh:round(tone*numpt)+spanh));
%har_bin=find(spectP(round(tone*numpt)-spanh:round(tone*numpt)+spanh)==har_peak);
%har_bin=har_bin+round(tone*numpt)-spanh-1;
%Ph=[Ph sum(spectP(har_bin-1:har_bin+1))];
Ph=[Ph sum(spectP(har_num*(fin-1):har_num*(fin-1)+2))];
end
%Determine the total distortion power
Pd=sum(Ph(2:10)); %Pd总失真功率 Ph(1) is fundamental harmonic谐波 power
%Determine the noise power
%Pn=sum(spectP(1:numpt/2))-Pdc-Ps-Pd; %Pn噪声功率 Ps信号功率
format;%设置输出格式
%SINAD=10*log10(Ps/(Pn+Pd)) % SINAD = 10*log10(Ps/(Pn+Pd)) 信号与噪声失真比
%SNR=10*log10(Ps/Pn) % SNR = 10*log10(Ps/Pn) 信噪比
%disp('THD is calculated from 2nd through 5th order harmonics');
%THD=10*log10(Pd/Ph(1)) % THD总谐波失真Pd总失真功率 Ph(1) 基波能量
SFDR=10*log10(Ph(1)/max(Ph(2:5))) %SFDR无杂散动态范围
disp('Signal & Harmonic Power Components:');
HD=10*log10(Ph(1:10)/Ph(1))
|
/2 
eetop公众号
创芯大讲堂
创芯人才网
发表于 2021-9-3 16:33:36
楼主
