Commit 944914f8 authored by Markus Mueller-Trapet's avatar Markus Mueller-Trapet

get correct sweep rate for exponential and linear sweeps

parent 6935219f
......@@ -368,16 +368,29 @@ classdef itaMSTF < itaMSPlaybackRecord
function sweepRate = sweepRate(this,value)
% get the sweep rate of the excitation
%% sweep rate from analytic calculation, only using sweep parameters / PDI
nSamples = ita_nSamples( this.fftDegree );
% MMT: use nSamples-1 here to be conform with sweep calculation
% based on timeVector and chirp function
finalExcitationLength = (nSamples-1)/this.samplingRate - this.stopMargin;
sweepRate(1) = log2(this.finalFreqRange(2)/this.finalFreqRange(1))/finalExcitationLength;
%% sweep rate of analysis of excitation signal
sweepRate(2) = ita_sweep_rate(this.raw_excitation,[2000 this.samplingRate/3]);
sweepRate = nan(2,1);
switch this.type
case 'exp'
% equation is: f(t) = f0*2^(sweeprate*t)
% sweep rate from analytic calculation, only using sweep parameters
nSamples = ita_nSamples( this.fftDegree );
% use nSamples-1 here to be conform with sweep calculation
% based on timeVector and chirp function
finalExcitationLength = (nSamples-1)/this.samplingRate - this.stopMargin;
sweepRate(1) = log2(this.finalFreqRange(2)/this.finalFreqRange(1))/finalExcitationLength;
case 'lin'
% equation is: f(t) = f0*(1 + sweeprate*t)
% sweep rate from analytic calculation, only using sweep parameters
nSamples = ita_nSamples( this.fftDegree );
% use nSamples-1 here to be conform with sweep calculation
% based on timeVector and chirp function
finalExcitationLength = (nSamples-1)/this.samplingRate - this.stopMargin;
sweepRate(1) = (this.finalFreqRange(2)/this.finalFreqRange(1) - 1)/finalExcitationLength;
otherwise
end
sweepRate(2) = ita_sweep_rate(this.raw_excitation,'freqRange',[max(this.finalFreqRange(1),2000) min(this.finalFreqRange(2),this.samplingRate/3)],'type',this.type,'f0',this.finalFreqRange(1));
if exist('value','var')
sweepRate = sweepRate(value);
end
......
......@@ -9,47 +9,46 @@ function varargout = ita_sweep_rate(varargin)
% Note: A too short sweep or a too small frequency range may be problematic.
% <ITA-Toolbox>
% This file is part of the ITA-Toolbox. Some rights reserved.
% You can find the license for this m-file in the license.txt file in the ITA-Toolbox folder.
% This file is part of the ITA-Toolbox. Some rights reserved.
% You can find the license for this m-file in the license.txt file in the ITA-Toolbox folder.
% </ITA-Toolbox>
% Author: Alexandre Bleus -- alexandre.bleus@akustik.rwth-aachen.de
% Created: 12-March-2010
%% Get ITA Toolbox preferences
thisFuncStr = [upper(mfilename) ':'];
%% Initialisation
narginchk(1,2);
if nargin == 1
varargin{2} = [1000 4000];
elseif nargin~=2
error([thisFuncStr 'Two input arguments are required.'])
end
if ~isa(varargin{1},'itaAudio')
error([thisFuncStr 'Oh Lord! An audio object (Exponential Sweep) is required.'])
else
inSweep = varargin{1};
freq_vec=varargin{2};
end
sArgs = struct('pos1_inSweep','itaAudio', 'freqRange',[1000 4000],'type','exp','f0',[]);
[inSweep,sArgs] = ita_parse_arguments(sArgs,varargin);
%% Calculate the sweep rate for some points
freqjump = round((freq_vec(2)-freq_vec(1))/10);
freq_gdel = freq_vec(1) + [(1:8).' (2:9).'].*freqjump;
freqjump = round((sArgs.freqRange(2)-sArgs.freqRange(1))/10);
freq_gdel = sArgs.freqRange(1) + [(1:8).' (2:9).'].*freqjump;
gdel = ita_groupdelay_ita(inSweep);
sweep_rate = log(freq_gdel(:,2)./freq_gdel(:,1)) ./ log(2) ./ (gdel.freq2value(freq_gdel(:,2)) - gdel.freq2value(freq_gdel(:,1)));
sweep_rate(sweep_rate < 0) = sweep_rate(end);
switch sArgs.type
case 'exp'
sweep_rate = log2(freq_gdel(:,2)./freq_gdel(:,1)) ./ (gdel.freq2value(freq_gdel(:,2)) - gdel.freq2value(freq_gdel(:,1)));
case 'lin'
sweep_rate = (freq_gdel(:,2) - freq_gdel(:,1)) ./ (gdel.freq2value(freq_gdel(:,2)) - gdel.freq2value(freq_gdel(:,1)));
if isempty(sArgs.f0)
error('Need f0 for linear sweep rate');
else
sweep_rate = sweep_rate./sArgs.f0;
end
otherwise
error('Type can only be exp or lin');
end
sweep_rate(sweep_rate < 0) = nan;
%% Verify that the sweep rates are similar
meanRate = round(10*mean(sweep_rate));
meanRate = round(10*mean(sweep_rate,'omitnan'));
if any(abs(round(sweep_rate*10) - meanRate) > 5)
ita_verbose_info('Oh Lord! Are you sure that you are using an Exponential Sweep?',0);
end
%% Set Output
inSweep.userData.sweep_rate = mean(sweep_rate);
inSweep.userData.sweep_rate = mean(sweep_rate,'omitnan');
varargout(1) = {inSweep.userData.sweep_rate};
end
\ No newline at end of file
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment