wolffd@0: function varargout = mirregularity(orig,varargin)
wolffd@0: %   i = mirregularity(x) calculates the irregularity of a spectrum, i.e.,
wolffd@0: %       the degree of variation of the successive peaks of the spectrum.
wolffd@0: %   Specification of the definition of irregularity:
wolffd@0: %       mirregularity(...,'Jensen') is based on (Jensen, 1999),
wolffd@0: %           where the irregularity is the sum of the square of the 
wolffd@0: %           difference in amplitude between adjoining partials.
wolffd@0: %           (Default approach)
wolffd@0: %       mirregularity(...,'Krimphoff') is based on (Krimphoff et al., 1994),
wolffd@0: %           where the irregularity is the sum of the amplitude minus the 
wolffd@0: %           mean of the preceding, same and next amplitude.
wolffd@0: %   If the input x is not already a spectrum with peak extracted, the peak
wolffd@0: %       picking is performed prior to the calculation of the irregularity.
wolffd@0: %       In this case the 'Contrast' parameter used in mirpeaks can be
wolffd@0: %       modified, and is set by default to .1.
wolffd@0: %
wolffd@0: %   [Krimphoff et al. 1994] J. Krimphoff, S. McAdams, S. Winsberg, 
wolffd@0: %       Caracterisation du timbre des sons complexes. II Analyses
wolffd@0: %       acoustiques et quantification psychophysique. Journal de Physique 
wolffd@0: %       IV, Colloque C5, Vol. 4. 1994. 
wolffd@0: %   [Jensen, 1999] K. Jensen, Timbre Models of Musical Sounds, Ph.D.
wolffd@0: %       dissertation, University of Copenhagen, Rapport Nr. 99/7.
wolffd@0: 
wolffd@0: 
wolffd@0:         meth.type = 'String';
wolffd@0:         meth.default = 'Jensen';
wolffd@0:         meth.choice = {'Jensen','Krimphoff'};
wolffd@0:     option.meth = meth;
wolffd@0:     
wolffd@0:         cthr.key = 'Contrast';
wolffd@0:         cthr.type = 'Integer';
wolffd@0:         cthr.default = .01;
wolffd@0:     option.cthr = cthr;
wolffd@0: 
wolffd@0: specif.option = option;
wolffd@0: 
wolffd@0: varargout = mirfunction(@mirregularity,orig,varargin,nargout,specif,@init,@main);
wolffd@0: 
wolffd@0: 
wolffd@0: function [x type] = init(x,option)
wolffd@0: if not(isamir(x,'mirdata')) || isamir(x,'miraudio')
wolffd@0:     x = mirspectrum(x);
wolffd@0: end
wolffd@0: if not(haspeaks(x))
wolffd@0:     x = mirpeaks(x,'Reso','SemiTone','Contrast',option.cthr);  %% FIND BETTER
wolffd@0: end
wolffd@0: type = 'mirscalar';
wolffd@0: 
wolffd@0: 
wolffd@0: function o = main(x,option,postoption)
wolffd@0: if iscell(x)
wolffd@0:     x = x{1};
wolffd@0: end
wolffd@0: m = get(x,'PeakVal');
wolffd@0: p = get(x,'PeakPos');
wolffd@0: y = cell(1,length(m));
wolffd@0: for h = 1:length(m)
wolffd@0:     y{h} = cell(1,length(m{h}));
wolffd@0:     for k = 1:length(m{h})
wolffd@0:         y{h}{k} = zeros(size(m{h}{k}));
wolffd@0:         for j = 1:size(m{h}{k},3)
wolffd@0:             for l = 1:size(m{h}{k},2)
wolffd@0:                 state = warning('query','MATLAB:divideByZero');
wolffd@0:                 warning('off','MATLAB:divideByZero');
wolffd@0:                 mm = m{h}{k}{1,l,j};
wolffd@0:                 pp = p{h}{k}{1,l,j};
wolffd@0:                 [pp oo] = sort(pp); % Sort peaks in ascending order of x abscissae. 
wolffd@0:                 mm = mm(oo);
wolffd@0:                 if strcmpi(option.meth,'Jensen')
wolffd@0:                     y{h}{k}(1,l,j) = sum((mm(2:end,:)-mm(1:end-1,:)).^2)...
wolffd@0:                                     ./sum(mm.^2);
wolffd@0:                 elseif strcmpi(option.meth,'Krimphoff')
wolffd@0:                     avrg = filter(ones(3,1),1,mm)/3;
wolffd@0:                     y{h}{k}(1,l,j) = log10(sum(abs(mm(2:end-1,:)-avrg(3:end))));
wolffd@0:                 end
wolffd@0:                 warning(state.state,'MATLAB:divideByZero');
wolffd@0:                 if isnan(y{h}{k}(1,l,j))
wolffd@0:                     y{h}{k}(1,l,j) = 0;
wolffd@0:                 end
wolffd@0:             end
wolffd@0:         end
wolffd@0:     end
wolffd@0: end
wolffd@0: if isa(x,'mirspectrum')
wolffd@0:     t = 'Spectral irregularity';
wolffd@0: else
wolffd@0:     t = ['Irregularity of ',get(x,'Title')];;
wolffd@0: end
wolffd@0: i = mirscalar(x,'Data',y,'Title',t);
wolffd@0: o = {i,x};