Fup Wavelet

Балдихина Вера Голубь Дима Как строить вейвлет на fup'ах Для начала нам нужно построить непосредственно сам fup, мы это делали по аналитической формуле: fup N − 1  x  =



4N



1 k 1  2 ∑ sinc N N k=1



N

∞

∏ sinc i =1



2−i  k 2k x cos N N



Код на матлабе который строит эту функцию для N −1 и точки x : function f=fup(N, x) summ=0; for k=1:20 proizv=1; for i=1:40 proizv = proizv*sin(pi*k*2^(-i)/N)/( pi*k*2^(-i)/N ); end summ = summ + proizv*cos(2*pi*k*x/N)*( sin(pi*k/N)/(pi*k/N) )^N; end

end

if abs(x)>(N+2)/2 f = 0; else f = 1/N*(1+2*summ); end

Сумма сдвигов fup обладает свойством:

∑ n ∈Z

fup N t  n = 1

[

Далее строим такую частичую сумму что бы supp = − r

fup N t =

∑

n = −r

fup N

]

4 4 ; , 3 3



 

a N2 r tn  2

Для атомарных функций выполняется: k

∑

f  t  n  = 0.5

n = −k

Используя это, находим коэффициент а , последнее уравнение для нашего случая будет выглядеть так (подставили в качестве t =  ):

r

∑

n = −r



fup N a

 

N 2 1 r n = 2 2

Получаем преобразование фурье искомого вейвлета  следующим производя формальную замену аргумента частичной суммы t  , не забывая извлечь квадратный корень:   =   fup N  =



r

∑

n = −r

fup N





a N2 r n  2



Приведу код который вычисляет заданное преобразование фурье отдавая массив координат [X ,Y ] : function [w, y] = vave() N = 22; a = 0.70833; freq=0.1; r = 8; w=-(N+2)/2:freq:(N+2)/2; for i=1:length(w) y(i) = 0; for n = -r:1:r y(i) = y(i) + fup(N+1, a/pi*((N+1)/2+1+r)*w(i)+n); end y(i) = y(i)^0.5; end end

Вот график который строит эта функция:

Далее все предельно просто, нам нужно вычислить преобразование фурье и получить искомую функцию  . Здесь ничего лучше нету, чем определения: ∞

f x =

1 f e i x  d  ∫  2  −∞

Воспользуемся самым простым алгоритмом численного интегрирования – методом прямоугольников, привожу алгоритм: function fourier_transform(x, y) h = x(2) - x(1); for k=1:length(x) sum = 0; for p=1:length(x) sum = sum + h*y(p)*exp(-j*x(p)*x(k)); end res(k) = real(sum); end plot(x, res); end

Что бы получить этот график надо выполнить команды в матлабе: [x, y] = vave(); fourier_transform(x, y)