 
function m0=m0(omega)
% this function is a mess! it does many things,
% depending what you uncomment 
% (1) gives m0(omega) if there is a nice expression
% or
% (2) takes $h$ filter and makes m0(omega);
%
%--------------------------------------------
% m0 - for haar <------------------
%m0 = 1/2 * (1+exp(-i * omega));
% m0 - for Daub2 <------------------
%m0 = ( 1/2 * (1+exp(-i * omega))).^2 * 1/2 .* ...
%   ( (1-sqrt(3))*exp(-i*omega) + (1+sqrt(3)) );
% some filters <------------------------------------
% h for Daub2
% sq3=sqrt(3); nor=4*sqrt(2); h=[(1+sq3)/nor ...
% (3+sq3)/nor (3-sq3)/nor (1-sq3)/nor];
% h for Daub4
 h = [0.2303778133 0.7148455706 0.63088076792 -0.0279837694 ...
    -0.1870348117 0.0308413818 0.03288301167 -0.01059740178];
% Given h this loop makes m0:
nn=length(h);
m0=linspace(0,0,length(omega));
for ii=1:nn
 m0 = m0 + h(ii) .* exp(-i * (ii-1) .* omega);
end
m0=1/sqrt(2) * m0;
%------------------------------------------------------------

