function dif_tone( idlPath, idlFile, strPath, strFile )
%--------------------------------------------------------------------------
% Error Spectrogram of Tone File
%
% Cooper Baker - 2014
%--------------------------------------------------------------------------
close all;
% Settings
%--------------------------------------------------------------------------
winSize = 512;
overlap = 4;
window = chebwin( winSize, 125 );
name = 'Error Spectrogram';
headLoc = 0.25;
midLoc = 0.5;
tailLoc = 0.75;
detail = 0.01;
% Initializations
%--------------------------------------------------------------------------
if any( exist( 'idlFile' ) ~= 1 )
[ idlFile, idlPath ] = uigetfile( '*.wav', 'Ideal Audio File' );
[ strFile, strPath ] = uigetfile( '*.wav', 'Stretched Audio File' );
end
[ idlBuf, sr ] = audioread( [ idlPath, idlFile ] );
[ strBuf, sr ] = audioread( [ strPath, strFile ] );
idlSize = length( idlBuf );
strSize = length( strBuf );
halfWinSize = winSize / 2;
overlapSamps = winSize - ( winSize / overlap );
hopSize = winSize / overlap;
hopMax = length( idlBuf ) - winSize;
idlMsec = ( idlSize / sr ) * 1000;
gridMsec = 500;
gridHops = gridMsec / ( ( ( winSize / overlap ) / sr ) * 1000 );
hopMsec = ( hopSize / sr ) * 1000;
frames = ( idlSize / winSize ) * overlap - ( overlap - 1 );
rangeMsec = detail * idlMsec;
headMsec = ( headLoc - ( detail / 2 ) ) * idlMsec;
midMsec = ( midLoc - ( detail / 2 ) ) * idlMsec;
tailMsec = ( tailLoc - ( detail / 2 ) ) * idlMsec;
headMsec = headMsec + rangeMsec * 0.25;
midMsec = midMsec + rangeMsec * 0.25;
tailMsec = tailMsec + rangeMsec * 0.25;
winSum = sum( window );
% trim or pad stretch buffer to match ideal buffer length
if idlSize > strSize
strBuf = [ strBuf ; zeros( idlSize - strSize, 1 ) ];
elseif strSize > idlSize
strBuf = strBuf( 1 : length( idlBuf ) );
end
% Normalization
%--------------------------------------------------------------------------
% get file metadata
idlInfo = audioinfo( [ idlPath, idlFile ] );
strInfo = audioinfo( [ strPath, strFile ] );
% get file normalization coefficient
idlFileNorm = str2num( idlInfo.Title );
strFileNorm = str2num( strInfo.Title );
% de-normalize normalized .wav file data
idlRaw = idlBuf / idlFileNorm;
strRaw = strBuf / strFileNorm;
% make hann window
normWin = hann( round( idlSize / 2 ) );
% copy and window middle of files
idlNormWin = idlRaw( idlSize * 0.25 : idlSize * 0.75 - 1 ) .* normWin;
strNormWin = strRaw( idlSize * 0.25 : idlSize * 0.75 - 1 ) .* normWin;
% calculate rms of windowed file middles
idlRms = rms( idlNormWin );
strRms = rms( strNormWin );
% calculate normalization scalars
strScale = 1 / ( strRms / idlRms );
% scale stretch buffer to match ideal buffer
idlNorm = idlRaw;
strNorm = strRaw * strScale;
% Analysis
%--------------------------------------------------------------------------
idlGram = spectrogram( idlNorm, window, overlapSamps, winSize, sr );
strGram = spectrogram( strNorm, window, overlapSamps, winSize, sr );
% calculate magnitude spectra
idlMag = abs( idlGram );
strMag = abs( strGram );
% calculate amplitude spectra
idlAmp = ( idlMag / winSum );
strAmp = ( strMag / winSum );
% calculate sones spectra
idlSones = idlAmp .^ 0.6;
strSones = strAmp .^ 0.6;
% copy spectra to plot arrays
idlPlot = idlSones;
strPlot = strSones;
% calculate ideal plot offset
offset = 0;
offset = abs( min( [ min( idlPlot ) min( strPlot ) ] ) );
% offset plots
idlPlot = idlPlot + offset;
strPlot = strPlot + offset;
% calculate difference between plots
diffPlot = abs( idlPlot - strPlot );
% calculate cropping values
side = round( frames * detail );
head = round( frames * headLoc );
mid = round( frames * midLoc );
tail = round( frames * tailLoc );
range = side * 2;
% copy portions of spectrogram
diffHead = diffPlot( :, [ ( head - side ) : ( head + side ) ] );
diffMid = diffPlot( :, [ ( mid - side ) : ( mid + side ) ] );
diffTail = diffPlot( :, [ ( tail - side ) : ( tail + side ) ] );
% find min and max values of difference plot
diffMax = max( diffPlot(:) );
diffMin = min( diffPlot(:) );
% Plot
%--------------------------------------------------------------------------
fontName = 'Times New Roman';
fontSize = 12;
fig = figure( 1 );
set( fig, 'Name', name );
set( fig, 'Position', [ 0 0 800 500 ] );
set( fig, 'defaultAxesFontName', fontName );
set( fig, 'defaultTextFontName', fontName );
colormap( [ [ 1 : -0.2/64 : 0.8 ]; ( [ 1 : -0.5/64 : 0.5 ] .^ 10 ) ; [ 1 : -0.5/64 : 0.5 ] .^ 20 ]' );
% Main Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 1, 1 ), 'Position', [ 0.1, 0.5, 0.8, 0.3666 ] );
imagesc( diffPlot );
caxis( 'auto' );
axis xy;
% grid and labels
axis( [ 0 length( diffPlot ) 1 halfWinSize ] );
set( gca, 'XTick', [ 0 : gridHops : length( diffPlot ) ] );
set( gca, 'XTickLabel', [ 0 : gridMsec : idlMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [ 0 : 2000 : ( sr / 2 ) ] );
% legend
bar = colorbar( 'location', 'East' );
set( get( bar, 'YLabel'), 'String', 'Sones' );
set( get( bar, 'YLabel'), 'Position', [ 0.75 ( diffMax * 0.11 ) ] );
% Head Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 4 ), 'Position', [ 0.1, 0.1, 0.26, 0.3666 ] );
imagesc( diffHead );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ headMsec : rangeMsec / 4 : headMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [ 0 : 2000 : ( sr / 2 ) ] );
ylabel( 'Frequency' );
yLab = ylabel( 'Frequency' );
text( range * 0.05, halfWinSize * 0.93, ' Head ', 'EdgeColor', 'k', 'BackgroundColor', 'w' );
% Mid Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 5 ), 'Position', [ 0.3666, 0.1, 0.2666, 0.3666 ] );
imagesc( diffMid );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ midMsec : rangeMsec / 4 : midMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [] );
xlabel( 'Milliseconds' );
text( range * 0.5, halfWinSize * 0.93, ' Middle ', 'EdgeColor', 'k', 'BackgroundColor', 'w', 'HorizontalAlignment', 'center' );
% Tail Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 6 ), 'Position', [ 0.6399, 0.1, 0.26, 0.3666 ] );
imagesc( diffTail );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ tailMsec : rangeMsec / 4 : tailMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [] );
text( range * 1.02, halfWinSize * 0.93, ' Tail ', 'EdgeColor', 'k', 'BackgroundColor', 'w', 'HorizontalAlignment', 'right' );
% tighten up figure borders
tightfig();
set( yLab, 'Position', [ range * -0.15, ( halfWinSize * 1.05 ) ] );
% Save Data to Files
%--------------------------------------------------------------------------
% get file name
[ a, fileName, b ] = fileparts( strFile );
% write plot to file
hgexport( fig, [ fileName, '.spect.eps' ] );
% write error data to files
csvFile = sprintf( '%s.spect.csv', fileName );
csvwrite( csvFile, diffPlot );
% EOF
%--------------------------------------------------------------------------
% Error Spectrogram of Tone File
%
% Cooper Baker - 2014
%--------------------------------------------------------------------------
close all;
% Settings
%--------------------------------------------------------------------------
winSize = 512;
overlap = 4;
window = chebwin( winSize, 125 );
name = 'Error Spectrogram';
headLoc = 0.25;
midLoc = 0.5;
tailLoc = 0.75;
detail = 0.01;
% Initializations
%--------------------------------------------------------------------------
if any( exist( 'idlFile' ) ~= 1 )
[ idlFile, idlPath ] = uigetfile( '*.wav', 'Ideal Audio File' );
[ strFile, strPath ] = uigetfile( '*.wav', 'Stretched Audio File' );
end
[ idlBuf, sr ] = audioread( [ idlPath, idlFile ] );
[ strBuf, sr ] = audioread( [ strPath, strFile ] );
idlSize = length( idlBuf );
strSize = length( strBuf );
halfWinSize = winSize / 2;
overlapSamps = winSize - ( winSize / overlap );
hopSize = winSize / overlap;
hopMax = length( idlBuf ) - winSize;
idlMsec = ( idlSize / sr ) * 1000;
gridMsec = 500;
gridHops = gridMsec / ( ( ( winSize / overlap ) / sr ) * 1000 );
hopMsec = ( hopSize / sr ) * 1000;
frames = ( idlSize / winSize ) * overlap - ( overlap - 1 );
rangeMsec = detail * idlMsec;
headMsec = ( headLoc - ( detail / 2 ) ) * idlMsec;
midMsec = ( midLoc - ( detail / 2 ) ) * idlMsec;
tailMsec = ( tailLoc - ( detail / 2 ) ) * idlMsec;
headMsec = headMsec + rangeMsec * 0.25;
midMsec = midMsec + rangeMsec * 0.25;
tailMsec = tailMsec + rangeMsec * 0.25;
winSum = sum( window );
% trim or pad stretch buffer to match ideal buffer length
if idlSize > strSize
strBuf = [ strBuf ; zeros( idlSize - strSize, 1 ) ];
elseif strSize > idlSize
strBuf = strBuf( 1 : length( idlBuf ) );
end
% Normalization
%--------------------------------------------------------------------------
% get file metadata
idlInfo = audioinfo( [ idlPath, idlFile ] );
strInfo = audioinfo( [ strPath, strFile ] );
% get file normalization coefficient
idlFileNorm = str2num( idlInfo.Title );
strFileNorm = str2num( strInfo.Title );
% de-normalize normalized .wav file data
idlRaw = idlBuf / idlFileNorm;
strRaw = strBuf / strFileNorm;
% make hann window
normWin = hann( round( idlSize / 2 ) );
% copy and window middle of files
idlNormWin = idlRaw( idlSize * 0.25 : idlSize * 0.75 - 1 ) .* normWin;
strNormWin = strRaw( idlSize * 0.25 : idlSize * 0.75 - 1 ) .* normWin;
% calculate rms of windowed file middles
idlRms = rms( idlNormWin );
strRms = rms( strNormWin );
% calculate normalization scalars
strScale = 1 / ( strRms / idlRms );
% scale stretch buffer to match ideal buffer
idlNorm = idlRaw;
strNorm = strRaw * strScale;
% Analysis
%--------------------------------------------------------------------------
idlGram = spectrogram( idlNorm, window, overlapSamps, winSize, sr );
strGram = spectrogram( strNorm, window, overlapSamps, winSize, sr );
% calculate magnitude spectra
idlMag = abs( idlGram );
strMag = abs( strGram );
% calculate amplitude spectra
idlAmp = ( idlMag / winSum );
strAmp = ( strMag / winSum );
% calculate sones spectra
idlSones = idlAmp .^ 0.6;
strSones = strAmp .^ 0.6;
% copy spectra to plot arrays
idlPlot = idlSones;
strPlot = strSones;
% calculate ideal plot offset
offset = 0;
offset = abs( min( [ min( idlPlot ) min( strPlot ) ] ) );
% offset plots
idlPlot = idlPlot + offset;
strPlot = strPlot + offset;
% calculate difference between plots
diffPlot = abs( idlPlot - strPlot );
% calculate cropping values
side = round( frames * detail );
head = round( frames * headLoc );
mid = round( frames * midLoc );
tail = round( frames * tailLoc );
range = side * 2;
% copy portions of spectrogram
diffHead = diffPlot( :, [ ( head - side ) : ( head + side ) ] );
diffMid = diffPlot( :, [ ( mid - side ) : ( mid + side ) ] );
diffTail = diffPlot( :, [ ( tail - side ) : ( tail + side ) ] );
% find min and max values of difference plot
diffMax = max( diffPlot(:) );
diffMin = min( diffPlot(:) );
% Plot
%--------------------------------------------------------------------------
fontName = 'Times New Roman';
fontSize = 12;
fig = figure( 1 );
set( fig, 'Name', name );
set( fig, 'Position', [ 0 0 800 500 ] );
set( fig, 'defaultAxesFontName', fontName );
set( fig, 'defaultTextFontName', fontName );
colormap( [ [ 1 : -0.2/64 : 0.8 ]; ( [ 1 : -0.5/64 : 0.5 ] .^ 10 ) ; [ 1 : -0.5/64 : 0.5 ] .^ 20 ]' );
% Main Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 1, 1 ), 'Position', [ 0.1, 0.5, 0.8, 0.3666 ] );
imagesc( diffPlot );
caxis( 'auto' );
axis xy;
% grid and labels
axis( [ 0 length( diffPlot ) 1 halfWinSize ] );
set( gca, 'XTick', [ 0 : gridHops : length( diffPlot ) ] );
set( gca, 'XTickLabel', [ 0 : gridMsec : idlMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [ 0 : 2000 : ( sr / 2 ) ] );
% legend
bar = colorbar( 'location', 'East' );
set( get( bar, 'YLabel'), 'String', 'Sones' );
set( get( bar, 'YLabel'), 'Position', [ 0.75 ( diffMax * 0.11 ) ] );
% Head Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 4 ), 'Position', [ 0.1, 0.1, 0.26, 0.3666 ] );
imagesc( diffHead );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ headMsec : rangeMsec / 4 : headMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [ 0 : 2000 : ( sr / 2 ) ] );
ylabel( 'Frequency' );
yLab = ylabel( 'Frequency' );
text( range * 0.05, halfWinSize * 0.93, ' Head ', 'EdgeColor', 'k', 'BackgroundColor', 'w' );
% Mid Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 5 ), 'Position', [ 0.3666, 0.1, 0.2666, 0.3666 ] );
imagesc( diffMid );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ midMsec : rangeMsec / 4 : midMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [] );
xlabel( 'Milliseconds' );
text( range * 0.5, halfWinSize * 0.93, ' Middle ', 'EdgeColor', 'k', 'BackgroundColor', 'w', 'HorizontalAlignment', 'center' );
% Tail Spectrogram
%--------------------------------------------------------------------------
set( subplot( 2, 3, 6 ), 'Position', [ 0.6399, 0.1, 0.26, 0.3666 ] );
imagesc( diffTail );
caxis( [ diffMin diffMax ] );
axis xy;
set( gca, 'XTick', [ range * 0.25 : ( range / 4 ) : range * 0.75 ] );
set( gca, 'XTickLabel', [ tailMsec : rangeMsec / 4 : tailMsec + rangeMsec ] );
set( gca, 'YTick', [ 1 : ( halfWinSize / ( sr / 4000 ) ) : halfWinSize ] );
set( gca, 'YTickLabel', [] );
text( range * 1.02, halfWinSize * 0.93, ' Tail ', 'EdgeColor', 'k', 'BackgroundColor', 'w', 'HorizontalAlignment', 'right' );
% tighten up figure borders
tightfig();
set( yLab, 'Position', [ range * -0.15, ( halfWinSize * 1.05 ) ] );
% Save Data to Files
%--------------------------------------------------------------------------
% get file name
[ a, fileName, b ] = fileparts( strFile );
% write plot to file
hgexport( fig, [ fileName, '.spect.eps' ] );
% write error data to files
csvFile = sprintf( '%s.spect.csv', fileName );
csvwrite( csvFile, diffPlot );
% EOF