Contents

Vibration Test Specification Demo

This DuraTools demo shows how field data from measurements on a truck are used to define a proper random test. The script uses the black list from the Data quality assessment demo to exclude data that were found to be bad.

% Copyright (c) 2003-2006, Axiom EduTech AB, Ljusterö, Sweden. All rights reserved.
% URL: http://www.vibratools.com Email: info@vibratools.com

Introduction

load blacklist
fs = 1000;
FontSize=10;

Calculate PSD Matrix and Turn Into 1/3 Octave Bands

CurrDir=pwd;
cd('DuraData')
m = 1;
for n = 1:102
    if isempty(find(n == blacklist)) == 1
        eval(['load test',int2str(n)]);
        [z,f] = psdnorm(Data,fs,16384);
        [zt,ft] = psd2ters(z,f);
        psdmatrix(:,m) = zt.*zt./ft/0.23;
        rmsvecn(n,1) = std(Data);
        rmsvec(m,1) = std(Data);
        m = m+1;
    end
end

Make a Scatterplot

figure(1)
loglog(ft,psdmatrix,'k*')
set(gca,'FontSize',FontSize)
xlabel('Frequency [Hz]','FontSize',FontSize)
title('PSD Scatter Plot','FontSize',FontSize)
axis([0.8 500 1.e-5 10])
xlogtics(0.8,500)

METHOD 1

We will show two different methods for producing the test spec. First, by creating a mean spectrum from the scatter, and then getting an RMS value for the test as mean(rms) + 1.96*std(rms). The statistics on rms are taken in log

Calculate a Mean Spectrum (log mean)

[M,N] = size(psdmatrix);
for m = 1:M
    meanpsd(m,1) = exp(mean(log(psdmatrix(m,:))));
end
meanrms = exp(mean(log(rmsvec)));
stdrms = exp(std(log(rmsvec)));
testval = meanrms+1.96*stdrms;
%

Calculate rms of Mean Spectrum

rmsmean = sqrt(sum(0.23*ft.*meanpsd));
%
% Correct the Amplitude
testspec1 = testval*testval/rmsmean/rmsmean*meanpsd;
figure(1)
loglog(ft,testspec1,'r','LineWidth',2)
set(gca,'FontSize',FontSize)
hold on
loglog(ft,meanpsd,'k','LineWidth',2)
loglog(ft,psdmatrix,'k*')
xlabel('Frequency [Hz]','FontSize',FontSize)
title('Test PSD #1','FontSize',FontSize)
legend('Test Spectrum #1','Mean PSD','Scatter')
axis([0.8 500 1.e-5 10])
xlogtics(0.8,500)

METHOD 2

Create a spectrum by taking mean and 1.96*std for each frequency. The statistics on rms are taken in log.

for m = 1:M
    testspec2(m,1) = exp(mean(log(psdmatrix(m,:)))+1.96*std(log(psdmatrix(m,:))));
end
figure(2)
loglog(ft,testspec2,'r','LineWidth',2)
set(gca,'FontSize',FontSize)
hold on
loglog(ft,psdmatrix,'k*')
xlabel('Frequency [Hz]','FontSize',FontSize)
title('Test PSD #2','FontSize',FontSize)
legend('Test Spectrum #2','Scatter')
axis([0.8 500 1.e-5 10])
xlogtics(0.8,500)

Compute Test Time

The truck was driving 30, 50, and 70 km/h during the test. To get the testing time, we calculate mean(rms) for different speeds. The speed is found in the first line of the Header.

speedvec = zeros(102,1);
for n = 1:102
    if isempty(find(n == blacklist)) == 1
        eval(['load test',int2str(n),' Header']);
        speed = str2num(Header(1,6:8));
        speedvec(n,1) = speed;
    end
end
k30 = find(speedvec == 30);
rms30 = mean(rmsvecn(k30));
k50 = find(speedvec == 50);
rms50 = mean(rmsvecn(k50));
k70 = find(speedvec == 70);
rms70 = mean(rmsvecn(k70));
%

Find rms of testspec2

rmstest = sqrt(sum(0.23*ft.*testspec2));
%
% Assume the following exposure for our equipment:
% 30 kmph 3000 h
% 50 kmph 2000 h
% 70 kmph 5000 h
%
% we use the formula:
% (test time)/(exposure time) = [(rms test)/(measured rms)]^(-8)
testtime30 = 3000*60*(rmstest/rms30)^(-8);
testtime50 = 2000*60*(rmstest/rms50)^(-8);
testtime70 = 5000*60*(rmstest/rms70)^(-8);
disp('Calculated test times in minutes')
disp(['Corresponding to 30 kmph  ',sprintf('%8.2f',testtime30)]);
disp(['Corresponding to 50 kmph  ',sprintf('%8.2f',testtime50)]);
disp(['Corresponding to 70 kmph  ',sprintf('%8.2f',testtime70)]);
disp(['Total Test Time           ',sprintf('%8.2f',testtime30+testtime50+testtime70)]);
disp(' ')
disp(['Total Test Time in hours  ',sprintf('%8.2f',ceil((testtime30+testtime50+testtime70)/60))]);
testtime = ceil((testtime30+testtime50+testtime70)/60);

Calculated test times in minutes
Corresponding to 30 kmph      0.03
Corresponding to 50 kmph      3.61
Corresponding to 70 kmph   2181.11
Total Test Time            2184.76
 
Total Test Time in hours     37.00

End Script

Go back to original directory and save results

cd(CurrDir)
save testspec2 testspec2 ft testtime
fprintf('VIBSPEC Done.\n')

VIBSPEC Done.

Published with MATLAB® 7.1