ivan@138: function MissingSampleTopologyExperiment(expParam) ivan@138: % For a total number of missing samples C in a frame, create several ivan@138: % configuration of B holes with length A, where A*B=C (i.e. the total ivan@138: % number of missing samples is constant). Test several values of C, several ivan@138: % solvers. For each C, test all possible combination of (A,B) such that ivan@138: % A*B=C. ivan@138: % Note that for each combination (A,B), a number of frames are tested at ivan@138: % random and SNR results are then averaged. ivan@138: % ivan@138: % Usage: MissingSampleTopologyExperiment(expParam) ivan@138: % ivan@138: % ivan@138: % Inputs: ivan@138: % - expParam is an optional structure where the user can define ivan@138: % the experiment parameters. ivan@138: % - expParam.soundDir: path to sound directory. All the .wav files ivan@138: % in this directory will be tested at random. ivan@138: % - expParam.destDir: path to store the results. ivan@138: % - expParam.N: frame length ivan@138: % - expParam.NFramesPerHoleSize: number of frames to use for each ivan@138: % testing configuration (A,B). Results are then averaged. ivan@138: % - expParam.totalMissSamplesList: list of all tested values C for ivan@138: % the total number of missing samples in a frame ivan@138: % - expParam.solvers: list of solvers with their parameters ivan@138: % ivan@138: % ------------------- ivan@138: % ivan@138: % Audio Inpainting toolbox ivan@138: % Date: June 28, 2011 ivan@138: % By Valentin Emiya, Amir Adler, Maria Jafari ivan@138: % This code is distributed under the terms of the GNU Public License version 3 (http://www.gnu.org/licenses/gpl.txt). ivan@138: ivan@138: if ~isdeployed ivan@138: addpath('../../Problems/'); ivan@138: addpath('../../Solvers/'); ivan@138: addpath('../../Utils/'); ivan@138: addpath('../../Utils/dictionaries/'); ivan@138: addpath('../../Utils/evaluation/'); ivan@138: % addpath('../../Utils/TCPIP_SocketCom/'); ivan@138: % javaaddpath('../../Utils/TCPIP_SocketCom'); ivan@138: dbstop if error ivan@138: close all ivan@138: end ivan@138: ivan@138: %% Set parameters ivan@138: if nargin<1 ivan@138: expParam = []; ivan@138: end ivan@138: % Path to audio files ivan@138: if ~isfield(expParam,'soundDir'), ivan@138: expParam.soundDir = '../../Data/testSpeech8kHz_from16kHz/'; ivan@138: end ivan@138: if ~isfield(expParam,'destDir'), ivan@138: expParam.destDir = '../../tmp/missSampTopoExp/'; ivan@138: end ivan@138: if ~exist(expParam.destDir,'dir') ivan@138: mkdir(expParam.destDir) ivan@138: end ivan@138: ivan@138: ivan@138: % frame length ivan@138: if ~isfield(expParam,'N'), ivan@138: expParam.N = 512; ivan@138: expParam.N = 256; ivan@138: warning('AITB:N','Frame length=256 is used to have faster computations. Recommended frame length is 512 at 8kHz.'); ivan@138: end ivan@138: ivan@138: % Number of random frames to test ivan@138: if ~isfield(expParam,'NFramesPerHoleSize'), ivan@138: expParam.NFramesPerHoleSize = 20; ivan@138: warning('AITB:NFrames','expParam.NFramesPerHoleSize = 20 is used to have faster computations. Recommended value: several hundreds.'); ivan@138: end ivan@138: ivan@138: % Number of missing samples: which numbers to test? ivan@138: if ~isfield(expParam,'totalMissSamplesList') ivan@138: expParam.totalMissSamplesList = [12,36,60,120,180,240]; ivan@138: expParam.totalMissSamplesList = [12,36]; ivan@138: warning('AITB:Miss','expParam.totalMissSamplesList = [12,36] is used to have faster computations. Recommended list: expParam.totalMissSamplesList = [12,36,60,120,180,240].'); ivan@138: end ivan@138: ivan@138: % Choose the solver methods you would like to test: OMP, L1, Janssen ivan@138: if ~isfield(expParam,'solvers'), ivan@138: nSolver = 0; ivan@138: nSolver = nSolver+1; ivan@138: expParam.solvers(nSolver).name = 'OMP-G'; ivan@138: expParam.solvers(nSolver).inpaintFrame = @inpaintFrame_OMP_Gabor; % solver function ivan@138: expParam.solvers(nSolver).param.N = expParam.N; % frame length ivan@138: expParam.solvers(nSolver).param.OMPerr = 0.001; ivan@138: expParam.solvers(nSolver).param.sparsityDegree = expParam.solvers(nSolver).param.N/4; ivan@138: expParam.solvers(nSolver).param.D_fun = @Gabor_Dictionary; % Dictionary (function handle) ivan@138: expParam.solvers(nSolver).param.redundancyFactor = 2; % Dictionary redundancy ivan@138: expParam.solvers(nSolver).param.wa = @wRect; % Analysis window ivan@138: ivan@138: nSolver = nSolver+1; ivan@138: expParam.solvers(nSolver).name = 'Janssen'; ivan@138: expParam.solvers(nSolver).inpaintFrame = @inpaintFrame_janssenInterpolation; % solver function ivan@138: expParam.solvers(nSolver).param.N = expParam.N; % frame length ivan@138: end ivan@138: ivan@138: ivan@138: ivan@138: ivan@138: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ivan@138: soundDir = expParam.soundDir; ivan@138: wavFiles = dir([soundDir '*.wav']); ivan@138: wavFiles = arrayfun(@(x)[soundDir x.name],wavFiles,'UniformOutput',false); ivan@138: ivan@138: %% Draw a list of random frames ivan@138: % Choose an audio file at random ivan@138: frameParam.kFrameFile = randi(length(wavFiles),expParam.NFramesPerHoleSize); ivan@138: ivan@138: % For each audio file, find maximum mean energy among all frames ivan@138: [dum fs] = wavread([soundDir wavFiles{1}],'size'); ivan@138: Ne = round(512/16000*fs); ivan@138: E2m = zeros(length(wavFiles),1); ivan@138: for kf = 1:length(wavFiles) ivan@138: x=wavread(wavFiles{kf}); ivan@138: xm = filter(ones(Ne,1)/Ne,1,abs(x.^2)); ivan@138: E2m(kf) = 10*log10(max(xm)); ivan@138: end ivan@138: ivan@138: % Choose the location of a frame at random, with a minimum energy ivan@138: maxDiffE2m = 10; ivan@138: frameParam.kFrameBegin = NaN(expParam.NFramesPerHoleSize,1); ivan@138: for kf = 1:expParam.NFramesPerHoleSize ivan@138: siz = wavread(wavFiles{frameParam.kFrameFile(kf)},'size'); ivan@138: while true ivan@138: frameParam.kFrameBegin(kf) = randi(siz(1)-expParam.N+1); ivan@138: x = wavread(wavFiles{frameParam.kFrameFile(kf)},[0,expParam.N-1]+frameParam.kFrameBegin(kf)); ivan@138: E2m0 = 10*log10(mean(abs(x.^2))); ivan@138: if E2m(frameParam.kFrameFile(kf))-E2m0 <= maxDiffE2m ivan@138: break ivan@138: end ivan@138: end ivan@138: end ivan@138: ivan@138: %% Test each number of missing samples ivan@138: PerfRes = cell(length(expParam.totalMissSamplesList),length(expParam.solvers)); ivan@138: factorsToTest = cell(length(expParam.totalMissSamplesList),length(expParam.solvers)); ivan@138: outputFile = [expParam.destDir 'missSampTopoExp.mat']; ivan@138: for kSolver = 1:length(expParam.solvers) ivan@138: fprintf('\n ------ Solver: %s ------\n\n',... ivan@138: expParam.solvers(kSolver).name); ivan@138: for kMiss = 1:length(expParam.totalMissSamplesList) ivan@138: NMissSamples = expParam.totalMissSamplesList(kMiss); ivan@138: factorsToTest{kMiss} = allFactors(NMissSamples); ivan@138: PerfRes{kMiss,kSolver} = zeros([length(factorsToTest{kMiss}),expParam.NFramesPerHoleSize]); ivan@138: for kFactor = 1:length(factorsToTest{kMiss}) ivan@138: holeSize = factorsToTest{kMiss}(kFactor); ivan@138: NHoles = NMissSamples/holeSize; ivan@138: fprintf('%d %d-length holes (%d missing samples = %.1f%%)\n',... ivan@138: NHoles,holeSize,NMissSamples,NMissSamples/expParam.N*100) ivan@138: problemParameters.holeSize = holeSize; ivan@138: problemParameters.NHoles = NHoles; ivan@138: for kFrame = 1:expParam.NFramesPerHoleSize ivan@138: %% load audio frame ivan@138: xFrame = wavread(... ivan@138: wavFiles{frameParam.kFrameFile(kFrame)},... ivan@138: frameParam.kFrameBegin(kFrame)+[0,expParam.N-1]); ivan@138: ivan@138: %% generate problem ivan@138: [problemData, solutionData] = ... ivan@138: generateMissingGroupsProblem(xFrame,problemParameters); ivan@138: ivan@138: %% solve problem ivan@138: xEst = ... ivan@138: expParam.solvers(kSolver).inpaintFrame(... ivan@138: problemData,... ivan@138: expParam.solvers(kSolver).param); ivan@138: ivan@138: %% compute and store performance ivan@138: [SNRAll, SNRmiss] = ... ivan@138: SNRInpaintingPerformance(... ivan@138: solutionData.xClean,... ivan@138: problemData.x,... ivan@138: xEst,... ivan@138: problemData.IMiss); ivan@138: PerfRes{kMiss,kSolver}(kFactor,kFrame) = SNRmiss(2); ivan@138: ivan@138: end ivan@138: end ivan@138: save(outputFile,'PerfRes','expParam'); ivan@138: end ivan@138: end ivan@138: ivan@138: figure ivan@138: Nrows = floor(sqrt(length(expParam.solvers))); ivan@138: Ncols = ceil(sqrt(length(expParam.solvers))/Nrows); ivan@138: cmap = lines; ivan@138: for kSolver = 1:length(expParam.solvers) ivan@138: subplot(Nrows,Ncols,kSolver) ivan@138: hold on,grid on ivan@138: for kMiss = 1:length(expParam.totalMissSamplesList) ivan@138: plot(factorsToTest{kMiss},mean(PerfRes{kMiss,kSolver},2),... ivan@138: 'color',cmap(kMiss,:)); ivan@138: end ivan@138: title(expParam.solvers(kSolver).name) ivan@138: end ivan@138: return ivan@138: ivan@138: function m = allFactors(n) ivan@138: % Find the list of all factors (not only prime factors) ivan@138: ivan@138: primeFactors = factor(n); ivan@138: ivan@138: degrees = zeros(size(primeFactors)); ivan@138: ivan@138: for k=1:length(degrees) ivan@138: degrees(k) = sum(primeFactors==primeFactors(k)); ivan@138: end ivan@138: ivan@138: [primeFactors, I] = unique(primeFactors); ivan@138: degrees = degrees(I); ivan@138: ivan@138: D = (0:degrees(1)).'; ivan@138: for k=2:length(degrees) ivan@138: Dk = ones(size(D,1),1)*(0:degrees(k)); ivan@138: D = [repmat(D,degrees(k)+1,1),Dk(:)]; ivan@138: end ivan@138: ivan@138: m = unique(sort(prod((ones(size(D,1),1)*primeFactors).^D,2))); ivan@138: ivan@138: return