matthiasm@0: function [f_audio_out,timepositions_afterDegr] = degradationUnit_applyWowResampling(f_audio, samplingFreq, timepositions_beforeDegr, parameter) matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: % Name: degradation_applyWowResampling matthiasm@0: % Date of Revision: 2013-03 matthiasm@0: % Programmer: Sebastian Ewert matthiasm@0: % matthiasm@0: % Description: matthiasm@0: % - This is useful for wow and flutter simulations. matthiasm@0: % - a playback speed of 1 is modulated by a cos signal via matthiasm@0: % f(x) = a_m * cos(2*pi*f_m*x)+1 matthiasm@0: % - the function mapping a position in the original recording matthiasm@0: % to a position in the generated recording: matthiasm@0: % F(x) = x + a_m * sin(2*pi*f_m*x) / (2*pi*f_m) matthiasm@0: % - an optimal solution now would employ a sinc reconstruction of f_audio, matthiasm@0: % sinc_fa, and sampling sinc_fa(F^-1(y)) equidistantly. matthiasm@0: % - This implementation employs only an approximation by first upsampling matthiasm@0: % f_audio and then sampling that using nearest neighbors matthiasm@0: % matthiasm@0: % Input: matthiasm@0: % f_audio - audio signal \in [-1,1]^{NxC} with C being the number of matthiasm@0: % channels matthiasm@0: % samplingFreq - sampling frequency of f_audio matthiasm@0: % timepositions_beforeDegr - some degradations delay the input signal. If matthiasm@0: % some points in time are given via this matthiasm@0: % parameter, timepositions_afterDegr will matthiasm@0: % return the corresponding positions in the matthiasm@0: % output. Set to [] if unavailable. Set f_audio matthiasm@0: % and samplingFreq to [] to compute only matthiasm@0: % timepositions_afterDegr. matthiasm@0: % matthiasm@0: % Input (optional): parameter matthiasm@0: % .intensityOfChange = 1.5 - a_m = parameter.intensityOfChange/100; stay below 100 matthiasm@0: % .frequencyOfChange = 0.5 - f_m matthiasm@0: % .upsamplingFactor = 5 - the higher the better the quality (and the more memory) matthiasm@0: % .normalizeOutputAudio = 0 - normalize audio matthiasm@0: % matthiasm@0: % Output: matthiasm@0: % f_audio_out - audio signal \in [-1,1]^{NxC} with C being the number matthiasm@0: % of channels matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: % Audio Degradation Toolbox matthiasm@0: % matthiasm@0: % Centre for Digital Music, Queen Mary University of London. matthiasm@0: % This file copyright 2013 Sebastian Ewert, Matthias Mauch and QMUL. matthiasm@0: % matthiasm@0: % This program is free software; you can redistribute it and/or matthiasm@0: % modify it under the terms of the GNU General Public License as matthiasm@0: % published by the Free Software Foundation; either version 2 of the matthiasm@0: % License, or (at your option) any later version. See the file matthiasm@0: % COPYING included with this distribution for more information. matthiasm@0: % matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: % Check parameters matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: if nargin<4 matthiasm@0: parameter=[]; matthiasm@0: end matthiasm@0: if nargin<3 matthiasm@0: timepositions_beforeDegr=[]; matthiasm@0: end matthiasm@0: if nargin<2 matthiasm@0: error('Please specify input data'); matthiasm@0: end matthiasm@0: matthiasm@0: if isfield(parameter,'intensityOfChange')==0 matthiasm@0: parameter.intensityOfChange = 1.5; % a_m = parameter.intensityOfChange/100; %stay below 100 matthiasm@0: end matthiasm@0: if isfield(parameter,'frequencyOfChange')==0 matthiasm@0: parameter.frequencyOfChange = 0.5; % f_m matthiasm@0: end matthiasm@0: if isfield(parameter,'upsamplingFactor')==0 matthiasm@0: parameter.upsamplingFactor = 5; % the higher the better the quality (and the more memory) matthiasm@0: end matthiasm@0: if isfield(parameter,'normalizeOutputAudio')==0 matthiasm@0: parameter.normalizeOutputAudio = 0; matthiasm@0: end matthiasm@0: matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: % Main program matthiasm@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% matthiasm@0: fs = samplingFreq; matthiasm@0: fsOverSampled = fs * parameter.upsamplingFactor; matthiasm@0: a_m = parameter.intensityOfChange / 100; matthiasm@0: f_m = parameter.frequencyOfChange; matthiasm@0: matthiasm@0: f_audio_out = f_audio; matthiasm@0: if ~isempty(f_audio) matthiasm@0: numSamples = size(f_audio,1); matthiasm@0: numChannels = size(f_audio,2); matthiasm@0: length_sec = numSamples/fs; matthiasm@0: numFullPeriods = floor(length_sec * f_m); matthiasm@0: numSamplesToWarp = round(numFullPeriods *fs / f_m); matthiasm@0: matthiasm@0: oldSamplePositions_to_newOversampledPositions = round(timeAssignment_newToOld([1:numSamplesToWarp]/fs,a_m,f_m)*fsOverSampled); matthiasm@0: matthiasm@0: for ch=1:numChannels matthiasm@0: audioUpsampled = resample(f_audio(:,ch),parameter.upsamplingFactor,1); matthiasm@0: matthiasm@0: f_audio_out(1:numSamplesToWarp,ch) = audioUpsampled(oldSamplePositions_to_newOversampledPositions); matthiasm@0: matthiasm@0: clear audioUpsampled matthiasm@0: end matthiasm@0: matthiasm@0: if parameter.normalizeOutputAudio matthiasm@0: f_audio_out = adthelper_normalizeAudio(f_audio_out, samplingFreq); matthiasm@0: end matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: %pitchShifted_semitones = getPitchShiftAtPositionsInNewFile([0:0.1:length_sec],a_m,f_m); % not used yet but might be useful later on matthiasm@0: matthiasm@0: %selftest(length_sec,a_m,f_m); matthiasm@0: matthiasm@0: % This degradation does impose a complex temporal distortion matthiasm@0: timepositions_afterDegr = []; matthiasm@0: if ~isempty(timepositions_beforeDegr) matthiasm@0: timepositions_afterDegr = timeAssignment_oldToNew(timepositions_beforeDegr,a_m,f_m); matthiasm@0: end matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: matthiasm@0: function timeAssigned = timeAssignment_oldToNew(x,a_m,f_m) matthiasm@0: % vectorized matthiasm@0: matthiasm@0: timeAssigned = x + a_m * sin(2*pi*f_m*x) / (2*pi*f_m); matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: function timeAssigned = timeAssignment_newToOld(y,a_m,f_m) matthiasm@0: % vectorized matthiasm@0: % SE: non linear equation: solved numerically (we are lucky: matthiasm@0: % As G(x) := x - (F(x)-y) is a contraction for a_m<1 we can simply employ the Banach fixpoint theorem) matthiasm@0: matthiasm@0: timeAssigned = y; matthiasm@0: for k=1:40 matthiasm@0: timeAssigned = y - a_m * sin(2*pi*f_m*timeAssigned) / (2*pi*f_m); matthiasm@0: end matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: function pitchShifted_semitones = getPitchShiftAtPositionsInNewFile(y,a_m,f_m) matthiasm@0: % vectorized matthiasm@0: matthiasm@0: % maps times from new to old and takes a look at the value of the derivative matthiasm@0: % there matthiasm@0: matthiasm@0: timeAssigned = timeAssignment_newToOld(y,a_m,f_m); matthiasm@0: matthiasm@0: derivative = 1 + a_m * cos(2*pi*f_m*timeAssigned); matthiasm@0: matthiasm@0: pitchShifted_semitones = log2(1./derivative)*12; matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: function selftest(length_sec,a_m,f_m) matthiasm@0: matthiasm@0: x0 = rand(10000,1) * length_sec; matthiasm@0: matthiasm@0: timeAssigned1 = timeAssignment_oldToNew(x0,a_m,f_m); matthiasm@0: timeAssigned2 = timeAssignment_newToOld(timeAssigned1,a_m,f_m); matthiasm@0: matthiasm@0: figure; matthiasm@0: plot(x0 - timeAssigned2); matthiasm@0: matthiasm@0: end matthiasm@0: matthiasm@0: matthiasm@0: matthiasm@0: matthiasm@0: matthiasm@0: