请教如何将功率谱转换为幅值谱

wdecai

请教各位采用功率谱模型函数生成的功率谱如何正确的转换为傅里叶幅值谱，我转换后的好像不正确，与实际幅值谱幅值之间相差很大。程序如下：
clc
clear all
close all
direc='Data';
records=fname(direc);
finame=char(records);
[xg,dt,NbPt]=peer2acc(direc,finame);%读取信号，dt为时间步长，xg为信号幅值，n为采样点数
EFz=1/dt;%采样频率
nfft=2^nextpow2(NbPt) ;
nbf=nfft-NbPt;
nzer=zeros(1,nbf);
xg=[xg,nzer];
n=nfft;
Fgen=EFz/nfft:EFz/nfft:EFz/2;%产生离散频率
SaCgen=Respec(Fgen,xg,dt);%计算反应谱
Wgen=2*pi*Fgen;
Z=0.05;%阻尼比
S=15;%强震持时
% 通过功率谱密度函数模型将反应谱转换为功率谱
GSUM=0.0;
for I=1:length(Fgen),
    POW=2*Z*Wgen(I)*S;
    TRANS=1-exp(-POW);
    ZS=Z/TRANS;
    QYS=sqrt(4*ZS/pi);
    XSP=Wgen(I)*S/2/pi/(-log(exp(-1)));
    RSTAR=sqrt(2*log(2*XSP));
    ET=-RSTAR*QYS*sqrt(pi/2);
    ARG=2*XSP*(1-exp(ET));
    RSP=sqrt(2*log(ARG));
    GW(I)=(4*ZS/(Wgen(I)*pi))*((SaCgen(I)/RSP)^2-GSUM);
    if GW(I)<=0.0, GW(I)=0.0;end
    if I==1,
        GSUM=0.5*GW(1)*Wgen(1);
    else GSUM=GSUM+GW(I)*(Wgen(I)-Wgen(I-1));
    end
end
%
% 将功率谱转换为傅里叶幅值谱
AA(1)=sqrt(GW(1)*Wgen(1)*1);
for I=2:length(Fgen)-1
    AA(I)=sqrt(GW(I)*(Wgen(I+1)-Wgen(I-1))*1);
end
I=length(Fgen);
AA(I)=sqrt(GW(I)*(Wgen(I)-Wgen(I-1))*1);

[ 本帖最后由 wdecai 于 2010-8-6 15:57 编辑 ]