Home > voicebox > v_ppmvu.m

v_ppmvu

PURPOSE ^

V_PPMVU calculate PPM, VU or EBU level of an audio signal [V,FX,FX1]=(X,FSX,M)

SYNOPSIS ^

function [v,fx,fx1]=v_ppmvu(x,fsx,m)

DESCRIPTION ^

V_PPMVU calculate PPM, VU or EBU level of an audio signal [V,FX,FX1]=(X,FSX,M)

 Usage: (1) v=v_ppmvu(x,fs,'a')      % calculate PPM level of signal x with sampling freq fs
        (2) [v,f]=v_ppmvu(x,fs,'aw') % calculate PPM + fast version as well
        (3) v=v_ppmvu(x,fs,'au')     % calculate VU level in linear units rather than dB
        (4) v=v_ppmvu(x,fs,'c')      % calculate VU level
        (5) v=v_ppmvu(x,fs,'e')      % calculate EBU loudness level
        (6) [v1,fx]=v_ppmvu(x1,fs,'a')   
                 v2=v_ppmvu(x2,fx)   % process in chunks: same as  v_ppmvu([x1; x2],fs,'a')

  Inputs: x = input signal, one **column** per channel
          fsx = sample frequency or fx output from a previous call
          m = either a structure with algorithm parameters (see below)
              or an attack/decay time constsant or a character string:
              'a' UK PPM characteristic [default]
              'b' DIN PPM
              'c' VU American
              'd' VU French
              'e' EBU short   [default toggle: 'q']
              'f' EBU medium  [default toggle: 'q']
            followed by any combination of modifier toggles:
              'z' remove mean (not yet implemented)
              'p' preemphasis (not yet implemented)
              's' set rise time to zero
              'o' oversample x 4
              'w' give fast output (with 0 decay time) as well as slow output
              'v' plot graph
              'q' average squared signal
              'u' output magnitude (or mean square value if 'q') instead of dB

 Outputs:  y = selected output (same size as x)
           fx = cell array holding algorithm state
                or, if 'w' option specified, the fast version of the output
           fx1 = cell array holding algorithm state (only if 'w' specified)

 Algorithm Parameters:
           mm   text string with options
           ta   attach time constant or zero if no attack smoothing (seconds)
           tm   averaging filter duration or zero (seconds)
           td   decay time constant or zero if no decay smoothing (seconds)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [v,fx,fx1]=v_ppmvu(x,fsx,m)
0002 %V_PPMVU calculate PPM, VU or EBU level of an audio signal [V,FX,FX1]=(X,FSX,M)
0003 %
0004 % Usage: (1) v=v_ppmvu(x,fs,'a')      % calculate PPM level of signal x with sampling freq fs
0005 %        (2) [v,f]=v_ppmvu(x,fs,'aw') % calculate PPM + fast version as well
0006 %        (3) v=v_ppmvu(x,fs,'au')     % calculate VU level in linear units rather than dB
0007 %        (4) v=v_ppmvu(x,fs,'c')      % calculate VU level
0008 %        (5) v=v_ppmvu(x,fs,'e')      % calculate EBU loudness level
0009 %        (6) [v1,fx]=v_ppmvu(x1,fs,'a')
0010 %                 v2=v_ppmvu(x2,fx)   % process in chunks: same as  v_ppmvu([x1; x2],fs,'a')
0011 %
0012 %  Inputs: x = input signal, one **column** per channel
0013 %          fsx = sample frequency or fx output from a previous call
0014 %          m = either a structure with algorithm parameters (see below)
0015 %              or an attack/decay time constsant or a character string:
0016 %              'a' UK PPM characteristic [default]
0017 %              'b' DIN PPM
0018 %              'c' VU American
0019 %              'd' VU French
0020 %              'e' EBU short   [default toggle: 'q']
0021 %              'f' EBU medium  [default toggle: 'q']
0022 %            followed by any combination of modifier toggles:
0023 %              'z' remove mean (not yet implemented)
0024 %              'p' preemphasis (not yet implemented)
0025 %              's' set rise time to zero
0026 %              'o' oversample x 4
0027 %              'w' give fast output (with 0 decay time) as well as slow output
0028 %              'v' plot graph
0029 %              'q' average squared signal
0030 %              'u' output magnitude (or mean square value if 'q') instead of dB
0031 %
0032 % Outputs:  y = selected output (same size as x)
0033 %           fx = cell array holding algorithm state
0034 %                or, if 'w' option specified, the fast version of the output
0035 %           fx1 = cell array holding algorithm state (only if 'w' specified)
0036 %
0037 % Algorithm Parameters:
0038 %           mm   text string with options
0039 %           ta   attach time constant or zero if no attack smoothing (seconds)
0040 %           tm   averaging filter duration or zero (seconds)
0041 %           td   decay time constant or zero if no decay smoothing (seconds)
0042 
0043 %      Copyright (C) Mike Brookes 2013
0044 %      Version: $Id: v_ppmvu.m 3387 2013-08-23 12:32:47Z dmb $
0045 %
0046 %   VOICEBOX is a MATLAB toolbox for speech processing.
0047 %   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
0048 %
0049 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0050 %   This program is free software; you can redistribute it and/or modify
0051 %   it under the terms of the GNU General Public License as published by
0052 %   the Free Software Foundation; either version 2 of the License, or
0053 %   (at your option) any later version.
0054 %
0055 %   This program is distributed in the hope that it will be useful,
0056 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0057 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0058 %   GNU General Public License for more details.
0059 %
0060 %   You can obtain a copy of the GNU General Public License from
0061 %   http://www.gnu.org/copyleft/gpl.html or by writing to
0062 %   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
0063 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0064 
0065 % define persistent constants
0066 persistent m4 k4 h1 g1 tamd deft p0
0067 if isempty(m4)
0068     % upsample by 4 using two cascaded half-band filters
0069     m4=5;   % one quarter of the first filter order
0070     h1f=kaiser(4*m4+1,5.3)'.*sinc(-m4:0.5:m4);
0071     h1=h1f(2:2:4*m4);
0072     k4=3;
0073     g1f=kaiser(4*k4+1,5.3)'.*sinc(-k4:0.5:k4);
0074     g1=g1f(2:2:4*k4);
0075     tamd=[0.00632 0 1.01335; % 'a' UK PPM characteristic [default]
0076         0.00316 0 0.73830; % 'b' DIN PPM
0077         0.10433 0 0; % 'c' VU American
0078         0.13089 0 0; % 'd' VU French
0079         0 0.4 0; % 'e' EBU momentary
0080         0 3 0]; % 'f' EBU short
0081     deft={'','','','','q','q'};
0082     p0.mm='a';
0083     p0.ta=0;
0084     p0.tm=0;
0085     p0.td=0;
0086 end
0087 % decode input arguments
0088 if iscell(fsx)
0089     fx=fsx;
0090     [nx,nc]=size(x);
0091     fsp=fx{1};
0092     fsx=fsp(1);
0093     fy=fsp(2);
0094     nt=fsp(3);
0095     nm=fsp(4);
0096     p=fx{2};
0097     mop=fx{3};
0098 else
0099     p=p0;  % default values
0100     if nargin<3 || isempty(m)
0101         m='a';
0102     end
0103     if ischar(m)
0104         p.mm=m;
0105         ix=m(1)-'a'+1;
0106         if ix<=0 || ix>6
0107             ix=1;
0108         end
0109         p.ta=tamd(ix,1);
0110         p.tm=tamd(ix,2);
0111         p.td=tamd(ix,3);
0112         p.mm=[p.mm deft{ix}]; % add in default toggles
0113     elseif isstruct(m)
0114         p=m;
0115     else
0116         p.mm=' ';
0117         p.ta=m;
0118         p.tm=0;
0119         p.td=m;
0120     end
0121     [nx,nc]=size(x);
0122     mop=rem(sum(repmat('zpsowvqu',length(p.mm),1)-p.mm(ones(8,1),:)'==0,1),2); % set option toggles
0123     fy=fsx*(1+3*mop(4)); % effective sample rate
0124     nt=round(fy*0.1);       % 0.1 s resolution in samples fro MA filter
0125     nm=round(fy*p.tm/nt);   % number of blocks for MA filter
0126     fsp=[fsx fy nt nm];
0127     fx={fsp p mop zeros(m4,nc),zeros(2*m4-1,nc),zeros(k4,nc),zeros(2*k4-1,nc), zeros(1,nc), 0, zeros(nm+1,nc), zeros(0,nc)};
0128 end
0129 
0130 % Stage 1: oversampling we use two cascaded half-band filters to upsample by 4
0131 
0132 if mop(4)
0133     nx2=2*nx;
0134     v=zeros(nx2,nc);
0135     v(1:2:2*m4,:)=fx{4};
0136     v(2*m4+1:2:nx2,:)=x(1:nx-m4,:);
0137     fx{4}=x(nx-m4+1:nx,:);
0138     [v(2:2:nx2,:),fx{5}]=filter(h1,1,x,fx{5});   % delayed by 2*m4 samples
0139     nx4=4*nx;
0140     y=zeros(nx4,nc);
0141     y(1:2:2*k4,:)=fx{6};
0142     y(2*k4+1:2:nx4-1,:)=v(1:nx2-k4,:);
0143     fx{6}=v(nx2-k4+1:nx2,:);
0144     [y(2:2:nx4,nc),fx{7}]=filter(g1,1,v,fx{7});   % delayed by 4*m4+2*k4 samples
0145 else
0146     y=x;
0147 end
0148 ty=1/fy;
0149 % Stage 2 make positive
0150 if mop(7)
0151     y=y.^2;
0152 else
0153     y=abs(y);
0154 end
0155 % Stage 2: attack filter
0156 if p.ta>0 && ~mop(3)
0157     za=exp(-ty/p.ta);
0158     [y,fx{8}]=filter(1-za,[1 -za],y,fx{8});
0159 end
0160 
0161 % Stage 3: moving average filter
0162 if nm>0
0163     ny=size(y,1);
0164     y=cumsum(y,1);
0165     nr=fx{9};
0166     jj=nt-nr:nt:ny; % end of frame indices
0167     kk=length(jj);
0168     if kk>0
0169         nmkk=nm+kk;
0170         v=zeros(nmkk,nc);
0171         v(nm+2:nmkk,:)=y(jj(2:kk),:)-y(jj(1:kk-1),:); % v(:,nc) contains the sum of 0.1 second blocks
0172         v(nm+1,:)=y(jj(1),:)+fx{10}(2,:);
0173         fx{10}(2,:)=y(ny,:)-y(jj(kk),:);
0174         v(2:nm,:)=fx{10}(3:nm+1,:);   % saved values
0175         fx{10}(3:nm+1,:)=v(2+kk:nmkk,:);
0176         v=cumsum(v,1);
0177         v(nm+1:nmkk,:)=v(nm+1:nmkk,:)-v(1:kk,:); % perform MA filter
0178         v(nm,:)=fx{10}(1,:);                % final MA output from previous chunk
0179         y=v(nm+floor((nr+1:nr+ny)/nt),:)/(nm*nt);   % copy MA ouptuts into output buffer
0180         fx{10}(1,:)=v(nmkk,:);                % save the final MA filter output for the next chunk
0181         fx{9}=ny-jj(kk);                    % update the length of the tail
0182     else                                    % no completed 100 ms blocks in this chunk
0183         fx{9}(1)=ny+nr;                     % update the length of the tail
0184         fx{10}(2,:)=y(ny,:)+fx{10}(2,:);    % update the sum of the tail
0185         y=fx{10}(ones(ny,1),:);             % set outputs to the final MA output from previous chunk
0186     end
0187 end
0188 % Stage 4: decay filter
0189 if p.td>0
0190     [v,vk,fx{11}]=maxfilt(y,exp(-ty/p.td),Inf,1,fx{11});
0191 else
0192     v=y;
0193 end
0194 
0195 % Stage 5: decimate and ouput
0196 if mop(5)    % if 'w' option specified, output y as well
0197     if mop(4) % if oversampling, we need to decimate
0198         y=y(1:4:end,:);
0199         v=v(1:4:end,:);
0200     end
0201     if ~mop(8)
0202         y=(20-10*mop(7))*log10(y);
0203         v=(20-10*mop(7))*log10(v);
0204     end
0205     fx1=fx;
0206     fx=y;
0207 else
0208     if mop(4) % if oversampling, we need to decimate
0209         v=v(1:4:end,:);
0210     end
0211     if ~mop(8)
0212         v=(20-10*mop(7))*log10(v);
0213     end
0214 end
0215 if mop(6) || nargout==0
0216     t=(1:length(x))/fsx;
0217     ax(1)=subplot(211);
0218     plot(t,x,'-b');
0219     if mop(7)>0 && mop(8)>0
0220         v=sqrt(v);
0221         if mop(5)
0222             y=sqrt(y);
0223         end
0224     end
0225     ax(2)=subplot(212);
0226     if mop(5)    % if 'w' option specified, output y as well
0227         plot(t,v,'-b',t,y,'-r');
0228     else
0229         plot(t,v,'-b');
0230     end
0231     if ~mop(8)
0232         set(gca,'ylim',[max(min(v)-1,max(v)-40) max(v)+1]);
0233         ylabel('dB');
0234     elseif mop(7)
0235         ylabel('x_{RMS}');
0236     else
0237         ylabel('|x|');
0238     end
0239     linkaxes(ax,'x');
0240     xlabel(['Time (' xticksi 's)']);
0241 end

Generated on Tue 10-Oct-2017 08:30:10 by m2html © 2003