Reply by divy...@gmail.com●January 29, 20072007-01-29

On Dec 15 2006, 5:48 pm, "divya_ratho...@gmail.com"
<divyarath...@gmail.com> wrote:
Sometime back I wrote a FFTW based plugin for Image Apprentice. You
may

A plugin alongwith source code has just been published. It does
Gaussian High Pass Filtering using convolution (utilizes FFTW)
http://www.adislindia.com/login/pluginlinks.php
One may use it as a reference.

Reply by tcharles●December 22, 20062006-12-22

>tcharles wrote:
>> double * fft::backward(fftw_complex * in) {
>> out = new double[nx*(ny/2+1)];
>> /* Backward transformation */
>> plan_backward = fftw_plan_dft_c2r_2d(nx, ny, in, out,

FFTW_ESTIMATE);

>> fftw_execute ( plan_backward );
>> return(out);
>> }
>
>> and I have a segmentation fault in the backward function, the debugger
>
>The real output array of an out-of-place c2r transform should be of
>size nx*ny, not nx*(ny/2+1). You are confusing the size of the complex
>array with the size of the real array. See the FFTW manual.
>
>You also have a memory leak because you don't call fftw_destroy_plan.
>
>Note also that using "new" instead of "fftw_malloc" may prevent FFTW
>from exploiting SIMD instructions (e.g. SSE/SSE2), since "new" is
>probably not 16-byte aligned.
>
>Steven G. Johnson

Thanks for your answer, I wrote a destructor to call fftw_destroy_plan und
fftw_free, und made the following changes:
fft::~fft() {
fftw_destroy_plan (plan_backward );
fftw_destroy_plan (plan_forward );
fftw_free ( out );
fftw_free ( R );
}
double * fft::backward(fftw_complex * in) {
out = (double *)fftw_malloc ( sizeof ( double ) *nx*ny );
/* Backward transformation */
plan_backward = fftw_plan_dft_c2r_2d(nx, ny, in, out, FFTW_ESTIMATE);
fftw_execute ( plan_backward );
return(out);
}
but I still have size fault. The debugger output is:
Address 0x464C028 is 0 bytes inside a block of size 20,808 alloc'd
==2095== at 0x401A4CE: malloc (vg_replace_malloc.c:149)
==2095== by 0x8096AF4: fft::backward(double (*) [2]) (fft_l.C:68)
the line 68 is out = (double *)fftw_malloc ( sizeof ( double ) *nx*ny );
Thanks for your help.
Charles Tchoutat

> and I have a segmentation fault in the backward function, the debugger

The real output array of an out-of-place c2r transform should be of
size nx*ny, not nx*(ny/2+1). You are confusing the size of the complex
array with the size of the real array. See the FFTW manual.
You also have a memory leak because you don't call fftw_destroy_plan.
Note also that using "new" instead of "fftw_malloc" may prevent FFTW
from exploiting SIMD instructions (e.g. SSE/SSE2), since "new" is
probably not 16-byte aligned.
Steven G. Johnson

Reply by tcharles●December 21, 20062006-12-21

>Sometime back I wrote a FFTW based plugin for Image Apprentice. You may
>consult that. Here's the link:
>http://www.adislindia.com/pdn.html
>or alternately (direct link)
>http://divyarathore.googlepages.com/fftwpluginforimageapprentice.htm
>
>I didn't do convolution there but it is self explainatory.
>In case you still need help, just write back.
>
>regards,
>Divya Rathore
>www.adislindia.com
>(remove underscores for email ID)
>

Thank you for your answer. I still have a problem with the fftw, to
perform the convolution I wrote the following functions:
fftw_complex * fft::forward(double * in) {
forward_result = (fftw_complex *) fftw_malloc ( sizeof ( fftw_complex )
* nx *(ny/2+1) );
plan_forward = fftw_plan_dft_r2c_2d( nx, ny, in, forward_result,
FFTW_ESTIMATE);
fftw_execute ( plan_forward );
return(forward_result);
}
double * fft::backward(fftw_complex * in) {
out = new double[nx*(ny/2+1)];
/* Backward transformation */
plan_backward = fftw_plan_dft_c2r_2d(nx, ny, in, out, FFTW_ESTIMATE);
fftw_execute ( plan_backward );
return(out);
}
//double * fft::multiply(double * A, double * B) {
fftw_complex * fft::multiply(fftw_complex * A, fftw_complex * B) {
double scale = 1.0 / (nx * ny);
R = (fftw_complex *) fftw_malloc ( sizeof ( fftw_complex ) * nx
*(ny/2+1) );
for(int i= 0; i < nx; i++)
for(int j= 0; j < ny/2+1; j++){
int ij = i*(ny/2+1) + j;
R[ij][0] = (A[ij][0] * B[ij][0]
- A[ij][1] * B[ij][1]) * scale;
R[ij][1] = (A[ij][0] * B[ij][0]
+ A[ij][1] * B[ij][0]) * scale;
}
return(R);
}
and I have a segmentation fault in the backward function, the debugger
output is
at 0x82175DA: apply_hc2r (in
/home/tchoutat/Work/projekt_tp1/nina/fdp/fdp)
==27739== by 0x821A687: apply (in
/home/tchoutat/Work/projekt_tp1/nina/fdp/fdp)
==27739== Address 0x459FD68 is 0 bytes after a block of size 10,608
alloc'd
==27739== at 0x401AD0D: operator new[](unsigned)
(vg_replace_malloc.c:195)
==27739== by 0x8096B14: fft::backward(double (*) [2]) (fft_l.C:66)
And I don't understand this size problem.
Thanks for your help.
Best regards.
Charles Tchoutat

Reply by divy...@gmail.com●December 15, 20062006-12-15

tcharles wrote:

> can anyone explain me how to perform the convolution of two dimentional
> array with fftw?

Sometime back I wrote a FFTW based plugin for Image Apprentice. You may
consult that. Here's the link:
http://www.adislindia.com/pdn.html
or alternately (direct link)
http://divyarathore.googlepages.com/fftwpluginforimageapprentice.htm
I didn't do convolution there but it is self explainatory.
In case you still need help, just write back.
regards,
Divya Rathore
www.adislindia.com
(remove underscores for email ID)

Reply by ●December 14, 20062006-12-14

tcharles wrote:

> Hello,
> can anyone explain me how to perform the convolution of two dimentional
> array with fftw?
> Thanks

This method of performing convolution is a relatively straightforward
application of the FFT; see these links:
http://www.dspguide.com/ch18/2.htmhttp://cnx.org/content/m10963/latest/
This would be no different with FFTW than any other environment; just
the mechanics of the source code used to perform the FFT would be
different.
Jason

Reply by tcharles●December 14, 20062006-12-14

Hello,
can anyone explain me how to perform the convolution of two dimentional
array with fftw?
Thanks