import string
import types

import Numeric
import LinearAlgebra
from UserArray import UserArray, asarray
from Numeric import dot, identity, multiply
from MLab import squeeze

# make translation table
_table = [None]*256
for k in range(256):
    _table[k] = chr(k)
_table = ''.join(_table)

_numchars = string.digits + ".-+jeEL"
_todelete = []
for k in _table:
    if k not in _numchars:
        _todelete.append(k)
_todelete = ''.join(_todelete)

def _eval(astr):
    return eval(astr.translate(_table,_todelete))

def _convert_from_string(data):
    data.find
    rows = data.split(';')
    newdata = []
    count = 0
    for row in rows:
        trow = row.split(',')
        newrow = []
        for col in trow:
            temp = col.split()
            newrow.extend(map(_eval,temp))
        if count == 0:
            Ncols = len(newrow)
        elif len(newrow) != Ncols:
            raise ValueError, "Rows not the same size."
        count += 1
        newdata.append(newrow)
    return newdata


_lkup = {'0':'000',
         '1':'001',
         '2':'010',
         '3':'011',
         '4':'100',
         '5':'101',
         '6':'110',
         '7':'111'}

def _binary(num):
    ostr = oct(num)
    bin = ''
    for ch in ostr[1:]:
        bin += _lkup[ch]
    ind = 0
    while bin[ind] == '0':
        ind += 1
    return bin[ind:]
        

class Matrix(UserArray):
    def __init__(self, data, typecode=None, copy=1, savespace=0):
        if type(data) is types.StringType:
            data = _convert_from_string(data)
        UserArray.__init__(self,data,typecode, copy, savespace)
        if len(self.array.shape) == 1:
            self.array.shape = (1,self.array.shape[0])
        self.shape = self.array.shape

    def __getitem__(self, index):        
        out = self._rc(self.array[index])
        try:
            n = len(index)
            if n > 1 and isinstance(index[0], types.SliceType) \
               and isinstance(index[1], types.IntType):
                sh = out.array.shape
                out.array.shape = (sh[1],sh[0])
                out.shape = out.array.shape
        except TypeError:  # Index not a sequence
            pass
        return out

    def __mul__(self, other):
        aother = asarray(other)
        if len(aother.shape) == 0:
            return self._rc(self.array*aother)
        else:
            return self._rc(dot(self.array, aother))

    def __rmul__(self, other):
        aother = asarray(other)
        if len(aother.shape) == 0:
            return self._rc(aother*self.array)
        else:
            return self._rc(dot(aother, self.array))

    def __imul__(self,other):
        aother = asarray(other)
        self.array = dot(self.array, aother)
        return self

    def __pow__(self, other):
        shape = self.array.shape
        if len(shape)!=2 or shape[0]!=shape[1]:
            raise TypeError, "matrix is not square"
        if type(other) in (type(1), type(1L)):
            if other==0:
                return Matrix(identity(shape[0]))
            if other<0:
                result=Matrix(LinearAlgebra.inverse(self.array))
                x=Matrix(result)
                other=-other
            else:
                result=self
                x=result
            if other <= 3:
                while(other>1):
                    result=result*x
                    other=other-1
                return result
            # binary decomposition to reduce the number of Matrix
            #  Multiplies for other > 3.
            beta = _binary(other)            
            t = len(beta)
            Z,q = x.copy(),0
            while beta[t-q-1] == '0':
                Z *= Z
                q += 1
            result = Z.copy()
            for k in range(q+1,t):
                Z *= Z
                if beta[t-k-1] == '1':
                    result *= Z
            return result
        else:
            raise TypeError, "exponent must be an integer"        

    def __rpow__(self, other):
        raise TypeError, "x**y not implemented for matrices y"

    def __invert__(self):
        return Matrix(Numeric.conjugate(self.array))
    
    def __setattr__(self,attr,value):
        if attr=='shape':
            if len(value) == 0:
                value = (1,1)
            if len(value) == 1:
                value = (1,) + value 
            if len(value) != 2:
                raise ValueError, "Can only reshape a Matrix as a 2-d array."
            self.array.shape=value
        self.__dict__[attr]=value

    def __getattr__(self, attr):
        if attr == 'A':
            return squeeze(self.array)
        elif attr == 'T':
            return Matrix(Numeric.transpose(self.array))
        elif attr == 'H':
            if len(self.array.shape) == 1:
                self.array.shape = (1,self.array.shape[0])
            return Matrix(Numeric.conjugate(Numeric.transpose(self.array)))
        elif attr == 'I':
            return Matrix(LinearAlgebra.inverse(self.array))
        elif attr == 'real':
            return Matrix(self.array.real)
        elif attr == 'imag':
            return Matrix(self.array.imag)
        elif attr == 'flat':
            return Matrix(self.array.flat)
        else:
            raise AttributeError, attr + " not found."

    def __setitem__(self, index, value):
        value = asarray(value, self._typecode)
        if multiply.reduce(value.shape) == 1:
            value = value[0,0]
        self.array[index] = value
    def __setslice__(self, i, j, value): self.array[i:j] = asarray(squeeze(value),self._typecode)

    def __float__(self):
        return float(squeeze(self.array))
    
    def m(self,b):
        return Matrix(self.array * asarray(b))
    
    def asarray(self,t=None):
        if t is None:
            return self.array
        else:
            return asarray(self.array,t)

if __name__ == '__main__':
        from Numeric import *
        m = Matrix( [[1,2,3],[11,12,13],[21,22,23]])
        print m*m
        print m.array*m.array
        print transpose(m)
        print m**-1
        m = Matrix([[1,1],[1,0]])
        print "Fibonacci numbers:",
        for i in range(10):
            mm=m**i
            print mm[0][0],
        print
        
            
            
