lebedov / scikit-cuda

Python interface to GPU-powered libraries
http://scikit-cuda.readthedocs.org/
Other
986 stars 179 forks source link

gpuarray.empty arguments #251

Closed saquibntt closed 6 years ago

saquibntt commented 6 years ago

Hi,

I am trying to do FFT on a wav source read using scipy like this. After following the demo provided I am a bit confused with the understanding of methods used in FFT computation. Please help me I didn't find the documentation for the methods.

The code:

import matplotlib.pyplot as plt
from scipy.io import wavfile as wav
#from scipy.fftpack import fft
import numpy as np
import pycuda.autoinit
import pycuda.gpuarray as gpuarray
import skcuda.fft as cu_fft

rate , data = wav.read('exam.wav')
x = np.asarray(data, np.float32)

print (data.shape[0]) #30045
print (data.shape[1]) # tuple index out of range
xf = np.fft.fft(data)
y = np.real(np.fft.ifft(xf))
x_gpu = gpuarray.to_gpu(data)
xf_gpu = gpuarray.empty(data.shape[0], data.shape[1]//2+1, np.complex64)
plan_forward = cu_fft.Plan(x_gpu.shape, np.float32, np.complex64)
cu_fft.fft(x_gpu, xf_gpu, plan_forward)

The error is

tuple index out of range

Now, even I know its because my data.shape[1] is null(as I printed the shape). I don't know what parameters I should pass.

Thanks.

lebedov commented 6 years ago

You are getting the exception because data is 1D and therefore data.shape[1] is undefined. The same problem would presumably occur on line xf_gpu = gpuarray.empty(data.shape[0], data.shape[1]//2+1, np.complex64); in any case, if data is 1D, then xf_gpu also only needs to be 1D.

saquibntt commented 6 years ago

I am trying the code below. I don't know why it gives false status. Is there any way to check if this is working properly. I am trying on Jetson TX2. python2 numpy 1.14.2 Cuda 9

from __future__ import print_function
import pycuda.autoinit
import pycuda.gpuarray as gpuarray
import numpy as np

from scipy.io import wavfile as wav
import skcuda.fft as cu_fft

print('Testing fft/ifft..')
#N = 4096 * 16
rate , data = wav.read('exam.wav')

x = np.asarray(data, np.float32)
xf = np.fft.fft(x)
y = np.real(np.fft.ifft(xf))

x_gpu = gpuarray.to_gpu(x)
xf_gpu = gpuarray.empty(data.shape[0], np.complex64)
plan_forward = cu_fft.Plan(x_gpu.shape, np.float32, np.complex64)
cu_fft.fft(x_gpu, xf_gpu, plan_forward)

y_gpu = gpuarray.empty_like(x_gpu)
plan_inverse = cu_fft.Plan(x_gpu.shape, np.complex64, np.float32)
cu_fft.ifft(xf_gpu, y_gpu, plan_inverse, True)

print('Success status: %r' % np.allclose(y, y_gpu.get(), atol=1e-6)) 

Please help! Thanks

lebedov commented 6 years ago

One quick check you can perform is plot y and y_gpu (or perhaps the lowest 100 frequency components) to see whether they actually are similar. Alternately, try gradually raising atol to see whether np.allclose returns True for a higher tolerance (I tried both of the above for some random 16 bit WAV file and had to raise atol to 1e-4 for the test to succeed). Different FFT implementations (i.e., numpy and CUFFT) will almost always produce slightly different results (especially when using single precision floating point) even if both implementations are mathematically correct.

saquibntt commented 6 years ago

Thank you @lebedov for your suggestions. I had to lower atol to 1e-2 for the status to be True. However, I could not compare the plots as plotting y_gpu gave a ValueError: setting an array element with a sequence. I am looking into the reason here.

Here is the tracelog:

Testing fft/ifft..
16000
Success status: True
Traceback (most recent call last):
  File "fft_demo.py", line 37, in <module>
    plt.plot(y_gpu)
  File "/usr/lib/python2.7/dist-packages/matplotlib/pyplot.py", line 3154, in plot
    ret = ax.plot(*args, **kwargs)
  File "/usr/lib/python2.7/dist-packages/matplotlib/__init__.py", line 1814, in inner
    return func(ax, *args, **kwargs)
  File "/usr/lib/python2.7/dist-packages/matplotlib/axes/_axes.py", line 1425, in plot
    self.add_line(line)
  File "/usr/lib/python2.7/dist-packages/matplotlib/axes/_base.py", line 1708, in add_line
    self._update_line_limits(line)
  File "/usr/lib/python2.7/dist-packages/matplotlib/axes/_base.py", line 1730, in _update_line_limits
    path = line.get_path()
  File "/usr/lib/python2.7/dist-packages/matplotlib/lines.py", line 925, in get_path
    self.recache()
  File "/usr/lib/python2.7/dist-packages/matplotlib/lines.py", line 621, in recache
    y = np.asarray(yconv, np.float_)
  File "/usr/local/lib/python2.7/dist-packages/numpy/core/numeric.py", line 492, in asarray
    return array(a, dtype, copy=False, order=order)
ValueError: setting an array element with a sequence.

Any suggestions?

lebedov commented 6 years ago

You need to explicitly transfer the contents of y_gpu back to system memory before you can plot them: plt.plot(y_gpu.get())

saquibntt commented 6 years ago

Thank you! It works