%《數字信號處理教程——MATLAB釋義與實現》第四章例4.7.1程序hc471
% 求矩形窗頻譜能量的程序hc471
% 電子工業出版社出版  陳懷琛編著 2004年9月
%
T=0.01;N=16384; a=2;    % 輸入原始數據
D=2*pi/(N*T);          % 頻率分辨率
M=floor(a/T);           % 求邊界點下標
x=[ones(1,2*M+1)];      % 給出時間序列
X=fft(x,N);            % 求傅立葉變換FFT
k1=input('k1=（書上取零） '),       % 求能量的起點頻率下標
k2=input('k2=（書上取floor(0.5*pi/D)）'),  % 求能量的終點頻率下標
k=k1:k2;                % 求能量頻段的下標序列
E=T/N*sum(abs(X(k+1)).^2)   % 求正頻率區間能量
E1=2*T/N*sum(abs(X(k+1)).^2) % 求全頻率區間能量
E1/4.01*100                 % 占總能量百分比