nikcleju@15
|
1 # -*- coding: utf-8 -*-
|
nikcleju@15
|
2 """
|
nikcleju@15
|
3 Created on Sat Nov 05 21:34:49 2011
|
nikcleju@15
|
4
|
nikcleju@15
|
5 @author: Nic
|
nikcleju@15
|
6 """
|
nikcleju@15
|
7 import numpy as np
|
nikcleju@15
|
8 import numpy.linalg
|
nikcleju@15
|
9 import scipy.io
|
nikcleju@15
|
10 import unittest
|
nikcleju@15
|
11 from pyCSalgos.SL0.SL0_approx import SL0_approx
|
nikcleju@15
|
12
|
nikcleju@15
|
13 class SL0results(unittest.TestCase):
|
nikcleju@15
|
14 def testResults(self):
|
nikcleju@15
|
15 mdict = scipy.io.loadmat('SL0approxtestdata.mat')
|
nikcleju@15
|
16
|
nikcleju@15
|
17 # A = system matrix
|
nikcleju@15
|
18 # Y = matrix with measurements (on columns)
|
nikcleju@15
|
19 # sigmamin = vector with sigma_min
|
nikcleju@15
|
20 for A,Y,eps,sigmamin,Xr in zip(mdict['cellA'].squeeze(),mdict['cellY'].squeeze(),mdict['cellEps'].squeeze(),mdict['sigmamin'].squeeze(),mdict['cellXr'].squeeze()):
|
nikcleju@15
|
21 for i in np.arange(Y.shape[1]):
|
nikcleju@15
|
22
|
nikcleju@15
|
23 # Fix numpy error "LapackError: Parameter a has non-native byte order in lapack_lite.dgesdd"
|
nikcleju@15
|
24 A = A.newbyteorder('=')
|
nikcleju@15
|
25 Y = Y.newbyteorder('=')
|
nikcleju@15
|
26 eps = eps.newbyteorder('=')
|
nikcleju@15
|
27 sigmamin = sigmamin.newbyteorder('=')
|
nikcleju@15
|
28 Xr = Xr.newbyteorder('=')
|
nikcleju@15
|
29
|
nikcleju@15
|
30 xr = SL0_approx(A, Y[:,i], eps.squeeze()[i], sigmamin)
|
nikcleju@15
|
31
|
nikcleju@15
|
32 # check if found solution is the same as the correct cslution
|
nikcleju@15
|
33 diff = numpy.linalg.norm(xr - Xr[:,i])
|
nikcleju@15
|
34 self.assertTrue(diff < 1e-12)
|
nikcleju@15
|
35 # err1 = numpy.linalg.norm(Y[:,i] - np.dot(A,xr))
|
nikcleju@15
|
36 # err2 = numpy.linalg.norm(Y[:,i] - np.dot(A,Xr[:,i]))
|
nikcleju@15
|
37 # norm1 = numpy.linalg.norm(xr,1)
|
nikcleju@15
|
38 # norm2 = numpy.linalg.norm(Xr[:,i],1)
|
nikcleju@15
|
39 #
|
nikcleju@15
|
40 # # Make a more robust condition:
|
nikcleju@15
|
41 # # OK; if solutions are close enough (diff < 1e-6)
|
nikcleju@15
|
42 # # or
|
nikcleju@15
|
43 # # (
|
nikcleju@15
|
44 # # Python solution fulfills the constraint better (or up to 1e-6 worse)
|
nikcleju@15
|
45 # # and
|
nikcleju@15
|
46 # # Python solution has l1 norm no more than 1e-6 larger as the reference solution
|
nikcleju@15
|
47 # # (i.e. either norm1 < norm2 or norm1>norm2 not by more than 1e-6)
|
nikcleju@15
|
48 # # )
|
nikcleju@15
|
49 # #
|
nikcleju@15
|
50 # # ERROR: else
|
nikcleju@15
|
51 # differr = err1 - err2 # intentionately no abs(), since err1` < err2 is good
|
nikcleju@15
|
52 # diffnorm = norm1 - norm2 # intentionately no abs(), since norm1 < norm2 is good
|
nikcleju@15
|
53 # if diff < 1e-6 or (differr < 1e-6 and (diffnorm < 1e-6)):
|
nikcleju@15
|
54 # isok = True
|
nikcleju@15
|
55 # else:
|
nikcleju@15
|
56 # isok = False
|
nikcleju@15
|
57 # self.assertTrue(isok)
|
nikcleju@15
|
58
|
nikcleju@15
|
59 #diff = numpy.linalg.norm(xr - Xr[:,i])
|
nikcleju@15
|
60 #if diff > 1e-6:
|
nikcleju@15
|
61 # self.assertTrue(diff < 1e-6)
|
nikcleju@15
|
62
|
nikcleju@15
|
63
|
nikcleju@15
|
64 if __name__ == "__main__":
|
nikcleju@15
|
65 #import cProfile
|
nikcleju@15
|
66 #cProfile.run('unittest.main()', 'profres')
|
nikcleju@15
|
67 unittest.main()
|
nikcleju@15
|
68 #suite = unittest.TestLoader().loadTestsFromTestCase(CompareResults)
|
nikcleju@15
|
69 #unittest.TextTestRunner(verbosity=2).run(suite)
|