hcrisp

11-25-2009, 12:32 PM

A user alerted me to the fact that use of FFT breaks in one case, whereas use of FFTCOMP breaks in another case. Looks like FFT has numerical limitations but FFTCOMP has memory limitations. Neither is fully-functional, and therefore PV-WAVE does not have a single FFT function that works for all cases. Let me illustrate.

Example where FFTCOMP works but FFT does not:

kernel_length = 4001L

h = FLTARR(kernel_length ) + 1.

h(0:2) = 0.

h(3) = 0.0300344

h(4) = 0.412832

h(5) = 0.882798

h(38) = 0.754520

h(39) = 0.254548

h(40) = 2.74063e-005

h(41) = 0.245479

h(42) = 0.745451

h(43) = 0.999973

h(64) = 0.972759

h(65) = 0.595404

h(66) = 0.122645

h(67:80) = 0.

h(81) = 0.240987

h(82) = 0.740877

h(83) = 0.999890

h(3918) = 0.999890

h(3919) = 0.740877

h(3920) = 0.240987

h(3921:3934) = 0.

h(3935) = 0.122645

h(3936) = 0.595404

h(3937) = 0.972759

h(3958) = 0.999973

h(3959) = 0.745451

h(3960) = 0.245479

h(3961) = 2.74063e-005

h(3962) = 0.254548

h(3963) = 0.754520

h(3996) = 0.882798

h(3997) = 0.412832

h(3998) = 0.0300344

h(3999:4000) = 0.

hinv_1 = FFT(h, 1) / kernel_length

hout_1 = FFT(hinv_1, -1)

hinv_2 = FFTCOMP(h, /BACKWARD, /COMPLEX) / kernel_length

hout_2 = FFTCOMP(hinv_2, /COMPLEX) / kernel_length

PLOT, ABS(hout_1)

OPLOT, ABS(hout_2), COLOR=230

; Notice the the FFTCOMP result is numerically correct, while the FFT result is not!

Example where FFT works but FFTCOMP does not:

kernel_length = 12000000L

h = DOUBLE(RANDOMU(s, kernel_length ))

hinv_1 = FFT(h, 1) / kernel_length

hout_1 = FFT(hinv_1, -1)

; Works!

kernel_length = 12000000L

h = DOUBLE(RANDOMU(s, kernel_length ))

hinv_2 = FFTCOMP(h, /BACKWARD, /COMPLEX) / kernel_length

;%%%FFTCOMP: Terminal error: MATH_OUT_OF_MEMORY_1

; Not enough memory, with "n" = 12000000.

;%%%

;%%%Execution halted at FFTCOMP <fftcomp.pro( 134)> (MATHSTAT_140).

;%%%Called from $MAIN$ .

Can anything be done to make both the FFT and/or FFTCOMP functions work for all cases? (I am using WAVE 8.00 on Windows XP.)

Example where FFTCOMP works but FFT does not:

kernel_length = 4001L

h = FLTARR(kernel_length ) + 1.

h(0:2) = 0.

h(3) = 0.0300344

h(4) = 0.412832

h(5) = 0.882798

h(38) = 0.754520

h(39) = 0.254548

h(40) = 2.74063e-005

h(41) = 0.245479

h(42) = 0.745451

h(43) = 0.999973

h(64) = 0.972759

h(65) = 0.595404

h(66) = 0.122645

h(67:80) = 0.

h(81) = 0.240987

h(82) = 0.740877

h(83) = 0.999890

h(3918) = 0.999890

h(3919) = 0.740877

h(3920) = 0.240987

h(3921:3934) = 0.

h(3935) = 0.122645

h(3936) = 0.595404

h(3937) = 0.972759

h(3958) = 0.999973

h(3959) = 0.745451

h(3960) = 0.245479

h(3961) = 2.74063e-005

h(3962) = 0.254548

h(3963) = 0.754520

h(3996) = 0.882798

h(3997) = 0.412832

h(3998) = 0.0300344

h(3999:4000) = 0.

hinv_1 = FFT(h, 1) / kernel_length

hout_1 = FFT(hinv_1, -1)

hinv_2 = FFTCOMP(h, /BACKWARD, /COMPLEX) / kernel_length

hout_2 = FFTCOMP(hinv_2, /COMPLEX) / kernel_length

PLOT, ABS(hout_1)

OPLOT, ABS(hout_2), COLOR=230

; Notice the the FFTCOMP result is numerically correct, while the FFT result is not!

Example where FFT works but FFTCOMP does not:

kernel_length = 12000000L

h = DOUBLE(RANDOMU(s, kernel_length ))

hinv_1 = FFT(h, 1) / kernel_length

hout_1 = FFT(hinv_1, -1)

; Works!

kernel_length = 12000000L

h = DOUBLE(RANDOMU(s, kernel_length ))

hinv_2 = FFTCOMP(h, /BACKWARD, /COMPLEX) / kernel_length

;%%%FFTCOMP: Terminal error: MATH_OUT_OF_MEMORY_1

; Not enough memory, with "n" = 12000000.

;%%%

;%%%Execution halted at FFTCOMP <fftcomp.pro( 134)> (MATHSTAT_140).

;%%%Called from $MAIN$ .

Can anything be done to make both the FFT and/or FFTCOMP functions work for all cases? (I am using WAVE 8.00 on Windows XP.)