Generalized pencil-of-function method

Extraction of two sinusoids from a noisy data through the GPOF method

Generalized pencil-of-function method (GPOF), also known as matrix pencil method, is a signal processing technique for estimating a signal or extracting information with complex exponentials. Being similar to Prony and original pencil-of-function methods, it is generally preferred to those for its robustness and computational efficiency. ^[1]

The method was originally developed by Yingbo Hua and Tapan Sarkar for estimating the behaviour of electromagnetic systems by its transient response, building on Sarkar's past work on the original pencil-of-function method. ^{[1] [2]} The method has a plethora of applications in electrical engineering, particularly related to problems in computational electromagnetics, microwave engineering and antenna theory. ^[1]

Method

Mathematical basis

A transient electromagnetic signal can be represented as: ^[3]

${\displaystyle y(t)=x(t)+n(t)\approx \sum _{i=1}^{M}R_{i}e^{s_{i}t}+n(t);0\leq t\leq T,}$

where

${\displaystyle y(t)}$ is the observed time-domain signal,

${\displaystyle n(t)}$ is the signal noise,

${\displaystyle x(t)}$ is the actual signal,

${\displaystyle R_{i}}$ are the residues (${\displaystyle R_{i}}$),

${\displaystyle s_{i}}$ are the poles of the system, defined as ${\displaystyle s_{i}=-\alpha _{i}+j\omega _{i}}$,

${\displaystyle z_{i}=e^{(-\alpha _{i}+j\omega _{i})T_{s}}}$ by the identities of Z-transform,

${\displaystyle \alpha _{i}}$ are the damping factors and

${\displaystyle \omega _{i}}$ are the angular frequencies.

The same sequence, sampled by a period of ${\displaystyle T_{s}}$, can be written as the following:

${\displaystyle y[kT_{s}]=x[kT_{s}]+n[kT_{s}]\approx \sum _{i=1}^{M}R_{i}z_{i}^{k}+n[kT_{s}];k=0,...,N-1;i=1,2,...,M}$,

Generalized pencil-of-function estimates the optimal ${\displaystyle M}$ and ${\displaystyle z_{i}}$'s. ^[4]

Noise-free analysis

For the noiseless case, two ${\displaystyle (N-L)\times L}$ matrices, ${\displaystyle Y_{1}}$ and ${\displaystyle Y_{2}}$, are produced: ^[3]

${\displaystyle [Y_{1}]={\begin{bmatrix}x(0)&x(1)&\cdots &x(L-1)\\x(1)&x(2)&\cdots &x(L)\\\vdots &\vdots &\ddots &\vdots \\x(N-L-1)&x(N-L)&\cdots &x(N-2)\end{bmatrix}}_{(N-L)\times L};}$ ${\displaystyle [Y_{2}]={\begin{bmatrix}x(1)&x(2)&\cdots &x(L)\\x(2)&x(3)&\cdots &x(L+1)\\\vdots &\vdots &\ddots &\vdots \\x(N-L)&x(N-L+1)&\cdots &x(N-1)\end{bmatrix}}_{(N-L)\times L}}$

where ${\displaystyle L}$ is defined as the pencil parameter. ${\displaystyle Y_{1}}$ and ${\displaystyle Y_{2}}$ can be decomposed into the following matrices: ^[3]

${\displaystyle [Y_{1}]=[Z_{1}][B][Z_{2}]}$

${\displaystyle [Y_{2}]=[Z_{1}][B][Z_{0}][Z_{2}]}$

where

${\displaystyle [Z_{1}]={\begin{bmatrix}1&1&\cdots &1\\z_{1}&z_{2}&\cdots &z_{M}\\\vdots &\vdots &\ddots &\vdots \\z_{1}^{(N-L-1)}&z_{2}^{(N-L-1)}&\cdots &z_{M}^{(N-L-1)}\end{bmatrix}}_{(N-L)\times M};}$ ${\displaystyle [Z_{2}]={\begin{bmatrix}1&z_{1}&\cdots &z_{1}^{L-1}\\1&z_{2}&\cdots &z_{2}^{L-1}\\\vdots &\vdots &\ddots &\vdots \\1&z_{M}&\cdots &z_{M}^{L-1}\end{bmatrix}}_{M\times L}}$

${\textstyle [Z_{0}]}$ and ${\textstyle [B]}$ are ${\textstyle M\times M}$ diagonal matrices with sequentially-placed ${\textstyle z_{i}}$ and ${\textstyle R_{i}}$ values, respectively. ^[3]

If ${\textstyle M\leq L\leq N-M}$, the generalized eigenvalues of the matrix pencil

${\displaystyle [Y_{2}]-\lambda [Y_{1}]=[Z_{1}][B]([Z_{0}]-\lambda [I])[Z_{2}]}$

yield the poles of the system, which are ${\displaystyle \lambda =z_{i}}$. Then, the generalized eigenvectors ${\displaystyle p_{i}}$ can be obtained by the following identities: ^[3]

${\displaystyle [Y_{1}]^{+}[Y_{1}]p_{i}=p_{i};}$    ${\displaystyle i=1,...,M}$

${\displaystyle [Y_{1}]^{+}[Y_{2}]p_{i}=z_{i}p_{i};}$    ${\displaystyle i=1,...,M}$

where the ${\displaystyle ^{+}}$ denotes the Moore–Penrose inverse, also known as the pseudo-inverse. Singular value decomposition can be employed to compute the pseudo-inverse.

Noise filtering

If noise is present in the system, ${\textstyle [Y_{1}]}$ and ${\textstyle [Y_{2}]}$ are combined in a general data matrix, ${\textstyle [Y]}$: ^[3]

${\displaystyle [Y]={\begin{bmatrix}y(0)&y(1)&\cdots &y(L)\\y(1)&y(2)&\cdots &y(L+1)\\\vdots &\vdots &\ddots &\vdots \\y(N-L-1)&y(N-L)&\cdots &y(N-1)\end{bmatrix}}_{(N-L)\times (L+1)}}$

where ${\displaystyle y}$ is the noisy data. For efficient filtering, L is chosen between ${\textstyle {\frac {N}{3}}}$ and ${\textstyle {\frac {N}{2}}}$. A singular value decomposition on ${\textstyle [Y]}$ yields:

${\displaystyle [Y]=[U][\Sigma ][V]^{H}}$

In this decomposition, ${\textstyle [U]}$ and ${\textstyle [V]}$ are unitary matrices with respective eigenvectors ${\textstyle [Y][Y]^{H}}$ and ${\textstyle [Y]^{H}[Y]}$ and ${\textstyle [\Sigma ]}$ is a diagonal matrix with singular values of ${\textstyle [Y]}$. Superscript ${\textstyle H}$ denotes the conjugate transpose. ^{[3] [4]}

Then the parameter ${\textstyle M}$ is chosen for filtering. Singular values after ${\textstyle M}$, which are below the filtering threshold, are set to zero; for an arbitrary singular value ${\textstyle \sigma _{c}}$, the threshold is denoted by the following formula: ^[1]

${\displaystyle {\frac {\sigma _{c}}{\sigma _{max}}}=10^{-p}}$,

${\textstyle \sigma _{max}}$ and p are the maximum singular value and significant decimal digits, respectively. For a data with significant digits accurate up to p, singular values below ${\textstyle 10^{-p}}$ are considered noise. ^[4]

${\textstyle [V_{1}']}$ and ${\textstyle [V_{2}']}$ are obtained through removing the last and first row and column of the filtered matrix ${\textstyle [V']}$, respectively; ${\textstyle M}$ columns of ${\textstyle [\Sigma ]}$ represent ${\textstyle [\Sigma ']}$. Filtered ${\textstyle [Y_{1}]}$ and ${\textstyle [Y_{2}]}$ matrices are obtained as: ^[4]

${\displaystyle [Y_{1}]=[U][\Sigma '][V_{1}']^{H}}$

${\displaystyle [Y_{2}]=[U][\Sigma '][V_{2}']^{H}}$

Prefiltering can be used to combat noise and enhance signal-to-noise ratio (SNR). ^[1] Band-pass matrix pencil (BPMP) method is a modification of the GPOF method via FIR or IIR band-pass filters. ^{[1] [5]}

GPOF can handle up to 25 dB SNR. For GPOF, as well as for BPMP, variance of the estimates approximately reaches Cramér–Rao bound. ^{[3] [5] [4]}

Calculation of residues

Residues of the complex poles are obtained through the least squares problem: ^[1]

${\displaystyle {\begin{bmatrix}y(0)\\y(1)\\\vdots \\y(N-1)\end{bmatrix}}={\begin{bmatrix}1&1&\cdots &1\\z_{1}&z_{2}&\cdots &z_{M}\\\vdots &\vdots &\ddots &\vdots \\z_{1}^{N-1}&z_{2}^{N-1}&\cdots &z_{M}^{N-1}\end{bmatrix}}{\begin{bmatrix}R_{1}\\R_{2}\\\vdots \\R_{M}\end{bmatrix}}}$

Applications

The method is generally used for the closed-form evaluation of Sommerfeld integrals in discrete complex image method for method of moments applications, where the spectral Green's function is approximated as a sum of complex exponentials. ^{[1] [6]} Additionally, the method is used in antenna analysis, S-parameter-estimation in microwave integrated circuits, wave propagation analysis, moving target indication, radar signal processing, ^{[1] [7] [8]} and series acceleration in electromagnetic problems. ^[9]

See also

• Estimation of signal parameters via rotational invariance techniques
• Generalized eigenvalue problem
• Matrix pencil
• MUSIC (algorithm)
• Prony's method

References

1. Sarkar, T. K.; Pereira, O. (February 1995). "Using the matrix pencil method to estimate the parameters of a sum of complex exponentials". IEEE Antennas and Propagation Magazine. 37 (1): 48–55. Bibcode:1995IAPM...37...48S. doi:10.1109/74.370583.
2. Sarkar, T.; Nebat, J.; Weiner, D.; Jain, V. (November 1980). "Suboptimal approximation/identification of transient waveforms from electromagnetic systems by pencil-of-function method". IEEE Transactions on Antennas and Propagation . 28 (6): 928–933. Bibcode:1980ITAP...28..928S. doi:10.1109/TAP.1980.1142411.
3. Hua, Y.; Sarkar, T. K. (February 1989). "Generalized pencil-of-function method for extracting poles of an EM system from its transient response". IEEE Transactions on Antennas and Propagation . 37 (2): 229–234. Bibcode:1989ITAP...37..229H. doi:10.1109/8.18710.
4. Hua, Y.; Sarkar, T. K. (May 1990). "Matrix pencil method for estimating parameters of exponentially damped/undamped sinusoids in noise". IEEE Transactions on Acoustics, Speech, and Signal Processing . 38 (5): 814–824. doi:10.1109/29.56027.
5. Hu, Fengduo; Sarkar, T. K.; Hua, Yingbo (January 1993). "Utilization of Bandpass Filtering for the Matrix Pencil Method". IEEE Transactions on Signal Processing . 41 (1): 442–446. Bibcode:1993ITSP...41..442H. doi:10.1109/TSP.1993.193174.
6. Dural, G.; Aksun, M. I. (July 1995). "Closed-form Green's functions for general sources and stratified media". IEEE Transactions on Microwave Theory and Techniques . 43 (7): 1545–1552. Bibcode:1995ITMTT..43.1545D. doi:10.1109/22.392913. hdl: 11693/10756 .
7. Kahrizi, M.; Sarkar, T. K.; Maricevic, Z. A. (January 1994). "Analysis of a wide radiating slot in the ground plane of a microstrip line". IEEE Transactions on Microwave Theory and Techniques . 41 (1): 29–37. doi:10.1109/22.210226.
8. Hua, Y. (January 1994). "High resolution imaging of continuously moving object using stepped frequency radar". Signal Processing. 35 (1): 33–40. doi:10.1016/0165-1684(94)90188-0.
9. Karabulut, E. Pınar; Ertürk, Vakur B.; Alatan, Lale; Karan, S.; Alişan, Burak; Aksun, M. I. (2016). "A novel approach for the efficient computation of 1-D and 2-D summations". IEEE Transactions on Antennas and Propagation . 64 (3): 1014–1022. Bibcode:2016ITAP...64.1014K. doi:10.1109/TAP.2016.2521860.