lm_transformation_factory.pyx 4.01 KB
Newer Older
1

csongor's avatar
csongor committed
2 3 4
import numpy as np
cimport numpy as np

csongor's avatar
csongor committed
5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33
def buildLm(inp, **kwargs):
    if inp.dtype == np.dtype('float32'):
        return _buildLm_f(inp, **kwargs)
    else:
        return _buildLm(inp, **kwargs)

def buildIdx(inp, **kwargs):
    if inp.dtype == np.dtype('complex64'):
        return _buildIdx_f(inp, **kwargs)
    else:
        return _buildIdx(inp, **kwargs)

cpdef np.ndarray[np.float32_t, ndim=1]  _buildIdx_f(np.ndarray[np
.complex64_t, ndim=1] nr, np.int_t lmax):
    cdef np.int size = (lmax+1)*(lmax+1)

    cdef np.ndarray res=np.zeros([size], dtype=np.complex64)
    cdef np.ndarray final=np.zeros([size], dtype=np.float32)
    res[0:lmax+1] = nr[0:lmax+1]

    for i in xrange(len(nr)-lmax-1):
        res[i*2+lmax+1] = nr[i+lmax+1]
        res[i*2+lmax+2] = np.conjugate(nr[i+lmax+1])
    final = _realify_f(res, lmax)
    return final

cpdef np.ndarray[np.float32_t, ndim=1] _realify_f(np.ndarray[np.complex64_t, ndim=1] nr, np.int_t lmax):
    cdef np.ndarray resi=np.zeros([len(nr)], dtype=np.float32)

34
    resi[0:lmax+1] = nr[0:lmax+1].real
csongor's avatar
csongor committed
35 36 37

    for i in xrange(lmax+1,len(nr),2):
        mi =  int(np.ceil(((2*lmax+1)-np.sqrt((2*lmax+1)*(2*lmax+1)-4*(i-lmax)+1))/2))
38 39
        resi[i]=np.sqrt(2)*(nr[i]).real*(-1)**(mi*mi)
        resi[i+1]=np.sqrt(2)*(nr[i]).imag*(-1)**(mi*mi)
csongor's avatar
csongor committed
40 41 42 43
    return resi

cpdef np.ndarray[np.float64_t, ndim=1]  _buildIdx(np.ndarray[np.complex128_t,
 ndim=1] nr, np.int_t lmax):
csongor's avatar
csongor committed
44 45 46 47 48 49 50 51 52
    cdef np.int size = (lmax+1)*(lmax+1)

    cdef np.ndarray res=np.zeros([size], dtype=np.complex128)
    cdef np.ndarray final=np.zeros([size], dtype=np.float64)
    res[0:lmax+1] = nr[0:lmax+1]

    for i in xrange(len(nr)-lmax-1):
        res[i*2+lmax+1] = nr[i+lmax+1]
        res[i*2+lmax+2] = np.conjugate(nr[i+lmax+1])
csongor's avatar
csongor committed
53
    final = _realify(res, lmax)
54 55
    return final

csongor's avatar
csongor committed
56
cpdef np.ndarray[np.float64_t, ndim=1] _realify(np.ndarray[np.complex128_t, ndim=1] nr, np.int_t lmax):
57 58
    cdef np.ndarray resi=np.zeros([len(nr)], dtype=np.float64)

59
    resi[0:lmax+1] = (nr[0:lmax+1]).real
60 61 62

    for i in xrange(lmax+1,len(nr),2):
        mi =  int(np.ceil(((2*lmax+1)-np.sqrt((2*lmax+1)*(2*lmax+1)-4*(i-lmax)+1))/2))
63 64
        resi[i]=np.sqrt(2)*(nr[i]).real*(-1)**(mi*mi)
        resi[i+1]=np.sqrt(2)*(nr[i]).imag*(-1)**(mi*mi)
65
    return resi
csongor's avatar
csongor committed
66

csongor's avatar
csongor committed
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82
cpdef np.ndarray[np.complex64_t, ndim=1] _buildLm_f(np.ndarray[np.float32_t,
ndim=1] nr, np.int_t lmax):
    cdef np.int size = (len(nr)-lmax-1)/2+lmax+1

    cdef np.ndarray res=np.zeros([size], dtype=np.complex64)
    cdef np.ndarray temp=np.zeros([len(nr)], dtype=np.complex64)
    res[0:lmax+1] = nr[0:lmax+1]

    temp = _inverseRealify_f(nr, lmax)

    for i in xrange(0,len(temp)-lmax-1,2):
        res[i/2+lmax+1] = temp[i+lmax+1]
    return res

cpdef np.ndarray[np.complex64_t, ndim=1] _inverseRealify_f(np.ndarray[np.float32_t, ndim=1] nr, np.int_t lmax):
    cdef np.ndarray resi=np.zeros([len(nr)], dtype=np.complex64)
83
    resi[0:lmax+1] = (nr[0:lmax+1]).real
csongor's avatar
csongor committed
84 85 86 87 88 89 90 91 92 93

    for i in xrange(lmax+1,len(nr),2):
        mi =  int(np.ceil(((2*lmax+1)-np.sqrt((2*lmax+1)*(2*lmax+1)-4*(i-lmax)+1))/2))
        resi[i]=(-1)**mi/np.sqrt(2)*(nr[i]+1j*nr[i+1])
        resi[i+1]=1/np.sqrt(2)*(nr[i]-1j*nr[i+1])
    return resi


cpdef np.ndarray[np.complex128_t, ndim=1] _buildLm(np.ndarray[np.float64_t,
ndim=1] nr, np.int_t lmax):
csongor's avatar
csongor committed
94 95 96 97 98 99
    cdef np.int size = (len(nr)-lmax-1)/2+lmax+1

    cdef np.ndarray res=np.zeros([size], dtype=np.complex128)
    cdef np.ndarray temp=np.zeros([len(nr)], dtype=np.complex128)
    res[0:lmax+1] = nr[0:lmax+1]

csongor's avatar
csongor committed
100
    temp = _inverseRealify(nr, lmax)
csongor's avatar
csongor committed
101 102 103

    for i in xrange(0,len(temp)-lmax-1,2):
        res[i/2+lmax+1] = temp[i+lmax+1]
104
    return res
csongor's avatar
csongor committed
105

csongor's avatar
csongor committed
106
cpdef np.ndarray[np.complex128_t, ndim=1] _inverseRealify(np.ndarray[np.float64_t, ndim=1] nr, np.int_t lmax):
csongor's avatar
csongor committed
107
    cdef np.ndarray resi=np.zeros([len(nr)], dtype=np.complex128)
108
    resi[0:lmax+1] = (nr[0:lmax+1]).real
csongor's avatar
csongor committed
109

110 111 112 113
    for i in xrange(lmax+1,len(nr),2):
        mi =  int(np.ceil(((2*lmax+1)-np.sqrt((2*lmax+1)*(2*lmax+1)-4*(i-lmax)+1))/2))
        resi[i]=(-1)**mi/np.sqrt(2)*(nr[i]+1j*nr[i+1])
        resi[i+1]=1/np.sqrt(2)*(nr[i]-1j*nr[i+1])
csongor's avatar
csongor committed
114 115
    return resi