| [rt_array2, lv_array2, lv_array_BG2]=rt_script(filenames_in, Fc, pbs, pab, outfile, num_x_filter, t_high, t_low, Room_Name, plot_all_mics, smd, num_std, method)
|
function [rt_array2, lv_array2, lv_array_BG2]=rt_script(filenames_in, Fc, pbs, pab, outfile, num_x_filter, t_high, t_low, Room_Name, plot_all_mics, smd, num_std, method)
% % Syntax;
% % [rt_array2, lv_array2, lv_array_BG2, ptsa, nptsa]=rt_script(filenames_in, Fc, pbs, pab, outfile, num_x_filter, t_high, t_low, Room_Name, plot_all_mics, smd, num_std, method)
% %
% % This program calculates the reverberation time in one third octave
% % bands. This program supports two measuremnt methods: 1) the speaker
% % on speaker off method and 2) the balloon pop method. The first method
% % using loud speakers is considered to be mmore accurate since it excites
% % the room modes more effectively. The first method generally has a
% % better signal to noise ratio than the second method. However somtimes
% % it is difficult to transport and setup the equipment for the speaker on
% % speakre off method. It is usually much easier anf more fun to pop
% % balloons. The balloon pop method may be sufficiently accurate
% % and is usually much more convenient.
% %
% % The data is processed by filtering the data into a third octave bands.
% % The instantaneous sound pressure level is calculated using the
% % hilbert transform.
% %
% % A moving time averaged mean value with a 40 ms sample rate for the
% % speaker on-off method and 10 ms sample rate for the balloon pop method
% % is used to smooth the instantaneous sound pressure level.
% %
% % For the speaker on-off method, the signal high level is calculated.
% % For the balloon pop method the maximum data point is found.
% %
% % The background level is estimated from the last t_low seconds of data.
% %
% % Cutoff times from the signal high to the signal low are calculated.
% % The signal cutoff are after the speaker was shut off or after the
% % balloon popped. The background cutoff is before the time record
% % reaches the background level.
% %
% % The instantaneous sound pressure level data between the cutoffs is
% % the transient reverberant decay.
% %
% % A decay line is robustly fit to the transient reverberant decay.
% %
% % The slope of the transient instantaneous amplitude is used to calculate
% % the T60 reverberation time.
% %
% % The program outputs a pdf file for each center frequency and one tab-
% % delimited text file which serve as documentation of the reverberation
% % time calculations.
% %
% % The program automatically loads the matlab files from the current
% % working directory. Each data file must contain a data variable which
% % is the time record of the sound pressure. The units must be linear.
% % The program is setup for units of Pa and a size eqaual to
% % (num_mics, num_data_points).
% %
% % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% % See reverb_time.m program for a baloon pop example.
% % Here is a speaker on-off example and a description of the input values.
% %
% Example=''; % Simulated reverberation data
% % For the Speaker on Speaker off method
%
% % This example is for a four microphone with each microphone recording at
% % a different level. The reverberation time is set to a different for each
% % microphone to illustrate the decay line fitting.
%
% filenames_in={'speaker_on_example.mat'};
% % cell array of filenames
%
% Fc=400; % Center Frequency array default has 31 bands
% % from 20 Hz to 20 KHz.
% %
% % There is a one-to-one correspondence between
% % The array Fc and the filenames
% %
% % Each file is assigned the next
% % center frequency in the array Fc
%
% pbs=15; % Percentage of the total range in dB below the average signal
% % level to discard from the curve fit. Increase pbs to avoid
% % spurious high level sound data at the speaker cutoff time
% %
% pab=40; % Percentage of the range above the average background level
% % to discard from the curve fit. ncrease pab to avoid spurious background
% % noise, where the reveberant decay approaches the background level.
% % This is more of a problem in environments with very low background levels;
% % for example laboratory hemi-anechoic chambers;
%
% pab=0; % 0 is a good number for acoustic environments with high continuous background levels;
% % for example animal shelters with dogs barking in the surrounding rooms.
%
% outfile='Reverb_data'; % Filename for the data for each microphone. The data includes
% % the reverberation times, high sound level at the beginning of the recording,
% % and the background sound pressure level at the end of the recording.
% % If more than 1 microphone is used then
% % standard deviation, confidence interval and outlier data are also printed.
%
% num_x_filter=2; % 2 is a typical value. Number of times to filter the data.
% % Minimum value is 1.
% % Typically a value of 2 to 10 is good. At low
% % frequencies (Fc < 100), num_x_filter=10 has a significant phase shift.
%
% t_high=0.5; % (s) duration of signal on at the beginning of the time record.
% % t_high is used to calculate the average
% % signal level. The speaker cutoff time is
% % estimated using the average signal level.
% % A longer t_high time is better when there is filter phase delay
% % especially at low frequencies (Fc < 100) or if the speaker turns on slowly.
%
% t_low=0.5; % (s) duration of quiet time at the end of the time record
% % used to calculate the average background
% % level, which is used to estimate the background cutoff time.
% % Background cutoff occurs when the reverberant decay approaches the
% % background noise level. A longer t_low time helps to avoid
% % the dips and peaks in the background level. Dips and peaks
% % can be caused by sudden quiet and background dog barks etc...
%
% Room_Name='Auditorium'; % The room name appears on the pdf plot
%
% plot_all_mics=1; % 1 plot all microphone time records
% % 0 plots only the microphone time record for the
% % microphone with the reverberation time
% % closest to the mean value
%
% smd=0; % acronym for "surpress microphone data" from the pdf plot legend.
% % smd=0; Allow all mics to have a level, reverberation time, and
% % quality marker on the plot legend.
% % smd=1; Allow only 1 mic to have a level, reverberation time, and
% % quality marker on the plot legend
% %
% num_std=2; % Keep all reverberation times within
% % 2 standard deviations of the geometric mean
% % num_std is the number of standard deviations
% % to accept data. 1 is good for data with high
% % background noise, 2 is better for laboratory
% % measurements.
%
% method=1; % 1 is the speaker on speaker off method.
% % otherwise is the balloon pop method
%
% % Now create the file
%
% tau=[2, 2.3, 3.0, 2.4]; % time constants for exponential decay
%
% tb=5; % seconds, duration of decay signal
%
% Fs=50000; % sampling frequency rate
%
% t=1/Fs*(1:(tb*Fs)); % seconds time array for calculating
% % decay signal
%
% A=[100, 110, 115, 130]; % RMS Amplitude in (Pa)
% decay1=A(1).*exp(-1./tau(1).*t).*randn(1,length(t));
% decay2=A(2).*exp(-1./tau(2).*t).*randn(1,length(t));
% decay3=A(3).*exp(-1./tau(3).*t).*randn(1,length(t));
% decay4=A(4).*exp(-1./tau(4).*t).*randn(1,length(t));
%
% % The test signal has
% % 1) nominal one second of signal high
% % 2) reverberant decay
% % 3) nominal one second of background noise
% %
% SP=zeros(4, 2*Fs+length(t));
% SP(1, :)=[A(1)*randn(1,Fs) decay1 A(1).*exp(-1./tau(1).*tb)*randn(1, Fs);];
% SP(2, :)=[A(2)*randn(1,1.5*Fs) decay2 A(2).*exp(-1./tau(2).*tb)*randn(1, 0.5*Fs);];
% SP(3, :)=[A(3)*randn(1,1.2*Fs) decay2 A(3).*exp(-1./tau(3).*tb)*randn(1, 0.8*Fs);];
% SP(4, :)=[A(4)*randn(1,0.8*Fs) decay2 A(4).*exp(-1./tau(4).*tb)*randn(1, 1.2*Fs);];
%
% filenames_in={'speaker_on_example.mat'};
% t_SP=1/Fs*(1:length(SP));
% save(filenames_in{1}, 'SP', 't_SP');
%
% [rt_array2, lv_array2, lv_array_BG2]=rt_script(filenames_in, Fc, pbs, ...
% pab, outfile, num_x_filter, t_high, t_low, Room_Name, plot_all_mics, smd, num_std, method);
%
% % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% % Output Variable Description
% %
% % rt_array2 cell array of reverberation times
% %
% % lv_array2 cell array of A-weighted average signal levels during t_high
% %
% % lv_array_BG2 cell array of A-weighted average background levels during t_low
% %
% % The size of the ouput cell arrays are {num_files, 1}
% % each cell has an array size (num_mics, 1)
% %
% %
% % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% %
% % Uses the Signal Processing Toolbox
% % Uses the following progrmas which can be found on matlab central
% % Uses the octave program set from Christophe Couvreur
% % Uses programs from Alexandros Leontitsis
% % Uses the modified programs from Alexandros Leontitsis
% %
% % List of Subprograms
% %
% % data_loader2
% % file_extension
% % convert_double
% % reverb_time
% % oct3dsgn From Christophe Couvreur found on matlab central
% % sub_mean
% % LMTSreg A modified program original from Alexandros Leontitsis
% % found on matlab central
% % LMSloc From Alexandros Leontitsis found on matlab central
% % LMS_trim
% % rand_int
% % LMTSregor A modified program original from Alexandros Leontitsis
% % found on matlab central
% % adsgn From Christophe Couvreur found on matlab central
% % cdsgn From Christophe Couvreur found on matlab central
% % Leq_all_calc
% % Aweight_time_filter
% % Cweight_time_filter
% % call_ci
% % ci from Raymond Reynolds found on matlab central
% %
% % save_a_plot_reverb_time
% % file_extension
% %
% %
% % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% %
% % Program was written by Edward L. Zechmann
% %
% % date 16 October 2007
% %
% % modified 9 November 2007 revision 1 revised comments
% % (grammar and clarity)
% %
% % modified 2 January 2008 updated comments
% %
% % modified 13 January 2008 Added the data_loader2
% %
% % modified 10 February 2008 Revised the balloon pop method
% %
% % modified 13 February 2008 Further revised the balloon pop method
% % Added an example. Updated comments
% %
% %
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% % Please feel free to modify this code as you please.
% %
figure(1);
delete(1);
% program assumes all matlab files in directory are recordings
num_files=length(filenames_in);
if isempty(Fc)
Fc=[20 25 31 40 50 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 10000 12500 16000 20000 ];
end
if length(Fc) < num_files
fprintf(1, '%s\r', 'Not enough Center Frequencies for input Fc');
fprintf(1, '%s\r', 'Adding more Center Frequencies for input Fc');
num1=ceil( (num_files-length(Fc))/31);
for e1=1:num1;
Fc=[Fc [20 25 31 40 50 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 10000 12500 16000 20000] ];
end
fprintf(1, '%s\r', 'New List of Center frequencies');
fprintf(1, '%d\t', Fc);
fprintf(1, '\r');
end
% The program processes the first file with the center frequency of Fc(1)
count=0;
rt_array2=cell(num_files, 1); % Reverberation time cell array
lv_array2=cell(num_files, 1); % Reverberation Level cell array
lv_array_BG2=cell(num_files, 1); % Background Level cell array
meanlevel=0;
mic_num_array=zeros(1,num_files);
ptsa=cell(1,num_files);
nptsa=cell(1,num_files);
default_wav_settings_in={};
default_mat_config_in={};
if num_x_filter < 1
num_x_filter=1;
end
for e1=1:num_files;
%reset filter coefficients
% B and A are for A-weighted time record
Ba=1;
Aa=1;
% Bc and Ac are for C-weighted time record
Bc=1;
Ac=1;
% load the data file
[SP_var, vibs_var, Fs_SP_var, Fs_vibs_var, default_wav_settings_out, default_mat_config_out]=data_loader2(filenames_in{e1}, default_wav_settings_in, default_mat_config_in);
% set the default wav file setttings
% for configuring the channels
default_wav_settings_in=default_wav_settings_out;
default_mat_config_in=default_mat_config_out;
num_vars=length(SP_var);
for e3=1:num_vars;
% Initialze the time data variable
SP=[];
if iscell(SP_var)
if length(SP_var) >= e3
SP=SP_var{e3};
end
if length(Fs_SP_var) >= e3
Fs_SP=Fs_SP_var{e3};
end
% currently there is not any code for reverberation time of
% vibrations
%if length(vibs_var) >= e7
% vibs=vibs_var{e7};
%end
%
%if length(Fs_vibs_var) >= e7
% Fs_vibs=Fs_vibs_var{e7};
%end
else
SP=SP_var;
Fs_SP=Fs_SP_var;
% currently there is not any code for reverberation time of
% vibrations
%vibs=vibs_var;
%Fs_vibs=Fs_vibs_var;
end
[num_mics, num_pts]=size(SP);
if num_mics > num_pts
SP=SP';
[num_mics, num_pts]=size(SP);
end
% Use the Fc array to input the center frequency of each recording file.
fc=Fc(e1);
dt=1/Fs_SP;
t_SP=dt*(1:num_pts);
title_str={Room_Name, 'For each Microphone', 'Level (dBA), T_6_0 (s), Quality', };
rt_array=zeros(num_mics,1); % buffer variable for the reverberation times
lv_array=zeros(num_mics,1); % buffer variable for the A-weighted signal level
lv_array_BG=zeros(num_mics,1); % buffer variable for the A-weighted background level
if logical(num_mics > 0) && logical(num_pts > 0)
count=count+1;
mic_num_array(count)=num_mics; % keep track of the number of mics in each file
figure(1);
delete(1);
color_a={[0 0 1], [0 1 0], [1 0 0], [1 0 1], [0 1 1], [0.75 0.5 1], [0.4 0.2 0.3], [0.1 0.5 0.1], [1 0.5 1], 0*[1 1 1], 0.4*[1 1 1], 0.7*[1 1 1], 0.9*[1 1 1]};
Lya=zeros(num_mics, num_pts); %
rtpa=zeros(num_mics, num_pts); %
legend_str=cell(num_mics, 1); % all the strings for the legend
t_r=max(t_SP)-min(t_SP);
for e2=1:num_mics;
[rt, Ly, rtp]=reverb_time(SP(e2, :), Fs_SP, fc, t_high, t_low, num_x_filter, pbs, pab, method);
Lya(e2, :)=Ly;
rtpa(e2, :)=rtp;
rt_array(e2)=rt;
if method == 1
% Calculate the A-weighted signal level
if isequal(Ba, 1) && isequal(Aa, 1)
[LeqA, LeqA8, LeqC, LeqC8, Leq, Leq8, Ba, Aa, Bc, Ac]=Leq_all_calc(SP(e2, 1:floor(Fs_SP*t_high)), Fs_SP, 1);
else
[LeqA, LeqA8, LeqC, LeqC8, Leq, Leq8, Ba, Aa, Bc, Ac]=Leq_all_calc(SP(e2, 1:floor(Fs_SP*t_high)), Fs_SP, 1, Ba, Aa, Bc, Ac);
end
else
LeqA=20*log10(1./0.00002.*max(abs(SP(e2, :))));
end
lv_array(e2)=LeqA;
% Calculate the A-weighted background level
if isequal(Ba, 1) && isequal(Aa, 1)
[LeqA_BG, LeqA8_BG, LeqC_BG, LeqC8_BG, Leq_BG, Leq8_BG, Ba, Aa, Bc, Ac]=Leq_all_calc(SP(e2, [(num_pts-floor(Fs_SP*t_low)):end]), Fs_SP, 1);
else
[LeqA_BG, LeqA8_BG, LeqC_BG, LeqC8_BG, Leq_BG, Leq8_BG, Ba, Aa, Bc, Ac]=Leq_all_calc(SP(e2, [(num_pts-floor(Fs_SP*t_low)):end]), Fs_SP, 1, Ba, Aa, Bc, Ac);
end
lv_array_BG(e2)=LeqA_BG;
% Concatentate the legend string parts
if rt < 0.001
legend_str{e2}=[ sprintf('NAN')];
else
legend_str{e2}=[ sprintf('%5.1f', 0.1*round(10*LeqA)), ' ', sprintf('%-6.3f', 1/1000*round(1000*rt))];
end
if e2 == num_mics
% removal of outliers
rt_array_cal=rt_array(rt_array > 0);
gm=geomean2(rt_array_cal); % calculate the geometric mean
stdrt=std(rt_array_cal); % calculate standard deviation
% find all values within num_std standard deviations of the geometric mean
if length(rt_array) > 3
min_lim=max([gm - num_std*stdrt 0]);
pts=intersect( find(rt_array > min_lim), find(rt_array < gm + num_std*stdrt));
npts=setdiff( 1:num_mics, pts); % find the indices of the outlying data
else
pts=1:num_mics; % small data set keep all data
npts=[]; % small data no outliers
end
outliers=zeros(size(rt_array));
if num_mics > 1
for e4=1:num_mics;
outliers(e4)=ismember(e4, npts);
if smd ~= 1 || num_mics > 1
if isequal(outliers(e4), 1)
legend_str{e4}=[legend_str{e4}, ' \oslash '];
else
legend_str{e4}=[legend_str{e4}, ' \surd '];
end
end
end
end
ptsa{count}=pts; % Cell array keeps track of all the indices of acceptable reverberation times
nptsa{count}=npts; % Cell array keeps track of all the indices of outlying reverberation times
rts=rt_array(pts); % The set of accepted reverberation times
% Calculate the geometric mean, standard deviation, and
% 95% confidence interval for the accepted reverberation times
mn_rt=geomean2(rts); % calculate the geometric mean
stdrt=std(rts); % calculate standard deviation
ci_int=call_ci(rts,95); % calculate 95% confidence interval of the standard error of the t-distribution with a two-sided test
% calculate the mean level of the signal during t_high
meanlevel=round(10*log10(mean(10.^(lv_array(pts, 1)./10))));
% Find the index of the median reverberation time for plotting
% purposes
[mbuf ix]=min(abs(rt_array-mn_rt));
if smd == 1
legend_str={legend_str{ix}};
end
figure(1);
% round to the nearest factor of 5 for plotting
% purposes
if plot_all_mics == 1
ymax=5*ceil(max(max(0.2*Lya)));
ymin=5*floor(min(min(0.2*Lya)));
else
ymax=5*ceil(1+max(max(0.2*Lya(ix, :))));
ymin=5*floor(-1+min(min(0.2*Lya(ix, :))));
end
ylim([ymin ymax]);
hv=0.05*(ymax-ymin);
if plot_all_mics == 1
% Plot the processed time averaged level of the time
% record.
for e5=1:num_mics;
e4=mod(e5-1, 12)+1;
plot(t_SP, Lya(e5, :), 'color', color_a{e4}, 'Linewidth', 1);
hold on;
end
% plot the solution of the reverberation time calculation
for e5=1:num_mics;
e4=mod(e5-1, 12)+1;
plot(t_SP, rtpa(e5, :), 'color', color_a{e4}, 'Linewidth', 2);
end
if smd == 1
% plot the first row after the header
plot(min(t_SP)+t_r*[0.62 0.65], (ymax-0.20*(ymax-ymin)-hv*(1-1))*[1 1], 'color', color_a{1}, 'Linewidth', 2);
else
for e5=1:num_mics; % plot a whie background for the colors of the lines
fill(min(t_SP)+t_r*[0.61 0.66 0.66 0.61 0.61], (ymax-0.20*(ymax-ymin)-hv*(e5-1))+[1-hv/2 1-hv/2 1+hv/2 1+hv/2 1-hv/2], [1 1 1], 'EdgeColor', 'none');
end
for e5=1:num_mics; % plot the colors of the lines
plot(min(t_SP)+t_r*[0.62 0.65], (ymax-0.20*(ymax-ymin)-hv*(e5-1))*[1 1], 'color', color_a{e5}, 'Linewidth', 2);
end
end
else
if smd == 1
% plot the median reverberation time
% add only the median mic data
plot(t_SP, Lya(ix, :), 'color', color_a{1}, 'Linewidth', 1);
hold on;
plot(t_SP, rtpa(ix, :), 'color', color_a{1}, 'Linewidth', 2);
fill(min(t_SP)+t_r*[0.61 0.66 0.66 0.61 0.61], (ymax-0.20*(ymax-ymin)-hv*(1-1))+[1-hv/2 1-hv/2 1+hv/2 1+hv/2 1-hv/2],[1 1 1], 'EdgeColor', 'none');
plot(min(t_SP)+t_r*[0.62 0.65], (ymax-0.20*(ymax-ymin)-hv*(1-1))*[1 1], 'color', color_a{1}, 'Linewidth', 2);
else
% plot the median reverberation time only add all of
% the mic data
plot(t_SP, Lya(ix, :), 'color', color_a{1}, 'Linewidth', 1);
hold on;
% plot the solution of the reverberation time calculations
plot(t_SP, rtpa(ix, :), 'color', color_a{1}, 'Linewidth', 2);
for e5=1:num_mics;
fill(min(t_SP)+t_r*[0.61 0.66 0.66 0.61 0.61], (ymax-0.20*(ymax-ymin)-hv*(e5-1))+[1-hv/2 1-hv/2 1+hv/2 1+hv/2 1-hv/2], [1 1 1], 'EdgeColor', 'none');
end
plot(min(t_SP)+t_r*[0.62 0.65], (ymax-0.20*(ymax-ymin)-hv*(ix-1))*[1 1], 'color', color_a{1}, 'Linewidth', 2);
end
end
% Give the pdf plot a title
% the time includes an estiamte of the reverberation time and
% if there are two or more mics then there will also be a two-sided 95%
% confidence interval
text(min(t_SP)+0.78*t_r, ymax-0.17*(ymax-ymin), title_str, 'Fontsize', 20, 'HorizontalAlignment', 'center', 'VerticalAlignment', 'bottom', 'BackgroundColor', 'w' );
% Give the pdf plot a legend
if smd == 1
text(min(t_SP)+0.66*t_r, ymax-0.20*(ymax-ymin)-hv*(1-1), legend_str{1}, 'Interpreter','Tex', 'Fontsize', 19, 'FontName','FixedWidth', 'HorizontalAlignment', 'left', 'VerticalAlignment', 'middle', 'BackgroundColor', 'w');
else
for e4=length(legend_str):(-1):1;
text(min(t_SP)+0.66*t_r, ymax-0.20*(ymax-ymin)-hv*(e4-1), legend_str{e4}, 'Interpreter','Tex', 'Fontsize', 19, 'FontName','FixedWidth', 'HorizontalAlignment', 'left', 'VerticalAlignment', 'middle', 'BackgroundColor', 'w');
end
end
set(gca, 'FontSize', 20);
xlabel('Time (Seconds)', 'Fontsize', 20);
ylabel('Instantaneous Amplitude ( dB ref. 20 \mu Pa)', 'Fontsize', 20);
if method == 1
title1=['Mean Level ', num2str(meanlevel), ' dBA, Freq = '];
else
title1=['Peak Level ', num2str(meanlevel), ' dB, Freq = '];
end
if num_mics > 1
title([title1, num2str(Fc(e1)), ' Hz, Reverb. Time = ', sprintf('%5.3f +- %5.3f (s)', [mn_rt ci_int])], 'Fontsize', 22);
else
title([title1, num2str(Fc(e1)), ' Hz, Reverb. Time = ', sprintf('%5.3f (s)', mn_rt)], 'Fontsize', 22);
end
rt_array2{count}=rt_array; % Collect the reverberation times
lv_array2{count}=lv_array; % Collect the A-weighted signal levels
lv_array_BG2{count}=lv_array_BG; % Collect the A-weighted background levels
cum_var_counter(count)=e1; % Cumulative variable Counter
hold on;
ylim([ymin ymax]);
% Give the PDF plot a file name with the format
% Room_Name_frequency_Hz_Var#
if e3 > 1
filename1=[ Room_Name, '_', num2str(floor(Fc(e1))), '_Hz', '_', num2str(e3) ];
else
filename1=[ Room_Name, '_', num2str(floor(Fc(e1))), '_Hz' ];
end
% Save the PDF plot
save_a_plot_reverb_time(1, filename1);
end
end
end
end
end
% Delete the PDF Plot from the screen.
% Speeds up processing and cleans up the screen.
delete(1);
% Open Output file to save the Reverberation Times, A-weighted signal
% Levels, and A-weighted background levels.
% Make sure outfile has a .txt extension
[filename_base, ext]=file_extension(outfile);
fid=fopen([filename_base '.txt'], 'w');
% % *********************************************************************
%Print Reverberation Times;
% % *********************************************************************
fprintf(fid, 'Reverberation Times\r\n');
fprintf(fid, 'Center Frequency\tMicrophone Number');
max_num_mics=max(mic_num_array);
tabs1='';
for e1=1:(max_num_mics);
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, 'Mean\tSTD\t95%%CI\tOutliers\r\n');
else
fprintf(fid, 'Mean\r\n');
end
fprintf(fid, 'Hz\t');
nums=1:(max_num_mics);
fprintf(fid, '%i\t', nums);
if num_mics > 1
fprintf(fid, '(s)\t(s)\t(s)\r\n');
else
fprintf(fid, '(s)\r\n');
end
for e3=1:count;
e1=cum_var_counter(count);
pts=ptsa{e3};
npts=nptsa{e3};
rts=rt_array2{e3}(pts);
rts_cal=rts(rts>0);
mn_rt=geomean2(rts_cal);
std_rt=std(rts_cal);
ci_int=call_ci(rts_cal,95);
fprintf(fid, '%i\t', floor(Fc(e3)));
fprintf(fid, '%s', sprintf('%6.3f\t', 0.001*round(1000*[rt_array2{e3}])) );
tabs1='';
for e2=1:(max_num_mics-mic_num_array(e3));
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, '%s', sprintf('%6.3f\t', 0.001*round(1000*[mn_rt, std_rt, ci_int])) );
for e2=1:length(npts);
fprintf(fid, '%d\t', npts(e2) );
end
else
fprintf(fid, '%s', sprintf('%6.3f\t', 0.001*round(1000*mn_rt)) );
end
fprintf(fid, '\r\n');
end
fprintf(fid, '\r\n\r\n\r\n');
% % *********************************************************************
%Print Signal High Levels;
% % *********************************************************************
if method == 1
fprintf(fid, 'Time Averaged Speaker on Sound Pressure Levels (dBA)\r\n');
else
fprintf(fid, 'Peak Sound Pressure Levels (dB)\r\n');
end
fprintf(fid, 'Center Frequency\tMicrophone Number');
max_num_mics=max(mic_num_array);
tabs1='';
for e1=1:(max_num_mics);
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, 'Mean\tSTD\t95%%CI\tOutliers\r\n');
else
fprintf(fid, 'Mean\r\n');
end
fprintf(fid, 'Hz\t');
nums=1:(max_num_mics);
fprintf(fid, '%i\t', nums);
if method == 1
if num_mics > 1
fprintf(fid, '(dBA)\t(dBA)\t(dBA)\r\n');
else
fprintf(fid, '(dBA)\r\n');
end
else
if num_mics > 1
fprintf(fid, '(dB)\t(dB)\t(dB)\r\n');
else
fprintf(fid, '(dB)\r\n');
end
end
for e3=1:count;
e1=cum_var_counter(count);
pts=ptsa{e3};
npts=nptsa{e3};
mn_rt=10*log10(mean(10.^(lv_array2{e3}(pts)./10)));
std_rt=std(lv_array2{e3}(pts));
ci_int=call_ci(lv_array2{e3}(pts),95);
fprintf(fid, '%i\t', floor(Fc(e3)));
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*[lv_array2{e3}])) );
tabs1='';
for e2=1:(max_num_mics-mic_num_array(e3));
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*[mn_rt, std_rt, ci_int])) );
for e2=1:length(npts);
fprintf(fid, '%d\t', npts(e2) );
end
else
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*mn_rt)) );
end
fprintf(fid, '\r\n');
end
fprintf(fid, '\r\n\r\n\r\n');
% % *********************************************************************
%Print Background Levels;
% % *********************************************************************
fprintf(fid, 'Time Averaged Background Sound Pressure Levels dBA\r\n');
fprintf(fid, 'Center Frequency\tMicrophone Number');
max_num_mics=max(mic_num_array);
tabs1='';
for e1=1:(max_num_mics);
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, 'Mean\tSTD\t95%%CI\tOutliers\r\n');
else
fprintf(fid, 'Mean\r\n');
end
fprintf(fid, 'Hz\t');
nums=1:(max_num_mics);
fprintf(fid, '%i\t', nums);
if num_mics > 1
fprintf(fid, '(dBA)\t(dBA)\t(dBA)\r\n');
else
fprintf(fid, '(dBA)\r\n');
end
for e3=1:count;
e1=cum_var_counter(count);
pts=ptsa{e3};
npts=nptsa{e3};
mn_rt=10*log10(mean(10.^(lv_array_BG2{e3}(pts)./10)));
std_rt=std(lv_array_BG2{e3}(pts));
ci_int=call_ci(lv_array_BG2{e3}(pts),95);
fprintf(fid, '%i\t', floor(Fc(e3)));
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*[lv_array_BG2{e3}])) );
tabs1='';
for e2=1:(max_num_mics-mic_num_array(e3));
tabs1=[tabs1, '\t'];
end
fprintf(fid, tabs1);
if num_mics > 1
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*[mn_rt, std_rt, ci_int])) );
for e2=1:length(npts);
fprintf(fid, '%d\t', npts(e2) );
end
else
fprintf(fid, '%s', sprintf('%6.1f\t', 0.1*round(10*mn_rt)) );
end
fprintf(fid, '\r\n');
end
% % *********************************************************************
% % Close the text file
% % *********************************************************************
fclose(fid);
fclose('all');
|
|