# Wiener deconvolution

Hi everyone.

I’ve been trying to develop a wiener deconvolution algorithm using CUDA. I’m following the form of it’s equation where to get the G(f) you divide the conjugate of the H(f) (the PSF after FFT) by the magnitude of H(f) plus a noise constant.
My PSF is zero padded to the size of the image, the FFT is performed as a 2D C2C cufft transform (I converse the real numbers to complex numbers first).
My functions look as follow:

global void conjugate(cufftComplex* signal_in) {

``````int id = blockDim.x*blockIdx.x + threadIdx.x;
float temp_imag[256];
float temp_real[256];

temp_imag[idx] = (-1)*(signal_in[id].y);
temp_real[idx] = signal_in[id].x;

signal_in[id].y = temp_imag[idx];
signal_in[id].x = temp_real[idx];
``````

};

global void mul(cufftComplex* signal_in, cufftComplex* signal_out, int* width, int* length) {

``````int id = blockDim.x*blockIdx.x + threadIdx.x;

cuComplex temp_real2[256];

float scale = 1.0f / (float)(*width*(*length));
temp_real2[idx] = cuCmulf(signal_in[id],signal_out[id]);

signal_in[id] = make_cuFloatComplex(scale*cuCrealf(temp_real2[idx]), scale*cuCimagf(temp_real2[idx]));

return;
``````

}

global void addnoise(cufftComplex* signal_in, cufftComplex* signal_out) {

``````int id = threadIdx.x;
int idx = threadIdx.x + blockDim.x*blockIdx.x;
float temp_real[256];
temp_real[id] = signal_in[idx].x;

signal_out[idx].x = (temp_real[id] + NOISECONTS) ;
``````

}

global void divide(cufftComplex* signal_in, cufftComplex* signal_in_noisy, int* width, int* length) {

``````int id = blockDim.x*blockIdx.x + threadIdx.x;
cuComplex temp[256];

temp[idx] = cuCdivf(signal_in[id], signal_in_noisy[id]);

signal_in[id] = temp[idx];
``````

};

global void magnitude(cufftComplex* signal_in, double* magnitude) {

``````int id = blockDim.x*blockIdx.x + threadIdx.x;

double temp_magRE[256];
double temp_magIM[256];

temp_magRE[idx] = signal_in[id].x;
temp_magIM[idx] = signal_in[id].y;