Quantcast

袁伟时:中国的奴性和戾气从哪里来?

“芯片大学”虚晃一枪,人才断层问题不能跑步解决

两大中国首富双双被重挫-释放信号强烈

民间帝王赖小民和性感女星舒淇与许晴

女性高潮时为什么会“喷水”?

Facebook Twitter

分享到微信朋友圈

点击图标下载本文截图到手机
即可分享到朋友圈。如何使用?

自由微信安卓APP发布,立即下载!
查看原文

振动台试验数据处理程序,内含超详细的注释

2016-07-01 声振之家 声振之家

% STDATA.m

% 振动台试验数据整理和处理软件

% CopyrightConcrete

% 江苏省土木工程与防灾减灾重点实验室

%本程序部分代码使用了王济、胡晓编著的《MATLAB在振动信号处理中的应用》一书中的代码,在此表示感谢;


%注意:

% 1. 因编程时间仓促,本软件为非通用软件,内部参数与具体试验相关;

% 2. 使用者在使用时应确保明了程序含义,并根据试验内容作相应修改;

% 3. 本软件为自由软件,使用者在注明出处后可以自由使用;

% 4. 对使用本程序造成的任何后果,本人不承担任何责任;

% 5. 欢迎指出程序错误之处,联系方法:concretelu@gmail.com


% 数据文件示例:

% 2 %需要处理的文件的数量

% AE1a.C02 %以下为需要处理的文件名

% AE3a.C02

% 128 %采样频率

% 0 %是否需要进行峰值调整;0:不进行;

% -1 %是否需要进行基线调整:-1:进行但不输出文件

% 3 %基线调整拟合多项式阶数

% 0 %是否需要进行滤波0:不进行;

% 1 %是否需要进行积分:1:进行且输出文件;

% 0.5 %最小截止频率

% 40 %最大截止频率

% 9800 %单位转换系数

% 2 %积分次数(加速度求位移)

% 1 %是否需要进行通道间计算

% 0 %是否需要进行绘图

clear

%清除内存中所有变量和函数

clc

%清除工作窗口中所显示的内容

close all hidden

%关闭所有隐藏的窗口

echo off

%关闭回显

fnc=input('数据处理控制文件名:''s')

fid=fopen(fnc'r')

nf=fscanf(fid'%f'1)

%需处理的文件数量

for i=1nf

fname(i,:)=fscanf(fid'%s'1)

%读入各数据文件名

end


sf=fscanf(fid'%f'1)

%读入采样频率

flag_adj=fscanf(fid'%d'1)

if flag_adj~=0

%是否需要进行峰值调整;0:不进行;1:进行且输出文件;-1:进行但不输出文件

adj=fscanf(fid'%f'[1nf])

%读入各文件调整系数

end


flag_bas=fscanf(fid'%d'1)

if flag_bas~=0

%是否需要进行基线调整0:不进行;1:进行且输出文件;-1:进行但不输出文件

basm=fscanf(fid'%d'1)

%读入拟合多项式阶数

end


flag_fil=fscanf(fid'%d'1)

if flag_fil~=0

%是否需要进行滤波0:不进行;1:进行且输出文件;-1:进行但不输出文件

fmin=fscanf(fid'%f'1) 

%最小截止频率

fmax=fscanf(fid'%f'1) 

%最大截止频率

end


flag_int=fscanf(fid'%d'1)

if flag_int~=0

%是否需要进行积分:0:不进行;1:进行且输出文件;-1:进行但不输出文件

fmin=fscanf(fid'%f'1)

%最小截止频率

fmax=fscanf(fid'%f'1)

%最大截止频率

c=fscanf(fid'%f'1)

%单位转换系数

it=fscanf(fid'%f'1)

%积分次数

end


flag_cmp=fscanf(fid'%d'1)

%是否需要进行通道间计算

flag_fig=fscanf(fid'%d'1)

%是否需要进行绘图

status=fclose(fid)

%关闭控制文件

delete('Acc.xls')

%删除存在的结果文件

delete('Disp.xls')

delete('Drift.xls')

delete('Model.xls')

delete('TotalAcc.xls')

delete('TotalDisp.xls')

delete('TotalDrift.xls')

Acc_string={'项目''台面加速度''台面位移''底层''二层''三层'...

'四层''五层''六层''七层''八层''九层''屋面''一层南''一层北'...

'六层南''六层北''屋面南''屋面北'} %建立加速度表头字符串

Disp_string={'项目',,'底层''二层''三层''四层''五层''六层''七层''八层'...

'九层''屋面''一层南''一层北''六层南''六层北''屋面南''屋面北'}

%建立位移表头字符串

Drift_string={'项目''底层''二层''三层''四层''五层''六层''七层'...

'八层''九层''一层南北''六层南北''屋面南北'} 

%建立层间位移表头字符串

for fi=1nf

x=dlmread(fname(fi,:)'\t'11) 

%从数据文件读入数据,去除标题行和时间列

[nch]=size(x) 

%确定数据长度n和通道数ch

t=01/sf(n-1)/sf 

%建立信号离散时间向量t

TableDisp=x(:,2) 

%读取第二通道:台面位移

if flag_adj~=0 

%确定是否需要进行加速度峰值调整

x=adj(fi).*x

%进行峰值调整

if flag_adj==1

filename=strcat('Adj_'fname(fi,:))

dlmwrite(filenamex)

%输出峰值调整后的数据

end


figure(1)

%绘制调整后的加速度信号

subplot(gl4[1 4])

%采用两栏显示台面加速度

plot(tx(:,1))

xlabel('时间 (s)')

ylabel('加速度 (g)')

title('台面加速度信号')

grid on

for i=3ch

%绘制各通道信号调整峰值后随时间变化的图形

subplot(gl4i+2)

plot(tx(:,i))

grid on

end

end


case_name=fname(fi,:)

%建立工况名称向量

x2=x.^2

%计算加速度有效值

xsum=cumsum(x2)

rmsx=sqrt(xsum(n,:)/n)

xlswrite('Acc.xls'Acc_stringcase_name'A1')

%输出加速度时程及统计值

xlswrite('Acc.xls'{'Max Acc.''Min Acc.''MaxAbs Acc.''Val Acc.'}...

case_name'A2')

acc_out=[max(x)min(x)max(-min(x)max(x))rmsx]

xlswrite('Acc.xls'acc_outcase_name'B2')

xlswrite('Acc.xls'[t'x]case_name'A6')

MaxAcc=[MaxAccmax(x)]

%建立每个工况的统计值,后期统一输出

MinAcc=[MinAccmin(x)]

AbsMaxAcc=[AbsMaxAccmax(-min(x)max(x))]

ValAcc=[ValAccrmsx]

if flag_bas~=0

%确定是否需要进行基线调整

for i=1ch

tt=(01/sf(n-1)/sf)'

%建立信号离散时间列向量t

a=polyfit(ttx(:,i)basm)

%计算趋势项的多项式待定系数向量a

x(:,i)=x(:,i)-polyval(att)

%形成消除趋势项的数据

end


if flag_bas==1

filename=['Bas_'fname(fi,:)]

dlmwrite(filenamex)

%输出基线调整后的数据

end

end


if flag_fil~=0 

%确定是否需要进行单独滤波

fmin=fscanf(fid'%f'1)

%最小截止频率

fmax=fscanf(fid'%f'1)

%最大截止频率

end


if flag_int~=0 

%确定是否需要进行积分

tt=(01/sf(n-1)/sf)'

%建立信号离散时间列向量t

nfft=2^nextpow2(n)

%取大于并接近n2的幂次方为FFT长度

df=sf/nfft

%计算频率间隔(Hz/s)

ni=round(fmin/df+1)

%计算指定频带对应数组的下标

na=round(fmax/df+1)

dw=2*pi*df

%计算圆频率间隔(rad/s)

w1=0dw2*pi*(0.5*sf-df)

%建立正的离散圆频率向量

w2=2*pi*(0.5*sf-df)-dw0

%建立负的离散圆频率向量

w=[w1w2]

%将正负圆频率向量合成单一向量

w=w.^it

%以积分次数为指数,建立圆频率变量向量

disp=zeros(nch)

%为避免不同向量长度影响,变量清零

adisp=zeros(nch)

Drift=zeros(nch-6)

for i=1ch

tx=x(:,i)'

y=fft(txnfft)

%FFT变换

a=zeros(1nfft)

%进行积分的频域变换

a(2nfft-1)=y(2nfft-1)./w(2nfft-1)

if it==2

%进行二次积分的相位变换

y=-a

else

real(y)=imag(a)

%进行一次积分的相位变换

imag(y)=-real(a)

end


a=zeros(1nfft)

a(nina)=y(nina)

%消除指定正频带外的频率成分

a(nfft-na+1nfft-ni+1)=y(nfft-na+1nfft-ni+1)

%消除指定负频带外的频率成分

y=ifft(anfft)

%IFFT变换

y=real(y(1n))*c

%取逆变换的实部n个元素并乘以单位变换系数作为积分结果

ab=polyfit(tty'2)

%计算趋势项的多项式待定系数向量a

adisp(:,i)=y'-polyval(abtt)

%形成消除趋势项的数据(绝对位移)

disp(:,i)=adisp(:,i)-TableDisp

%计算各点的相对位移

ac=polyfit(ttdisp(:,i)2)

%计算趋势项的多项式待定系数向量ac

disp(:,i)=disp(:,i)-polyval(actt)

%形成消除趋势项的数据(相对位移)

end

end


if flag_int==1

%如果需要输出位移文件

y1=min(disp(:,3ch))

%搜寻每层相对位移最小值(负向)

y2=max(disp(:,3ch))

%搜寻每层相对位移最大值(正向)

d2=disp(:,3ch).^2

dsum=cumsum(d2)

rmsd=sqrt(dsum(n,:)/n)

%搜寻每层相对位移有效值

xlswrite('Disp.xls'Disp_stringcase_name'A1')

%输出位移时程及统计值

xlswrite('Disp.xls'{'Max Disp''Min Disp''MaxAbs Disp''Val Disp'}...

case_name'A2')

dispout=[y2y1max(-y1y2)rmsd]

xlswrite('Disp.xls'dispoutcase_name'B2')

xlswrite('Disp.xls'[t' disp(:,3ch)]case_name'A6')

%输出基线调整后的相对位移数据

MaxDisp=[MaxDispy2]

MinDisp=[MinDispy1]

AbsMaxDisp=[AbsMaxDispmax(-y1y2)]

ValDisp=[ValDisprmsd]

end


if flag_cmp~=0

%确定是否需要进行通道间计算

for i=19

Drift(:,i)=disp(:,i+3)-disp(:,i+2)

%计算层间位移

end


%计算两角点位移差

Drift(:,10)=disp(:,14)-disp(:,13)

Drift(:,11)=disp(:,16)-disp(:,15)

Drift(:,12)=disp(:,18)-disp(:,17)

end


if flag_cmp==1

y3=min(Drift)

%搜寻各层间位移最小值

y4=max(Drift)

%搜寻各层间位移最大值

dr2=Drift.^2

drsum=cumsum(dr2)

rmsdr=sqrt(drsum(n,:)/n)

%计算层间位移有效值

xlswrite('Drift.xls'Drift_stringcase_name'A1')

xlswrite('Drift.xls'{'Max Drift''Min Drift''MaxAbs Drift''Val Drift'}...

case_name'A2')

driftout=[y4y3max(-y3y4)rmsdr]

xlswrite('Drift.xls'driftoutcase_name'B2')

xlswrite('Drift.xls'[t' Drift]case_name'A6')

%输出层间位移数据

MaxDrift=[MaxDrifty4]

MinDrift=[MinDrifty3]

AbsMaxDrift=[AbsMaxDriftmax(-y3y4)]

ValDrift=[ValDriftrmsdr]

end


[TxyF]=tfestimate(x(:,3)x(:,12)nfft/2[][]40)

%求台面信号和顶层信号间的传递函数

%aTxy=abs(Txy)

%输出模态幅值

%xlswrite('Model.xls'{'频率''幅值''相位角'}case_name'A1')

%xlswrite('Model.xls'[F aTxy angle(Txy)]case_name'A2')

if flag_fig~=0 

%确定是否需要进行绘图

gl=ceil((ch-2)/4)+1

%计算图形显示所需要的行数,列数为4

figure(2)

subplot(gl4[1 4])

%采用两栏显示台面加速度

plot(tadisp(:,1)tTableDisp'r')

xlabel('时间 (s)')

ylabel('位移 (mm)')

title([case_name' 积分位移信号'])

grid on

ymin=floor(min(y1'))

%搜寻全部相对位移最小值

ymax=ceil(max(y2'))

%搜寻全部相对位移最大值

xmax=ceil((n-1)/sf)

for i=3ch

%绘制各通道相对位移随时间变化的图形

subplot(gl4i+2)

plot(tdisp(:,i))

axis([0 xmax ymin ymax])

%xlabel('时间 (s)')

%ylabel('加速度 (g)')

%title('通道'+int2str(ch))

grid on

end


figure(3)

ymin=floor(min(y3'))

%搜寻全部相对位移最小值

ymax=ceil(max(y4'))

%搜寻全部相对位移最大值

xmax=ceil((n-1)/sf)

for i=112

%绘制各通道信号x随时间变化的图形

subplot(43i)

plot(tDrift(:,i))

axis([0 xmax ymin ymax])

grid on

end


figure(4)

subplot(211)

plot(FaTxy)

grid on

subplot(212)

plot(Fangle(Txy)*180/pi())

%输出相位角

grid on

end

end


%建立数据表名称

xlswrite('TotalAcc.xls'Acc_string'MaxAcc''A1')

xlswrite('TotalAcc.xls'Acc_string'MinAcc''A1')

xlswrite('TotalAcc.xls'Acc_string'AbsMaxAcc''A1')

xlswrite('TotalAcc.xls'Acc_string'ValAcc''A1')

xlswrite('TotalDisp.xls'Disp_string'MaxDisp''A1')

xlswrite('TotalDisp.xls'Disp_string'MinDisp''A1')

xlswrite('TotalDisp.xls'Disp_string'AbsMaxDisp''A1')

xlswrite('TotalDisp.xls'Disp_string'ValDisp''A1')

xlswrite('TotalDrift.xls'Drift_string'MaxDrift''A1')

xlswrite('TotalDrift.xls'Drift_string'MinDrift''A1')

xlswrite('TotalDrift.xls'Drift_string'AbsMaxDrift''A1')

xlswrite('TotalDrift.xls'Drift_string'ValDrift''A1')

xlswrite('TotalAcc.xls'MaxAcc'MaxAcc''B2') %输出数据

xlswrite('TotalAcc.xls'MinAcc'MinAcc''B2')

xlswrite('TotalAcc.xls'AbsMaxAcc'AbsMaxAcc''B2')

xlswrite('TotalAcc.xls'ValAcc'ValAcc''B2')

xlswrite('TotalDisp.xls'MaxDisp'MaxDisp''B2')

xlswrite('TotalDisp.xls'MinDisp'MinDisp''B2')

xlswrite('TotalDisp.xls'AbsMaxDisp'AbsMaxDisp''B2')

xlswrite('TotalDisp.xls'ValDisp'ValDisp''B2')

xlswrite('TotalDrift.xls'MaxDrift'MaxDrift''B2')

xlswrite('TotalDrift.xls'MinDrift'MinDrift''B2')

xlswrite('TotalDrift.xls'AbsMaxDrift'AbsMaxDrift''B2')

xlswrite('TotalDrift.xls'ValDrift'ValDrift''B2')

close all


本文转载自新浪Glulam的博客,程序版权见程序注释

关联阅读

数据处理名词术语中英文对照,一网打尽

信号处理中各种噪声定义、分类和性质,种类多多!

实例解析FFT在非平稳信号分析中的局限性

数字滤波的优点及10种常用数字滤波方法比较



声明:本微信转载文章出于非商业性的教育和科研目的,并不意味着赞同其观点或证实其内容的真实性。版权归原作者所有,如转载稿涉及版权等问题,请立即联系我们,我们会予以更改或删除相关文章,保证您的权利!

文章有问题?点此查看未经处理的缓存