Whitening in Frequency Domain

On 3/9/13 11:22 AM, radams2000@gmail.com wrote:
>
> But the question that the OP wants to answer is, lets assume that through some process you end up modifying each fft frame in a way that would cause aliasing in the reconstruction, for example by multiplying bins by 1 or 0. What's the best way to fix up the spectrum such that, after inverse transform and overlap-add,  the aliases are completely eliminated?

whatever it is that you're multiplying the FFT frame with (let's call 
that a "frequency response"), that's convolution in the time domain. 
but there is no guarantee that the impulse response that corresponds 
with that frequency response is short enough to prevent time-aliasing.

so you *force* that impulse response to be short enough so that you get 
*some* number of samples that are not affected by time aliasing.  you 
impose shortening on that impulse response by multiplying some of the 
samples by zero.  that is a form of a window.  then you may as well make 
the window a "nice" window.

multiplying the original impulse response by any window is the same as 
convolving the frequency response with *something* (the DFT of the 
window).  if it's a "nice" window, that convolving can legitimately be 
called "smoothing".  applying the window in the time domain to the 
impulse response is the same as convolution with the DFT of the window 
to the frequency response.  the longer the window, the tighter the 
convolution is.  the shorter the window, the more smearing in the 
smoothing occurs in the frequency domain.

so i dunno how to fix it after the whole overlap-add is done (if you 
have glitches from time-aliasing, i think you're already screwed), but i 
*think* i know how to fix it *inside* the overlap-add method.

> I think Dale's answer is to smooth the spectrum directly in the frequency domain by convolving with a simple 3-term weighting vector.  Can this eliminate ALL the aliasing?

only if it corresponds to an impulse response that is short enough to 
prevent time aliasing for the output samples that you keep.  i don't see 
that the simple 3-term vector would correspond to *any* particular 
window, whereas we need a window that shortens the impulse response to 
be guaranteed to avoid time-aliasing for at least *some* output samples.

>
> I don't know.
>
> Bob


-- 

r b-j                  rbj@audioimagination.com

"Imagination is more important than knowledge."

On Saturday, March 9, 2013 5:06:59 AM UTC-8, radam...@gmail.com wrote:
> ...
> 
> Actually I have a question about the freq-domain convolution method; how to you handle the end points?
> 
...
> Bob

Real world signals have had some kind of anti-aliasing filter applied. Usually a low pass filter. That leaves a region at Fs/2 where there is no signal and the convolution doesn't need to be calculated because there is no valid information to be smoothed or whitened.

Dale B. Dalrymple

On 3/9/13 1:42 PM, dbd wrote:
> On Saturday, March 9, 2013 5:06:59 AM UTC-8, radam...@gmail.com wrote:
>> ...
>>
>> Actually I have a question about the freq-domain convolution method; how to you handle the end points?
>>

i was thinking about time-domain endpoints.  looks like i missed 
something in Bob's question.


>
> Real world signals have had some kind of anti-aliasing filter applied. Usually a low pass filter. That leaves a region at Fs/2 where there is no signal and the convolution doesn't need to be calculated because there is no valid information to be smoothed or whitened.

but about the endpoints in the freq domain, i said this before (without 
referring to "endpoints"): "(... we expect reflection about DC anyway), 
but if the convolution is short (or "small") then smoothing everywhere 
other than around Nyquist should behave exactly as we would expect."

now Dale, rather than asking Bob, i'll ask you:

how is convolution in the frequency domain (for the purpose of 
smoothing) effectively different than multiplication in the time domain. 
  and how is multiplication in the time domain different than windowing? 
  and how are you going to guarantee the avoidance of time-aliasing 
without making certain that the impulse response is short enough?

just curious.

-- 

r b-j                  rbj@audioimagination.com

"Imagination is more important than knowledge."

Robert, Dale

This all gets to the heart of my question.

If in order to completely eliminate aliasing, the frequency-domain vector that I need to use to convolve with the fft frame is very long, then I might as well go back to the time domain, window so that I get zeros in all the right places, and then go back to the frequency domain to get the properly smoothed spectrum. If on the other hand a 3-term freq-domain vector can have an inverse fft that is half zeros then there's no point in doing it in the time domain. 

But in thinking about this I am getting even more confused because it seems like there are some redundant steps.

1) original signal, window, take fft
2) modify fft with non-smoothed frequency response
3) inv fft to time domain
4) window to force 0's for N/2 samples
5) fft to freq domain. Now I have smoothed response.
6) inv fft and overlap-add to to time domain

But don't steps 5 and 6 just cancel each other out? Couldn't I just stay in the time domain after applying my window and use this in the overlap-add?


Bob

On 3/9/13 3:33 PM, radams2000@gmail.com wrote:
> Robert, Dale
>
> This all gets to the heart of my question.
>
> If in order to completely eliminate aliasing, the frequency-domain vector that I need to use to convolve with the fft frame is very long, then I might as well go back to the time domain, window so that I get zeros in all the right places, and then go back to the frequency domain to get the properly smoothed spectrum. If on the other hand a 3-term freq-domain vector can have an inverse fft that is half zeros then there's no point in doing it in the time domain.
>
> But in thinking about this I am getting even more confused because it seems like there are some redundant steps.
>
> 1) original signal, window, take fft

for regular overlap-add, the only windowing here would be a rect().  but 
i think it can be a complementary window (like the Hann), if you cut the 
hop length by 1/2.

> 2) modify fft with non-smoothed frequency response

i didn't realize that.  i thought it was different.

> 3) inv fft to time domain
> 4) window to force 0's for N/2 samples
> 5) fft to freq domain. Now I have smoothed response.
> 6) inv fft and overlap-add to to time domain
>
> But don't steps 5 and 6 just cancel each other out? Couldn't I just stay in the time domain after applying my window and use this in the overlap-add?
>


i thought it might be

0) FFT length = N.  decide on a fixed FIR length, L, so that the 
overlap-add will be the same each frame.  the frame hop, N-L+1, will 
remain constant.

1) input one frame of N-L+1 samples.  zero-pad (rect() window). FFT.

2) compute    G[k] =  ( (Xref)^2 / |X[k]|^2 )^(alpha/2)

    alpha is the compression amount.  alpha=0 has no effect alpha=1
    means that it fully on and alpha = 1/2 is what the OP suggests.
    G[k] >= 0 has even symmetry:  G[N-k] = G[k]

3) smooth G[k] by inv FFT and windowing g[n] with even symmetry window

          h[n]  =  g[n] w[n]     where w[N-n] = w[n]

          w[n] = 0  for all |n - N/2| <  (N-L+1)/2

4) FFT h[n] to get H[k],   again even symmetry and real.

5) multiply  Y[k] = X[k] H[k]

6) inv FFT to get y[n]

7) overlap first L-1 samples of current frame y[n] with last L-1 samples 
of previous frame and add.

8) advance input and output pointers by N-L+1 samples, go to 1)

now they have a word for this sorta thing: "multiband compressor". 
okay, 2 words.  lotta bands (unless h[n] is very short).

if the impulse response h[n] is long, the "smoothing" is kinda tight 
(not a lot of smoothing between bins in G[k]).  if the impulse response 
is short (L small), then the smoothing between bins is kinda sloppy and 
there are effectively fewer bands.

i think the window could be a Kaiser and that would limit the reach of 
the smoothing for a fixed length L (of h[n]).  (bins of G[k] far away 
would be weighted by sidelobe gains of the F.T. of the Kaiser which are 
low in magnitude.)  but it could be whatever window with even symmetry 
you like.

perhaps H[k] for the current frame be "smoothed" a little with the H[k] 
(same k) from the previous frame (low pass filtering, in time, on each 
H[k] over the frames), if you didn't want to hear this multiband 
compressor pump real fast.

the cool thing about overlap-add w.r.t. overlap-save is that with OLA, 
over the L-1 samples that overlap, there is a sorta linear crossfade 
from the filter of the previous frame to the filter of the current 
frame.  with overlap-scrap, it's more like a butt splice (razor blade 
that is right angle with the recording tape).  so use OLA.

if you want a non-rectangular window on the input, i think you can do it 
if it's complementary like a Hann, but there are some overlap-adding 
details that i am not sure about.  i think maybe the number of samples 
in the frame hop which was N-L+1 would be cut in half to (N-L+1)/2 which 
is doable because i think that L is an odd number so that w[n] can be 
symmetrical about w[0].


-- 

r b-j                  rbj@audioimagination.com

"Imagination is more important than knowledge."

On 3/9/13 10:16 PM, robert bristow-johnson wrote:
...
> 3) smooth G[k] by inv FFT and windowing g[n] with even symmetry window
>
>           h[n]  =  g[n] w[n]     where w[N-n] = w[n]
>
>           w[n] = 0  for all |n - N/2| <  (N-L+1)/2

i think the h[n] should be delayed to make it "causal" so that the 
overlapping regions of y[n] of the current frame and previous frame line 
up for overlap-add.  so i think it should be

       h[n + (L-1)/2]  =  g[n] w[n]      for  0 <= n < N - (L-1)/2

       h[n + (L-1)/2]  =  g[N+n] w[N+n]  for -(L-1)/2 <= n < 0




-- 

r b-j                  rbj@audioimagination.com

"Imagination is more important than knowledge."

On Saturday, March 9, 2013 10:53:48 AM UTC-8, robert bristow-johnson wrote:
> ...
> now Dale, rather than asking Bob, i'll ask you:

As is traditional on comp.dsp, neither the OP or Bob have explained what all things their envisioned complete system would do in which domain. I have no idea what Bob's system is or how it compares to what the OP (thinks) his system may be (now).

> how is convolution in the frequency domain (for the purpose of 
> smoothing) effectively different than multiplication in the time domain. 
>   and how is multiplication in the time domain different than windowing? 

In theory they are the same, but Bob suggested an fft, applying a window and an ifft. I didn't know why he chose that, but given that he did, for the windows I stated, the frequency domain convolution would cost less than the two transforms. Does Bob still want to do that? Ask him. I'm waiting (at least) until the OP provides some information before engaging in speculative system design.

>   and how are you going to guarantee the avoidance of time-aliasing 
> without making certain that the impulse response is short enough?

So far, the OP has given no system details adequate to determine what effect time domain aliasing would have on the intended purpose or operation of his system. How much and how long matters is yet to be determined. I've had customers specify frequency domain operations, but they also knew how to specify aliasing requirements for their application. (Maybe because they knew what the function of their system was.)

Radar and sonar systems often perform frequency domain whitening, interference rejection and continue to perform the detection functions in the frequency domain. In that case there is no time aliasing generated and smoothers like the boxcar, split boxcar and median can be used to save calculations and/or provide robust interference rejection. Those tasks can be much easier in the frequency domain. (No Tim this isn't an attempt to use the fft for magic, it's a matter of choosing the domain that makes signals, noises and interferences of interest easy to observe and manipulate for the system's intended functions.)

Dale B. Dalrymple

Robert

Thanks. Does this also work if you window the input, assuming the window is properly constrained?

Dale; my field is audio so the tolerance for aliasing is low (at least in the case of music reproduction) . I agree there are plenty of situations where you just need to reduce the aliasing and where a short convolution to smmoth the spectrum would do the trick. 

Bob

On 3/10/13 7:29 AM, radams2000@gmail.com wrote:
> Robert
>
> Thanks. Does this also work if you window the input, assuming the window is properly constrained?
>

okay, i have to say the obvious.  i hadn't ever done this "whitening" 
thing.  but on multiple occasions i have done frequency-domain 
processing, both the basic fast convolution (where the input window is 
rectangular for overlap-add and nonexistent for overlap-scrap) and for 
phase vocoder.  other than "analysis" (not "processing", no audio out), 
i don't do much frequency domain processing anymore.

that said, i am more confident of what i speculated last night:

(A)  if you want a non-rectangular window on the input, i think you can 
do it if it's complementary like a Hann.  i think the number of samples 
in the frame hop which was N-L+1 (for the rect) would be cut in half to 
(N-L+1)/2 which is doable because L is an odd number so that w[n] can be 
symmetrical about w[0].  the overlap would be greater than L-1 samples. 
  i think it would be  L-1 + (N-L+1)/2 = (N+L-1)/2  samples of overlap 
adding.

so you get more smoothness between frames because you're doing twice as 
many frames, for the same size of FFT N.  the non-zero width of your 
complementary window would be the same length as the rect which is
N-L+1 so that there are still L-1 zeros padded (to prevent time-aliasing 
in the circular convolution).  but the hop length is cut in half to 
(N-L+1)/2 so the last half of the N-L+1 samples in one frame are reused 
as the first half of the following frame.  the window should be 
complementary so that

     v[n]  +  v[n + (N-L+1)/2]  =  1    for 0 <= n < (N-L+1)/2

this window, v[n], is not the same window w[n] that is mentioned in my 
previous post.  that w[n] window need not be complementary, just even.


(B)  the other thing that occurs to me now, is that your original input 
x[n] can be used twice, one is rectangularly windowed for the 
overlap-add fast convolution previously referred and then you can window 
the same data and FFT it for the purposes of analysis, to get that G[k]. 
  it might smooth out G[k] a little to begin with.  that's just an idea 
offa the top of my head.


so are you experimenting with this, Bob?


-- 

r b-j                  rbj@audioimagination.com

"Imagination is more important than knowledge."

Robert

Not really, I just like to educate myself on the basics of frequency-domain manipulation without aliasing, and I've always been slightly confused on the topic so this thread has been helpful. 

Bob