V_HORIZDIFF - Estimates the horizontal difference between two functions of x Usage: (1) z=v_horizdiff(y,v,x); % Approximately: y(x) = v(x+z) (2) SNRgain=v_horizdiff(PESQenh,PESQraw,SNRs); % SNRs are the test SNRs, PESQraw and PESQenh are the % quatlity metrics of the raw and enhanced signal, SNRgain % is the effective SNR improvement of the enhancer. (3) x=(1:10)'; v_horizdiff(1.4*x,x,x); % plots illustrative example Inputs: y(n,m) each column is a function of x v(k,1) reference function of u (this will be extrapolated if necessary) x(n,1) x values for y [default: x=(1:n)'] u(k,1) x values for v [default: v=x] q interpolation mode 'l' linear [default] 'p' pchip (n>=2) 's' spline (n>=4) Outputs: z(n,m) horizontal difference between v and y zm(1,m) MMSE horizontal difference that minimizes sum((y(x)-v(x+z)).^2)
0001 function [z,zm]=v_horizdiff(y,v,x,u,q) 0002 %V_HORIZDIFF - Estimates the horizontal difference between two functions of x 0003 % 0004 % Usage: (1) z=v_horizdiff(y,v,x); % Approximately: y(x) = v(x+z) 0005 % 0006 % (2) SNRgain=v_horizdiff(PESQenh,PESQraw,SNRs); % SNRs are the test SNRs, PESQraw and PESQenh are the 0007 % % quatlity metrics of the raw and enhanced signal, SNRgain 0008 % % is the effective SNR improvement of the enhancer. 0009 % 0010 % (3) x=(1:10)'; v_horizdiff(1.4*x,x,x); % plots illustrative example 0011 % 0012 % Inputs: y(n,m) each column is a function of x 0013 % v(k,1) reference function of u (this will be extrapolated if necessary) 0014 % x(n,1) x values for y [default: x=(1:n)'] 0015 % u(k,1) x values for v [default: v=x] 0016 % q interpolation mode 0017 % 'l' linear [default] 0018 % 'p' pchip (n>=2) 0019 % 's' spline (n>=4) 0020 % 0021 % Outputs: z(n,m) horizontal difference between v and y 0022 % zm(1,m) MMSE horizontal difference that minimizes sum((y(x)-v(x+z)).^2) 0023 % 0024 0025 % Copyright (C) Mike Brookes 2012-2024 0026 % Version: $Id: v_axisenlarge.m 10865 2018-09-21 17:22:45Z dmb $ 0027 % 0028 % VOICEBOX is a MATLAB toolbox for speech processing. 0029 % Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html 0030 % 0031 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0032 % This program is free software; you can redistribute it and/or modify 0033 % it under the terms of the GNU General Public License as published by 0034 % the Free Software Foundation; either version 2 of the License, or 0035 % (at your option) any later version. 0036 % 0037 % This program is distributed in the hope that it will be useful, 0038 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0039 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0040 % GNU General Public License for more details. 0041 % 0042 % You can obtain a copy of the GNU General Public License from 0043 % http://www.gnu.org/copyleft/gpl.html or by writing to 0044 % Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA. 0045 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0046 % 0047 % 0048 % choose interpolation method 0049 % 0050 [n,m]=size(y); 0051 if nargin<5 || isempty(m) 0052 q=''; 0053 end 0054 if n>=4 && any(q=='s') 0055 im='spline'; 0056 elseif n>=2 && (any(q=='s') || any(q=='p')) 0057 im='pchip'; 0058 else 0059 im='linear'; 0060 end 0061 if nargin<3 || isempty(x) 0062 x=(1:n)'; 0063 end 0064 if nargin<4 || isempty(u) 0065 u=x; 0066 end 0067 % 0068 % Interpolate inverse function, u(v) 0069 % 0070 z=interp1(v,u,y,im,'extrap')-repmat(x,1,m); 0071 if nargout>1 0072 zm=zeros(1,m); 0073 ff=@(zm,c)sum((c-interp1(u,v,x+zm,im,'extrap')).^2); 0074 for i=1:m 0075 zm(i)=fminbnd(@(zm) ff(zm,y(:,i)),u(1)-x(end),u(end)-x(1)); 0076 end 0077 end 0078 if ~nargout 0079 ax=subplot(212); 0080 if m==1 0081 plot(x,z,'-r'); 0082 else 0083 plot(x,z); 0084 end 0085 xlabel('x'); 0086 ylabel('\Delta{x} Gain'); 0087 if all(abs(z-z(1))<=z(1)*1e-4) 0088 set(gca,'ylim',2*[min(z(1),0) max(z(1),0)]); 0089 end 0090 ax(2)=subplot(211); 0091 plot(x,y); 0092 hold on 0093 legax=plot(u,v,':k'); 0094 if m==1 0095 plot([x x+z]',[y y]','-r',x+z,y,'>r'); 0096 end 0097 hold off 0098 linkaxes(ax,'x'); 0099 ylabel('y(x) and v(x)'); 0100 legend(legax,'Reference, v(x)','location','best'); 0101 end