Here's the mentioned file:
crd.lbl.gov/~dhbailey/dhbpapers/convops.pdf
Martin
Reply by Ron N.●September 6, 20072007-09-06
On Sep 3, 11:36 pm, "jravi" <ravi.jand...@nxp.com> wrote:
> I want to compute an M point DFT of a sequence x(n)where M is not a power
> of 2. One way to do this would be:
> 1. Pad (N - M) zeros and compute an N point DFT where N = power of 2.
> 2. Find the (continuous) DTFT (discrete time fourier transform) of x(n)
> from the N DFT samples by interpolating the DFT with the sinc
> function.
>
> X(w) = Summation {Y(k)* sinc[N/2 * (w - 2*pi*k/N)] }
> k=0 to N-1
> where
> X(w) = DTFT of the sequence x(n)
> Y(k) = N point DFT samples of x(n)
> N = power of 2
> w = omega (radians/s)
>
> 3. Sample the DTFT at M equally spaced points (2*pi/M apart) which would
> give the M point DFT.
>
> X(p) = Summation {Y(k)* sinc[pi*N * (p/M -k/N)] } 0<= p <= M-1
> k=0 to N-1
> where
> X(p) = M point DFT of x(n)
> Y(k) = N point DFT of x(n)
>
> This method would require computing the next higher power of 2 FFT and
> additional M * N multiplications after that. Also, I need to store the M*N
> sinc values in a table.
>
> Is there any simpler method to achieve the same?? I want to compute the
> exact M point DFT ...
The simplest method to code would be a plain explicit DFT.
If you want an exact M point DFT, then your above method won't
work as stated. Non-windowed Sinc interpolation requires a infinite
width Sinc kernel. You would need to circularly pre-wrap a Sinc an
infinite number of times around your N point window (or at least until
it decayed below your definition of "exact"), and put that summation
in your N*M filter tap table.
Are you sure that a windowed Sinc or an interpolated approximation
to a windowed Sinc won't do? The table would be smaller, or you
could do it without any table at a higher computational cost.
IMHO. YMMV.
--
rhn A.T nicholson d.0.t C-o-M
Reply by ●September 6, 20072007-09-06
On 4 Sep., 08:36, "jravi" <ravi.jand...@nxp.com> wrote:
> I want to compute an M point DFT of a sequence x(n)where M is
> not a power of 2.
You can do an FFT of the next *lower* power-of-2 length, compute the
values missing from one end by time-domain convolution, and subtract
the latter from the other FFT end to remove the time aliasing. I once
found a write-up online investigating the running-time properties of
this approach, but unfortunately I can't remember the author or locate
it without that info right now. I'll look it up later unless someone
beats me to it.
Martin
Reply by John E. Hadstate●September 5, 20072007-09-05
"John Hadstate" <jh113355@hotmail.com> wrote in message
news:1188999271.501266.164990@k79g2000hse.googlegroups.com...
> On Sep 4, 2:36 am, "jravi" <ravi.jand...@nxp.com> wrote:
>> I want to compute an M point DFT of a sequence x(n)where
>> M is not a power
>> of 2.
>
>
> Rader and Gold has a section on how to do exactly what
> you're asking
> for.
Rader, Charles M. and Gold, Bernard, Digital Processing of
Signals, (C) McGraw-Hill, 1969, ISBN 07-023642-9
Chapter 7, "High-speed Convolution and Correlation with
Applications to Digital Filtering, was contributed by Thomas
G. Stockham, Jr. The relevant text appears in Section 7.8,
"Discrete Fourier Transforms Expressed as Convolutions or
Correlations", on pages 213-215. Quoting a piece:
"The major step, computing the correlation, can be
implemented by high-speed techniques, so that we have
arrived at the seemingly futile result that we can use an
FFT to compute a DFT in a somewhat roundabout manner.
However, this result is not really useless, because we have
bypassed several of the restrictions inherent in the
computation of a DFT.... For example, there is no
restriction that the number of points ... be highly
composite since an FFT based on a power of 2 can be used to
compute lagged products with any number of points."
I doubt that you can find a copy of Rader and Gold, but you
might be able to find Stockham's chapter as a paper
somewhere else. Good luck!
Reply by John Hadstate●September 5, 20072007-09-05
On Sep 4, 2:36 am, "jravi" <ravi.jand...@nxp.com> wrote:
> I want to compute an M point DFT of a sequence x(n)where M is not a power
> of 2.
> I want to compute an M point DFT of a sequence x(n)where M is not a power
> of 2. One way to do this would be:
> 1. Pad (N - M) zeros and compute an N point DFT where N = power of 2.
> 2. Find the (continuous) DTFT (discrete time fourier transform) of x(n)
> from the N DFT samples by interpolating the DFT with the sinc
> function.
Continuous? How. Digital computers don't do continuous anything.
> X(w) = Summation {Y(k)* sinc[N/2 * (w - 2*pi*k/N)] }
> k=0 to N-1
> where
> X(w) = DTFT of the sequence x(n)
> Y(k) = N point DFT samples of x(n)
> N = power of 2
> w = omega (radians/s)
>
> 3. Sample the DTFT at M equally spaced points (2*pi/M apart) which would
> give the M point DFT.
>
> X(p) = Summation {Y(k)* sinc[pi*N * (p/M -k/N)] } 0<= p <= M-1
> k=0 to N-1
> where
> X(p) = M point DFT of x(n)
> Y(k) = N point DFT of x(n)
>
> This method would require computing the next higher power of 2 FFT and
> additional M * N multiplications after that. Also, I need to store the M*N
> sinc values in a table.
>
> Is there any simpler method to achieve the same?? I want to compute the
> exact M point DFT and not just increase the frequency resolution through
> zero padding. Please let me know.
Jerry
--
Engineering is the art of making what you want from things you can get.
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Reply by jravi●September 4, 20072007-09-04
I want to compute an M point DFT of a sequence x(n)where M is not a power
of 2. One way to do this would be:
1. Pad (N - M) zeros and compute an N point DFT where N = power of 2.
2. Find the (continuous) DTFT (discrete time fourier transform) of x(n)
from the N DFT samples by interpolating the DFT with the sinc
function.
X(w) = Summation {Y(k)* sinc[N/2 * (w - 2*pi*k/N)] }
k=0 to N-1
where
X(w) = DTFT of the sequence x(n)
Y(k) = N point DFT samples of x(n)
N = power of 2
w = omega (radians/s)
3. Sample the DTFT at M equally spaced points (2*pi/M apart) which would
give the M point DFT.
X(p) = Summation {Y(k)* sinc[pi*N * (p/M -k/N)] } 0<= p <= M-1
k=0 to N-1
where
X(p) = M point DFT of x(n)
Y(k) = N point DFT of x(n)
This method would require computing the next higher power of 2 FFT and
additional M * N multiplications after that. Also, I need to store the M*N
sinc values in a table.
Is there any simpler method to achieve the same?? I want to compute the
exact M point DFT and not just increase the frequency resolution through
zero padding. Please let me know.