Lecture notes, lecture 4 PDF

Preview

Summary

Download Lecture notes, lecture 4 PDF

Description

DEPARTMENT OF ELECTRICAL AND COMPUTER ENGINEERING COLLEGE OF ENGINEERING AND APPLIED SCIENCES

Generalized Analysis: Mathematical derivation and implementation notes for a generalized multirate digital filterbank channelizer.

Draft 27 March 2012 Modified October 17 2013

Dr. Bradley J. Bazuin Associate Professor

Technical Report: TN03272012BJB

Abstract The following provides time sample based mix-filter-decimation derivation with MATLAB vector considerations. In additions the difference between an inverse-transform analysis and forward transform analysis are derived and the “modulo-time” or constant phase difference shown between the two derivations. Both FFT/IFFT and OFFT/IOFFT algorithms are presented.

©2013 Bradley J. Bazuin Disclaimer This is a technical report generated by the author as a record of personal research interest and activity. Western Michigan University makes no representation that the material contained in this report is correct, error free or complete in any or all respects. Thus, Western Michigan University, it’s faculty, it’s administration, or students make no recommendation for the use of said material and take no responsibility for such usage. Therefore, persons or organizations who choose to use this material do so at their own risk.

Analysis Derivation 1: using the IDFT The analysis processing operations is formed as shown in Fig. 1. w n  x n hn M

y m 

W Knk Figure 1: Generalized complex Mixing filter decimation The algorithm to process this structure can be defined for a particular frequency, k, as y k  m 



 hl   xmM  l   W 

mM l k K

(1)

l  

where



WKnk  exp  i2  nk K



(2)

For this equation, define a causal filter of length K  where

0  hn  h n  0 

n 0 0  n  K  K   n

(3)

Describe the structure of the filter in terms of K rows of  columns for, n → (rK+ ). Here, r goes from 0 to  and  goes from 0 to K-1. h K   h 0   h 1 h K  1 h rK           h K  1 h2 K  1 The equation becomes

h  1K    h   1K  1      h K  1   

(4)

 1 K 1

mM rK  k y k m   hrN     x mM  rK     W K   

(5)

r  0  0

The input signal is causal, the past history and current sample is known, but the future samples are unknown. Describe the structure of the incoming signal in terms of K rows of  columns to appropriately align with the filter coefficients.

Page 1 of 13

©2013 Bradley J. Bazuin

xmM  K  0   xmM    1K  0   x mM  0   x mM  1 x mM  K  1  x mM    1K  1       x mM rK          x mM   K  1   xmM  K  1 x mM  2 K  1  This describes the appropriate filtering relationship. Rearranging the equation, K 1

 1

 0

r 0

y k m  W KmMk   W Kk  hrK     x mM  rK     W K rKk

(6)

(7)

The generalized form of the equation results after adjustments for the exponential factors as K 1

 1

 0

r 0

y k m   WKmk   W Kk  h rK    x mM  rK    M

(8)

The processing elements then consist of a polyphase filter of dimensions K x λ, an inverse discrete Fourier transform, and a post multiplication related to the desired frequency of the output signal for k=0 to K-1. As will be discussed later, the post multiplication can be absorbed into the IDFT by circularly shifting the polyphase filter outputs prior to performing the IDFT.

SpecificCases For Nyquist sampling considerations, the number of frequencies K and decimation rate M are related. This can be described in a ratio as a  KM In interpretation, the frequencies for the demodulation are considered frequency bins based on the input sample rate of the data. If the data is complex, the sample rate and frequency content may be equivalent based on complex Nyquist sampling. If the original frequency band is evenly divided into K segments, the decimation rate and filter bandwidths must be appropriately defined to avoid aliasing in the resulting frequency bins. For example, if M=K, the bin width and decimated signal bandwidth are the same, whether the filter allows aliasing or not. If the bins frequency content (passband bandwidth) is to be equal to the bin width, then the output sample rate would want to allow a bandwidth greater than the bin width. As an example 2M=K would provide a complex sampling rate of fs/M = 2 (fs/K) or twice the fs/K bin spacing. This would facilitate an fs/K bandwidth filter passband, transition bands and stop bands within the 2 fs/K output signal frequency band. This concept readily supports the processing of closely spaced RF channels such as radio or television.

If the K/M substitution is performed, the result reduces to K 1

 1

 0

r 0

mk k y k  m  Wa   WK  h rK    x mM  rK  

Page 2 of 13

(9)

©2013 Bradley J. Bazuin

As a note, maximal rate decimation for an input signal would define a=1. The resulting equation would be K 1

 1

 0

r0

y k m    W K k h rK    x mK  rK   

(10)

which becomes a more familiar compact notation form of polyphase filters K 1

 1

0

r 0

y k  m   W K k  h r   x   m  r 

(11)

If x(n) is real, you would compute frequency channels for k=0 to K/2. For x(n) complex, all frequency channels may be of interest. As a concern, the decimated frequency channels are very narrow. If each channel is to contain frequency content across the entire decimated frequency bin, it may be desirable to select a=2 or greater. For a=2, the decimated frequency band would contain frequency elements from the band center to the center of the next adjacent frequency band as compared to the edge between the frequency bands (+/- fdecimate instead of +/- fdecimate/2). This would allow the applied filter to be 3dB at +/- fdecimate/2 , have a reasonable transition band and an excellent stopband attenuation. If signal reconstruction is required, this is a desirable characteristic! For a=2, the derived equation becomes K 1  1 K  yk  m  W2mk   WK k  h rK     x m  rK     2  r 0  0 which can be simplified in terms if the complex multiplication to K 1 1  K  y k m    1mk   W Kk  h rK     x  m  rK     2  r 0 0

Page 3 of 13

(12)

(13)

©2013 Bradley J. Bazuin

IDFT Implementation for arbitrary K and M Repeating the developed equations: xmM  K  0   x mM  0   x mM  1 x mM  K  1 xmM  rK           x mM  K  1 x mM  2 K  1 K 1

 1

 0

r 0

 xmM    1K  0   x mM    1K  1      x mM   K  1  

mMk k y k m  W K   WK  h rK     x mM  rK   

(6)

(8)

The first interest is in the polyphase elements, which can be defined as  1

PP  m   h rK     x mM  rK    ,

for k  0 : K  1

(41)

r 0

Based on the definition of the data matrix shown in Equ. 6, when K and M are not equal the matrix must be reformed for each computation. In reforming the matrix, the oldest M samples are removed and then newest M samples are added. This removal and addition can be seen as a “raster shift” of the data, where data shift off the bottom of columns and appear at the top of the next column. The oldest data is shifted off the bottom of the last column, while new data is shifted into the start of the first column. The process described can be easily performed using a software based circular buffer, where a starting point for data is identified by a pointer and offsets of K from the pointer are identified by “modulo- K  ” addition (add K to the previous pointer, if the sum is greater than or equal to K  than subtract K ). The result forms a K length vector that must be further processed. For a direct implementation following Equ. 7,

y k  m  W mMk  W K k  PP  m K K 1

(42)

 0

Here, the inverse discrete Fourier transform of the polyphase elements is computed followed by a phase rotation based on the frequency, k, and decimated time sample, m, of the desired output.

Page 4 of 13

©2013 Bradley J. Bazuin

PhaseShiftusingaCircularShift Studying the Fourier matrix form of the transformation,

WK 0  WK 0 WK 0 WK 0    0 K 1 1 2  WK WK WK    WK 2 4 2 K 1 WKnk  WK0 (43)  WK     WK WK          W 0 W  K 1  W 2  K 1   W  K 1  K 1  K K K   K each column provide a different “k-based” multiplication or shifting of the previous column. Instead of following this form, the twiddle factor may be combined with the complex shift as K 1





(44)

 W PP

m 

(45)

y k m   WK  mM  k  PP  m  0

y k  m 

K 1  mM

rk K

r mM

r  mM

Using the form of Equ. 45, it can be seen that the new summation is a –mM circularly shifted version of the previous Fourier matrix. As the order of summation in a Fourier transform is arbitrary and the “twiddle factors” form a circular symmetry (-mM:-1) is the same as K-mM:K1), this equation could also be solved as K 1





yk m   WK rk  PPr mM K m 

(46)

r 0

This implies that the data is circular shifted by mM prior to performing the transformation! If this is done, than no “post-multiplication” is required. A block diagram of the implementation sequence is shown in Fig. 2.

Figure 2: Analysis filter processing sequence for M and K integers. Note that his analysis was performed for the FFT Analysis (1) which results in the use of an inverse-discrete-Fourier transform. For the other derivations, a similar analysis can and should be performed in order to get the buffering and circular shifting correct. Another way to visualize the processing structure is shown in Fig. 3.

Page 5 of 13

©2013 Bradley J. Bazuin

xn ym 

Figure 3: Visualization of the Analysis Processing. With this structure, the commutated input of an M sample vector is placed in a K by  matrix after appropriately shifting all previous M vectors forward in the matrix. The polyphase row computations are performed and output to the modulo mM circular shifter. After shifting the IDFT is applied to produce K distinct y k m  outputs for further processing.

Page 6 of 13

©2013 Bradley J. Bazuin

Analysis Derivation 2: using the DFT There is an alternate ordering of filter and data that can be made. In this ordering, the filter goes from 1 to K  and the data is collected as a sequence from 0 to K-1, with K-1 representing the most recent time sample. For a filter of length K  where

n 0 0  n  K K  n

0  hn  hn 0 

(21)

Describe the structure of the filter in terms of K rows of  columns, n → (rK- ). The r goes from 1 to  and  goes from 0 to K-1. h K  h 2K     h rK       h 2  h K  2    h 1  h K  1  The equation becomes 

h K        h   1K  2    h   1 K  1 

(22)

K

z k m    hrK     xmM  rK     W kmM rk  k

(23)

r 1  1

The input signal is causal, the past history and current sample is known, but the future samples are unknown. Describe the structure of the incoming signal in terms of K rows of  columns.  xmM  K  1  x mM  2 K  1  x mM   K  1   xmM K 2  xmM 2 K 2   xmM K 2          x mM  rK               xmM     1K  xmM  xmM  K   This describes the appropriate filtering relationship. Rearranging the equation, K



 1

r 1

z k m   WKmMk   WKk h rK     x mM  rK     WK rKk

(24)

(25)

Which reduces to K



1

r1

zk  m  W Kmk   W Kk  h rK     x mM  rK    M

(26)

For critically sampled analysis, a=1, we have K



 1

r 1

z k m    W Kk  hrK     xmK  rK   

Page 7 of 13

(27)

©2013 Bradley J. Bazuin

While for 2x oversampling of a=2, 

K

K z k m    1mk   WKk  hrK    x m  rK    (28)  2   1 r 1 Of note is that the summation in  goes from 1 to K. To form a “proper” FFT or DFT, Equ. 27 should be changed to K



 1

r 1

z k m  W Kk   W K 1 k  h rK    x mK  rK    Note that the vector and matrix data ordering remains as shown, but when compared to the previous derivation the index on the filter matrix has changed and the values have flipped from top to bottom. In addition, the data matrix has flipped top to bottom.

Comparingderivations1IDFTand2DFT The first derivation uses data numbered from the current sample backward in time and requires an IFFT or IDFT to produce a full channelizer. In contrast, the second derivation reverses the time order of the data and uses a forward transform. Also of note, the filter and data from one derivation to the next are “flipped” in the matrix forms from top to bottom. While either derivation will produce results, it is important when performing a “re-synthesis” of the analyzed channels to perform the matching “inverse” operation (IDFT follows by DFT or DFT followed by IDFT). Comparing the IDFT and the DFT of the flipped column. K 1

X  k   x K  1  n  WKnk

(29)

n 0

Performing a variable substitution, m=K-1-n for m=0 to K-1 K 1

X k    x m   WKK

1 m k

(30)

m0

Isolating the exponential terms and cancelling where possible K 1

X  k  WKk   x m  WKmk

(31)

m0

By flipping the sequence, there is a residual phase element or complex constant that multiplies each of the frequency terms. When performed across multiple analysis outputs, the term for each frequency becomes a complex constant. As a result, the two formulations will not produce identical output results as “time waveforms”, but each is mathematically valid with a fixed multiplicative phase term when they are compared. As a final note, the “magnitude” of the outputs should be identical.

Page 8 of 13

©2013 Bradley J. Bazuin

Analysis Derivation using the Odd-Frequency Fourier Transform, Derivation 1 IOFFT An OFFT analysis can be performed using an odd frequency mixing 

y k m  

 h l   x mM

mM  l k 1 2 

 l   WK

(14)

l  

Applying similar processes to the FFT, the intermediate value becomes mM rK   k  12 

1 K 1

y k m    h rN     x mM  rK     W K

(15)

r 0  0

Expanding the exponential and allocating to appropriate summations,



y k  m  W K

mM k 1 2



K 1



  k 12

  WK  0



 1

  hrN    xmM  rK    WK



 rK k  1 2



(16)

r 0

This becomes









 1 r  K    hrN    x m  rK      1  0 r 0  a  The sign change based in r can be applied to the filter coefficients as

y k  m  Wa

m k 1

2





K 1

 k 1 2

 WK





1





 1





(18)





(19)

  K     1r  hrN     x m  rK     0 r 0  a  For critically sampled analysis, a=1, we have

yk m  W a

m k 1

2

K1

  k  12

  WK

K 1

y k m   1 m   WK

  k  12

0

(17)

   1r  h rN     x mK  rK   r 0

For a=2,

y k m    i  m

 2 k1 

K 1

  WK



 k 1 2

 0







 1  K    1r h rN     x m  rK    r 0  2 

(20a)

or









 1  K     1r  h rN     x  m  rK    (20b) r 0  2   0 Special note: multiplication by a fixed phase is indicative of a time-delayed signal when taking the Fourier transform. The equivalent for a IDFT or IFFT is to perform a circular shift on the data prior to performing the IDFT or IFFT. K 1

y k m    i  m   1 mk   WK 

  k  12

Page 9 of 13

©2013 Bradley J. Bazuin

Analysis Derivation using the Odd-Frequency Fourier Transform, Derivation 2 OFFT An OFFT analysis can be performed using an odd frequency mixing 

y k m  

 h l   x mM

mM  l k 1 2 

 l   WK

(31)

l  

Applying similar processes to the FFT, the intermediate value becomes mM rK  k  12 

 K 1

y k m   h rN     x mM  rK     W K

(32)

r 1  0

Expanding the exponential and allocating to appropriate summations,



y k m   WK

mM k 1 2



K 1





 k 1 2

 WK





h rN    x mM  rK    W



 0





rK k 1 2 K

(33)

r 1

This becomes









 r  K    h rN    x m  rK      1   0 r 1  a  The sign change based in r can be applied to the filter coefficients as

y k  m  W a

m k1



2

K 1

 k 1

 WK



K 1

y k m    1m 

K 1



2







  K     1r  h rN     x m  rK     0 r 1  a  For critically sampled analysis, a=1, we have

y k m   W a

m k1 2

...