Code covered by the BSD License  

Highlights from
deOliveira wavelets

deOliveira wavelets

by

Hélio de Oliveira

 

cdeo is the complex wavelet (type 5), deo is an orthogonal infinitely regular wavelet (type 3)

cdeowavf(varargin);
function varargout = cdeowavf(varargin);
% cdeo  wavelet -- TYPE 5 WAVELET according to waveinfo
%   [PSI,T] = cdeo(LOWB,UPPB,N,alpha) returns deO
%   complex wavelet function evaluated on an
%   N point regular grid on the interval [LOWB,UPPB].
%   Output arguments are the wavelet function PSI
%   computed on the grid T (defined by LOWB, UPPB and N)
%
%   These functions have [-12 12] as effective support.
%
%   N must be a power of two.
%   See also CDEOINFO, WAVEFUN, WAVEINFO.
%
%   Luciana R. Soares, H.M. de Oliveira.
%   Last Revision: 10-Dec-2002

% Check arguments
if errargn(mfilename,nargin,[4],nargout,[0:2]), error('*'), end

% Translate input arguments...

    [lowb,uppb,N,label] = deal(varargin{:});
    ind = strncmpi('cdeo',label,4);
    if isequal(ind,1) label([1:4]) = []; end
    alpha = wstr2num(deblank(label));
 
% Check alpha factor (roll-off)
if isempty(alpha) error('** deO: Invalid wavelet number!'); end
if (alpha < 0) | (alpha > 1/3)
    error('** deO: Invalid value for alpha **')
end
if (alpha == 0)
   alpha = 10.^-5; % alpha ~ 0
   disp('deO: Shannon square root Wavelet');
end

% Check N
if errargt(mfilename,log(N)/log(2),'int')
    error('** deO: Invalid value for N **')
end

% Check Interval
if errargt(mfilename,uppb-lowb,'re0')
   error('** deO: Verify upper/lower bound... **')
end

% Transform interval bounds to grid
lint = (uppb-lowb);
step = lint/N;
if step > 0.3
   disp('** deO: the number of grid points is not enough for that support... Be carefull **');
end
t  = [lowb:step:uppb];
td = t - 1/2;

% Initial Variables
va = 2*(1 - alpha);
vb = 1 + alpha;
vc = 1 - alpha;
vd = 2*(1 + alpha);
epslon = 0.0001;

% Computing critical values...

crit_cb1 = 1/(2*(vc-vb));
crit_cb2 = -1/(2*(vc-vb));

y = find(td==0);
ycb1 = find(abs(td - crit_cb1) < epslon);
ycb2 = find(abs(td - crit_cb2) < epslon);

crit_da1 = 1/(2*(vd-va));
yda1 = find(abs(td - crit_da1) < epslon);
crit_da2 = -1/(2*(vd-va));
yda2 = find(abs(td - crit_da2) < epslon);

% Avoiding critical values...

td(y) = epslon;
td(ycb1) = crit_cb1 + epslon;
td(ycb2) = crit_cb2 + epslon;
td(yda1) = crit_da1 + epslon;
td(yda2) = crit_da2 + epslon;

% Wavelet function (psi)

   % Real part
   hra = sin(va*pi*td)./(pi*td);
   hrb = sin(vb*pi*td)./(pi*td);
   mrbc = (2*abs(vc-vb)*(cos(vc*pi*td) + 2*(vc-vb).*td.*sin(vb*pi*td)))./(pi.*(1-(2*td*(vc-vb)).^2));
   mrad = (2*abs(vd-va)*(cos(vd*pi*td) + 2*(vd-va).*td.*sin(va*pi*td)))./(pi.*(1-(2*td*(vd-va)).^2));
   RdeO = (hra - hrb + mrbc + mrad)/(2*sqrt(2*pi));
   % Imaginary part
   denom = pi * td;
   hia = cos(va*pi*td)./denom;
   hib = cos(vb*pi*td)./denom;
   den1 = pi.*(1-(2*td*(vc-vb)).^2);
   den2 = pi.*(1-(2*td*(vd-va)).^2);
   mibc = (2*abs(vc-vb)*(sin(vc*pi*td) - 2*(vc-vb).*td.*cos(vb*pi*td)))./(den1);
   miad = (2*abs(vd-va)*(sin(vd*pi*td) - 2*(vd-va).*td.*cos(va*pi*td)))./(den2);
   IdeO = (hia - hib + mibc + miad)/(2*sqrt(2*pi));
   % deO Wavelet
   psi = 2*sqrt(pi)*(RdeO + j*IdeO);
end

% Set output arguments
if nargout ~= 2 
   error('**deO: too few output parameters: should be [psi, t]=...'); 
end
varargout = {psi,t};
   

%% typical deo computation
%% load a signal (e.g. cuspamax, nelec, lelecum etc.)
% load cuspamax;
% dados = cuspamax(1:1024);
% lv = length(dados);
% subplot(311), plot(dados); title('signal...');
%% performing a cwt by deo0.33333
%% all integer scales from 1 to 64
% subplot(313)
% wt = cwt(dados,1:64,'deo0.33333','plot');
% title('CWT, absolute coefficients')
% colormap(bone(64));
% ylabel('scale');

Contact us