Efficient algorithm for calculating of the small part of DFT

[dora]

fft fixed-point dsp 16 winograd

Hello friends,

I have 24 long real data vector and I would like to calculate the
1,2,3,4 Fourier coefficients only. (the zero coeff is the DC term)

I would like to do this as efficiently as possible using 16 bits fixed point DSP.

What I am having in mind is the Goertzel algorithm for calculation of the single FFT bin.
It seems to me that the 4 Goertzel filters can be calculate more efficiently together then
separately but still I don't have the formula. Currently I am thinking in this respect.

In addition I need only the magnitude value of the 1,2,3,4 Fourier coeffs so probably this
can simplify the things even more.

So can someone propose solution of my problem?

Thank you in advance!
dora

 

[zorro]

calculating dft

Hi Dora,

I would take the flow graph of a FFT algorithm and mark back the branches from the desired results towards the data. In this way I would have the subset of operations that have to be performed for a reduced FFT-like algorithm. Note that there are several possible decompositions of the data length. Maybe this needs less computation that several Goertzel filters or the direct computation of DFT bins.
Regards

Z

 

[berraco]

Re: Efficient algorithm for calculating of the small part of

Hi, Dora.
Perhaps the modified Goertzel algorithm proposed here

**broken link removed**

would suit you. I hope it. Anyway, I will continue thinking on it.

 

[dora]

Re: Efficient algorithm for calculating of the small part of

Hi zorro,

yes i also thought about something like this.
Unfortunately my data are 24 points long (it is part image projection)
so I can not use the classical FFT. Probably I can pad to 32 points but in this case willget the interpolated spectrum in other frequency grid. Probably this modifiead grid will be also acceptable for me however.

Zorro do you know fast algoritham for FFT calculation of the 24 real points? I know that such algorithams exists but most probably they are quite complex for implementation and less efficient than FFT.

Best Regards
dora

 

[zorro]

Re: Efficient algorithm for calculating of the small part of

Hi Dora,

See the section “FFT Algorithms for composite N” in “Discrete Time Signal Processing” by Oppenheim and Schafer (section 9.6 in the 1st edition). The general algorithm is explained, and in fig. 9.30 - 9.31 it is shown the flow graph of a 12-point FFT. References for prime factor algorithms are given.
In addition, if you have purely real data, you can use a N-point algorithm for calculating a 2N-point transform, so this simple 12-point algorithm could be useful for you. See section 6.2.2 in “Digital Signal Processing (2nd. Ed)” by Proakis and Manolakis.
Let me know your progress.
Regards

Z

 

[GeorgeM]

Re: Efficient algorithm for calculating of the small part of

As zorro said, to calculate 4 bins (spectral lines) of 24 use DFT (Direct Fourier Transform) that's easiest way to get them in terms of the processing time and your efforts.

Regards,
GeorgeM

 

[Aser]

Re: Efficient algorithm for calculating of the small part of

There are a lot of fast DFT algorithms which have the minimized number of multiplications and additions like Winograd algorithm, algorithms based on the fast convolution.
But all of them suffer from their unregularity.
Therefore they have to be programmed in DSP microprocessor using inlining.
As a result, they are very hard to progam and debug, they are not hardware efficient.
Therefore when the vector length is less than 32-64 the best way is to calculate the sraightforward DFT.
DSP microprocessor provides the loop support for such an algorithm, memory addressing, and MAC operations.
And it will turn that it is faster than some "optimized algorithm" especially when only 4 bins are needed.

 

[dora]

Re: Efficient algorithm for calculating of the small part of

Hi Zorro,

>Let me know your progress.
Currently i selected the modified Goertzel approach. Basically for 2 reasons:

1. I have 4 Goertzel filters, one for each bins and its appears that the filter for the 4-th bin uses 2*cos(2*pi*4/24) ==1 so it is very simple (no multiplications).

2. Using filter approach it is simpler for me to implement the algorithm using fix point arithmetic for the precision which I need.
I have clean understanding of the filter gain for each freq , step response, DC gain ... so I can estimate the sailing factors needed to avoid overflows.

I managed to do the job and now it seems that the time bottleneck is in the other part. But let me describe more in details.

I need to estimate the direction field of the fingerprint image. So I take 24x24 window, put it at the certain place on the image , rotate it at certain angle, sum along the rows of the rotated window (made image projection at the selected direction) so I got 24 long vector. Calculate the DFT of this vector to see if there is big spectral energy in one of the first 4 bins (covering the ridge spatial frequency of interest). Then for the next direction ... like this

Currently I made the DFT much faster then making of the projection sums. Unfortunately the speed which I got (mainly due to projections sum)is not acceptable for my fingerprint application about 4 sec using 144Mhz
tms320c5509.

So currently I reconsider using of the spectral methods for fingerprint direction fields estimation. Probably some other methods will be robust for this application?

All suggestions will be highly appreciated
Dora

 

[zorro]

Re: Efficient algorithm for calculating of the small part of

Hi Dora,

An idea: what about performing a double (2-D) DFT over the image? If I’m not confused, an energy peak should appear in the spectral domain at a frequency pair (Fx, Fx); the pitch between lines would be 1/sqrt(Fx^2 + Fy^2), and the orientation angle would be atan(Fx/Fy).
If N bins out of M in the one-dimensional DFT is good in your case, approximately pi/4*N^2 out of M^2 of the 2-D DFT values are needed with this method. The algorithm for calculate DFT’s could be FFT or Goertzel. No projection has to be calculated.
What do you think?
Best regards

Z

 

[dora]

Re: Efficient algorithm for calculating of the small part of

Hello Zorro,

You gave me very good idea (as usually)
Since I am not so familiar with the 2D Fourier transform I decided to test you approach with the example.
It seems that the spectral peak indeed is located at (Fx, Fy) as you described.

Several days ago I Implemented the gradient approach. Let me describe it
The image is divided on let say 12x12 blocks. For each block the direction is estimated as follows:

For each pixel in the block the Sobel operator is applied


dI | -1 -2 -1 | dI | -1 0 1 |
--- = | 0 0 0 | --- = | -2 0 2 |
dy | 1 2 1 | dx | -1 0 2 |


So we have the approximation of the gradient (dI_dx, dI_dy), rotate it 90deg and get normal vector at the point
(-dI_dy, dI_dx). Then in the literature they advice LMS averaging to be used (to speed up the procedure I use
simple averaging) So I got speed several times bigger then the spectral approach (with projections)

Now I am going to implement the spectral algorithm you proposed and will compare times and results.
About your idea: If I need good accuracy (In fact I don't need good accuracy, let say 10 deg) the
spectral peak should be not so close to the origin (DC is at the origin (middle) of the 2D spectrum image).
This will impose some image block size limitations (if I don't use padding). I will check !


Thank you Zorro!
Dora

 

[zorro]

Re: Efficient algorithm for calculating of the small part of

Hello Dora,

Very interesting.

About your idea: If I need good accuracy (In fact I don't need good accuracy, let say 10 deg) the
spectral peak should be not so close to the origin (DC is at the origin (middle) of the 2D spectrum image).
This will impose some image block size limitations (if I don't use padding). I will check !

Right. If the peak is at the origin, it is a uniform image and no direction can be defined. If the peak is at a distance of 1 bin from the origin, that means that there is only one (relatively wide) line, and the direction is not as "clear" as in the case there were several (narrow) lines.
I think that interpolation can be used for estimating the location of the peak between bins, e.g., the gravity center ot the 2-D power spectrum.
If necessary, windowing the image previous to transform could improve the estimation.
Regards

Z