nikcleju@10: # -*- coding: utf-8 -*-
nikcleju@10: """
nikcleju@10: Created on Sat Nov 05 20:29:02 2011
nikcleju@10: Test SL0 algorithm
nikcleju@10: @author: Nic
nikcleju@10: """
nikcleju@10: 
nikcleju@10: 
nikcleju@10: import numpy as np
nikcleju@10: import numpy.linalg
nikcleju@10: import scipy.io
nikcleju@10: import unittest
nikcleju@10: from pyCSalgos.SL0.SL0 import SL0
nikcleju@10: 
nikcleju@10: class SL0results(unittest.TestCase):
nikcleju@10:   def testResults(self):
nikcleju@10:     mdict = scipy.io.loadmat('SL0testdata.mat')
nikcleju@10:     
nikcleju@10:     # A = system matrix
nikcleju@10:     # Y = matrix with measurements (on columns)
nikcleju@10:     # sigmamin = vector with sigma_min
nikcleju@10:     for A,Y,sigmamin,Xr in zip(mdict['cellA'].squeeze(),mdict['cellY'].squeeze(),mdict['sigmamin'].squeeze(),mdict['cellXr'].squeeze()):
nikcleju@10:       for i in np.arange(Y.shape[1]):
nikcleju@10:         
nikcleju@10:         # Fix numpy error "LapackError: Parameter a has non-native byte order in lapack_lite.dgesdd"
nikcleju@10:         A = A.newbyteorder('=')
nikcleju@10:         Y = Y.newbyteorder('=')
nikcleju@10:         sigmamin = sigmamin.newbyteorder('=')
nikcleju@10:         Xr = Xr.newbyteorder('=')
nikcleju@10:         
nikcleju@10:         xr = SL0(A, Y[:,i], sigmamin)
nikcleju@10:         
nikcleju@10:         # check if found solution is the same as the correct cslution
nikcleju@10:         diff = numpy.linalg.norm(xr - Xr[:,i])
nikcleju@10:         self.assertTrue(diff < 1e-12)
nikcleju@10:     #        err1 = numpy.linalg.norm(Y[:,i] - np.dot(A,xr))
nikcleju@10:     #        err2 = numpy.linalg.norm(Y[:,i] - np.dot(A,Xr[:,i]))
nikcleju@10:     #        norm1 = numpy.linalg.norm(xr,1)
nikcleju@10:     #        norm2 = numpy.linalg.norm(Xr[:,i],1)
nikcleju@10:     #                
nikcleju@10:     #        # Make a more robust condition:
nikcleju@10:     #        #  OK;    if   solutions are close enough (diff < 1e-6)
nikcleju@10:     #        #              or
nikcleju@10:     #        #              (
nikcleju@10:     #        #               Python solution fulfills the constraint better (or up to 1e-6 worse)
nikcleju@10:     #        #                 and
nikcleju@10:     #        #               Python solution has l1 norm no more than 1e-6 larger as the reference solution
nikcleju@10:     #        #                 (i.e. either norm1 < norm2   or   norm1>norm2 not by more than 1e-6)
nikcleju@10:     #        #              )
nikcleju@10:     #        #        
nikcleju@10:     #        #  ERROR: else        
nikcleju@10:     #        differr  = err1 - err2    # intentionately no abs(), since err1` < err2 is good
nikcleju@10:     #        diffnorm = norm1 - norm2  # intentionately no abs(), since norm1 < norm2 is good
nikcleju@10:     #        if diff < 1e-6 or (differr < 1e-6 and (diffnorm < 1e-6)):
nikcleju@10:     #          isok = True
nikcleju@10:     #        else:
nikcleju@10:     #          isok = False
nikcleju@10:     #        self.assertTrue(isok)
nikcleju@10:         
nikcleju@10:         #diff = numpy.linalg.norm(xr - Xr[:,i])
nikcleju@10:         #if diff > 1e-6:
nikcleju@10:         #    self.assertTrue(diff < 1e-6)
nikcleju@10: 
nikcleju@10:   
nikcleju@10: if __name__ == "__main__":
nikcleju@10:     #import cProfile
nikcleju@10:     #cProfile.run('unittest.main()', 'profres')
nikcleju@10:     unittest.main()    
nikcleju@10:     #suite = unittest.TestLoader().loadTestsFromTestCase(CompareResults)
nikcleju@10:     #unittest.TextTestRunner(verbosity=2).run(suite)