On Jan 18, 8:12�am, Greg Heath <he...@alumni.brown.edu> wrote: ...> More simply, > > 1. fftshift shifts and reflects the positive time points to the left > hand > � � side of the window causing the original middle point at t=0 > � � to become the right most point. > whereas > 2. ifftshift shifts and reflects the negative time components to the > � � right hand side of the window causing the original middle point > � � at t=0 to become the left most point.okay, Greg, i went to the MW site to read about ifftshift() http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ifftshift.html since i didn't know about it before. i've known and used fftshift() for even-lengthed arrays for 15+ years. i know precisely what fftshift () does to even-length arrays. okay. the word "reflect" and "reflects" is in error and not in the MATLAB doc. for even-lengthed arrays, fftshift() *swaps* the first half data with the second half along every dimension. equivalent to: Y = fftshift(X); or N = length(X); Y = [X(N/2+1:N) X(1:N/2)]; it just swaps the two halves and using it again will swap the halves back. now the MATLAB lit says that fftshift() does not undo itself for odd- lengthed arrays and that's why you need ifftshift(). that's very well the case and i don't know what fftshift (or ifftshift) does for odd- length arrays anyway. but, assuming N even, the purpose of fftshift() in the input to fft() is precisely what i said. let's say you have a large audio file and you want to make a sliding spectrum analyzer for it with MATLAB. try reading a big long wav file by frames (say your frames are 1024 long). pick out a 1024 snippet of data and window it with hann or hamming or kaiser or the tapered window of your choice. in your frame array (of length 1024) the data is bigger in the middle than at the edges, right? now fft that frame and do it again but with the two half-frames swapped (use fftshift, if you want, or do it "manually") going into the fft. look at the two fft outputs. the magnitudes should be identical, but comparing the phases might be a learning opportunity. one of them will look better than the other. and there's a reason for that. r b-j
MATLABs FFT and fftshifting the input
Started by ●January 11, 2010
Reply by ●January 19, 20102010-01-19
Reply by ●January 20, 20102010-01-20
On Jan 19, 12:47 am, robert bristow-johnson <r...@audioimagination.com> wrote:> On Jan 18, 8:12 am, Greg Heath <he...@alumni.brown.edu> wrote: > ... > > > More simply, > > > 1. fftshift shifts and reflects the positive time points to the leftWHOOPS! delete "and reflects"> > hand side of the window causing the original middle point at t=0 > > to become the right most point. > > whereas > > 2. ifftshift shifts and reflects the negative time components to theWHOOPS! delete "and reflects"> > right hand side of the window causing the original middle point > > at t=0 to become the left most point. > > okay, Greg, i went to the MW site to read about ifftshift() > > http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ifftshift.html > > since i didn't know about it before. i've known and used fftshift() > for even-lengthed arrays for 15+ years. i know precisely what fftshift > () does to even-length arrays.For even-length arrays ifftshift and fftshift are equivalent. The source code is exactly the same except fftshift uses ceil(N/2) and ifftshift uses floor(N/2). See below.> okay. the word "reflect" and "reflects" is in error and not in the > MATLAB doc.I agree. Thanks for the correction.> for even-lengthed arrays, fftshift() *swaps* the first half data with > the second half along every dimension. equivalent to: > > Y = fftshift(X); > or > N = length(X); > Y = [X(N/2+1:N) X(1:N/2)]; > > it just swaps the two halves and using it again will swap the halves > back.Yes. ifftshift and ffthift are always inverses. However, when N is even they are equivalent. Of course, when N is even, you have to remember that there are more negative time points than positive time points.> now the MATLAB lit says that fftshift() does not undo itself for odd- > lengthed arrays and that's why you need ifftshift(). that's very well > the case and i don't know what fftshift (or ifftshift) does for odd- > length arrays anyway.ifftshift(x(1:N)) = [ x(1+floor(N/2):N), x(1:floor(N/2))] fftshift(X(1:N)) = [ X(1+ceil(N/2):N), X(1:ceil(N/2))] Again, I recommend that MATLAB newbies use ifftshift in the ifft (time) domain and fftshift in the fft (frequency) domain. Then most of the worry about whether N is even or odd is mitigated .> but, assuming N even, the purpose of fftshift() in the input to fft() > is precisely what i said. > > let's say you have a large audio file and you want to make a sliding > spectrum analyzer for it with MATLAB. > > try reading a big long wav file by frames (say your frames are 1024 > long). pick out a 1024 snippet of data and window it with hann or > hamming or kaiser or the tapered window of your choice. in your frame > array (of length 1024) the data is bigger in the middle than at the > edges, right? > > now fft that frame and do it again but with the two half-frames > swapped (use fftshift, if you want, or do it "manually") going into > the fft. look at the two fft outputs. the magnitudes should be > identical, but comparing the phases might be a learning opportunity. > one of them will look better than the other. and there's a reason for > that.There have been examples in comp.soft-sys.matlab where an odd number of points were used so that the symmetry axis t = 0 could be EXACTLY in the middle of the waveform. When fftshift was used in the time domain to attempt to reference the phase to the symmetry point, the phase plot did not have the expected behavior. After I resolved the problem, I found the best advice was to only use fftshift in the fft domain and only use ifftshift in the ifft domain. I have enclosed one of the codes used to investigate these issues. The time waveform is Gaussian. When the phase reference is chosen "correctly" the transform is also Gaussian ( real and positive with zero phase). However, when the phase is not chosen "correctly" the imaginary part can be nonzero and/or the real part can be negative. In a few cases when the real part should be positive, roundoff error yields small ( ~1e-5) negative values. As a result, there will be phase spikes from zero to pi. Hope this helps. Greg close all, clear all, clc, k=0 format short g for N = 14:15 T = 5, dt = T/N if ceil(N/2)==floor(N/2) t = dt*[ -N/2 : N/2-1 ]; %[ -T/2 : dt : T/2-dt ] else t = dt*[-(N-1)/2:(N-1)/2 ]; % [-(T-dt)/2:dt:(T-dt)/2] end minmaxt = minmax(t) z1 = exp(-t.^2); [z1max i1max] = max(z1) z2 = ifftshift(z1); [z2max i2max] = max(z2) z3 = fftshift(z1); [z3max i3max] = max(z3) maxabsz2z3 = maxabs(z2-z3) k=k+1, figure(k) subplot(3,1,1) plot(z1) ylabel('x(t)') subplot(3,1,2) plot(z2,'r') ylabel('ifftshift(x(t))') subplot(3,1,3) plot(z3,'g') ylabel('fftshift(x(t))') Z1 = fftshift(fft(z1))/N; X1 = real(Z1); Y1 = imag(Z1); A1 = abs(Z1); P1 = angle(Z1); minmaxX1 = minmax(X1) minmaxY1 = minmax(Y1) minmaxP1 = minmax(P1) k=k+1,figure(k) subplot(2,2,1) plot(X1) ylabel('Real') title('FOURIER TRANSFORM ') subplot(2,2,2) plot(Y1) ylabel('Imaginary') title('OF exp(-(t-T/2).^2') subplot(2,2,3) plot(A1) ylabel('Magnitude') subplot(2,2,4) plot(P1) ylabel('Phase') Z2 = fftshift(fft(z2))/N; X2 = real(Z2); Y2 = imag(Z2); A2 = abs(Z2); P2 = angle(Z2); minmaxX2 = minmax(X2) minmaxY2 = minmax(Y2) minmaxP2 = minmax(P2) k=k+1,figure(k) subplot(2,2,1) plot(X2) ylabel('Real') title('FOURIER TRANSFORM ') subplot(2,2,2) plot(Y2) ylabel('Imaginary') title(' OF IFFTSHIFT(exp(-(t-T/2).^2)') subplot(2,2,3) plot(A2) ylabel('Magnitude') subplot(2,2,4) plot(P2) ylabel('Phase') Z3 = fftshift(fft(z3))/N; X3 = real(Z3); Y3 = imag(Z3); A3 = abs(Z3); P3 = angle(Z3); minmaxX3 = minmax(X3) minmaxY3 = minmax(Y3) minmaxP3 = minmax(P3) k=k+1,figure(k) subplot(2,2,1) plot(X3) ylabel('Real') title('FOURIER TRANSFORM ') subplot(2,2,2) plot(Y3) ylabel('Imaginary') title(' OF FFTSHIFT(exp(-(t-T/2).^2)') subplot(2,2,3) plot(A3) ylabel('Magnitude') subplot(2,2,4) plot(P3) ylabel('Phase') end for j = k:-1:1 figure(j) end
Reply by ●January 20, 20102010-01-20
On Jan 20, 4:07�am, Greg Heath <he...@alumni.brown.edu> wrote:> On Jan 19, 12:47 am, robert bristow-johnson > > <r...@audioimagination.com> wrote: > > On Jan 18, 8:12 am, Greg Heath <he...@alumni.brown.edu> wrote: > > ... > > > > More simply, > > > > 1. fftshift shifts and reflects the positive time points to the left > > WHOOPS! delete "and reflects" > > > > � �hand side of the window causing the original middle point at t=0 > > > � �to become the right most point. > > > whereas > > > 2. ifftshift shifts and reflects the negative time components to the > > WHOOPS! delete "and reflects" > > > > � � right hand side of the window causing the original middle point > > > � � at t=0 to become the left most point. > > > okay, Greg, i went to the MW site to read about ifftshift() > > >http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ifftshift.html > > > since i didn't know about it before. �i've known and used fftshift() > > for even-lengthed arrays for 15+ years. �i know precisely what fftshift > > () does to even-length arrays. > > For even-length arrays ifftshift and fftshift are equivalent. > The source code is exactly the same except fftshift uses > ceil(N/2) and ifftshift uses floor(N/2). See below. > > > okay. �the word "reflect" and "reflects" is in error and not in the > > MATLAB doc. > > I agree. Thanks for the correction. > > > for even-lengthed arrays, fftshift() *swaps* the first half data with > > the second half along every dimension. �equivalent to: > > > � Y = fftshift(X); > > or > > � N = length(X); > > � Y = [X(N/2+1:N) X(1:N/2)]; > > > it just swaps the two halves and using it again will swap the halves > > back. > > Yes. ifftshift and ffthift are always inverses. However, > when N is even they are equivalent. Of course, when N is > even, you have to remember that there are more negative > time points than positive time points.really, the point x[N/2] (crap! this is stupid dumb-ass MATLAB that makes the DFT off by 1, so the point at x(N/2+1)) is negative time? the DFT will treat that one as ambiguous. in frequency, X[N/2] (or in crappy MATLAB, it's X(N/2+1)) is called the "Nyquist point" and it is also equally ambiguous. but the main thing you're missing now, Greg, is that the DFT really only transforms a periodic discrete sequence in the "time" domain to another periodic discrete sequence in the "frequency" domain. those points (using horrible MATLAB indexing convention), X(N/2+2:N) are just as plausible being positive frequency above Nyquist as being negative frequency. same for the time-domain sequence, the points x(N/ 2+2:N) can be considered the "negative-time" portion or just the latter half of the periodic waveform starting at x(1) (which is really x[0] since MATLAB is so stupid). now, Greg, i am not going through your code at all. take a *big* chunk of signal (say you got it out of a wav file, a million samples): x_input(1:1000000) from that extract, by use of windowing, a smaller chunk to send to the FFT. n = 100000; % or whatever big number of your choice. N = 4096; % FFT size x = x_input(n+1:n+N); % extract it x = x .* hanning(size(x)); % window it i don't care what the data is, as long as it's not noise (there has to be some sense of smoothness or predictability of the signal). now look at the results of these two operations: X = fft(x); Y = fft( fftshift(x) ); or maybe we would say... Y = fft( ifftshift(x) ); which should not make a difference for N even. so Greg, can you tell us what the difference is between X and Y? why might it be useful to swap the two halves of x going into the FFT?> > now the MATLAB lit says that fftshift() does not undo itself for odd- > > lengthed arrays and that's why you need ifftshift(). �that's very well > > the case and i don't know what fftshift (or ifftshift) does for odd- > > length arrays anyway. > > ifftshift(x(1:N)) = [ x(1+floor(N/2):N), x(1:floor(N/2))] > fftshift(X(1:N)) �= [ X(1+ceil(N/2):N), X(1:ceil(N/2))]that's useful to know. so, if it's odd length, ifftshift will move the point precisely in the middle to the "t=0" location whereas fftshift does not.> Again, I recommend that MATLAB newbies use ifftshift in > the ifft (time) domain and fftshift in the fft (frequency) > domain. Then most of the worry about whether N is even or > odd is mitigated .i just would not recommend using odd-length FFT for the most part. then there is no "Nyquist point" to swap with the DC point. so, Greg, i am still curious if you understand what use one might derive from swapping the two halves of the windowed input going into the FFT? if you don't get it, that's fine, i'll leave it to your future. r b-j
Reply by ●January 21, 20102010-01-21
On Jan 20, 11:40 am, robert bristow-johnson <r...@audioimagination.com> wrote:> On Jan 20, 4:07 am, Greg Heath <he...@alumni.brown.edu> wrote: > > On Jan 19, 12:47 am, robert bristow-johnson > > > <r...@audioimagination.com> wrote: > > > On Jan 18, 8:12 am, Greg Heath <he...@alumni.brown.edu> wrote: > > > ... > > > > > More simply, > > > > > 1. fftshift shifts and reflects the positive time points to the left > > > WHOOPS! delete "and reflects" > > > > > hand side of the window causing the original middle point at t=0 > > > > to become the right most point. > > > > whereas > > > > 2. ifftshift shifts and reflects the negative time components to the > > > WHOOPS! delete "and reflects" > > > > > right hand side of the window causing the original middle point > > > > at t=0 to become the left most point. > > > > okay, Greg, i went to the MW site to read about ifftshift() > > > >http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ifftshift.html > > > > since i didn't know about it before. i've known and used fftshift() > > > for even-lengthed arrays for 15+ years. i know precisely what fftshift > > > () does to even-length arrays. > > > For even-length arrays ifftshift and fftshift are equivalent. > > The source code is exactly the same except fftshift uses > > ceil(N/2) and ifftshift uses floor(N/2). See below. > > > > okay. the word "reflect" and "reflects" is in error and not in the > > > MATLAB doc. > > > I agree. Thanks for the correction. > > > > for even-lengthed arrays, fftshift() *swaps* the first half data with > > > the second half along every dimension. equivalent to: > > > > Y = fftshift(X); > > > or > > > N = length(X); > > > Y = [X(N/2+1:N) X(1:N/2)]; > > > > it just swaps the two halves and using it again will swap the halves > > > back. > > > Yes. ifftshift and ffthift are always inverses. However, > > when N is even they are equivalent. Of course, when N is > > even, you have to remember that there are more negative > > time points than positive time points. > > really, the point x[N/2] (crap! this is stupid dumb-ass MATLAB that > makes the DFT off by 1, so the point at x(N/2+1)) is negative time? > the DFT will treat that one as ambiguous. in frequency, X[N/2] (or in > crappy MATLAB, it's X(N/2+1)) is called the "Nyquist point" and it is > also equally ambiguous.I agree that MATLAB not using nonpositive indexing results in wierd stuff... Don't shoot the messenger! For example, regardless of how the time sample is circularly shifted, when N is even, MATLAB's fftshift puts the Nyquist point at the leftmost point of the negative frequency interval. However, that is really irrelevant to the main point. Regardless of the indexing convention, when N is odd, MATLAB's fftshift is not self inverting. Therefore it's use should be confined to the frequency domain for centering spectra that are created over a nonnegative frequency interval by the fft. That answers the question posed by the OP who is trying to grasp how to use MATLAB's version of fft. So why is this thread dragging on?> but the main thing you're missing now, Greg, is that the DFT really > only transforms a periodic discrete sequence in the "time" domain to > another periodic discrete sequence in the "frequency" domain.No. The DFT really only transforms a finite discrete sequence in the "time" domain to a periodic discrete sequence in the "frequency" domain. The IDFT transforms a finite discrete sequence in the "frequency" domain to a periodic discrete sequence in the "time" domain. In other words, the representation obtained by the transform and its inverse creates periodicity even when it is not present in the original data.> those > points (using horrible MATLAB indexing convention), X(N/2+2:N) are > just as plausible being positive frequency above Nyquist as being > negative frequency.I have a problem with that statement. It's akin to saying that X(2:N/2) are just as plausible being negative frequency below negative Nyquist as being positive frequency.> same for the time-domain sequence, the points x(N/ > 2+2:N) can be considered the "negative-time" portion or just the > latter half of the periodic waveform starting at x(1) (which is really > x[0] since MATLAB is so stupid). > > now, Greg, i am not going through your code at all.That was not the purpose of the inclusion. It was to yield plots, for those with MATLAB, that visually demonstrate the difference between using ifftshift, fftshift or neither. After all, this thread has dragged on because I recommend using ifftshift, and not fftshift, in the time domain.> take a *big* > chunk of signal (say you got it out of a wav file, a million samples): > > x_input(1:1000000) > > from that extract, by use of windowing, a smaller chunk to send to the > FFT. > > n = 100000; % or whatever big number of your choice. > N = 4096; % FFT size > > x = x_input(n+1:n+N); % extract it > x = x .* hanning(size(x)); % window it > > i don't care what the data is, as long as it's not noise (there has to > be some sense of smoothness or predictability of the signal). now > look at the results of these two operations: > > X = fft(x); > > Y = fft( fftshift(x) ); > > or maybe we would say... > > Y = fft( ifftshift(x) ); > > which should not make a difference for N even. so Greg, can you tell > us what the difference is between X and Y? why might it be useful to > swap the two halves of x going into the FFT?If it is something that my examples didn't cover, just come out and tell me because I don't have a clue to what you are getting at.> > > now the MATLAB lit says that fftshift() does not undo itself for odd- > > > lengthed arrays and that's why you need ifftshift(). that's very well > > > the case and i don't know what fftshift (or ifftshift) does for odd- > > > length arrays anyway. > > > ifftshift(x(1:N)) = [ x(1+floor(N/2):N), x(1:floor(N/2))] > > fftshift(X(1:N)) = [ X(1+ceil(N/2):N), X(1:ceil(N/2))] > > that's useful to know. so, if it's odd length, ifftshift will move > the point precisely in the middle to the "t=0" location whereas > fftshift does not. > > > Again, I recommend that MATLAB newbies use ifftshift in > > the ifft (time) domain and fftshift in the fft (frequency) > > domain. Then most of the worry about whether N is even or > > odd is mitigated . > > i just would not recommend using odd-length FFT for the most part. > then there is no "Nyquist point" to swap with the DC point.In typical situations using a power of two or zero-padding to a power of two is the norm. The deal here is the atypical situation instigated by an inexperienced MATLAB user. Most of the time it is someone who is committed to N odd in order to make t = 0 exactly in the middle of the time waveform and/or someone who doesn't have much data and is reluctant to throw away one measurement point.> so, Greg, i am still curious if you understand what use one might > derive from swapping the two halves of the windowed input going into > the FFT? if you don't get it, that's fine, i'll leave it to your > future.I have only recommended the swap to take advantage of symmetry about the "middle" of the waveform. If there is another reason, you'd better tell me now. I'm too old to count on the future. Greg
Reply by ●January 21, 20102010-01-21
On Jan 21, 8:42�am, Greg Heath <he...@alumni.brown.edu> wrote: ...> I agree that MATLAB not using nonpositive indexing results > in wierd stuff... > > Don't shoot the messenger!not planning to (unless your name is Cleve Moler). this is an old and tired problem with MathWorks spouting their old and tired excuses for not fixing nor even recognizing the problem.> For example, regardless of how the time sample is circularly > shifted, when N is even, MATLAB's fftshift puts the Nyquist > point at the leftmost point of the negative frequency interval.yes. and there is no Nyquist point for N odd.> However, that is really irrelevant to the main point. Regardless > of the indexing convention, when N is odd, MATLAB's fftshift is > not self inverting.well, since i have never done an FFT in MATLAB with odd N, i never got around to wondering what would happen. i believe the main point from the OP on is "why would someone use fftshift() on the *input* to the fft()?" and there is an answer to that.> Therefore it's use should be confined to the > frequency domain for centering spectra that are created over a > nonnegative frequency interval by the fft.ya know. these are numbers, arrays of numbers. the DFT does not know nor care if the input to it is time or frequency domain. the iDFT is just like the DFT except for a scaling factor.> That answers the question posed by the OP who > is trying to grasp how to use MATLAB's version > of fft. > > So why is this thread dragging on? > > > but the main thing you're missing now, Greg, is that the DFT really > > only transforms a periodic discrete sequence in the "time" domain to > > another periodic discrete sequence in the "frequency" domain. > > No. > > The DFT really only transforms a finite discrete sequence in the > "time" domain to a periodic discrete sequence in the "frequency" > domain.that's an old debate we had here at comp.dsp periodically, that you are sure to lose. there is not one mathematical difference between what we call the "Discrete Fourier Transform" and what we call the "Discrete Fourier Series". i might recommend reading what Oppenheim and Schafer have to say about it.> The IDFT transforms a finite discrete sequence in the "frequency" > domain to a periodic discrete sequence in the "time" domain.you *better* put those quotes there because the iDFT doesn't know and doesn't care if it's "time" or "frequency". and there is *no* qualitative difference between the DFT and the iDFT. most definitions put the 1/N factor completely in one or the other (usually in the iDFT) but you can have a legit definition that puts 1/sqrt(N) in both the DFT and iDFT. there is no qualitative difference between -j and +j (they both have equal claim to squaring to be -1 and no other numbers can make that claim). so, if you're gonna build a case that ifftshift() should be used only with one and that fftshift() should be used only with the other, it's already a loser.> In other words, the representation obtained by the transform and > its inverse creates periodicity even when it is not present in the > original data.perhaps a shorter way to say that is that the DFT periodically extends the input presented to it. it "assumes" (to anthropomorphize a bit) that the data passed to it is exactly on period of a periodic sequence. and the DFT output is one period of a periodic sequence of the same period, N.> > > those > > points (using horrible MATLAB indexing convention), X(N/2+2:N) are > > just as plausible being positive frequency above Nyquist as being > > negative frequency. > > I have a problem with that statement. It's akin to saying that > X(2:N/2) are just as plausible being negative frequency below > negative Nyquist as being positive frequency.positive frequency above Nyquist. the data (going in or coming outa the DFT) is periodic. always. it never is not periodic. so (not using horrible MATLAB indexing) the data from X[N/2] to X[N] can just as well be from X[-N/2] tp X[0]. it's interpretation. if you look at x[0] to x[N-1] as being sampled from a real, continuous- time function x(t), then we expect *before* sampling that there was no energy at or above Nyquist (and at and below -Nyquist for negative continuous frequency). that makes it a little confusing for X[N/2] if there is any energy in it. if there is 0 in X[N/2], then we're okay and we can safely say that X[N/2+1] to X[N-1] meaningfully map only to the negative frequencies x[-(N/2-1)] to X[-1]. X[0] remains the DC component and X[1] to X[N/2-1] remain the positive frequency components. but after sampling, the continuous spectrum itself is periodic, that is another reason to understand X[k] as being periodic in k with period N.> > same for the time-domain sequence, the points x(N/2+2:N) > > can be considered the "negative-time" portion or just the > > latter half of the periodic waveform starting at x(1) (which is really > > x[0] since MATLAB is so stupid). > > > now, Greg, i am not going through your code at all. > > That was not the purpose of the inclusion. It was to > yield plots, for those with MATLAB, that visually > demonstrate the difference between using ifftshift, > fftshift or neither. > > After all, this thread has dragged on because I > recommend using ifftshift, and not fftshift, in > the time domain.what's "time domain" what's "frequency domain". what difference, in MATLAB or anywhere else, is there between the DFT and iDFT?> > > > take a *big* > > chunk of signal (say you got it out of a wav file, a million samples): > > > � � x_input(1:1000000) > > > from that extract, by use of windowing, a smaller chunk to send to the > > FFT. > > > � � n = 100000; �% or whatever big number of your choice. > > � � N = 4096; � �% FFT size > > > � � x = x_input(n+1:n+N); � � � % extract it > > � � x = x .* hanning(size(x)); �% window it > > > i don't care what the data is, as long as it's not noise (there has to > > be some sense of smoothness or predictability of the signal). �now > > look at the results of these two operations: > > > � � X = fft(x); > > > � � Y = fft( fftshift(x) ); > > > or maybe we would say... > > > � � Y = fft( ifftshift(x) ); > > > which should not make a difference for N even. �so Greg, can you tell > > us what the difference is between X and Y? �why might it be useful to > > swap the two halves of x going into the FFT? > > If it is something that my examples didn't cover, just > come out and tell me because I don't have a clue to > what you are getting at.i said it before. if you pass the above x as defined by x = x_input(n+1:n+N); % extract it x = x .* hanning(size(x)); % window it directly to the FFT and look at the results, and compare that to passing it to the FFT by way of fftshift(), the latter will have the center of the window positioned around n=0 whereas the former will have the center of the window positioned around n=N/2. that's a delay of N/2. that will effectively multiply every odd sample of X[k] by -1. that will add a phase of +pi or -pi (take your pick) to *only* the X[k] with k odd. a nice smooth phase response won't look so smooth because pulses centered at t=0 have less phasiness than pulses centered at t=N/2.> > > > > > now the MATLAB lit says that fftshift() does not undo itself for odd- > > > > lengthed arrays and that's why you need ifftshift(). �that's very well > > > > the case and i don't know what fftshift (or ifftshift) does for odd- > > > > length arrays anyway. > > > > ifftshift(x(1:N)) = [ x(1+floor(N/2):N), x(1:floor(N/2))] > > > fftshift(X(1:N)) �= [ X(1+ceil(N/2):N), X(1:ceil(N/2))] > > > that's useful to know. �so, if it's odd length, ifftshift will move > > the point precisely in the middle to the "t=0" location whereas > > fftshift does not. > > > > Again, I recommend that MATLAB newbies use ifftshift in > > > the ifft (time) domain and fftshift in the fft (frequency) > > > domain. Then most of the worry about whether N is even or > > > odd is mitigated . > > > i just would not recommend using odd-length FFT for the most part. > > then there is no "Nyquist point" to swap with the DC point. > > In typical situations using a power of two or > zero-padding to a power of two is the norm. > > The deal here is the atypical situation instigated by an > inexperienced MATLAB user. Most of the time it is someone > who is committed to N odd in order to make t = 0 exactly > in the middle of the time waveform and/or someone who > doesn't have much data and is reluctant to throw away > one measurement point. > > > so, Greg, i am still curious if you understand what use one might > > derive from swapping the two halves of the windowed input going into > > the FFT? �if you don't get it, that's fine, i'll leave it to your > > future. > > I have only recommended the swap to take advantage of > symmetry about the "middle" of the waveform.sometimes you want that "middle of the waveform" positioned around t=0.> If there is > another reason, you'd better tell me now. I'm too old > to count on the future.i can relate. r b-j
Reply by ●January 21, 20102010-01-21
robert bristow-johnson wrote: ...> so, if you're gonna build a case that ifftshift() should be used only > with one and that fftshift() should be used only with the other, it's > already a loser.Greg's claim is that for naive users who need one simple rule that always works properly, it's a good way to go. That seems plausible. ... Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
Reply by ●January 21, 20102010-01-21
On Jan 21, 3:48�pm, Jerry Avins <j...@ieee.org> wrote:> robert bristow-johnson wrote: > > � �... > > > so, if you're gonna build a case that ifftshift() should be used only > > with one and that fftshift() should be used only with the other, it's > > already a loser. > > Greg's claim is that for naive users who need one simple rule that > always works properly, it's a good way to go. That seems plausible.i just dunno if the difference between fftshift and ifftshift has anything to do with the difference between fft and ifft (of which the latter has no qualitative difference other than scaling). the issue is, for N odd, do we want the DC point X[0] to go to the center, X [(N-1)/2] or to the point right after it? the consequence is that it's the point right after the center, X[(N+1)/2] that gets mapped to DC, X[0]. in the "time domain", i dunno whether it's better, for mathematical reasons, to map x[(N+1)/2] to x[0] or to map X[(N-1)/2] to x[0]. i think maybe it's the latter. so i think it's not better to use ifftshift to fix the phasey issue than fftshift. but i am not entirely sure. hey Jerry, i am coming down to NYC this Sunday to put up posters for an event that i am promoting (FYI it's http://eve-fest.org ). i am staying Sunday overnight in NYC, but want to check up on our house in NJ on Monday. would you like me to pay you a visit? L8r, r b-j
Reply by ●January 21, 20102010-01-21
robert bristow-johnson wrote:> On Jan 21, 3:48 pm, Jerry Avins <j...@ieee.org> wrote: >> robert bristow-johnson wrote: >> >> ... >> >>> so, if you're gonna build a case that ifftshift() should be used only >>> with one and that fftshift() should be used only with the other, it's >>> already a loser. >> Greg's claim is that for naive users who need one simple rule that >> always works properly, it's a good way to go. That seems plausible. > > i just dunno if the difference between fftshift and ifftshift has > anything to do with the difference between fft and ifft (of which the > latter has no qualitative difference other than scaling). the issue > is, for N odd, do we want the DC point X[0] to go to the center, X > [(N-1)/2] or to the point right after it? the consequence is that > it's the point right after the center, X[(N+1)/2] that gets mapped to > DC, X[0]. > > in the "time domain", i dunno whether it's better, for mathematical > reasons, to map x[(N+1)/2] to x[0] or to map X[(N-1)/2] to x[0]. i > think maybe it's the latter. so i think it's not better to use > ifftshift to fix the phasey issue than fftshift. but i am not > entirely sure. > > hey Jerry, i am coming down to NYC this Sunday to put up posters for > an event that i am promoting (FYI it's http://eve-fest.org ). i am > staying Sunday overnight in NYC, but want to check up on our house in > NJ on Monday. would you like me to pay you a visit?Sure. There's a bed for you. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
Reply by ●January 23, 20102010-01-23
On Jan 21, 10:26 pm, robert bristow-johnson <r...@audioimagination.com> wrote:> On Jan 21, 3:48 pm, Jerry Avins <j...@ieee.org> wrote: > > robert bristow-johnson wrote: > > > so, if you're gonna build a case that ifftshift() should be used only > > > with one and that fftshift() should be used only with the other, it's > > > already a loser. > > > Greg's claim is that for naive users who need one simple rule that > > always works properly, it's a good way to go. That seems plausible. > > i just dunno if the difference between fftshift and ifftshift has > anything to do with the difference between fft and ifft (of which the > latter has no qualitative difference other than scaling).An important difference is that f = 0 has physical significance whereas t = 0 does not. No matter how x(t) is circularly shifted, X = fft(x) always yields X(1) = X(f=0).> the issue is, for N odd, do we want the DC point X[0] > to go to the center, X[(N-1)/2] or to the point right after it?In MATLAB, the DC point is X(1) and the center is X((N+1)/2) fftshift([1 2 3 4 5]) = [ 4 5 1 2 3 ]> the consequence is that it's the point right after the center, >X[(N+1)/2] that gets mapped to DC, X[0].Yes. fftshift([ 4 5 1 2 3 ]) = [ 2 3 4 5 1 ] However, ifftshift([ 4 5 1 2 3 ]) = [ 1 2 3 4 5 ] That is because the MATLAB functions fftshift and ifftshift are not self inverting when N is odd. I'm not familiar with the source code of fftshift for other languages. It is possible that fftshift could be coded so that it is always self inverting. If so, ... Your issue is with the MATLAB convention; not, given the MATLAB convention and noninformative documentation, what is the best way to use it.> in the "time domain", i dunno whether it's better, for mathematical > reasons, to map x[(N+1)/2] to x[0] or to map X[(N-1)/2] to x[0]. i > think maybe it's the latter. so i think it's not better to use > ifftshift to fix the phasey issue than fftshift. but i am not > entirely sure.The important point is, for symmetric time functions, to map the symmetry point to x(1). fftshift([ 2 1 0 1 2 ]) = [1 2 2 1 0] wheras fftshift([ 2 1 0 1 2 ]) = [0 1 2 2 1] Since it is desirable that the same function be used when N is even or odd and large or small, the choice is pretty much limited to the approach I have recommended. Of course, there is always the option of writing alternative code and naming it fftshiftrjb.m Hope this helps. Greg.
Reply by ●January 23, 20102010-01-23
On Jan 23, 4:55�am, Greg Heath <he...@alumni.brown.edu> wrote:> On Jan 21, 10:26 pm, robert bristow-johnson > > <r...@audioimagination.com> wrote: > > On Jan 21, 3:48 pm, Jerry Avins <j...@ieee.org> wrote: > > > robert bristow-johnson wrote: > > > > so, if you're gonna build a case that ifftshift() should be used only > > > > with one and that fftshift() should be used only with the other, it's > > > > already a loser. > > > > Greg's claim is that for naive users who need one simple rule that > > > always works properly, it's a good way to go. That seems plausible. > > > i just dunno if the difference between fftshift and ifftshift has > > anything to do with the difference between fft and ifft (of which the > > latter has no qualitative difference other than scaling). > > An important difference is that f = 0 has physical significance > whereas t = 0 does not.that's making an assumption of the physical meaning of x[] or X[]. i consider the DFT and iDFT to be qualitatively equivalant.> No matter how x(t) is circularly shifted, > X = fft(x) always yields X(1) = X(f=0).you can say the same about the inverse. no matter how X(f) is circularly shifted, x = ifft(X) always yields x(1) = x(t=0).> I'm not familiar with the source code of fftshift for > other languages.i've never seen a library with a function named "fftshift" other than MATLAB.> > in the "time domain", i dunno whether it's better, for mathematical > > reasons, to map x[(N+1)/2] to x[0] or to map X[(N-1)/2] to x[0]. �i > > think maybe it's the latter. �so i think it's not better to use > > ifftshift to fix the phasey issue than fftshift. �but i am not > > entirely sure. > > The important point is, for symmetric time functions, to map the > symmetry point to x(1). > > fftshift([ 2 1 0 1 2 ]) = �[1 � � 2 � � 2 � � 1 � � 0] > > wheras > > fftshift([ 2 1 0 1 2 ]) = �[0 � � 1 � � 2 � � 2 � � 1] >there's a missing "i", so let's see what you meant. quoting from before: ifftshift(x(1:N)) = [ x(1+floor(N/2):N), x(1:floor(N/2)) ] fftshift(X(1:N)) = [ X(1+ceil(N/2):N), X(1:ceil(N/2)) ] i think it's fftshift([ 2 1 0 1 2 ]) = [1 2 2 1 0] wheras ifftshift([ 2 1 0 1 2 ]) = [0 1 2 2 1] is that not correct?> Since it is desirable that the same function be used > when N is even or odd and large or small, the > choice is pretty much limited to the approach I > have recommended.for the purpose of mapping the point of symmetry to x(1) (what the rest of the world calls x[0]), then you should use ifftshift for odd N (or even N). so does this make sense?: y_output = zeros(size(x_input)); N = length(x); framehop = floor(N/2); n = 0; while (n < length(x_input)-N) x = x_input(n+1:n+N); x = x .* hanning(size(x)); X = fft(ifftshift(x)); % % do your processing of X here. % phase of X is nicely ironed out. % Y = X; % null processing % % y = fftshift(ifft(Y)); y_output(n+1:n+N) = y_output(n+1:n+N) + y; n = n + framehop; end; would you say that's correct, Greg, or should the first one be fftshift and the second be ifftshift? and do you see why one might want to use one or the other on the input to the fft? without considering if further, i am not so sure how the overlap- adding will work for odd N. but it should work transparently (except for the two ends of y_output) for even N. maybe it should be framehop = ceil(N/2), but that doesn't seem right either. r b-j
