##A428

'''
Draw the graph of D4(f) for your favorite signal f. Compare the result with D2(f). 
'''

#From L428
import sympy as s

r2 = s.sqrt(2)
r3 = s.sqrt(3)
d= 4*r2

a1 = (1+r3)/d; a2 = (3+r3)/d 
a3 = (3-r3)/d; a4 = (1-r3)/d

# D4trend
def D4trend(f) :
    N = len(f)
    if N%2: return 'D4trend: N is not divisible by 2'%N
    return [a1*f[2*j] +a2*f[2*j+1]+a3*f[(2*j+2)%N]+a4*f[(2*j+3)%N] \
    for j in range(N//2)]
    
#D4fluct
def D4fluct(f) :
    N = len(f)
    if N%2: return 'D4trend: N is not divisible by 2'%N
    return [a4*f[2*j]-a3*f[2*j+1]+a2*f[(2*j+2)%N]-a1*f[(2*j+3)%N] \
    for j in range(N//2)]

#D4 transform
def D4(f) :
    N = len(f)
    if N%2: return 'D4trend: N is not divisible by 2'%N
    return D4trend(f) + D4fluct(f)
    

from cdi_haar import *


#From cdi-graphics:
def write(X,T,ha='center', va = 'center', fs='12', rot=0.0):
    return plt.text(X[0],X[1], T, horizontalalignment=ha,verticalalignment=va, fontsize=fs, rotation=rot)
def hline(a,b,h,lw=1,dash='k-'):
    x = np.arange(a,b+0.0001,b-a)
    return plt.plot(x,h+0*x,dash, lw=lw)
    
    
x = np.arange(0,1,1/(2**9))
def func(x) : return (30*(x**2))*((1-x)**6)*np.sin(pi*5*x)
plt.figure("30*(x^2)*((1-x)^6)*sin(5*pi*x)")
plt.plot(x,func(x),'r-', lw=2)
plt.plot(x,haar(func(x),1),'g-', lw=2)
plt.plot(x,D4(func(x)),'b--', lw=2)
hline(0.5,0.65,-0.2,2,'r-')
write((0.8,-0.2),'Original function')
hline(0.5,0.65,-0.3,2,'g-')
write((0.8,-0.3),'Haar transform')
hline(0.5,0.65,-0.4,2,'b--')
write((0.82,-0.4),'Daubechies transform',)



