##// END OF EJS Templates
jro_soft » compressed_sensing
RSS
Partial implementation of wavelet transform. Implemented functions: dualfilt1, FSfarras, irls_dn (partial), irls_dn2 (partial), deb4_basis (partial).
yamaro -
r17:18
parent child
Show More
Show file before | Show file after | Hide comments
@@ -0,0 +1,65
1 '''
No newline at end of file
2 Created on May 26, 2014
No newline at end of file
3
No newline at end of file
4 @author: Yolian Amaro
No newline at end of file
5 '''
No newline at end of file
6
No newline at end of file
7 import pywt
No newline at end of file
8 import numpy as np
No newline at end of file
9
No newline at end of file
10 def FSfarras():
No newline at end of file
11 #function [af, sf] = FSfarras
No newline at end of file
12
No newline at end of file
13 # Farras filters organized for the dual-tree
No newline at end of file
14 # complex DWT.
No newline at end of file
15 #
No newline at end of file
16 # USAGE:
No newline at end of file
17 # [af, sf] = FSfarras
No newline at end of file
18 # OUTPUT:
No newline at end of file
19 # af{i}, i = 1,2 - analysis filters for tree i
No newline at end of file
20 # sf{i}, i = 1,2 - synthesis filters for tree i
No newline at end of file
21 # See farras, dualtree, dualfilt1.
No newline at end of file
22 #
No newline at end of file
23 # WAVELET SOFTWARE AT POLYTECHNIC UNIVERSITY, BROOKLYN, NY
No newline at end of file
24 # http://taco.poly.edu/WaveletSoftware/
No newline at end of file
25 #
No newline at end of file
26 # Translated to Python by Yolian Amaro
No newline at end of file
27
No newline at end of file
28
No newline at end of file
29
No newline at end of file
30 a1 = np.array( [
No newline at end of file
31 [ 0, 0],
No newline at end of file
32 [-0.08838834764832, -0.01122679215254],
No newline at end of file
33 [ 0.08838834764832, 0.01122679215254],
No newline at end of file
34 [ 0.69587998903400, 0.08838834764832],
No newline at end of file
35 [ 0.69587998903400, 0.08838834764832],
No newline at end of file
36 [ 0.08838834764832, -0.69587998903400],
No newline at end of file
37 [-0.08838834764832, 0.69587998903400],
No newline at end of file
38 [ 0.01122679215254, -0.08838834764832],
No newline at end of file
39 [ 0.01122679215254, -0.08838834764832],
No newline at end of file
40 [0, 0]
No newline at end of file
41 ] );
No newline at end of file
42
No newline at end of file
43 a2 = np.array([
No newline at end of file
44 [ 0.01122679215254, 0],
No newline at end of file
45 [ 0.01122679215254, 0],
No newline at end of file
46 [-0.08838834764832, -0.08838834764832],
No newline at end of file
47 [ 0.08838834764832, -0.08838834764832],
No newline at end of file
48 [ 0.69587998903400, 0.69587998903400],
No newline at end of file
49 [ 0.69587998903400, -0.69587998903400],
No newline at end of file
50 [ 0.08838834764832, 0.08838834764832],
No newline at end of file
51 [-0.08838834764832, 0.08838834764832],
No newline at end of file
52 [ 0, 0.01122679215254],
No newline at end of file
53 [ 0, -0.01122679215254]
No newline at end of file
54 ]);
No newline at end of file
55
No newline at end of file
56 #print a2.shape
No newline at end of file
57
No newline at end of file
58 af = np.array([ [a1,a2] ], dtype=object)
No newline at end of file
59
No newline at end of file
60 s1 = a1[::-1]
No newline at end of file
61 s2 = a2[::-1]
No newline at end of file
62
No newline at end of file
63 sf = np.array([ [s1,s2] ], dtype=object)
No newline at end of file
64
No newline at end of file
65 return af, sf No newline at end of file
Show file before | Show file after | Hide comments
@@ -0,0 +1,38
1 '''
No newline at end of file
2 Created on May 26, 2014
No newline at end of file
3
No newline at end of file
4 @author: Yolian Amaro
No newline at end of file
5 '''
No newline at end of file
6
No newline at end of file
7 import numpy as np
No newline at end of file
8 import FSfarras
No newline at end of file
9 import dualfilt1
No newline at end of file
10
No newline at end of file
11 def deb4_basis(N):
No newline at end of file
12
No newline at end of file
13 Psi = np.zeros(shape=(N,2*N+1));
No newline at end of file
14 idx = 1;
No newline at end of file
15
No newline at end of file
16 J = 4;
No newline at end of file
17 [Faf, Fsf] = FSfarras;
No newline at end of file
18 [af, sf] = dualfilt1;
No newline at end of file
19 # #
No newline at end of file
20 # # # compute transform of zero vector
No newline at end of file
21 # # x = zeros(1,N);
No newline at end of file
22 # # w = dualtree(x, J, Faf, af);
No newline at end of file
23 # #
No newline at end of file
24 # # # Uses both real and imaginary wavelets
No newline at end of file
25 # # for i in range (1, J+1):
No newline at end of file
26 # # for j in range (1, 2):
No newline at end of file
27 # # for k in range (1, (w[i][j]).size):
No newline at end of file
28 # # w[i][j](k) = 1;
No newline at end of file
29 # # y = idualtree(w, J, Fsf, sf);
No newline at end of file
30 # # w[i][j](k) = 0;
No newline at end of file
31 # # # store it
No newline at end of file
32 # # Psi(:,idx) = y.T.conj();
No newline at end of file
33 # # idx = idx + 1;
No newline at end of file
34 # #
No newline at end of file
35 # # # Add uniform vector (seems to be useful if there's a background
No newline at end of file
36 # # Psi(:,2*N+1) = 1/np.sqrt(N);
No newline at end of file
37 #
No newline at end of file
38 # return Psi No newline at end of file
Show file before | Show file after | Hide comments
@@ -0,0 +1,65
1 '''
No newline at end of file
2 Created on May 29, 2014
No newline at end of file
3
No newline at end of file
4 @author: Yolian Amaro
No newline at end of file
5 '''
No newline at end of file
6
No newline at end of file
7 import numpy as np
No newline at end of file
8
No newline at end of file
9 def dualfilt1():
No newline at end of file
10
No newline at end of file
11 # Kingsbury Q-filters for the dual-tree complex DWT
No newline at end of file
12 #
No newline at end of file
13 # USAGE:
No newline at end of file
14 # [af, sf] = dualfilt1
No newline at end of file
15 # OUTPUT:
No newline at end of file
16 # af{i}, i = 1,2 - analysis filters for tree i
No newline at end of file
17 # sf{i}, i = 1,2 - synthesis filters for tree i
No newline at end of file
18 # note: af{2} is the reverse of af{1}
No newline at end of file
19 # REFERENCE:
No newline at end of file
20 # N. G. Kingsbury, "A dual-tree complex wavelet
No newline at end of file
21 # transform with improved orthogonality and symmetry
No newline at end of file
22 # properties", Proceedings of the IEEE Int. Conf. on
No newline at end of file
23 # Image Proc. (ICIP), 2000
No newline at end of file
24 # See dualtree
No newline at end of file
25 #
No newline at end of file
26 # WAVELET SOFTWARE AT POLYTECHNIC UNIVERSITY, BROOKLYN, NY
No newline at end of file
27 # http://taco.poly.edu/WaveletSoftware/
No newline at end of file
28
No newline at end of file
29 # These cofficients are rounded to 8 decimal places.
No newline at end of file
30
No newline at end of file
31 a1 = np.array([
No newline at end of file
32 [ 0.03516384000000, 0],
No newline at end of file
33 [ 0, 0],
No newline at end of file
34 [-0.08832942000000, -0.11430184000000],
No newline at end of file
35 [ 0.23389032000000, 0],
No newline at end of file
36 [ 0.76027237000000, 0.58751830000000],
No newline at end of file
37 [ 0.58751830000000, -0.76027237000000],
No newline at end of file
38 [ 0, 0.23389032000000],
No newline at end of file
39 [-0.11430184000000, 0.08832942000000],
No newline at end of file
40 [ 0, 0],
No newline at end of file
41 [ 0, -0.03516384000000]
No newline at end of file
42 ]);
No newline at end of file
43
No newline at end of file
44 a2 = np.array([
No newline at end of file
45 [ 0, -0.03516384000000],
No newline at end of file
46 [ 0, 0],
No newline at end of file
47 [-0.11430184000000, 0.08832942000000],
No newline at end of file
48 [ 0, 0.23389032000000],
No newline at end of file
49 [ 0.58751830000000, -0.76027237000000],
No newline at end of file
50 [ 0.76027237000000, 0.58751830000000],
No newline at end of file
51 [ 0.23389032000000, 0],
No newline at end of file
52 [ -0.08832942000000, -0.11430184000000],
No newline at end of file
53 [ 0, 0],
No newline at end of file
54 [ 0.03516384000000, 0]
No newline at end of file
55 ]);
No newline at end of file
56
No newline at end of file
57 af = np.array([ [a1,a2] ], dtype=object)
No newline at end of file
58
No newline at end of file
59 s1 = a1[::-1]
No newline at end of file
60 s2 = a2[::-1]
No newline at end of file
61
No newline at end of file
62 sf = np.array([ [s1,s2] ], dtype=object)
No newline at end of file
63
No newline at end of file
64
No newline at end of file
65 return af, sf No newline at end of file
Show file before | Show file after | Hide comments
@@ -0,0 +1,80
1 '''
No newline at end of file
2 Created on May 27, 2014
No newline at end of file
3
No newline at end of file
4 @author: Yolian Amaro
No newline at end of file
5 '''
No newline at end of file
6
No newline at end of file
7 #from scipy.sparse import eye
No newline at end of file
8 from scipy import linalg
No newline at end of file
9 import scipy.sparse as sps
No newline at end of file
10 import numpy as np
No newline at end of file
11
No newline at end of file
12 def irls_dn(A,b,p,lambda1):
No newline at end of file
13
No newline at end of file
14
No newline at end of file
15 # Minimize lambda*||u||_p + ||A*u-b||_2, 0 < p <= 1
No newline at end of file
16 # using Iterative Reweighted Least Squares
No newline at end of file
17 # (see http://math.lanl.gov/Research/Publications/Docs/chartrand-2008-iteratively.pdf
No newline at end of file
18 # and http://web.eecs.umich.edu/~aey/sparse/sparse11.pdf)
No newline at end of file
19
No newline at end of file
20 # Note to self: I found that "warm-starting" didn't really help too much.
No newline at end of file
21
No newline at end of file
22 [M,N] = A.shape;
No newline at end of file
23 # Initialize and precompute:
No newline at end of file
24 eps = 1e-2; # damping parameter
No newline at end of file
25 [Q,R] = linalg.qr(A.T.conj(),0);
No newline at end of file
26 print A.shape
No newline at end of file
27 print R.shape
No newline at end of file
28 print b.shape
No newline at end of file
29 c = linalg.solve(R.T.conj(),b); # will be used later also
No newline at end of file
30 u = Q*c; # minimum 2-norm solution
No newline at end of file
31 I = sps.eye(M);
No newline at end of file
32
No newline at end of file
33 # Spacing of floating point numbers
No newline at end of file
34 eps = np.spacing(1)
No newline at end of file
35
No newline at end of file
36 # Loop until damping parameter is small enough
No newline at end of file
37 while (eps > 1e-7):
No newline at end of file
38 epschange = 0;
No newline at end of file
39 # Loop until it's time to change eps
No newline at end of file
40 while (~epschange):
No newline at end of file
41 # main loop
No newline at end of file
42 # u_n = W*A'*(A*W*A'+ lambda*I)^-1 * b
No newline at end of file
43 # where W = diag(1/w)
No newline at end of file
44 # where w = (u.^2 + eps).^(p/2-1)
No newline at end of file
45
No newline at end of file
46 # Update
No newline at end of file
47 w = (u**2 + eps)**(1-p/2);
No newline at end of file
48
No newline at end of file
49 # Empty temporary N x N matrix
No newline at end of file
50 temp = np.zeros(shape=(N,N))
No newline at end of file
51
No newline at end of file
52 # Sparse matrix
No newline at end of file
53 for i in range (1, N):
No newline at end of file
54 for j in range (1,N):
No newline at end of file
55 if(i==j):
No newline at end of file
56 temp[i,j] = w
No newline at end of file
57
No newline at end of file
58 # Compressed Sparse Matrix
No newline at end of file
59 W = sps.csr_matrix(temp); #Compressed Sparse Row matrix
No newline at end of file
60
No newline at end of file
61
No newline at end of file
62 WAT = W*A.T.conj();
No newline at end of file
63 u_new = WAT * ( linalg.solve (A*WAT + lambda1*I), b);
No newline at end of file
64
No newline at end of file
65 # See if this subproblem is converging
No newline at end of file
66 delu = np.linalg.norm(u_new-u)/np.linalg.norm(u);
No newline at end of file
67 epschange = delu < (np.sqrt(eps)/100);
No newline at end of file
68
No newline at end of file
69 # Make update
No newline at end of file
70 u = u_new;
No newline at end of file
71
No newline at end of file
72
No newline at end of file
73 eps = eps/10; # decrease eps
No newline at end of file
74 # Print info
No newline at end of file
75 print 'eps =',eps;
No newline at end of file
76
No newline at end of file
77 return u
No newline at end of file
78
No newline at end of file
79
No newline at end of file
80 No newline at end of file
Show file before | Show file after | Hide comments
@@ -0,0 +1,50
1 '''
No newline at end of file
2 Created on May 27, 2014
No newline at end of file
3
No newline at end of file
4 @author: Yolian Amaro
No newline at end of file
5 '''
No newline at end of file
6
No newline at end of file
7 from irls_dn import *
No newline at end of file
8
No newline at end of file
9 def irls_dn2(A,b,p,G):
No newline at end of file
10
No newline at end of file
11 # Minimize ||u||_p subject to ||A*u-b||_2^2 <= G (0 < p <= 1)
No newline at end of file
12
No newline at end of file
13 # What this function actually does is finds the lambda1 so that the solution
No newline at end of file
14 # to the following problem satisfies ||A*u-b||_2^2 <= G:
No newline at end of file
15 # Minimize lambda1*||u||_p + ||A*u-b||_2
No newline at end of file
16
No newline at end of file
17 # Start with a large lambda1, and do a line search until fidelity <= G.
No newline at end of file
18 # (Inversions with large lambda1 are really fast anyway).
No newline at end of file
19
No newline at end of file
20 # Then spin up fzero to localize the root even better
No newline at end of file
21
No newline at end of file
22 # Line Search
No newline at end of file
23
No newline at end of file
24 alpha = 2; # Line search parameter
No newline at end of file
25 lambda1 = 1e5; # What's a reasonable but safe initial guess?
No newline at end of file
26 u = irls_dn(A,b,p,lambda1);
No newline at end of file
27 # fid = np.norm(A*u-b)^2;
No newline at end of file
28 #
No newline at end of file
29 # print '----------------------------------\n';
No newline at end of file
30 #
No newline at end of file
31 # while (fid >= G)
No newline at end of file
32 # lambda1 = lambda1 / alpha; # Balance between speed and accuracy
No newline at end of file
33 # u = irls_dn(A,b,p,lambda1);
No newline at end of file
34 # fid = np.norm(A*u-b)^2;
No newline at end of file
35 # print 'lambda1 = #2e \t ||A*u-b||^2 = #.1f\n',lambda1,fid);
No newline at end of file
36 #
No newline at end of file
37 # # Refinement using fzero
No newline at end of file
38 # lambda10 = [lambda1 lambda1*alpha]; # interval with zero-crossing
No newline at end of file
39 # f = @(lambda1) np.norm(A*irls_dn(A,b,p,lambda1) - b)^2 - G;
No newline at end of file
40 # opts = optimset('fzero');
No newline at end of file
41 # # opts.Display = 'iter';
No newline at end of file
42 # opts.Display = 'none';
No newline at end of file
43 # opts.TolX = 0.01*lambda1;
No newline at end of file
44 # lambda1 = fzero(f,lambda10,opts);
No newline at end of file
45 #
No newline at end of file
46 # u = irls_dn(A,b,p,lambda1);
No newline at end of file
47 #
No newline at end of file
48 #
No newline at end of file
49 # return u;
No newline at end of file
50 No newline at end of file
Show file before | Show file after | Hide comments
@@ -9,20 +9,17
9 #---------------------------------------------------------- No newline at end of file
9 #----------------------------------------------------------
10 No newline at end of file
10
11 import math
11 import math
No newline at end of file
12 #import cmath
No newline at end of file
13 #import scipy
No newline at end of file
14 #import matplotlib No newline at end of file
15 import numpy as np
12 import numpy as np
No newline at end of file
16 #from numpy import linalg No newline at end of file
17 import matplotlib.pyplot as plt No newline at end of file
13 import matplotlib.pyplot as plt
18 from scipy import linalg No newline at end of file
14 from scipy import linalg
19 import time No newline at end of file
15 import time
20 from y_hysell96 import* No newline at end of file
16 from y_hysell96 import*
21 from deb4_basis import * No newline at end of file
17 from deb4_basis import *
22 from modelf import *
18 from modelf import *
No newline at end of file
19 #from scipy.optimize import fsolve No newline at end of file
23 from scipy.optimize import fsolve No newline at end of file
24 from scipy.optimize import root No newline at end of file
20 from scipy.optimize import root
25 import pywt No newline at end of file
21 import pywt
22 from irls_dn2 import * No newline at end of file
26 No newline at end of file
23
27 No newline at end of file
24
28 ## Calculate Forward Model No newline at end of file
25 ## Calculate Forward Model
General Comments 0
You need to be logged in to leave comments. Login now