[GRASS-dev] GRASS inefficiency and FFTW

[Glynn Clements]

Dylan Beaudette wrote:

> Hmm, it doesn't look like there are all that many files actually calling or 
> defining fft():
> 
> /include/grass/gmath.h
> imagery/i.fft/fftmain.c
> imagery/i.fft/local_proto.h
> imagery/i.ifft/ifftmain.c
> imagery/i.ifft/local_proto.h
> lib/gmath/del2g.c
> lib/gmath/fft.c
> raster/r.surf.fractal/frac.h
> raster/r.surf.fractal/spec_syn.c

For dependency tracking, you're usually better off using the databases
created by tools/sql.sh (assuming that you have PostgreSQL installed):

	grass=> SELECT * FROM obj_imp WHERE symbol = 'fft' ;
	                         object                         | symbol 
	--------------------------------------------------------+--------
	 imagery/i.fft/OBJ.i686-pc-linux-gnu/fftmain.o          | fft
	 imagery/i.ifft/OBJ.i686-pc-linux-gnu/ifftmain.o        | fft
	 lib/gmath/OBJ.i686-pc-linux-gnu/del2g.o                | fft
	 raster/r.surf.fractal/OBJ.i686-pc-linux-gnu/spec_syn.o | fft
	(4 rows)
	
	grass=> SELECT * FROM obj_imp WHERE symbol = 'del2g' ;
	                  object                   | symbol 
	-------------------------------------------+--------
	 imagery/i.zc/OBJ.i686-pc-linux-gnu/main.o | del2g
	(1 row)


> Does anyone gave a general idea of the tasks which would need to be completed 
> to remove the backwards-compatibility layer? Perhaps a couple people can 
> perform these tasks. I am not a great programmer, but if the tasks are fairly 
> simple I might be able to help, with some oversight.

First, if the code is padding or scaling the data so that the array
dimensions are powers of two, remove that; FFTW doesn't have this
restriction.

Second, rather than having separate arrays of real and imaginary
components, have one array of complex values, i.e.:

change:
	double *data[2];
to:
	double (*data)[2];

change:
        data[0] = G_malloc((rows*cols)*sizeof(double));
        data[1] = G_malloc((rows*cols)*sizeof(double));
to:
        data = G_malloc((rows*cols)*2*sizeof(double));

change:
	data[0][i] = value;
	data[1][i] = 0.0;
to:
	data[i][0] = value;
	data[i][1] = 0.0;

Finally, call fft2() (which I've just added) instead of fft().

-- 
Glynn Clements <glynn at gclements.plus.com>