holger@0: % This is an example script to illustrate the use of the CAdaptInstrSpec
holger@0: % class. The script loads the spectra extracted two recordings of two
holger@0: % different violins. It then applied the CAdaptInstrSpec class to adapt the
holger@0: % spectra of the first violin to those of the second.
holger@0: 
holger@0: clearvars;
holger@0: close all;
holger@0: 
holger@0: addpath('./functions/');
holger@0: 
holger@0: 
holger@0: %% constants
holger@0: 
holger@0: % file paths of instrument spectra
holger@0: noteSpectraFilePaths = {'./instrSpectra/Violin.CQT.mat'; ...
holger@0:                         './instrSpectra/Violin2.CQT.mat'};
holger@0: numInstr = length(noteSpectraFilePaths);
holger@0: 
holger@0: % analysis
holger@0: tuningFreqInHz = 440;
holger@0: costFctName = 'LS';
holger@0: beta = 2; % only needed for plotting. should match the cost function above
holger@0: numIter = 5;
holger@0: 
holger@0: 
holger@0: 
holger@0: 
holger@0: %% load instrument spectra
holger@0: instrSpec = cell(numInstr,1);
holger@0: f0Idcs    = cell(numInstr,1);
holger@0: 
holger@0: for instrIdx = 1:numInstr
holger@0:     load(noteSpectraFilePaths{instrIdx}); % loads variables: noteSpectra, midiPitches, freqsInHz
holger@0:     instrSpec{instrIdx} = noteSpectra;
holger@0:     f0Idcs{instrIdx} = midiPitch2Shift(midiPitches, tuningFreqInHz, freqsInHz)+1;
holger@0: end
holger@0: numBinsPerSemitone = getNumBinsPerSemitoneFromFreqVec(freqsInHz);
holger@0: 
holger@0: 
holger@0: %% estimate filter curve
holger@0: 
holger@0: % create instrSpecFilter object
holger@0: instrFiltObj = CAdaptInstrSpec(instrSpec{1}, instrSpec{2}, f0Idcs{1}, f0Idcs{2}, numBinsPerSemitone, costFctName);
holger@0: 
holger@0: % apply update function for h
holger@0: for iterIdx = 1:numIter;
holger@0:     instrFiltObj = instrFiltObj.updateH;
holger@0: end
holger@0: 
holger@0: h = instrFiltObj.getH;
holger@0: h_smooth = instrFiltObj.getSmoothedH;
holger@0: 
holger@0: 
holger@0: %% estimate spectra with given filter curve
holger@0: estSpectra = instrFiltObj.estimateSpectra(f0Idcs{1});
holger@0: 
holger@0: 
holger@0: 
holger@0: 
holger@0: %% plot results
holger@0: 
holger@0: % plot original spectra and estimated spectra
holger@0: yticks = [250 500 1000 2000 4000];
holger@0: yticklabel = {'250'; '500'; '1k'; '2k'; '4k'};
holger@0: 
holger@0: figure;
holger@0: subplot(311);
holger@0: imagesc(f0Idcs{1}, freqsInHz, instrSpec{1});
holger@0: axis xy;
holger@0: set(gca, 'YScale', 'log', 'YTick', yticks, 'YTickLabel', yticklabel);
holger@0: xlabel('pitch index');
holger@0: ylabel('frequency [Hz]');
holger@0: title('spectra of instrument 1');
holger@0: 
holger@0: 
holger@0: subplot(312);
holger@0: imagesc(f0Idcs{1}, freqsInHz, instrSpec{2});
holger@0: axis xy;
holger@0: set(gca, 'YScale', 'log', 'YTick', yticks, 'YTickLabel', yticklabel);
holger@0: xlabel('pitch index');
holger@0: ylabel('frequency [Hz]');
holger@0: title('spectra of instrument 2');
holger@0: 
holger@0: subplot(313);
holger@0: imagesc(f0Idcs{1}, freqsInHz, estSpectra);
holger@0: axis xy;
holger@0: set(gca, 'YScale', 'log', 'YTick', yticks, 'YTickLabel', yticklabel);
holger@0: xlabel('pitch index');
holger@0: ylabel('frequency [Hz]');
holger@0: title('spectra of instrument 2 adapted to instrument 1')
holger@0: 
holger@0: 
holger@0: % plot filter curve and errors
holger@0: xticks = [125 250 500 1000 2000 4000];
holger@0: xticklabel = {'125'; '250'; '500'; '1k'; '2k'; '4k'};
holger@0: 
holger@0: % compute beta divergence for each element
holger@0: betaDivsOriginal = betaDivergencePerElement(instrSpec{1}, instrSpec{2}, beta);
holger@0: betaDivsEstimate = betaDivergencePerElement(instrSpec{1}, estSpectra,   beta);
holger@0: 
holger@0: % scale per-element beta divergence for image plot
holger@0: numColors = size(colormap,1);
holger@0: maxDist     = max([betaDivsOriginal(:); betaDivsEstimate(:)]);
holger@0: betaDivsOriginal = betaDivsOriginal / maxDist * numColors;
holger@0: betaDivsEstimate = betaDivsEstimate / maxDist * numColors;
holger@0: 
holger@0: 
holger@0: 
holger@0: figure;
holger@0: subplot(311)
holger@0: %stem(freqsInHz(h ~= 0), h(h ~= 0), 'k', 'filled');
holger@0: plot(freqsInHz(h ~= 0), h(h ~= 0), 'k.');
holger@0: hold on;
holger@0: plot(freqsInHz, h_smooth, 'k');
holger@0: hold off;
holger@0: axis tight;
holger@0: set(gca, 'XScale', 'log', 'XTick', xticks, 'XTickLabel', xticklabel);
holger@0: title('estimated filter curve ''h''');
holger@0: legend('original', 'smoothed', 'Location', 'NorthWest');
holger@0: 
holger@0: subplot(312);
holger@0: image(f0Idcs{1}, freqsInHz, betaDivsOriginal);
holger@0: axis xy;
holger@0: set(gca, 'YScale', 'log', 'YTick', yticks, 'YTickLabel', yticklabel);
holger@0: xlabel('pitch index');
holger@0: ylabel('frequency [Hz]');
holger@0: title('elementwise differences between original sets of spectra');
holger@0: 
holger@0: subplot(313);
holger@0: image(f0Idcs{1}, freqsInHz, betaDivsEstimate);
holger@0: axis xy;
holger@0: set(gca, 'YScale', 'log', 'YTick', yticks, 'YTickLabel', yticklabel);
holger@0: xlabel('pitch index');
holger@0: ylabel('frequency [Hz]');
holger@0: title('differences between original and adapted set of spectra');