CMSIS-DSP: Update to Python wrapper

New functions added.

.gitignore

@@ -1,6 +1,7 @@
 DSP_Lib_TestSuite/build/
 PythonWrapper/build/
 PythonWrapper/cmsisdsp.cp36-win_amd64.pyd
+PythonWrapper/*.so
 PythonWrapper/rec_2.dat
 *.pickle
 build_*/

PythonWrapper/cmsisdsp_pkg/src/cmsismodule.c

@@ -3,13 +3,13 @@
  * Title:        cmsismodule.c
  * Description:  C code for the CMSIS-DSP Python wrapper
  *
- * $Date:        25. March 2019
+ * $Date:        27 April 2021
- * $Revision:    V0.0.1
+ * $Revision:    VV1.0
  *
  * Target Processor: Cortex-M cores
  * -------------------------------------------------------------------- */
 /*
- * Copyright (C) 2010-2019 ARM Limited or its affiliates. All rights reserved.
+ * Copyright (C) 2010-2021 ARM Limited or its affiliates. All rights reserved.
  *
  * SPDX-License-Identifier: Apache-2.0
  *
@@ -186,10 +186,10 @@ Method_##NAME##_##FIELD(ml_##NAME##Object *self, PyObject *ignored)\
     PyObject *OBJ=PyArray_SimpleNewFromData(ND, dims, NPY_FLOAT, DATA);
 
 #define FLOATARRAY1(OBJ,NB1,DATA)                                          \
-    npy_intp dims[1];                                                  \
+    npy_intp dims##OBJ[1];                                                  \
-    dims[0]=NB1;                                                       \
+    dims##OBJ[0]=NB1;                                                       \
-    const int ND=1;                                                    \
+    const int ND##OBJ=1;                                                    \
-    PyObject *OBJ=PyArray_SimpleNewFromData(ND, dims, NPY_FLOAT, DATA);
+    PyObject *OBJ=PyArray_SimpleNewFromData(ND##OBJ, dims##OBJ, NPY_FLOAT, DATA);
 
 #define FLOAT64ARRAY1(OBJ,NB1,DATA)                                          \
     npy_intp dims[1];                                                  \
@@ -221,6 +221,11 @@ Method_##NAME##_##FIELD(ml_##NAME##Object *self, PyObject *ignored)\
     const int ND=1;                                                     \
     PyObject *OBJ=PyArray_SimpleNewFromData(ND, dims, NPY_BYTE, DATA);
 
+#define TYP_ARRAY1(OBJ,NB1,DATA,NPYTYPE)                                          \
+    npy_intp dims[1];                                                   \
+    dims[0]=NB1;                                                        \
+    const int ND=1;                                                     \
+    PyObject *OBJ=PyArray_SimpleNewFromData(ND, dims, NPYTYPE, DATA);
 
 #define MATRIXFROMNUMPY(EXT,TYP,SRCTYPE,NUMPYTYPE)                                   \
 arm_matrix_instance_##EXT *EXT##MatrixFromNumpy(PyObject *o)                   \
@@ -257,6 +262,7 @@ MATRIXFROMNUMPY(f32,float32_t,double,NPY_DOUBLE);
 MATRIXFROMNUMPY(f64,float64_t,double,NPY_DOUBLE);
 MATRIXFROMNUMPY(q31,q31_t,int32_t,NPY_INT32);
 MATRIXFROMNUMPY(q15,q15_t,int16_t,NPY_INT16);
+MATRIXFROMNUMPY(q7,q7_t,int8_t,NPY_BYTE);
 
 #define CREATEMATRIX(EXT,TYP)                                        \
 arm_matrix_instance_##EXT *create##EXT##Matrix(uint32_t r,uint32_t c)\
@@ -274,6 +280,7 @@ CREATEMATRIX(f32,float32_t);
 CREATEMATRIX(f64,float64_t);
 CREATEMATRIX(q31,q31_t);
 CREATEMATRIX(q15,q15_t);
+CREATEMATRIX(q7,q7_t);
 
 #define NUMPYVECTORFROMBUFFER(EXT,CTYPE,NUMPYTYPE_FROMC)                     \
 PyObject *NumpyVectorFrom##EXT##Buffer(CTYPE *ptr,int nb)                    \
@@ -310,6 +317,7 @@ NUMPYARRAYFROMMATRIX(f32,NPY_FLOAT);
 NUMPYARRAYFROMMATRIX(f64,NPY_DOUBLE);
 NUMPYARRAYFROMMATRIX(q31,NPY_INT32);
 NUMPYARRAYFROMMATRIX(q15,NPY_INT16);
+NUMPYARRAYFROMMATRIX(q7,NPY_BYTE);
 
 //#include "specific.h"
 #include "cmsismodule.h"

PythonWrapper/fixedpoint.py

@@ -0,0 +1,49 @@
+import numpy as np
+
+def q31sat(x):
+     if x > 0x7FFFFFFF:
+          return(np.int32(0x7FFFFFFF))
+     elif x < -0x80000000:
+          return(np.int32(0x80000000))
+     else:
+          return(np.int32(x))
+
+q31satV=np.vectorize(q31sat)
+
+def toQ31(x):
+     return(q31satV(np.round(x * (1<<31))))
+
+def q15sat(x):
+     if x > 0x7FFF:
+          return(np.int16(0x7FFF))
+     elif x < -0x8000:
+          return(np.int16(0x8000))
+     else:
+          return(np.int16(x))
+
+q15satV=np.vectorize(q15sat)
+
+def toQ15(x):
+     return(q15satV(np.round(x * (1<<15))))
+
+def q7sat(x):
+     if x > 0x7F:
+          return(np.int8(0x7F))
+     elif x < -0x80:
+          return(np.int8(0x80))
+     else:
+          return(np.int8(x))
+
+q7satV=np.vectorize(q7sat)
+
+def toQ7(x):
+     return(q7satV(np.round(x * (1<<7))))
+
+def Q31toF32(x):
+     return(1.0*x / 2**31)
+
+def Q15toF32(x):
+     return(1.0*x / 2**15)
+
+def Q7toF32(x):
+     return(1.0*x / 2**7)

PythonWrapper/setup.py

@@ -58,12 +58,20 @@ interpolation = glob.glob(os.path.join(ROOT,"Source","InterpolationFunctions","*
 interpolation.remove(os.path.join(ROOT,"Source","InterpolationFunctions","InterpolationFunctions.c"))
 interpolation.remove(os.path.join(ROOT,"Source","InterpolationFunctions","InterpolationFunctionsF16.c"))
 
+quaternion = glob.glob(os.path.join(ROOT,"Source","QuaternionMathFunctions","*.c"))
+quaternion.remove(os.path.join(ROOT,"Source","QuaternionMathFunctions","QuaternionMathFunctions.c"))
+
+#distance = glob.glob(os.path.join(ROOT,"Source","DistanceFunctions","*.c"))
+#distance.remove(os.path.join(ROOT,"Source","DistanceFunctions","DistanceFunctions.c"))
+
+
 #modulesrc = glob.glob(os.path.join("cmsisdsp_pkg","src","*.c"))
 modulesrc = []
 modulesrc.append(os.path.join("cmsisdsp_pkg","src","cmsismodule.c"))
 
 allsrcs = support + fastmath + filtering + matrix + statistics + complexf + basic
 allsrcs = allsrcs + controller + transform + modulesrc + common+ interpolation
+allsrcs = allsrcs + quaternion
 
 def notf16(number):
   if re.search(r'f16',number):
@@ -88,10 +96,10 @@ module1 = Extension(config.extensionName,
                               )
 
 setup (name = config.setupName,
-       version = '0.0.1',
+       version = '1.0.0',
        description = config.setupDescription,
        ext_modules = [module1],
-       author = 'Copyright (C) 2010-2020 ARM Limited or its affiliates. All rights reserved.',
+       author = 'Copyright (C) 2010-2021 ARM Limited or its affiliates. All rights reserved.',
        url="https://github.com/ARM-software/CMSIS_5",
        classifiers=[
         "Programming Language :: Python",

PythonWrapper/testdsp2.py

@@ -0,0 +1,459 @@
+import cmsisdsp as dsp
+import numpy as np
+from scipy import signal
+from scipy.fftpack import dct 
+import fixedpoint as f
+from pyquaternion import Quaternion
+
+import colorama
+from colorama import init,Fore, Back, Style
+import statsmodels.tsa.stattools
+
+import scipy.spatial
+
+
+init()
+
+def printTitle(s):
+    print("\n" + Fore.GREEN + Style.BRIGHT +  s + Style.RESET_ALL)
+
+def printSubTitle(s):
+    print("\n" + Style.BRIGHT + s + Style.RESET_ALL)
+
+
+def imToReal2D(a):
+    ar=np.zeros(np.array(a.shape) * [1,2])
+    ar[::,0::2]=a.real
+    ar[::,1::2]=a.imag
+    return(ar)
+
+def realToIm2D(ar):
+    return(ar[::,0::2] + 1j * ar[::,1::2])
+
+def normalize(a):
+  return(a/np.max(np.abs(a)))
+
+def autocorr(x):
+    result = np.correlate(x, x, mode='full')
+    return result[result.size//2:]
+
+#################### MAX AND ABSMAX ##################################
+printTitle("Max and AbsMax")
+a=np.array([1.,-3.,4.,0.,-10.,8.])
+
+printSubTitle("Float tests")
+i=dsp.arm_max_f32(a)
+print(i)
+
+i=dsp.arm_absmax_f32(a)
+print(i)
+
+printSubTitle("Fixed point tests")
+
+# Normalize for fixed point tests
+a = a / i[0]
+
+a31 = f.toQ31(a)
+i=dsp.arm_absmax_q31(a31)
+print(f.Q31toF32(i[0]),i[1])
+
+a8 = f.toQ15(a)
+i=dsp.arm_absmax_q15(a8)
+print(f.Q15toF32(i[0]),i[1])
+
+a7 = f.toQ7(a)
+i=dsp.arm_absmax_q7(a7)
+print(f.Q7toF32(i[0]),i[1])
+
+################### MIN AND ABSMIN ################################
+
+printTitle("Min and AbsMin")
+a=np.array([1.,-3.,4.,0.5,-10.,8.])
+
+printSubTitle("Float tests")
+i=dsp.arm_min_f32(a)
+print(i)
+
+i=dsp.arm_absmin_f32(a)
+print(i)
+
+printSubTitle("Fixed point tests")
+
+# Normalize for fixed point tests
+idx=i[1]
+i=dsp.arm_absmax_f32(a)
+a = a / i[0]
+print(a)
+print(a[idx])
+
+a31 = f.toQ31(a)
+i=dsp.arm_absmin_q31(a31)
+print(f.Q31toF32(i[0]),i[1])
+
+a8 = f.toQ15(a)
+i=dsp.arm_absmin_q15(a8)
+print(f.Q15toF32(i[0]),i[1])
+
+a7 = f.toQ7(a)
+i=dsp.arm_absmin_q7(a7)
+print(f.Q7toF32(i[0]),i[1])
+
+##################### CLIPPING ###################
+printTitle("Clipping tests tests")
+a=np.array([1.,-3.,4.,0.5,-10.,8.])
+i=dsp.arm_absmax_f32(a)
+
+minBound =-5.0 
+maxBound =6.0
+b=dsp.arm_clip_f32(a,minBound,maxBound)
+print(a)
+print(b)
+
+a = a / i[0]
+print(a)
+minBound = minBound / i[0]
+maxBound = maxBound / i[0]
+print(minBound,maxBound)
+
+b=dsp.arm_clip_q31(f.toQ31(a),f.toQ31(minBound),f.toQ31(maxBound))
+print(f.Q31toF32(b))
+
+b=dsp.arm_clip_q15(f.toQ15(a),f.toQ15(minBound),f.toQ15(maxBound))
+print(f.Q15toF32(b))
+
+b=dsp.arm_clip_q7(f.toQ7(a),f.toQ7(minBound),f.toQ7(maxBound))
+print(f.Q7toF32(b))
+
+############### MAT VECTOR MULT
+
+printTitle("Matrix x Vector")
+a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
+b=np.array([-2,-1,3,4])
+
+c = np.dot(a,b)
+print(c)
+c = dsp.arm_mat_vec_mult_f32(a,b)
+print(c)
+
+printSubTitle("Fixed point")
+normalizationFactor=2.0*np.sqrt(np.max(np.abs(c)))
+a=a/normalizationFactor
+b=b/normalizationFactor
+print(np.dot(a,b))
+
+c=dsp.arm_mat_vec_mult_q31(f.toQ31(a),f.toQ31(b))
+print(f.Q31toF32(c))
+
+c=dsp.arm_mat_vec_mult_q15(f.toQ15(a),f.toQ15(b))
+print(f.Q15toF32(c))
+
+c=dsp.arm_mat_vec_mult_q7(f.toQ7(a),f.toQ7(b))
+print(f.Q7toF32(c))
+
+############### MATRIX MULTIPLY
+
+printTitle("Matrix x Matrix")
+
+a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
+b=np.array([[1.,2,3],[5.1,6,7],[9.1,10,11],[5,8,4]])
+print(np.dot(a , b))
+c=dsp.arm_mat_mult_f32(a,b)
+print(c[1])
+
+printSubTitle("Fixed point")
+
+normalizationFactor=2.0*np.sqrt(np.max(np.abs(c[1])))
+a = a / normalizationFactor
+b = b / normalizationFactor
+c=dsp.arm_mat_mult_f32(a,b)
+print(c[1])
+
+print("")
+af = f.toQ31(a)
+bf = f.toQ31(b)
+c = dsp.arm_mat_mult_q31(af,bf)
+print(f.Q31toF32(c[1]))
+
+print("")
+af = f.toQ15(a)
+bf = f.toQ15(b)
+s=bf.shape 
+nb=s[0]*s[1]
+tmp=np.zeros(nb)
+c = dsp.arm_mat_mult_q15(af,bf,tmp)
+print(f.Q15toF32(c[1]))
+
+print("")
+af = f.toQ7(a)
+bf = f.toQ7(b)
+s=bf.shape 
+nb=s[0]*s[1]
+tmp=np.zeros(nb)
+c = dsp.arm_mat_mult_q7(af,bf,tmp)
+print(f.Q7toF32(c[1]))
+
+################# MAT TRANSPOSE #################
+
+printTitle("Transposition")
+a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
+normalizationFactor=np.max(np.abs(c[1]))
+a = a / normalizationFactor
+
+print(np.transpose(a))
+print("")
+r=dsp.arm_mat_trans_f32(a)
+print(r[1])
+print("")
+
+r=dsp.arm_mat_trans_q31(f.toQ31(a))
+print(f.Q31toF32(r[1]))
+print("")
+
+r=dsp.arm_mat_trans_q15(f.toQ15(a))
+print(f.Q15toF32(r[1]))
+print("")
+
+r=dsp.arm_mat_trans_q7(f.toQ7(a))
+print(f.Q7toF32(r[1]))
+print("")
+
+################## FILL FUNCTIONS #################
+
+v=0.22 
+nb=10 
+a=np.full((nb,),v)
+print(a)
+
+a=dsp.arm_fill_f32(v,nb)
+print(a)
+
+a=f.Q31toF32(dsp.arm_fill_q31(f.toQ31(v),nb))
+print(a)
+
+a=f.Q15toF32(dsp.arm_fill_q15(f.toQ15(v),nb))
+print(a)
+
+a=f.Q7toF32(dsp.arm_fill_q7(f.toQ7(v),nb))
+print(a)
+
+################# COMPLEX MAT TRANSPOSE #################
+
+printTitle("Complex Transposition")
+a=np.array([[1. + 0.0j ,2 + 1.0j,3 + 0.0j,4 + 2.0j],
+            [5 + 1.0j,6 + 2.0j,7 + 3.0j,8 + 1.0j],
+            [9 - 2.0j,10 + 1.0j,11 - 4.0j,12 + 1.0j]])
+normalizationFactor=np.max(np.abs(c[1]))
+a = a / normalizationFactor
+
+print(np.transpose(a))
+print("")
+r=dsp.arm_mat_cmplx_trans_f32(imToReal2D(a))
+print(realToIm2D(r[1]))
+print("")
+
+r=dsp.arm_mat_cmplx_trans_q31(f.toQ31(imToReal2D(a)))
+print(realToIm2D(f.Q31toF32(r[1])))
+print("")
+
+r=dsp.arm_mat_cmplx_trans_q15(f.toQ15(imToReal2D(a)))
+print(realToIm2D(f.Q15toF32(r[1])))
+print("")
+
+################ Levinson ##################
+
+printTitle("Levinson Durbin")
+na=5
+s = np.random.randn(na+1)
+s = normalize(s)
+phi = autocorr(s)
+phi = normalize(phi)
+
+sigmav,arcoef,pacf,sigma,phi1=statsmodels.tsa.stattools.levinson_durbin(phi,nlags=na,isacov=True)
+      
+print(arcoef)
+print(sigmav)
+
+(a,err)=dsp.arm_levinson_durbin_f32(phi,na)
+print(a)
+print(err)
+
+phiQ31 = f.toQ31(phi)
+(aQ31,errQ31)=dsp.arm_levinson_durbin_q31(phiQ31,na)
+print(f.Q31toF32(aQ31))
+print(f.Q31toF32(errQ31))
+
+################## Bitwise operations #################
+
+printTitle("Bitwise operations")
+def genBitvectors(nb,format):
+    if format == 31:
+       maxVal = 0x7fffffff
+    if format == 15:
+       maxVal = 0x7fff
+    if format == 7:
+       maxVal = 0x7f 
+
+    minVal = -maxVal-1
+    
+    return(np.random.randint(minVal, maxVal, size=nb))
+
+NBSAMPLES=10
+
+
+
+printSubTitle("u32")
+su32A=genBitvectors(NBSAMPLES,31)
+su32B=genBitvectors(NBSAMPLES,31)
+ffff = (np.ones(NBSAMPLES)*(-1)).astype(np.int)
+
+
+ref=np.bitwise_and(su32A, su32B)
+#print(ref)
+result=dsp.arm_and_u32(su32A, su32B).astype(int)
+print(result-ref)
+
+ref=np.bitwise_or(su32A, su32B)
+#print(ref)
+result=dsp.arm_or_u32(su32A, su32B).astype(int)
+print(result-ref)
+
+ref=np.bitwise_xor(su32A, su32B)
+#print(ref)
+result=dsp.arm_xor_u32(su32A, su32B).astype(int)
+print(result-ref)
+
+ref=np.bitwise_xor(ffff, su32A)
+#print(ref)
+result=dsp.arm_not_u32(su32A).astype(int)
+print(result-ref)
+
+printSubTitle("u16")
+su16A=genBitvectors(NBSAMPLES,15)
+su16B=genBitvectors(NBSAMPLES,15)
+
+ffff = (np.ones(NBSAMPLES)*(-1)).astype(np.int)
+
+
+ref=np.bitwise_and(su16A, su16B)
+#print(ref)
+result=dsp.arm_and_u16(su16A, su16B).astype(np.short)
+print(result-ref)
+
+ref=np.bitwise_or(su16A, su16B)
+#print(ref)
+result=dsp.arm_or_u16(su16A, su16B).astype(np.short)
+print(result-ref)
+
+ref=np.bitwise_xor(su16A, su16B)
+#print(ref)
+result=dsp.arm_xor_u16(su16A, su16B).astype(np.short)
+print(result-ref)
+
+ref=np.bitwise_xor(ffff, su16A)
+#print(ref)
+result=dsp.arm_not_u16(su16A).astype(np.short)
+print(result-ref)
+
+printSubTitle("u8")
+
+su8A=genBitvectors(NBSAMPLES,7)
+su8B=genBitvectors(NBSAMPLES,7)
+
+ref=np.bitwise_and(su8A, su8B)
+#print(ref)
+result=dsp.arm_and_u8(su8A, su8B).astype(np.byte)
+print(result-ref)
+
+ref=np.bitwise_or(su8A, su8B)
+#print(ref)
+result=dsp.arm_or_u8(su8A, su8B).astype(np.byte)
+print(result-ref)
+
+ref=np.bitwise_xor(su8A, su8B)
+#print(ref)
+result=dsp.arm_xor_u8(su8A, su8B).astype(np.byte)
+print(result-ref)
+
+ref=np.bitwise_xor(ffff, su8A)
+#print(ref)
+result=dsp.arm_not_u8(su8A).astype(np.byte)
+print(result-ref)
+
+#################### Quaternion tests ##################
+NBSAMPLES=3
+
+def flattenQuat(l):
+    return(np.array([list(x) for x in l]).reshape(4*len(l)))
+
+def flattenRot(l):
+    return(np.array([list(x) for x in l]).reshape(9*len(l)))
+
+# q and -q are representing the same rotation.
+# So there is an ambiguity for the tests.
+# We force the real part of be positive.
+def mkQuaternion(mat):
+    q=Quaternion(matrix=mat)
+    if q.scalar < 0:
+        return(-q)
+    else:
+        return(q)
+
+a=[2.0*Quaternion.random() for x in range(NBSAMPLES)]
+src=flattenQuat(a)
+
+
+res=flattenQuat([x.normalised for x in a])
+print(res)
+output=dsp.arm_quaternion_normalize_f32(src)
+print(output)
+print("")
+
+res=flattenQuat([x.conjugate for x in a])
+print(res)
+output=dsp.arm_quaternion_conjugate_f32(src)
+print(output)
+print("")
+
+res=flattenQuat([x.inverse for x in a])
+print(res)
+output=dsp.arm_quaternion_inverse_f32(src)
+print(output)
+print("")
+
+res=[x.norm for x in a]
+print(res)
+output=dsp.arm_quaternion_norm_f32(src)
+print(output)
+print("")
+
+a=[x.normalised for x in a]
+ra=[x.rotation_matrix for x in a]
+rb=[mkQuaternion(x) for x in ra]
+
+srca=flattenQuat(a)
+resa=dsp.arm_quaternion2rotation_f32(srca)
+resb=dsp.arm_rotation2quaternion_f32(resa)
+
+
+print(ra)
+print(resa)
+print("")
+print(rb)
+print(resb)#
+
+a=[2.0*Quaternion.random() for x in range(NBSAMPLES)]
+b=[2.0*Quaternion.random() for x in range(NBSAMPLES)]
+
+c = np.array(a) * np.array(b)
+print(c)
+
+srca=flattenQuat(a)
+srcb=flattenQuat(b)
+resc=dsp.arm_quaternion_product_f32(srca,srcb)
+
+print(resc)
+
+print(a[0]*b[0])
+res=dsp.arm_quaternion_product_single_f32(srca[0:4],srcb[0:4])
+print(res)
+