testdsp.py 8.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370
  1. import cmsisdsp as dsp
  2. import numpy as np
  3. from scipy import signal
  4. import matplotlib.pyplot as plt
  5. from scipy.fftpack import dct
  6. r = dsp.arm_add_f32(np.array([1.,2,3]),np.array([4.,5,7]))
  7. print(r)
  8. r = dsp.arm_add_q31([1,2,3],[4,5,7])
  9. print(r)
  10. r = dsp.arm_add_q15([1,2,3],[4,5,7])
  11. print(r)
  12. r = dsp.arm_add_q7([-1,2,3],[4,127,7])
  13. print(r)
  14. r = dsp.arm_scale_f32([1.,2,3],2)
  15. print(r)
  16. r = dsp.arm_scale_q31([0x7FFF,0x3FFF,0x1FFF],1 << 20,2)
  17. print(r)
  18. r = dsp.arm_scale_q15([0x7FFF,0x3FFF,0x1FFF],1 << 10,2)
  19. print(r)
  20. r = dsp.arm_scale_q7([0x7F,0x3F,0x1F],1 << 5,2)
  21. print(r)
  22. r = dsp.arm_negate_f32([1.,2,3])
  23. print(r)
  24. r = dsp.arm_negate_q31([1,2,3])
  25. print(r)
  26. r = dsp.arm_negate_q15([1,2,3])
  27. print(r)
  28. r = dsp.arm_negate_q7(np.array([0x80,0x81,0x82]))
  29. print(r)
  30. r = dsp.arm_cmplx_conj_f32([1.,2,3,4])
  31. print(r)
  32. r = dsp.arm_cmplx_conj_q31([1,2,3,4])
  33. print(r)
  34. r = dsp.arm_cmplx_conj_q15([1,2,3,4])
  35. print(r)
  36. r = dsp.arm_cmplx_dot_prod_f32([1.,2,3,4],[1.,2,3,4])
  37. print(r)
  38. r = dsp.arm_cmplx_dot_prod_q31([0x1FFF,0x3FFF,0x1FFF,0x3FFF],[0x1FFF,0x3FFF,0x1FFF,0x3FFF])
  39. print(r)
  40. r = dsp.arm_cmplx_mult_real_f32([1.0,2,3,4],[5.,5.,5.,5.])
  41. print(r)
  42. pidf32 = dsp.arm_pid_instance_f32(Kp=1.0,Ki=1.2,Kd=0.4)
  43. print(pidf32.Kp())
  44. print(pidf32.Ki())
  45. print(pidf32.Kd())
  46. print(pidf32.A0())
  47. dsp.arm_pid_init_f32(pidf32,0)
  48. print(pidf32.A0())
  49. print(dsp.arm_cos_f32(3.14/4.))
  50. print(dsp.arm_sqrt_q31(0x7FFF))
  51. firf32 = dsp.arm_fir_instance_f32()
  52. dsp.arm_fir_init_f32(firf32,3,[1.,2,3],[0,0,0,0,0,0,0])
  53. print(firf32.numTaps())
  54. filtered_x = signal.lfilter([3,2,1.], 1.0, [1,2,3,4,5,1,2,3,4,5])
  55. print(filtered_x)
  56. print(dsp.arm_fir_f32(firf32,[1,2,3,4,5]))
  57. print(dsp.arm_fir_f32(firf32,[1,2,3,4,5]))
  58. def q31sat(x):
  59. if x > 0x7FFFFFFF:
  60. return(np.int32(0x7FFFFFFF))
  61. elif x < -0x80000000:
  62. return(np.int32(0x80000000))
  63. else:
  64. return(np.int32(x))
  65. q31satV=np.vectorize(q31sat)
  66. def toQ31(x):
  67. return(q31satV(np.round(x * (1<<31))))
  68. def q15sat(x):
  69. if x > 0x7FFF:
  70. return(np.int16(0x7FFF))
  71. elif x < -0x8000:
  72. return(np.int16(0x8000))
  73. else:
  74. return(np.int16(x))
  75. q15satV=np.vectorize(q15sat)
  76. def toQ15(x):
  77. return(q15satV(np.round(x * (1<<15))))
  78. def q7sat(x):
  79. if x > 0x7F:
  80. return(np.int8(0x7F))
  81. elif x < -0x80:
  82. return(np.int8(0x80))
  83. else:
  84. return(np.int8(x))
  85. q7satV=np.vectorize(q7sat)
  86. def toQ7(x):
  87. return(q7satV(np.round(x * (1<<7))))
  88. def Q31toF32(x):
  89. return(1.0*x / 2**31)
  90. def Q15toF32(x):
  91. return(1.0*x / 2**15)
  92. def Q7toF32(x):
  93. return(1.0*x / 2**7)
  94. firq31 = dsp.arm_fir_instance_q31()
  95. x=np.array([1,2,3,4,5])/10.0
  96. taps=np.array([1,2,3])/10.0
  97. xQ31=toQ31(x)
  98. tapsQ31=toQ31(taps)
  99. dsp.arm_fir_init_q31(firq31,3,tapsQ31,[0,0,0,0,0,0,0])
  100. print(firq31.numTaps())
  101. resultQ31=dsp.arm_fir_q31(firq31,xQ31)
  102. result=Q31toF32(resultQ31)
  103. print(result)
  104. a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
  105. b=np.array([[1.,2,3,4],[5.1,6,7,8],[9.1,10,11,12]])
  106. print(a+b)
  107. print("OK")
  108. v=dsp.arm_mat_add_f32(a,b)
  109. print(v)
  110. a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
  111. b=np.array([[1.,2,3],[5.1,6,7],[9.1,10,11],[5,8,4]])
  112. print(np.dot(a , b))
  113. v=dsp.arm_mat_mult_f32(a,b)
  114. print(v)
  115. def imToReal2D(a):
  116. ar=np.zeros(np.array(a.shape) * [1,2])
  117. ar[::,0::2]=a.real
  118. ar[::,1::2]=a.imag
  119. return(ar)
  120. def realToIm2D(ar):
  121. return(ar[::,0::2] + 1j * ar[::,1::2])
  122. a=np.array([[1. + 2j,3 + 4j],[5 + 6j,7 + 8j],[9 + 10j,11 + 12j]])
  123. b=np.array([[1. + 2j, 3 + 5.1j ,6 + 7j],[9.1 + 10j,11 + 5j,8 +4j]])
  124. print(np.dot(a , b))
  125. # Convert complex array to real array for use in CMSIS DSP
  126. ar = imToReal2D(a)
  127. br = imToReal2D(b)
  128. v=dsp.arm_mat_cmplx_mult_f32(ar,br)
  129. print(v)
  130. a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]]) / 30.0
  131. b=np.array([[1.,2,3,4],[5.1,6,7,8],[9.1,10,11,12]]) / 30.0
  132. print(a+b)
  133. aQ31=toQ31(a)
  134. bQ31=toQ31(b)
  135. v=dsp.arm_mat_add_q31(aQ31,bQ31)
  136. rQ31=v[1]
  137. r=Q31toF32(rQ31)
  138. print(r)
  139. a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
  140. print(np.transpose(a))
  141. print(dsp.arm_mat_trans_f32(a))
  142. a = np.array([[1., 2.], [3., 4.]])
  143. print(np.linalg.inv(a))
  144. print(dsp.arm_mat_inverse_f32(a))
  145. a = np.array([[1., 2.], [3., 4.]])
  146. print(np.linalg.inv(a))
  147. print(dsp.arm_mat_inverse_f64(a))
  148. a=np.array([[1.,2,3,4],[5,6,7,8],[9,10,11,12]])
  149. print(2.5*a)
  150. print(dsp.arm_mat_scale_f32(a,2.5))
  151. a=np.array([1.,2,3,4,5,6,7,8,9,10,11,12])
  152. print(np.max(a))
  153. print(np.argmax(a))
  154. print(dsp.arm_max_f32(a))
  155. print(np.mean(a))
  156. print(dsp.arm_mean_f32(a))
  157. print(np.dot(a,a))
  158. print(dsp.arm_power_f32(a))
  159. def imToReal1D(a):
  160. ar=np.zeros(np.array(a.shape) * 2)
  161. ar[0::2]=a.real
  162. ar[1::2]=a.imag
  163. return(ar)
  164. def realToIm1D(ar):
  165. return(ar[0::2] + 1j * ar[1::2])
  166. nb = 16
  167. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  168. result=np.fft.fft(signal)
  169. print(result)
  170. signalR = imToReal1D(signal)
  171. cfftf32=dsp.arm_cfft_instance_f32()
  172. status=dsp.arm_cfft_init_f32(cfftf32,nb)
  173. print(status)
  174. resultR = dsp.arm_cfft_f32(cfftf32,signalR,0,1)
  175. resultI = realToIm1D(resultR)
  176. print(resultI)
  177. signal = signal / 10.0
  178. result=np.fft.fft(signal)
  179. print(result)
  180. signalR = imToReal1D(signal)
  181. signalRQ31=toQ31(signalR)
  182. cfftq31=dsp.arm_cfft_instance_q31()
  183. status=dsp.arm_cfft_init_q31(cfftq31,nb)
  184. print(status)
  185. resultR = dsp.arm_cfft_q31(cfftq31,signalRQ31,0,1)
  186. resultI = realToIm1D(Q31toF32(resultR))*16
  187. print(resultI)
  188. signal = signal / 10.0
  189. result=np.fft.fft(signal)
  190. print(result)
  191. signalR = imToReal1D(signal)
  192. signalRQ15=toQ15(signalR)
  193. cfftq15=dsp.arm_cfft_instance_q15()
  194. status=dsp.arm_cfft_init_q15(cfftq15,nb)
  195. print(status)
  196. resultR = dsp.arm_cfft_q15(cfftq15,signalRQ15,0,1)
  197. resultR=Q15toF32(resultR)
  198. resultI = realToIm1D(resultR)*16
  199. print(resultI)
  200. nb = 128
  201. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  202. result=np.fft.fft(signal)
  203. print(result)
  204. cfftradix4f32=dsp.arm_cfft_radix4_instance_f32()
  205. rfftf32=dsp.arm_rfft_instance_f32()
  206. status=dsp.arm_rfft_init_f32(rfftf32,cfftradix4f32,nb,0,1)
  207. print(status)
  208. resultI = dsp.arm_rfft_f32(rfftf32,signal)
  209. print(result)
  210. nb = 128
  211. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  212. signalRQ31=toQ31(signal)
  213. result=np.fft.fft(signal)
  214. print(result)
  215. rfftq31=dsp.arm_rfft_instance_q31()
  216. status=dsp.arm_rfft_init_q31(rfftq31,nb,0,1)
  217. print(status)
  218. resultI = dsp.arm_rfft_q31(rfftq31,signalRQ31)
  219. resultI=Q31toF32(resultI)*(1 << 7)
  220. print(result)
  221. nb = 128
  222. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  223. signalRQ15=toQ15(signal)
  224. result=np.fft.fft(signal)
  225. print(result)
  226. rfftq15=dsp.arm_rfft_instance_q15()
  227. status=dsp.arm_rfft_init_q15(rfftq15,nb,0,1)
  228. print(status)
  229. resultI = dsp.arm_rfft_q15(rfftq15,signalRQ15)
  230. resultI=Q15toF32(resultI)*(1 << 7)
  231. print(result)
  232. nb = 128
  233. nb2=64
  234. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  235. result=dct(signal,4,norm='ortho')
  236. print(result)
  237. cfftradix4f32=dsp.arm_cfft_radix4_instance_f32()
  238. rfftf32=dsp.arm_rfft_instance_f32()
  239. dct4f32=dsp.arm_dct4_instance_f32()
  240. status=dsp.arm_dct4_init_f32(dct4f32,rfftf32,cfftradix4f32,nb,nb2,0.125)
  241. print(status)
  242. state=np.zeros(2*nb)
  243. resultI = dsp.arm_dct4_f32(dct4f32,state,signal)
  244. print(resultI)
  245. signal = signal / 10.0
  246. result=dct(signal,4,norm='ortho')
  247. signalQ31=toQ31(signal)
  248. cfftradix4q31=dsp.arm_cfft_radix4_instance_q31()
  249. rfftq31=dsp.arm_rfft_instance_q31()
  250. dct4q31=dsp.arm_dct4_instance_q31()
  251. status=dsp.arm_dct4_init_q31(dct4q31,rfftq31,cfftradix4q31,nb,nb2,0x10000000)
  252. print(status)
  253. state=np.zeros(2*nb)
  254. resultI = dsp.arm_dct4_q31(dct4q31,state,signalQ31)
  255. resultI=Q31toF32(resultI)*(1 << 7)
  256. nb = 128
  257. nb2=64
  258. signal = np.cos(2 * np.pi * np.arange(nb) / nb)
  259. signal = signal / 10.0
  260. result=dct(signal,4,norm='ortho')
  261. signalQ15=toQ15(signal)
  262. cfftradix4q15=dsp.arm_cfft_radix4_instance_q15()
  263. rfftq15=dsp.arm_rfft_instance_q15()
  264. dct4q15=dsp.arm_dct4_instance_q15()
  265. status=dsp.arm_dct4_init_q15(dct4q15,rfftq15,cfftradix4q15,nb,nb2,0x1000)
  266. print(status)
  267. state=np.zeros(2*nb)
  268. resultI = dsp.arm_dct4_q15(dct4q15,state,signalQ15)
  269. resultI=Q15toF32(resultI)*(1 << 7)
  270. from pylab import figure, clf, plot, xlabel, ylabel, xlim, ylim, title, grid, axes, show
  271. figure(1)
  272. plot(np.absolute(signal))
  273. t = np.arange(nb)
  274. freq = np.fft.fftfreq(t.shape[-1])
  275. resultmag=np.absolute(result)
  276. figure(2)
  277. plot(resultmag)
  278. figure(3)
  279. cmsigmag=np.absolute(resultI)
  280. plot(cmsigmag)
  281. show()
  282. biquadf32 = dsp.arm_biquad_casd_df1_inst_f32()
  283. numStages=1
  284. state=np.zeros(numStages*4)
  285. coefs=[1.,2,3,4,5]
  286. dsp.arm_biquad_cascade_df1_init_f32(biquadf32,1,coefs,state)
  287. print(dsp.arm_biquad_cascade_df1_f32(biquadf32,[1,2,3,4,5]))