1/6 Excited to announce the McGill Radio Lab’s latest paper on digital signal processing for real time astronomical radio interferometers (radio telescopes) was submitted for publication and just hit the Arxiv. https://t.co/w3AqqDQ21G

Inverting and Up-channelizing Critically-Sampled Polyphase Filter Banks

Source: https://arxiv.org/html/2608.23441

Inverting and Up-channelizing Critically-Sampled Polyphase Filter BanksNumpy(11), Scipy(31), Matplotlib(14).

Stephen FayAffiliation:Department of Physics, McGill University, Montreal, Quebec H3A 2T8, CanadaEmail:[email protected]Mohan AgrawalAffiliation:Department of Physics, McGill University, Montreal, Quebec H3A 2T8, CanadaAffiliation:Trottier Space Institute, McGill University, Montreal, Quebec H3A 2A7, CanadaEmail:[email protected]Email:[email protected]Corresponding author:Mohan AgrawalH. Cynthia ChiangAffiliation:Department of Physics, McGill University, Montreal, Quebec H3A 2T8, CanadaAffiliation:Trottier Space Institute, McGill University, Montreal, Quebec H3A 2A7, CanadaEmail:[email protected]Aman ChokshiAffiliation:Department of Physics, McGill University, Montreal, Quebec H3A 2T8, CanadaAffiliation:Trottier Space Institute, McGill University, Montreal, Quebec H3A 2A7, CanadaEmail:[email protected]Jonathan SieversAffiliation:Department of Physics, McGill University, Montreal, Quebec H3A 2T8, CanadaAffiliation:Trottier Space Institute, McGill University, Montreal, Quebec H3A 2A7, CanadaEmail:[email protected]Simon TartakovskyAffiliation:Department of Physics, Princeton University, Princeton, New Jersey 08540, USAEmail:[email protected]

Abstract

Polyphase filter banks (PFBs) are widely used to channelize digitized time-domain data in real time. However, applications such as radio astronomy often require higher spectral resolution than real-time systems can provide. We introduce a circulant inverse-PFB formalism and present a new, scalable, fast Fourier transform (FFT)-based algorithm for inverting critically sampled PFBs. This formalism explicitly identifies the origin of quantization noise amplification and enables two complementary mitigation strategies: (1) a fast, practical Wiener filter and (2) a rigorous maximum-likelihood reconstruction with time-domain priors. We evaluate both approaches using simulated 4-bit quantized PFB data. For autospectra, Wiener filtering confines reconstruction errors to less than 10% of the channel width, with a peak error below 10%. The maximum-likelihood approach further reduces reconstruction errors as additional time-domain prior information is incorporated, with priors spanning 10% of the time-domain samples reducing worst-case errors to below 2%. This framework enables high-resolution spectral analysis of archival PFB data, facilitating applications including long-baseline interferometry, pulsar and fast radio burst searches, and ultra-narrow-band physics experiments.

IIntroduction

A polyphase filter bank(12;7;29;3;21, PFB;)is a signal processing technique, common in fields like radio astronomy and telecommunications, for splitting a broadband digital signal into many narrow frequency channels. PFB channelizers allow custom sculpting of the frequency response of the output channels with user-defined filters while also being much faster than traditional channelizing schemes. Their distinct advantage is manifest in radio astronomy where minimizing spectral leakage(20)between adjacent channels is necessary to prevent strong human-generated interference from contaminating frequencies of scientific interest. A PFB channelizer can provide any desired level of passband size and sidelobe suppression with an appropriate length of filter used in its construction. Typical PFBs used in modern radio telescopes provide excellent stop-band attenuation (Fig.1), and for the same level attenuation, do so at a fraction of the computational cost associated with naive channelization techniques.

Radio telescopes typically process large instantaneous radio-frequency (RF) bandwidths, often ranging up to several GHz. Processing large RF bandwidths requires fast sampling rates, which results in high throughput data rates. For this reason, hardware constraints at the input, such as buffer size limitations on field programmable gate arrays (FPGA), or at the output, such as data transmission or storage capacity, limit the frequency resolution achievable with real-time channelization. Therefore, many instruments record only coarsely-channelized PFB data, even though finer frequency resolution could be scientifically valuable.

Recovering the original timestream from archival PFB data lifts these restrictions: once inverted, the data can be re-channelized to arbitrarily fine frequency resolution to suit the science case at hand. This capability is desirable across a range of applications. In long baseline interferometry, for example, coarse channelization limits achievable frequency resolution, which in turn limits the maximum usable baseline length and field-of-view because of the bandwidth smearing effect(4). Re-channelization of coarse PFB data mitigates this effect, thereby enabling ultra-wide-field long-baseline interferometry(5). An accurate PFB inversion scheme could also enhance the sensitivity of pipelines designed for pulsar searches, and fast radio burst search and localization(6;15). Lastly, re-channelization capability could potentially enable searches for new physics, such as ultra-narrowband (few-Hz) signals from axion-like particles in the ionosphere(2), with existing radio telescopes.

The most common implementation of a PFB channelizer iscritically sampled, meaning that the frequency spacing between channels is the same as a channel’s Nyquist bandwidth. Passing a real-valued signal of bandwidthBB, sampled at the Nyquist rate of2​B2B, through a PFB that uses anNN-point FFT at the channelization stage producesnc​h​a​n=N/2n_{chan}=N/2positive frequency channels that are spacedB/nc​h​a​nB/n_{chan}apart. If the PFB is critically-sampled, then for each channel, a new complex-valued sample also arrives at the rate ofB/nc​h​a​nB/n_{chan}, which makes the input rate and the output rate (summed across all channels) equal. Critical-sampling causes any signal leaking outside a channel’s Nyquist width ofB/nc​h​a​nB/{n_{chan}}to alias back into the channel, corrupting the channel edges (see Fig.1). This aliasing makes critically-sampled PFBs difficult to invert. As will be discussed later, aliasing alone is not sufficient to make inversion impossible: while information has been scrambled, it is not lost.

Traditionally, invertible PFBs have been designed from the ground-up using so-calledperfect-reconstructionfilters. These filters exactly cancel aliasing artifacts from neighbouring channels, but this property imposes strict constraints on their length and shape(17;16;29;30). Such restrictions limit the filter’s versatility in providing the desired level of stop-band attenuation and pass-band flatness with the available hardware resources. Certain classes of these filters also lack a linear phase response(26), which makes them further unsuitable for phase-sensitive applications where linearity is required.

Another way to sidestep the difficulties posed by critical sampling is to use oversampled PFBs(19;1;28, OS-PFBs;), where the output timestream is sampled faster thanB/nc​h​a​nB/n_{chan}. A higher output sampling rate makes a channel’s Nyquist bandwidth larger than the channel width. The out-of-channel ripples thus fall within a channel’s Nyquist band, reducing the aliasing. However, the benefits of OS-PFBs come at the cost of additional hardware resources that are required to process the now higher output data rate, and correct for systematic phase rotations introduced by the oversampling(27).

In this paper, we focus on the challenges associated with inverting critically-sampled PFBs. First we pinpoint the reasons for invertibility issues that fundamentally arise from symmetry properties of the filter weights (windows) typically used in PFBs. We then present a computationally efficient inversion algorithm. The organization of the rest of the paper is as follows. SectionIIestablishes the foundations of PFBs, from ideal channelizers to their interpretation as independent correlations. SectionIIIpresents an FFT-based algorithm for fast PFB inversion and addresses the effects of quantization noise on the inverted PFB. Finally, sectionIVaugments the inversion algorithm with two de-noising methods: a Wiener filter approach, and a time-domain maximum-likelihood approach, and compares their performance using simulated data.

II.1Ideal Channelizer

The process of channelizing a timestream signal refers to splitting the information content into narrow frequency bands, or “channels.” A band-limited, real timestreamx⁡(t)x(t)can be split intonc​h​a​nn_{chan}channels, each of which has its own corresponding timestream. The sum of these per-channel timestreams yields the originalx⁡(t)x(t). Ideally, any individual channelkkhaving widthΔ​f\Delta fhas information that is strictly contained within the frequency range|f−fk|≤Δ​f/2|f-f_{k}|\leq\Delta f/2, with no information content from other frequencies. The ideal channelizer is straightforward. To obtain the timestream for channelkk, the Fourier transformX⁡(f)X(f)of the original timestream is calculated and then multiplied by a boxcar filter, defined asHk​(f)=1H_{k}(f)=1for|f−fk|≤Δ​f/2|f-f_{k}|\leq\Delta f/2, and zero otherwise. This boxcar simply selects frequencies within the vicinity offkf_{k}(thekk’th channel) and rejects frequencies outside this channel. By the convolution theorem, multiplyingX⁡(f)X(f)with a boxcar centred atfkf_{k}is equivalent to convolving the timestream with the shifted boxcar’s Fourier transform: a sinc function times a phaseej​2​π​fk​te^{j2\pi f_{k}t}.

If the bandwidth ofx⁡(t)x(t)isBB, such thatB=nc​h​a​n​Δ​fB=n_{chan}\Delta f, then the timestream can be digitized by sampling at the Nyquist rate of2​B2B. The convolution for obtaining the timestream of channelkkcan then be written as a sum over the samples ofx⁡(t)x(t):

yk​(t)=∑mx⁡[tm]​sinc⁡(tm−t)​e−i​2​π​fk​(tm−t)⏟Convolution with the filter,y_{k}(t)=\sum_{m}x[t_{m}]\underbrace{{\rm sinc}(t_{m}-t)e^{-i2\pi f_{k}(t_{m}-t)}}_{\text{Convolution with the filter}},(1)wheretm=m/2​Bt_{m}=m/2Bare the timestream samples andyk​(t)y_{k}(t)are the timestream samples for channelkk. As the frequency content of the resultingnc​h​a​nn_{chan}channelized timestreams is limited toΔ​f=B/nc​h​a​n\Delta f=B/n_{chan}, each channelized timestream can be sampled more sparsely thanx⁡(t)x(t). The lowest possible sampling rate is set by the Nyquist limit to beB/nc​h​a​n=Δ​fB/n_{chan}=\Delta f, which reduces the number of output samples ofyk​(t)y_{k}(t)by a factor ofnc​h​a​nn_{chan}relative to the length ofx⁡(t)x(t). Now, if the channelizer output for channelkkis evaluated at times that are integer multiples for1/Δ​f1/\Delta f, then (1) reduces to an infinitely-long discrete Fourier transform (DFT)111This formulation simplifies the formal boundary mapping between continuous spectra and discrete DFT bins. The underlying mechanism is that sampling the channelized output at the downsampled rate ofB/nchanB/n_{\text{chan}}frequency-translates each channel to baseband, causing the complex phase factore−j​2​π​νk​te^{-j2\pi\nu_{k}t}to evaluate to unity at these sampling instants.. Moreover, as the output of the channelizer is now complex, the total data rate across all channels is twicenc​h​a​nn_{chan}timesB/nn​c​h​a​nB/n_{nchan}, equal to2​B2B, which is the input data rate ofx⁡(t)x(t). This property makes the channelizer critically sampled.

II.2Introduction to Polyphase Filter Banks

Channelizing a real timestream with finite lengthMMrequires some modifications to the ideal channelizer. First, the sinc filter is truncated to the same lengthMMas the timestream. To obtainnc​h​a​nn_{chan}complex-valued output channels, a minimum ofN=2​nc​h​a​nN=2n_{chan}real-valued input timestream points are required, such thatM≥NM\geq N. We refer to both a unit ofNNtime-domain samples and a set ofN/2N/2complex channelized values as aframeof data. A sequence ofppconsecutive frames is then defined as asegmentofMMsamples. As M increases, the Fourier transform of the truncated sinc progressively approaches the ideal boxcar filter. The ratiop=M/Np=M/Nis known as the number of PFBtaps. In practice,ppis typically determined by available computational memory, withp=4p=4being common in radio astronomy.

WithM=p​NM=pNset as the number of timestream data points to channelize, (1) can be rewritten for one segment of data of sizeMMas

y[k]=∑m=0M−1x[m]sinc[m−M/2]e−i2πmk/N.y[k]=\sum_{m=0}^{M-1}x[m]{\rm sinc}[m-M/2]e^{-i2\pi mk/N}.(2)This form can be obtained by using the following substitutions:fk=k​B/nn​c​h​a​nf_{k}=kB/n_{nchan},tm=m/2​Bt_{m}=m/2B,N=2​nc​h​a​nN=2n_{chan}, along with the fact thatttis an integer multiple ofnc​h​a​n/Bn_{chan}/B. Moreover, as the summation starts fromm=0m=0, the sinc is shifted byM/2M/2samples so that it is symmetric across one block of data222A filter having an even frequency response (top-hat), must necessarily be even in time..

IfM>NM>N, then (2) is equivalent to taking a longerMM-sized DFT, but with only everypp’th frequency channel retained:

X[k′=pk]=∑m=0M−1x[m]sinc[m−M/2]e−i2πmpk/M⏟e−i2πmk/N.X[k^{\prime}=pk]=\sum_{m=0}^{M-1}x[m]{\rm sinc}[m-M/2]\underbrace{e^{-i2\pi mpk/M}}_{e^{-i2\pi mk/N}}.(3)The efficiency of this processing step can be improved by noting that the exponential in (2) repeats with a period ofNN. Thus, a summation overM=p​NM=pNhasppidentical exponentials. This redundancy can be exploited by dividing the sum overMMsamples intoppdifferent sums, using the substitutionm=r​N+nm=rN+n, wherer<pr<pandn<Nn<N. Lettingh⁡[m]≡sinc⁡[m−M/2]h[m]\equiv{\rm sinc}[m-M/2]for notational brevity, it follows that

y[k]=∑n=0N−1e−i2πmk/N∑r=0p−1x[rN+n]h[rN+M].y[k]=\sum_{n=0}^{N-1}e^{-i2\pi mk/N}\sum_{r=0}^{p-1}x[rN+n]h[rN+M].(4)Equation (4) represents channelization with a polyphase filter bank. The sequence of operations involved in producing one PFB spectrum, as illustrated in Fig.2, are as follows. A timestream segment of lengthp​NpNis multiplied by a sinc and divided intoppsequential frames of lengthNN, starting at samples00,NN, …(p−1)​N(p-1)Nrespectively. Corresponding samples of thesespptimestreams are summed, and the resultingNNsamples are Fourier transformed. This process is repeated for the adjacent timestream segment, starting a frame (NNsamples) later. Each data point is therefore usedpptimes in the computation. The power of the PFB lies in achieving the filtering fidelity of a Fourier transformpptimes longer than the transform actually computed. For a given number of frequency channels, increasing the number of tapsppenhances both the flatness of the response as well as out-of-band rejection.

In practice, the filterh⁡[m]h[m]is multiplied by a smooth taper (or windowing) function(22)to avoid abruptly truncating the sinc, allowing the shapes of the pass and side-bands to be fine tuned. Because the length of the filter is necessarily finite, its response outside the nominal channel width ofB/nc​h​a​nB/n_{chan}cannot be zero. When the sample rate of a channel timestream is reduced toB/nc​h​a​nB/n_{chan}, this spectral leakage aliases back into the channel. Fig.1compares the spectral response of various PFB filters varying in their length and taper functions. In the rest of this text,windowandfilterare used interchangeably to mean the sinc function (optionally) multiplied by a conventional taper function.

Figure 1:Comparison of frequency responses of a channel of PFB and FFT filter banks. All responses are normalized to unity at the center of the reference channel. The gray shaded area indicates the ideal pass band of a digital channel. The output of a channelizer at any given instant is the input spectrum multiplied by the frequency response. Thus, the plots indicate the amount of power contributed to the reference channel by its neighbors (referred to as spectral leakage). Compared to a simple FFT, a windowed FFT reduces the level of sidelobes at the expense of a wider main-lobe. In contrast, a PFB achieves both a narrow main-lobe and low side-lobes.

II.3The PFB as Independent Correlations

The first step towards PFB inversion is rewriting the channelizer in (4) as a matrix operation:

𝐲=ℱ​𝒮​𝒲​𝐱.\mathbf{y}=\mathcal{FSW}\mathbf{x}.(5)Here, the PFB operatorℱ​𝒮​𝒲\mathcal{FSW}acts on an input timestream vector,𝐱\mathbf{x}, of lengthp​NpN. First, the square, diagonal window matrix𝒲\mathcal{W}element-wise multiplies𝐱\mathbf{x}with the filter weights. Next,𝒮\mathcal{S}is a summing matrix of shape[N,p​N][N,pN]that comprisesppidentity matrices of sizeNNstacked horizontally. This matrix pointwise-sums samples ofppconsecutive time-domain frames, outputting anNN-long vector. Finally,ℱ\mathcal{F}is aNN-point DFT operator withℱn,k=exp(−j2πnk/N)\mathcal{F}_{n,k}=\exp(-j2\pi nk/N). To compute the PFB of a longer timestream, the operatorℱ​𝒮​𝒲\mathcal{FSW}is repeatedly applied to successive timestream segments, each shifted by a frame ofNNsamples. This section provides a detailed breakdown ofℱ​𝒮​𝒲\mathcal{FSW}, which encodes a limited set of independent operations. By reshaping pieces of the constituent matrices, these operations can be isolated from each other, yielding a framework for PFB inversion (and its subsequent implementation with parallel computations).

The only non-zero entries in rownnof𝒮​𝒲\mathcal{SW}are:h⁡[n],h⁡[n+N],…,h⁡[n+(p−1)​N]h[n],h[n+N],\ldots,h[n+(p-1)N]in columnsn,n+N,…,n+(p−1)​Nn,n+N,\ldots,n+(p-1)N, respectively. Thus, each𝒮​𝒲\mathcal{SW}row combines only the correspondingpptimestream samples. To streamline the application of𝒮​𝒲\mathcal{SW}to an arbitrarily long timestream, two quantities are reshaped. First, the input timestream is written as a 2-D array withNNcolumns and an arbitrary number of rows, and second, thep​NpN-long window is reshaped into an array of dimensions[p,N][p,N]. With this restructuring, a single application of𝒮​𝒲\mathcal{SW}encompasses two steps. First,pprows (frames) of the data array are elementwise multiplied by the window array. Second, allpprows of the resulting product are summed along columns, yielding aNN-long output which is passed toℱ\mathcal{F}, resulting in one PFB spectrum. To obtain the next spectrum, the multiply (𝒲\mathcal{W}) and add (𝒮\mathcal{S}) operations are repeated again, but with the data array segment shifted by one frame (row). Repeated application of𝒮​𝒲\mathcal{SW}can be envisioned as the window array sliding down the data array, as illustrated in Fig.2.

Figure 2:Sequence of steps involved in a PFB, shown without the final Fourier transform. The illustrations demonstrate a 4-tap,NN-point PFB. The top panels show 6 frames of timestream data, and the4​N4N-long filter window composed of sinc function pointwise-multiplied with a Hamming taper function. The three lower panels demonstrate the operations performed by the PFB on three successive overlapping segments of four frames: first, multiply the four frames with a window weights, then sum the transformed frames pointwise in the vertical direction, producing one frame of filtered data. The operation is repeated on an input shifted by a frame. Though not illustrated here, each of the final output frames are Fourier transformed to complete the PFB channelization.The reshaped representation reduces the𝒮​𝒲\mathcal{SW}operation to independent column-wise multiply-and-sum operations that together produce a single output frame (NNvalues). This output corresponds to the “polyphase” step (the inner sum in (4)), which can be understood by considering the operations acting ononeof theNNindependent columns. Consider columnnnof the data array, which contains timestream samplesxn,xn+N,xn+2​N,…x_{n},x_{n+N},x_{n+2N},\dotsetc. These samples are drawn from the input timestream starting at samplennand sub-sampled by a factor ofNN. Similarly, columnnnof the window array is theNN-fold sub-sampled window function. The output corresponding to columnnnat every𝒮​𝒲\mathcal{SW}application (once everyNNsamples) is then simply the correlation333Correlation with the sub-sampled window function is equivalent to convolution with the full-rate filter described by (1). The window coefficients are the time-reversed filter coefficients; see (2) for example.of the sub-sampled timestream with the corresponding sub-sampled filter. Filtered outputs from allNNsub-filters are then passed to the Fourier transform block. In summary, in a PFB,NNsub-sampled timestreams are filtered by their respectiveNN-distinctpp-sized filters in parallel. In the remainder of this text, a sub-sampled filter is referred to as asub-filter, and a sub-sampled timestream as asub-stream.

As a concrete example, consider aNN-point PFB withp=4p=4taps acting on a long timestream of lengthb​NbN(b≫pb\gg p). The output of sub-filter pathnncan be written as

𝐝n=Wn​𝐱n,\mathbf{d}_{n}=W_{n}\mathbf{x}_{n},(6)where𝐱n\mathbf{x}_{n}is the sub-sampled input timestream corresponding to columnnnand the operatorWnW_{n}is rectangular matrix of shape[b−3,b][b-3,b]given by

Wn=[wnwn+Nwn+2​Nwn+3​N⋯00wnwn+Nwn+2​N⋯0⋱0⋯wnwn+Nwn+2​Nwn+3​N]W_{n}=\begin{bmatrix}w_{n}&w_{n+N}&w_{n+2N}&w_{n+3N}&\cdots&0\\ 0&w_{n}&w_{n+N}&w_{n+2N}&\cdots&0\\ \vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\ 0&\cdots&w_{n}&w_{n+N}&w_{n+2N}&w_{n+3N}\\ \end{bmatrix}(7)The Toeplitz structure of operatorWnW_{n}clearly demonstrates the correlation of the sub-sampled timestream of pathnnwith a filter of sizep=4p=4. By decomposing a PFB intoNNindependent and identical sub-problems, the task of inverting the PFB is reduced to that of inverting one sub-problem. We exploit this decomposition to develop a fast inversion strategy discussed in the next section.

IIIInverse Polyphase Filter Banks

PFB inversion begins with applying the inverse FFT (IFFT) to the recorded PFB spectra. This step undoes the action ofℱ\mathcal{F}, leaving the task of decorrelating the output of the stacked𝒮​𝒲\mathcal{SW}operator, i.e., solving for𝐱n\mathbf{x}_{n}in (6) for allNNindependent sub-filters. The system of equations given by (6) is underdetermined since there are more input values than output values. One way to solve this system is to use the Moore-Penrose pseudo-inverse, e.g., as explored by24. However, in practice, the PFB output often contains additional noise introduced by quantization, and this noise is amplified when a naive pseudo-inverse is applied without any additional regularization (details in §III.2). This section describes a method for solving (6) using circulant boundary conditions. This method enables a fast inverse PFB (IPFB) implementation and elucidates the conditions under which quantization noise dominates. (Methods for suppressing this noise will be addressed in §IV.)

III.1Circulant inversion

Supposeb−p+1b-p+1PFB frames were recorded usingbbframes of input data according to (6). Asb≫pb\gg p, the matrixWnW_{n}is nearly circulant, or equivalently, the correlation in (6) is nearly a circular correlation. By appendingp−1p-1rows toWnW_{n}(three additional rows forp=4p=4), it can be promoted to a circulant matrixWncircW_{n}^{\mathrm{circ}}of sizebb. Equation (6) then becomes a circulant correlation,

𝐝ncirc=Wncirc​𝐱n=𝐰n⊛𝐱n,\mathbf{d}^{\mathrm{circ}}_{n}=W_{n}^{\mathrm{circ}}\mathbf{x}_{n}=\mathbf{w}_{n}\circledast\mathbf{x}_{n},(8)where⊛\circledastdenotes a circular correlation,𝐝ncirc\mathbf{d}^{\mathrm{circ}}_{n}is the original data𝐝n\mathbf{d}_{n}extended byp−1p-1entries corresponding to the circulant output, and the kernel𝐰n\mathbf{w}_{n}is the first row ofWnW_{n}. The advantage of adopting this circulant model is thatWncircW_{n}^{\mathrm{circ}}is diagonalized by the DFT operator, thus allowing FFTs to be used for decorrelation. To implement circulant decorrelation in practice, we extend the size of𝐝n\mathbf{d}_{n}such that it containsbbframes of real data, i.e.,dim(𝐝n)=dim(𝐝ncirc)=b\dim(\mathbf{d}_{n})=\dim(\mathbf{d}^{\mathrm{circ}}_{n})=b. An estimate of the original timestream is then formed with the following procedure. LetFFdenote a DFT of sizebb, andF†F^{\dagger}its inverse. By the correlation theorem,

𝐱^n=F†​[F​𝐝n(F​𝐰n)∗]=F†​[1(F​𝐰n)∗]⊛𝐝n,\hat{\mathbf{x}}_{n}=F^{\dagger}\!\left[\frac{F\mathbf{d}_{n}}{(F\mathbf{w}_{n})^{*}}\right]=F^{\dagger}\!\left[\frac{1}{(F\mathbf{w}_{n})^{*}}\right]\circledast\mathbf{d}_{n},(9)where𝐱^n\hat{\mathbf{x}}_{n}denotes an estimate of𝐱n\mathbf{x}_{n}, the asterisk (∗*) denotes an element-wise conjugate, and all multiplication and division operations in (9) are also applied element-wise.

As real data are not circulant, the lastp−1p-1entries of the extended𝐝n\mathbf{d}_{n}differ from𝐝ncirc\mathbf{d}^{\mathrm{circ}}_{n}. By enforcing circularity during decorrelation, the samples of the reconstructed𝐱^n\hat{\mathbf{x}}_{n}are corrupted near the array boundaries. However, the error in the circulant model approaches zero for samples that are sufficiently far from the edges (explained in more detail in Fig.4). In use cases with high data rates, as is typical for radio astronomy, the loss from edge corruption is negligibly small even for modest lengths of data.

In (9), the Fourier transformF​𝐰nF\mathbf{w}_{n}of the sub-sampled window is the frequency response of thennth sub-filter, with Nyquist frequencyB/NB/N. If the response of a sub-filter reaches zero, the reconstruction of the corresponding sub-stream𝐱n\mathbf{x}_{n}is ill defined at that frequency since inverting the correlation requires dividing by zero. In situations whereF​𝐰nF\mathbf{w}_{n}is sufficiently small, round-off error and quantization noise in𝐝n\mathbf{d}_{n}dominate the recovered𝐱n\mathbf{x}_{n}.

The entries ofF​𝐰nF\mathbf{w}_{n}are also the eigenvalues of the circulant matrixWncircW_{n}^{\text{circ}}. SinceWnW_{n}is obtained fromWncircW_{n}^{\text{circ}}by deleting its lastp−1p-1rows, the singular values ofWnW_{n}are well approximated by the eigenvalue magnitudes|F​𝐰n||F\mathbf{w}_{n}|. A vanishing eigenvalue therefore corresponds to a vanishing singular value, and any decorrelation algorithm suffers noise amplification along that mode. We accordingly use the termseigenvalue spectrumandfrequency spectrumof a sub-filter interchangeably.

Invertibility of the PFB is assessed by examining the spectra of allNNsub-filters, shown in Fig.3for a PFB withN=2048N=2048,p=4p=4, and a sinc–Hamming window. Each of the 2048 columns of Fig.3is the FFT of the zero-padded kernel𝐰n\mathbf{w}_{n}(positive frequencies only), i.e., the frequency spectrum of that sub-filter. The spectra of sub-filter paths nearn=N/2n=N/2(center of the plot) approach zero at the band edge, suppressing those frequency modes of the corresponding sub-sampled timestream. Because each channel of the final PFB spectrum depends on the summed output ofallsub-filters, information lost in any one sub-filter prevents perfect reconstruction of the underlying spectrum when the inverted data are re-channelized at finer resolution. AppendixAgives details of how this degeneracy is unavoidable for critically sampled PFBs with any symmetric window function.

Figure 3:Magnitudes of the sub-filter eigenvalues for ap=4p=4tap,N=2048N=2048point PFB with a sinc–Hamming window. Columnnnshows|F​𝐰n||F\mathbf{w}_{n}|, the frequency spectrum of sub-filternn(positive frequencies only, increasing upward); these are the eigenvalue magnitudes of the circulant operatorsWncircW_{n}^{\text{circ}}of Eq. (8). The normalized Nyquist frequency corresponds to a physical frequency ofB/NB/N, the width of one PFB channel. For sub-filters nearn=N/2n=N/2, the spectrum approaches zero toward the band edge (top center). Reconstruction in this region is degenerate: quantization noise dominates the decorrelation, limiting invertibility.

III.2Quantization Effects

Raw PFB data from an FPGA are often quantized using a limited number of bits because of constraints on data volume or bandwidth. Quantization introduces noise to the recorded spectra, and this noise is amplified during the PFB inversion process because of singularities in the convolution kernel. Circular decorrelation enables straightforward quantitative characterization of this corruption.

If the PFB output in (5) is quantized using a high number of bits (≳4\gtrsim 4), then under the Additive Quantization Noise Model(10;8), the recorded spectrum isℱ​𝒮​𝒲​𝐱+𝐧𝐪\mathcal{FSW}\mathbf{x}+\mathbf{n_{q}}, where𝐧𝐪\mathbf{n_{q}}is the quantization noise (difference between each sample inℱ​𝒮​𝒲​𝐱\mathcal{FSW}\mathbf{x}and the nearest quantized level.) As different values of𝐧𝐪\mathbf{n_{q}}are uncorrelated with each other and are approximately uniformly distributed,𝐧𝐪\mathbf{n_{q}}can be treated as white noise when assessing its impact on the IPFB.

After invertingℱ\mathcal{F}, the decorrelation of𝒮​𝒲\mathcal{SW}is performed by dividing the DFTs of the sub-filters by the conjugate DFT ofWn\mathrm{W}_{n}, but in the presence of noise:

𝐱^n=F†​{F​𝐝n(F​𝐰n)∗+F​ϵn(F​𝐰n)∗}.\hat{\mathbf{x}}_{n}=F^{\dagger}\left\{\frac{F\mathbf{d}_{n}}{(F\mathbf{w}_{n})^{*}}+\frac{F\boldsymbol{\epsilon}_{n}}{(F\mathbf{w}_{n})^{*}}\right\}.(10)Hereϵn\boldsymbol{\epsilon}_{n}is the time-domain version of the original quantization noise𝐧𝐪\mathbf{n_{q}}, obtained after application ofNN-point IFFT to the PFB spectra, i.e.ϵn=(ℱ−1​𝐧𝐪)n\boldsymbol{\epsilon}_{n}=(\mathcal{F}^{-1}\mathbf{n_{q}})_{n}. For sub-filter paths close toN/2N/2,F​𝐰nF\mathbf{w}_{n}is nearly singular, causing the quantization noise term in (10) to be strongly amplified during reconstruction. Consequently, the reconstruction error,𝐱^n−𝐱n\hat{\mathbf{x}}_{n}-\mathbf{x}_{n}, in samples corresponding to these sub-filter paths is much higher than would be naively expected from the quantization noise in the PFB data. For example, when a PFB withN=2048N=2048-point DFTs is inverted, the RMS reconstruction error in samples corresponding to sub-filter pathN/2=1024N/2=1024is many times higher than in samples from other filter paths. The reconstruction error therefore appears in the timestream as periodic noise spikes concentrated near the center of each frame. When the IPFB timestream is channelized again at higher frequency resolution, these periodic timestream spikes produce a corresponding picket-fence pattern of noise spikes in the frequency spectrum.

This decorrelation operation can be equivalently viewed as correlation with a kernel having a frequency response of1/F​𝐰n1/F\mathbf{w}_{n}. Examining the time-domain response of this “inverse” kernel for various sub-filter paths provides an alternate approach for understanding the noise amplification in the reconstructed timestream. Decorrelation kernels for two sub-filter paths are contrasted in Fig.4. Plots labeled “unfiltered” correspond to inverse kernels used for decorrelation in (10), and are discussed in this section. Since the inverse kernel for path 1024 has a wide region of support, one sample of reconstructed sub-stream𝐱^N/2\hat{\mathbf{x}}_{N/2}receives noise contributions from several thousand samples around it, thereby boosting the noise floor. Inverse kernels for sub-filter paths far away fromN/2N/2, however, approach a delta function. Therefore, a sample of sub-stream𝐱^0\hat{\mathbf{x}}_{0}, for example, contains quantization noise from that sample alone.

The qualitative behaviour of reconstruction error in the presence of quantization noise described here will be shown to agree with detailed numerical simulations in the next section, which describes various methods to mitigate the effects of quantization noise.

IVMitigating Quantization Effects

This section describes and compares two noise mitigation methods. The first is a Wiener filter approach, which employs prior knowledge of the noise and signal statistical properties and can be applied because the noise is additive. The second is a maximum-likelihood optimization formulation of the IPFB, in which a time-domain prior containing information on the degenerate regions of the filter response is imposed to improve the conditioning of the correlation matrices.

IV.1Wiener Filter

Figure 4:Time-domain representation of IPFB decorrelation kernels for a 2048-point, 4-tap PFB that uses a sinc-Hamming window function. Decorrelation kernels are obtained as the FFT of the inverse filter response,Hn/(F​𝐰n)∗H_{n}/(F\mathbf{w}_{n})^{*}. In this example,HnH_{n}is either 1 (naive decorrelation) or the Wiener filter, discussed in sectionIV. The two columns of panels show two different sub-filter paths. The bottom panels are zoomed-in versions of the top panels. Samples on the horizontal axis correspond to a sampling rate of2​B/N2B/N(raw timestream sub-sampled by factor ofNN). The kernels show that sub-filters near the frame center (1024) have a broad but finite domain of support: the worst-conditioned filter path decays to 0 in approximately 20,000 samples. As the sub-filter path moves away from the center, the corresponding kernel approaches a delta function. Therefore, in the absence of quantization noise, circulant decorrelation approaches floating-point precision away from the edges of the recovered timestream. Upon decorrelation, the 3-sample error mentioned in sectionIIIcorrupts only as many samples as the width of the worst-conditioned filter path. When inverting large volumes of PFB data, the fraction lost to edge-corruption is insignificant.When the covariance, or equivalently the power spectral density (PSD) of the signal and noise are known, the Wiener filter yields the optimal444in a least-squares sense.estimate of the underlying signal from a noisy observation(32;23). The general form of the filter isS/(S+N)S/(S+N), whereSSandNNare, respectively, the PSD of signal and noise. For the case of inverting a PFB, the appropriate form of the filter for decorrelating the effect of the sub-sampled window function from sub-filter pathnnis

Hn​(f)=|Wn​(f)|2|Wn​(f)|2+ϕ2​(1+ϕ2),H_{n}(f)=\frac{|W_{n}(f)|^{2}}{|W_{n}(f)|^{2}+\phi^{2}}(1+\phi^{2}),(11)whereWn​(f)=F​𝐰nW_{n}(f)=F\mathbf{w}_{n}is the frequency spectrum of sub-filternn(vertical slices in Fig.3), andϕ\phiis a threshold parameter that controls the strength of noise suppression. The multiplicative factor1+ϕ21+\phi^{2}is a normalization term that keeps power unchanged whereWn=1W_{n}=1. Fig.5shows the effects of Wiener filtering on reducing the percentage RMS error on the recovered time-stream samples, and on the auto-spectrum of the re-channelized PFB data. The fractional error on auto-spectrum is calculated as

δ⁡(f)=⟨|X⁡(f)−X^​(f)|2⟩⟨|X⁡(f)|2⟩,\delta(f)=\frac{\langle|X(f)-\hat{X}(f)|^{2}\rangle}{\langle|X(f)|^{2}\rangle},(12)whereX⁡(f)X(f)andX^​(f)\hat{X}(f)are the high-frequency-resolution spectra of the original timestream and the IPFB timestream, respectively. The ensemble average is taken by averaging many spectra along the time axis (described in more detail in sectionIV.3). Physically,ϕ\phiis related to the quantization SNR555defined as the ratio of signal variance to quantization noise varianceasϕ2=1/SNR\phi^{2}=1/\mathrm{SNR}. Numerical simulations show that, as a rule of thumb,ϕ∼0.1\phi\sim 0.1generally provides sufficient noise suppression for 4-bit-quantized PFB data.

The periodic nature of time- and frequency-domain noise spikes described in the previous section is apparent from the simulation results in Fig.5. Both the time- and frequency-domain spikes are broadened because inversion problems arise from a range of sub-filter eigenvalues that lie near zero. From (10), invertibility issues are expected when the response of any sub-filter drops (the region near top-center in Fig.3) below the threshold set by the quantization noise. Consequently, corruption appears in small regions near a frame center and near the boundaries of the original channels in the re-channelized spectrum. By Wiener filtering, power in the low SNR modes is suppressed, preventing noise amplification. Equivalently, from the time-domain decorrelation viewpoint, a Wiener filter shrinks and restricts the region of support of inverse kernels (plots labelled “filtered” kernels in Fig.4). The exact form of the Wiener filter, the optimal value ofϕ\phi, and an analytical expression forδ⁡(f)\delta(f)are derived in appendixB. The computational cost of the Wiener-filtered IPFB algorithm is compared with that of other filtering-based IPFB approaches in appendixC.

The Wiener filter is one statistical approach to down-weighting noisy eigenmodes. More generally, because the circulant formalism explicitly provides the eigenspectrum of the PFB transform matrix, arbitrary eigenvalue filtering schemes can be applied at negligible additional computational cost.

Figure 5:Numerically simulated comparison of percentage error in time-stream reconstruction and re-channelized power spectrum. The simulation PFB’d a white-noise timestream, quantized the spectra to 4-bits, and inverted the spectra using the two shown methods. The reconstructed timestream re-channelized using 64-times longer frame sizes for finer frequency resolution. Top plot shows percentage error in the frequency spectrum of re-channelized spectrum compared to channelization of the original timestream, averaged over many spectra. Middle plot shows reconstruction error in the timestream as compared to the original timestream. It can be seen that samples near frame centre are corrupted. Bottom plot zooms-in on a frame of filtered and unfiltered IPFB timestreams to highlight the different levels of reconstruction noise. The channels with 10% error can either be discarded, or the number of quantization bits can be increased, which will rapidly improve the re-construction error. (See appendixB.)

IV.2Maximum-Likelihood Optimization

IV.2.1Time Domain Prior

An alternate approach to spectral filtering is formulating timestream reconstruction as a maximum likelihood (ML) estimation problem. While more computationally intensive, ML-method can improve upon the Wiener filter results. Given some datadd, a linear modelAA, model parametersxx, and a covariant, positive-definite noise matrixNN, we write

χ2=(d−A​x)T​N−1​(d−A​x).\chi^{2}=(d-Ax)^{T}N^{-1}(d-Ax).(13)Minimizingχ2\chi^{2}yields the model parameters most likely to have produced the observed data. For the case of computing an IPFB,ddis the quantized, channelized data after the PFB is applied,xxis the unknown time-domain data, andAAis the PFB operator given by (5). Quantization noise is assumed to be uniform and uncorrelated such that the data noise covariance matrix,NN, is proportional to the identity matrix. This least-squares system can be solved using standard techniques like conjugate gradient (CG)(25)to obtain an estimate of the original timestream. In practice, an inverse-FFT (IFFT) of the PFB data is taken (undoingℱ\mathcal{F}) before solving the system, such thatA=𝒮​𝒲A=\mathcal{SW}, which makes theχ2\chi^{2}evaluation faster.

A direct least-squares solution to (13) inverts the near-zero singular values of the𝒮​𝒲\mathcal{SW}matrix, which, again, leads to the undesirable noise amplification. To better constrain the problem, a prior/regularization term must be included to fill in the attenuated singular modes. In §IV.1, prior information was used by specifying the PSD of the signal and the noise (assumed to be known a-priori), which led to the Wiener filter formalism. Alternately, the requirement of known PSDs can be sidestepped by specifying a time-domain prior. Since reconstruction noise is highest near the centre of a PFB frame, if a subset of the original timestream samples corresponding to the corrupted sub-filter paths are retained, they can be incorporated as a prior in the ML optimization as follows:

χ2=(d−A​x)T​N−1​(d−A​x)+(x−p)T​Q−1​(x−p).\chi^{2}=(d-Ax)^{T}N^{-1}(d-Ax)+(x-p)^{T}Q^{-1}(x-p).(14)The pre-saved time-stream samples are captured by the vectorpp, and the inverse noise matrix in the prior term,Q−1Q^{-1}, is 1 for time indices saved and 0 for indices discarded. The prior term then evaluates to∑i(xi−pi)2\sum_{i}(x_{i}-p_{i})^{2}, whereiiiterates over the saved timestream samples. Practical implementation of this approach requires modifying the digital logic for conditionally re-routing a fraction of raw timestream samples. Consequently, the choice of which samples to retain from each sub-filter path directly determines the additional digital bandwidth required by the hardware for enabling this algorithm.

IV.2.2Implementing the Prior

Consider a toy model of anNN-point,pp-tap PFB that uses a symmetric window and has no output quantization. For this PFB, the sub-sampled window array for filter pathN/2N/2,𝐰N/2\mathbf{w}_{N/2}, is of the form[a,b,b,a,0,0,…]T[a,b,b,a,0,0,\ldots]^{T}. Due to symmetry, the eigenvalue spectrum of this filter path is exactly zero at the Nyquist frequency (see appendixA). The correlation operation of the PFB, given by (6), therefore loses all information about the Nyquist mode of the sub-sampled timestream𝐱N/2\mathbf{x}_{N/2}. Decorrelating the Nyquist mode𝐱N/2\mathbf{x}_{N/2}is now ill-defined because it involves division by a zero eigenvalue. All other eigenmodes, however, can be successfully recovered.

Suppose we now take a pseudo-inverse, where the offending zero eigenvalue is truncated. The reconstructed𝐱^N/2\hat{\mathbf{x}}_{N/2}then differs from the original𝐱N/2\mathbf{x}_{N/2}only in the Nyquist frequency. In particular, the residual is𝐫=α​[−1,1,−1,1,…]T\mathbf{r}=\alpha[-1,1,-1,1,\ldots]^{T}, for someα∈ℝ.\alpha\in\mathbb{R}.Moreover,𝐫\mathbf{r}lies in the null-space ofWN/2W_{N/2}, so that moving in that direction of makes no difference inχ2\chi^{2}given by (14). This missing information (or degenerate direction) can be constrained with a prior. Including prior information in this example is straightforward. If we save just a single sample of the original timestream, we knowα\alpha, which in turn means we know𝐫\mathbf{r}. Subtracting𝐫\mathbf{r}from our reconstructed𝐱N/2\mathbf{x}_{N/2}then gives us the true𝐱N/2\mathbf{x}_{N/2}.

Figure 6:SNR (in dB) relative to the maximum possible SNR that can be obtained upon a naive “no filter” IPFB. For a given quantization scheme, the relative SNR is simply by1/|F​𝐰n|21/|F\mathbf{w}_{n}|^{2}. The abscissa spans the width of one original frequency channel. The two curves show the SNR behaviour for two different sub-filter paths for a 2048-point PFB. The width of the frequency band where SNR drops more than 6 dB is marked as corrupt, and used to guide the sample saving strategy for time-domain prior method described in sectionIV.2. Gray highlighted region illustrates corrupted region for the filter pathN/2N/2. Filter paths away fromN/2N/2have a smaller corrupted bandwidth.In a real PFB, the introduction of quantization noise corrupts a small but finite band of frequencies as shown in Fig.5. In this case, we need to save more than just one sample, but still only a tiny fraction of total samples. If the width of the corrupted band isΔ​ν\Delta\nu, the correlation length of the resulting residual timestream is of the order1/Δ​ν1/\Delta\nu. Thus, samples of𝐱N/2\mathbf{x}_{N/2}saved at a rate of roughlyΔ​ν\Delta\nucontain sufficient information to enable reconstruction of the corrupted modes. Instead of taking a pseudo-inverse and then solving for the missing modes, as in the toy example, the saved samples can now be used as a principled prior in the full-least squares solution to (14).

If one out of everyRRsamples is saved fromLLout ofNNsub-filter paths, the fraction of the original timestream saved is

β=LR​N\beta=\frac{L}{RN}(15)In our numerical experiments (described in the following section), we fixR=4R=4, i.e. save every 4th time-sample of a channel, and varyLLto obtain the desired level of saved fraction,β\beta. For guiding the eye, the choice ofR=4R=4corresponds to a corrupted bandwidth of approximately 20% for the sub-filter pathN/2N/2, as depicted in Fig.6. We emphasize, however, thatRRandLLare tunable parameters, and other, more efficient sample saving strategies are possible. For example, filter paths away fromN/2N/2have lower levels of noise corruption (smallerΔ​ν\Delta\nu), such that their samples can be saved at a slower rate.

IV.3Quantization Mitigation Results

Figure 7:Percentage error in re-channelized autospectrum for varying levels of time-domain samples saved. The simulation used here is identical to the one used for the Wiener-filtering results shown in Fig.5, and described in detail in sectionIV.3. For this plot, simulated 4-bit-quantized PFB data with 2048 channels was inverted with a conjugate gradient (CG)-based maximum likelihood optimizer, where the fraction of samples retained from the original timestream was varied from 1% to 10%. The IPFB timestream was re-channelized with 64-times finer channel resolution. The displayed bandwidth is one original channel or 64 fine channels wide. The plots are centered at the edge of two original channels to highlight the reconstruction error. With only 3% of timestream samples saved, the worst-case error drops to less than half of Wiener-filter error. Saving 10% of timestream samples further brings the error down to 2%. For comparison, the noise floor set by 4-bit quantization is approximately 1%. The secondary y-axis shows the effective number of bits (ENOB) obtained from effective quantization SNR (1/error) using the relation in (B1).The effects of quantization noise on the reconstructed timestream and the re-channelized spectrum, obtained using the two IPFB methods, are assessed using an end-to-end simulation consisting of the following steps. A digital timestream of length 227is generated from the standard normal distribution to simulate a signal incident upon an antenna. This timestream is passed through a 2048-point, 4-tap PFB, and the real/imaginary components of the resulting complex spectra are quantized using a 15-level (4 bits), mid-tread quantizer (consisting of uniformly-spaced levels from -7 to 7, with a 0 level.) The level-spacing for the quantizer was chosen to be 0.353 of the RMS of the channelized spectra as this value minimizes the mean-squared quantization error for a 4-bit uniform quantizer(18). The PFB is then inverted in three ways: circulant inversion without filtering, circulant inversion with Wiener filtering, and CG inversion with varying fractions of raw timestream samples saved. Finally, the reconstructed timestream is channelized again using a PFB with 64-times longer frame size, to produce 64-times finer frequency resolution. The fractional error in the rechannelized autospectrum, denoted byδ⁡(f)\delta(f), is calculated using (12), where the true spectra is obtained by directly channelizing the original timestream to the same, 64-times higher frequency resolution. Asδ⁡(f)\delta(f)is the ratio of jitter in power relative to the true power, the metric1/δ⁡(f)1/\delta(f)indicates the effective SNR achieved in a channel after rechannelization. To assess the quality of raw timestream reconstruction, fractional error in timestream is calculated asVar⁡[x−x^]/σx\sqrt{\mathrm{Var}[x-\hat{x}]}/\sigma_{x}, whereσx\sigma_{x}is the standard deviation of the original timestream. Both error metrics are estimated by averaging several thousand PFB frames.

Error in rechannelized autospectrum obtained after naive and Wiener filter IPFB are compared in the top plot of Fig.5. The result shows that near the center of an original channel, where all correlation kernels are well-behaved,δ⁡(f)\delta(f)is approximately 1.1%. In units of decibels, the relative error is10​log10​[δ⁡(f)]≈−19​dB10\mathrm{log_{10}}[\delta(f)]\approx-19\;\mathrm{dB}, which is consistent with the SNR ceiling of approximately 19 dB set by our 4-bit quantization scheme(13). Near the channel edges, whereF​𝐰nF\mathbf{w}_{n}drops precipitously for sub-filter paths close toN/2N/2, the Wiener filter is seen to suppress noise amplification. With Wiener filtering, the maximum value ofδ⁡(f)\delta(f)near the channel edges is approximately 10%, with the worst-affected region occupying less than 10% of the channel width. This error level also corresponds to salvaging 2.5 bits of information near the channel edges (see Fig.7). In comparison, the worst-case error is more than 100% without filtering, resulting in near-total loss of information for some sub-filter paths.

The fractional RMS error in timestream obtained using naive and Wiener filter IPFB is compared in the middle plot of Fig.5. In agreement with the autospectra results, Wiener filter decorrelation limits the noise amplification for samples near the center of the original PFB frame (highlighted in gray), which correspond to near-singular sub-filter paths. Most of the PFB frame samples have errors of about of 10%, which is the 4-bit quantization noise floor in amplitude.

Results obtained with CGχ2\chi^{2}minimization are displayed in Fig.7. The autospectrum error,δ⁡(f)\delta(f), is shown for four levels of saved fraction,β\beta, along with the Wiener filter and no-filter results for comparison. For eachβ\beta, the CG solver was provided the Wiener filter solution as the initial guess and was run for 25 iterations. The rechannelization error is seen to reduce smoothly, and approaches the quantization noise floor, as samples from more sub-filter paths are included in the prior term. ForL=4L=4andN=2048N=2048, an increase inβ\betafrom 1% to 10% corresponds to an increase in the number of saved paths from 82 to 820 (symmetric around center of the frame, 1024), where every 4th sample from each path is retained. Near the edges of the channel, the CG error is seen to drop to less than half of the Wiener filter solution with retention of only 3% of raw timestream samples. Equivalently, 3 out of 4 bits of information can be recovered in the worst-corrupted region by retaining 3% of timestream samples. With 10% of samples saved, the error reaches approximately 2%, or within 3 dB of the 4-bit quantization floor of approx. 1%.

To put the additional resource demand of ML method in perspective, the floating-point oversampled-PFB reconstruction reported by19usex an oversampling factor of 4/3, which requires 33% higher digital bandwidth than an equivalent critically-sampled PFB. For quantized PFBs, our sparse time-domain prior method achieves reconstruction to within 6 dB of the maximum achievable SNR with 4-bit quantization with only 3–5% of samples saved, which is a substantially lower requirement on additional bandwidth. We expect that more sophisticated sample saving schemes might enable comparable results with an even smaller fraction of samples saved, but we defer a full exploration of these schemes to future work.

VConclusion

This paper presents an efficient, Fourier-transform-based algorithm to invert critically sampled PFB measurements, enabling offline “up-channelization” of data to finer frequency resolution. By imposing circulant boundary conditions, the PFB inversion can be expressed as a set of independent decorrelations that are diagonal in the Fourier basis, providing a direct path to an efficient FFT-based implementation. This formalism also exposes a fundamental limitation: the problem is intrinsically ill-conditioned at frequencies where a polyphase sub-filter’s response approaches zero, and quantization noise is subsequently amplified strongly during inversion. The circulant formalism makes it straightforward to characterize the exact behaviour of this noise amplification in the time and frequency domains .

We explore two approaches for mitigating this noise amplification. The first applies a Wiener filter directly to the sub-filter eigenmodes in Fourier space, suppressing modes in which quantization noise dominates the reconstruction. For 4-bit quantized data, Wiener filtering reduces the worst re-channelized spectral errors from greater than 100% to approximately 10%, with significant errors confined to less than 10% of the original channel width. The reconstruction fidelity improves rapidly with increasing quantization depth. The Fourier-domain implementation allows Wiener filtering to be applied with negligible additional computational cost.

The second approach casts PFB inversion as a maximum-likelihood optimization problem, using a sparse prior of retained time-domain samples to constrain poorly measured modes. In simulations of 4-bit quantized data, retaining only 3% of the original timestream reduces the worst-case spectral error to less than half that obtained with Wiener filtering, and retaining 5–10% of samples further brings down to within 3 dB of the 4-bit noise floor. Although this approach requires additional data acquisition and greater computational effort compared to the FFT-based Wiener filter method, only a small amount of complementary time-domain information is needed to substantially improve reconstruction fidelity.

For conventionally recorded PFB data, the Wiener-filtered inverse provides a robust and computationally efficient approach to inversion/up-channelization while controlling the noise amplification inherent to the inverse problem. In a future paper, we will demonstrate the performance of a scalable GPU implementation of the Wiener filter IPFB algorithm, which can invert many days of data collected with the ALBATROS radio telescope in a continuous, streaming fashion. In applications where greater reconstruction fidelity is required, our maximum-likelihood results demonstrate that retaining only a small fraction of the original timestream can recover much of the information lost to poorly constrained modes. Together, these results provide a practical framework for understanding and mitigating the limitations of critically sampled PFB inversion, enabling finer spectral resolution to be recovered from channelized radio-astronomical data.

We acknowledge the support of the Natural Sciences and Engineering Research Council of Canada (NSERC), RGPIN-2019-04506, RGPNS 534549-19; Canada Foundation for Innovation John R. Evans Leaders Fund #40824; National Geographic Society Explorer Grant NGS-94983T-22. We acknowledge the support of the Government of Canada’s New Frontiers in Research Fund (NFRF) NFRFE-2023-00069. This research was undertaken, in part, thanks to funding from the Canada 150 Research Chairs Program. MA acknowledges support from the Vanier Canada Graduate Scholarship Program.

{contribution}

SF wrote the initial version of Wiener filter and CG solvers for inverse PFB, investigated their performance and behavior, created the visualizations, and wrote an initial draft of the manuscript. MA conducted the literature review, performed formal mathematical analysis of the IPFB methods, analyzed and validated simulation results, and co-authored this manuscript (including the appendices). He also contributed to the development of IPFB methods and ideas for visualization. HCC reviewed and edited the manuscript, provided feedback on the visualizations, and supervised the presentation of technical aspects of the IPFB methods. AC contributed in creating and refining the visualizations, participated in technical discussions on IPFB methods, and collaborated on writing the manuscript. JS conceived the initial research idea, contributed to the development of IPFB methods and interpretation of results, reviewed and edited the manuscript, and supervised the project from its conception to completion. ST conducted initial investigation into possible approaches for inverting PFBs, reviewed the manuscript, and provided feedback.

Appendix AEigenvalues of a symmetric PFB window

The singularities encountered in sub-filter responses shown in Fig.3arise, fundamentally, due to symmetry properties of the window function used in the PFB. Here, we provide a heuristic proof of why any reasonable choice of window function will always lead to singular values in the critically-sampled PFB operator. We require a sensible window arrayWWof lengthMMto be: (1) symmetric about its center, such thatWk=WM−kW_{k}=W_{M-k}, and (2) go to zero on the boundary, so thatW0=WM−1=0W_{0}=W_{M-1}=0. The latter requirement is reasonable because windows are designed to suppress spectral leakage by gently approaching zero towards the boundaries of a finite timestream.

Consider aNN-point,pp-tap PFB withM=p​NM=pNconstructed usingWW, and letppbe an even number. From (8), a PFB correlates a sub-sampled input timestream with the sub-sampled window array,𝐰n\mathbf{w}_{n}, which comprises onlyppnon-zero entries. For the sub-filter pathN/2N/2, theseppentries are exactly symmetric. For example, withp=4p=4taps,𝐰N/2\mathbf{w}_{N/2}is of the form[a,b,b,a,0,0,…][a,b,b,a,0,0,\ldots]due to the symmetry requirement. To evaluate the the eigenvalues (or frequency response) of this sub-filter, aMM-point Fourier transform,F​𝐰nF\mathbf{w}_{n}, is performed. The eigenvalue corresponding to the Nyquist frequency is then seen to be a dot product of[a,b,b,a,0,0,…][a,b,b,a,0,0,\ldots]with the alternating sequence[1,−1,1,−1,…][1,-1,1,-1,\ldots]; the latter being the DFT phasorej​π​me^{j\pi m}for the Nyquist frequency. This dot-product evidently vanishes, corresponding to the point of singularity at the top-center of the plot in Fig.3. For an odd number of taps, it is the eigenvalue of path00, instead of theN/2N/2, that vanishes. BecauseW0=0W_{0}=0, the dot product of remainingp−1p-1elements of𝐰0\mathbf{w}_{0}evaluates to zero due to symmetry666Consider a 5-tap PFB as a concrete example. The window will be symmetric around sample5​N/25N/2whereNNis the number of samples per frame. Path 0 will consist of samples at[0,N,2​N,3​N,4​N][0,N,2N,3N,4N]. If the window goes to 0 at the edge then sample 0 will be zero. The symmetry around5​N/25N/2means the window at samples2​N2Nand3​N3Nare the same, and that the window at samplesNNand4​N4Nwill be the same. That leaves𝐰𝟎=[0,a,b,b,a,0​…]\mathbf{w_{0}}=[0,a,b,b,a,0...], which will Fourier transform to zero at the Nyquist frequency..

The more closely a window adheres to the two requirements of symmetry and zero at the boundaries, the closer the smallest PFB eigenvalue is to zero. Consequently, inverse PFB suffers a correspondingly higher noise amplification. Common windows implemented in theNumpypackage, such as Hamming, Hanning, and Blackman are not exactly symmetric inMMdue to the normalization choice ofc​o​s​[2​π​m/(M−1)]cos[2\pi m/(M-1)], but they approach perfect symmetry with increasingMM. The IPFB noise amplification problem is, therefore, worse for PFBs that use high FFT frame sizes,NN, or too many taps,pp.

Appendix BOptimal Wiener Threshold

A radio telescope’s signal processing chain quantizes the real and imaginary parts of the output of a PFB. In what follows, letssbe the real (or imaginary) component of a PFB channel and𝒬⁡(s)\mathcal{Q}(s)its quantized version. In an optimally-tuned quantizer, all quantization levels are well-exercised and the clipped fraction is kept low. For such a quantizer, the RMS of quantization noise,σq\sigma_{q}, defined asσq=Var⁡[s−𝒬⁡(s)]\sigma_{q}=\sqrt{\mathrm{Var}[s-\mathcal{Q}(s)]}is approximatelyΔ/12\Delta/\sqrt{12}, whereΔ\Deltais the level spacing or one least-signifcant bit (LSB) of the quantizer. For the example of 4-bit (15 level, symmetric) quantization used in sectionIV.2, the signal777by “signal” we mean antenna timestream. Noise refers specifically to quantization noise, and all references to SNR refer to quantization SNR in this section.RMSσs\sigma_{s}andΔ\Deltaare related asΔ=0.353​σs\Delta=0.353\sigma_{s}in the optimal case. The signal-to-noise ratio (SNR), defined asσs2/σq2\sigma_{s}^{2}/\sigma_{q}^{2}, achieved with this tuning is approximately 80. In decibels,SNRdB=10​log10​[δ⁡(f)]\mathrm{SNR_{dB}}=10\mathrm{log_{10}}[\delta(f)], the optimal SNR is 19 dB, and for a given signal level, it increases by approximately 6 dB per additional bit. Another useful metric assessing the SNR achieved by a given signal processing system is effective number of bits (ENOB), which quantifies SNR in terms of bits of information retained or lost. For example, an ENOB of 3 implies that the SNR of the system is equivalent to that achieved by an optimally-tuned 3 bit digitizer. To convert between quantization SNR and ENOB, we use the following relation

SNRdB=5.2​R−1.8,\mathrm{SNR_{dB}}=5.2\mathrm{R}-1.8,(B1)where,R=l​o​g2​Nb​i​t​sR=log_{2}N_{bits}. This approximate scaling relation was obtained fitting a line to the results obtained for uniform quantization of a gaussian signal by13, displayed in their figure 2. In §IVin the main text, ENOB is used to asses the amount of information successfully recovered by various IPFB methods in the worst-corrupted portions of the rechannelized spectrum.

After taking an N-point IFFT of the PFB spectra, the decorrelation problem for a sub-filter pathnn(n<Nn<N) is𝐝n=Wncirc​𝐱n+ϵn\mathbf{d}_{n}=W_{n}^{\mathrm{circ}}\mathbf{x}_{n}+\boldsymbol{\epsilon}_{n}, whereϵn\boldsymbol{\epsilon}_{n}is the original quantization noise scaled by the IFFT and window normalization factors. Since the scaling factors are common to both the quantization noise,ϵn\boldsymbol{\epsilon}_{n}, and the signal,𝐱n\mathbf{x}_{n}, the SNR in time-domain does not change. At this stage, Wiener filtering𝐝n\mathbf{d}_{n}would suppress power in the noisy eigenmodes, preventing noise amplification upon decorrelation. To obtain the exact expression in (11), let the signal power spectral density (PSD) beS⁡(f)=|Wn​(f)|2​⟨|Xn​(f)|2⟩S(f)=|W_{n}(f)|^{2}\langle|X_{n}(f)|^{2}\rangle, and the noise PSD beN⁡(f)=⟨|Qn​(f)|2⟩N(f)=\langle|Q_{n}(f)|^{2}\rangle. HereWn​(f)W_{n}(f),Xn​(f)X_{n}(f), andQn​(f)Q_{n}(f)are the Fourier transforms of thenn’th sub-filter, thenn’th sub-sampled timestream, and the quantization noise in that timestream, respectively. The frequency variableffthroughout this section is in units of cycles/sample (ork/Mk/M, whereMMis the number of samples being Fourier transformed).

We assume that quantization noise in time samples is uncorrelated, such thatN⁡(f)N(f)is flat. Moreover, if the PFB channels are narrow enough that the PSD of the underlying signal is approximately flat within one channel, then⟨|Xn​(f)|2⟩\langle|X_{n}(f)|^{2}\ranglewill also be approximately flat888Sub-sampling the original timestreamxxby a factorNNaliasesNNNyquist zones, eachB/nc​h​a​nB/n_{chan}-wide, on top of each other, realizing an approximately flat spectrum.. Flat PSDs of signal and noise are proportional to their variances, which allows the Wiener filter,S/(S+N)S/(S+N), for sub-filter pathnnto be written exactly as

Hn​(f)=|Wn​(f)|2​σs2|Wn​(f)|2​σs2+σn2=|Wn​(f)|2|Wn​(f)|2+1/SNR\displaystyle\begin{split}H_{n}(f)&=\frac{|W_{n}(f)|^{2}\sigma_{s}^{2}}{|W_{n}(f)|^{2}\sigma_{s}^{2}+\sigma_{n}^{2}}\\ &=\frac{|W_{n}(f)|^{2}}{|W_{n}(f)|^{2}+1/\mathrm{SNR}}\\ \end{split}(B2)whereσn\sigma_{n}is RMS of quantization noise. The dimensionless parameterϕ2\phi^{2}from sectionIV.1can now be recognized as1/SNR1/\text{SNR}. For the 4-bit example, the optimal SNR of 19 dB corresponds toϕ≈0.1\phi\approx 0.1. A high value ofϕ\phiforces the Wiener filter to use an artificially low SNR, causing signal suppression near channel edges. However, we have found in our numerical experiments that the filtered spectrum is not overly sensitive to the exact choice ofϕ\phi. Optionally, in (11), we have approximately normalized the filter response by noting thatWn​(f)∼1W_{n}(f)\sim 1for most sub-filters, such that the integral ofHn​(f)H_{n}(f)over[−0.5,0.5][-0.5,0.5]is roughly1/(1+ϕ2)1/(1+\phi^{2}).

Wiener filtering𝐝n\mathbf{d}_{n}followed by decorrelation is equivalent to decorrelation with an effective kernelHn​(f)/Wn∗​(f)H_{n}(f)/W_{n}^{*}(f). With the Wiener filter specified by (B2), the decorrelated estimate (denoted with a hat) ofnn’th timestream’s spectrum can be written as

X^n​(f)=Hn​(f)​Xn​(f)+Hn​(f)Wn∗​(f)​Qn​(f).\hat{X}_{n}(f)=H_{n}(f)X_{n}(f)+\frac{H_{n}(f)}{W_{n}^{*}(f)}Q_{n}(f).(B3)However, of greater practical interest is the spectrum of the reconstituted full-rate timestream, and the bias introduced in its PSD due to quantization noise and Wiener filtering. Conceptually, the full-rate timestream can be obtained by upsampling all𝐱n\mathbf{x}_{n}’s (appendingN−1N-1zeros between their samples), delaying each upsampled𝐱n\mathbf{x}_{n}bynnsamples, and then adding them together. This convertsx0,0,…x_{0},0,...andx1,0,..x_{1},0,..intox0,x1,…x_{0},x_{1},.... Following this process, the spectrum of the full-rate timestream can be written in terms of its sub-sampled components as

X⁡(f)=∑n=0N−1Xn​(N​f)​e−j​2​π​N​f​n.X(f)=\sum_{n=0}^{N-1}X_{n}(Nf)e^{-j2\pi Nfn}.(B4)The process of upsampling compresses the original spectrum and creates images, so that full-rate spectrum containsNNreplicas of the low-rate spectrum in[−0.5,0.5][-0.5,0.5]. A proof of this identity can be found in standard texts on multi-rate signal processing(29;22). The relative bias or error in the power spectrum of filtered data is then defined as:

δ⁡(f)=⟨|X⁡(f)−X^​(f)|2⟩⟨|X⁡(f)|2⟩.\delta(f)=\frac{\langle|X(f)-\hat{X}(f)|^{2}\rangle}{\langle|X(f)|^{2}\rangle}.(B5)Using (B2), (B3), and (B4), the bias can be evaluated easily for the case of white input signal. In this case, the spectrum of two sub-filter paths are uncorrelated, such that⟨Xn​(f)​Xm​(f)⟩∝σs2​δm​n\langle X_{n}(f)X_{m}(f)\rangle\propto\sigma_{s}^{2}\delta_{mn}. For the case of naive decorrelation, the relative error is

δ⁡(f)=1SNR​1N​∑n=0N−11|Wn​(N​f)|2,\delta(f)=\frac{1}{\mathrm{SNR}}\frac{1}{N}\sum_{n=0}^{N-1}\frac{1}{|W_{n}(Nf)|^{2}},(B6)which shows that the error will increase sharply if one of the sub-filters’ response approaches zero. This is illustrated in Fig.6. Similarly, relative error in the case of Wiener filtered decorrelation is

δ⁡(f)=1−1N​∑n=0N−1Hn​(N​f).\delta(f)=1-\frac{1}{N}\sum_{n=0}^{N-1}H_{n}(Nf).(B7)Evaluation ofδ⁡(f)\delta(f)in practice involves summation of allNNfilter responses, followed byNNrepeats of the summed spectrum, which spans the sampling rate2​B2Bof the original timestream. The expression in (B7) is the analytical estimate in Fig.5, which lines up well with the Monte-Carlo estimate. It is worth mentioning here that the task of calculating the full-rate spectrum from spectra ofNNpolyphase branches is mathematically identical to calculating the spectrum of interleaved analog-to-digital converters(9). Finally, near the center of any PFB channel, which corresponds to the region aroundHn​(f=0)H_{n}(f=0)or the bottom edge of Fig.3,Hn∼1H_{n}\sim 1, such that (B7) can be simplified to

δ⁡(f)|N​fmodN=0=11+SNR.\delta(f)|_{Nf\bmod{N}=0}=\frac{1}{1+\mathrm{SNR}}.(B8)This simplified expression agrees with the expectation that the relative error in the well-behaved region of is set by maximum achievable quantization SNR, provided the assumptions about signal and noise PSD continue to hold. SinceSNR≫1\mathrm{SNR}\gg 1, the error in majority of the re-channelized spectra goes as∼SNR−1\sim\mathrm{SNR}^{-1}, whereSNR\mathrm{SNR}increases approximately exponentially with the number of bits. Near the edges of original channels, the exact level ofSNR\mathrm{SNR}degradation can be computed from (B7).

Appendix CCirculant Inverse PFB Time Complexity

Using the circulant algorithm to invert a chunk ofpp-tapNN-point PFB spectra of shape (M spectra and N/2 unique channels) requires three stages of FFTs. First stage appliesMMindependent N-point FFTs to invertℱ\mathcal{F}, costing𝒪⁡(M​N​l​o​g​N)\mathcal{O}(MNlogN). This step producesNNcorrelated timestreams𝐝n\mathbf{d}_{n}’s, each of sizeMM. The subsequent decorrelation of these timestreams in Fourier space, according to equation (9), requires two stages ofNNindependentMM-point FFTs, with an optional point-wise multiplication step for Wiener filtering, costing roughly𝒪⁡(N​M​l​o​g​M)\mathcal{O}(NMlogM). Thus, the total computational cost for the circulant inversion scales as𝒪⁡(M​N​l​o​g​M​N)\mathcal{O}(MNlogMN).

The Wiener-filter decorrelation step can also be implemented by passing each of the𝐝n\mathbf{d}_{n}’s through a finite impulse response (FIR) filter, where the impulse response of the FIR filter is given by the “decorrelation kernels” shown in Fig.4). Cost of FIR filtering one of the sub-sampled timestreams scales as𝒪⁡(M​L)\mathcal{O}(ML)whereLLis the length of the decorrelation kernel. The total cost of FIR filtering based IPFB is, therefore,𝒪⁡(M​N​max​(l​o​g​N,L))\mathcal{O}(MN\mathrm{max}(logN,L)). From Fig.4, the length of the longest required FIR kernel (corresponding to sub-filter pathN/2N/2) is approximately 50 samples for the example case of a 4-tap PFB with 4-bit output quantization (corresponding toϕ∼0.1\phi\sim 0.1). This length increases steadily with higher taps or higher bit-depth. In a typical radio astronomy use-cases, PFB frame sizes,NN, are in the range210−2142^{10}-2^{14}, and the number PFB frames being inverted in a single function call,MM, might be in the range104−10610^{4}-10^{6}, such thatl​o​g​(M​N)<Llog(MN)<L. Thus, FIR filtering will be slower than the FFT-method by a factor of a few, unlessLLis artificially reduced by choosing a high value ofϕ\phi.

Finally, the Wiener-filter machinery can be replicated using a pseudo-inverse solver with Tikhonov regularization for noise suppression. A well-optimized solver can, in principle, solve equation (6) in𝒪⁡(M​p2)\mathcal{O}(Mp^{2})time using banded Cholesky decomposition, for a total cost of𝒪⁡(M​N​p2)\mathcal{O}(MNp^{2}). Withp=4p=4taps,log2​M∼p2\mathrm{log_{2}}M\sim p^{2}, such that the circulant algorithm and the Tikhonov solver have similar computational cost. However, the circulant algorithm offers two important benefits. First, its cost does not scale with the number of taps. Due to itsp2p^{2}scaling, the Tikhonov solver becomes quickly unwieldy for other common configurations likep=8p=8or1616, and prohibitively expensive for PFB configurations found in telecommunication sector, where the number of taps can be as high as 128. Second, the circulant algorithm allowsanyeigenvalue filtering scheme, not just Wiener filtering. In contrast, custom filtering of singular values of equation (6) requires a full SVD that costs𝒪⁡(M3)\mathcal{O}(M^{3}), making it intractable for inverting large amounts of PFB data.

References

• L. H. Arnaldi and H. DellavaleOversampled filter bank channelizer for cryogenic detectors.Review of Scientific Instruments92(2).Cited by:§I.
• Beadleet al.(2024)C. Beadle, A. Caputo, and S. A. EllisResonant conversion of wave dark matter in the ionosphere.Physical Review Letters133(25),pp. 251001.Cited by:§I.
• Bellangeret al.(1976)M. Bellanger, G. Bonnerot, and M. CoudreuseDigital filtering by polyphase network: application to sample-rate alteration and filter banks.IEEE Transactions on acoustics, speech, and signal processing24(2),pp. 109–114.Cited by:§I.
• Bridle and Schwab (1999)A. H. Bridle and F. SchwabBandwidth and time-average smearing.InSynthesis Imaging in Radio Astronomy II,Vol.180,pp. 371.Cited by:§I.
• Chianget al.(2020)H. C. Chiang, T. Dyson, E. Egan, S. Eyono, N. Ghazi, J. Hickish, J. M. Jauregui-Garcia, V. Manukha, T. Menard, T. Moso, J. Peterson, L. Philip, J. L. Sievers, and S. TartakovskyThe array of long baseline antennas for taking radio observations from the sub-antarctic.External Links:2008.12208,LinkCited by:§I.
• Choet al.(2020)H. Cho, J. Macquart, R. M. Shannon, A. T. Deller, I. S. Morrison, R. D. Ekers, K. W. Bannister, W. Farah, H. Qiu, M. W. Sammons,et al.Spectropolarimetric analysis of frb 181112 at microsecond resolution: implications for fast radio burst emission mechanism.The Astrophysical Journal Letters891(2),pp. L38.Cited by:§I.
• Crochiere and Rabiner (2005)R. E. Crochiere and L. R. RabinerInterpolation and decimation of digital signals—a tutorial review.Proceedings of the IEEE69(3),pp. 300–331.Cited by:§I.
• Demir and Bjornson (2020)O. T. Demir and E. BjornsonThe bussgang decomposition of nonlinear systems: basic theory and mimo extensions [lecture notes].IEEE Signal Processing Magazine38(1),pp. 131–136.Cited by:§III.2.
• El-Chammas (2010)M. El-ChammasBackground calibration of timing skew in time-interleaved a/d converters.Stanford University.Cited by:Appendix B.
• Fletcheret al.(2007)A. K. Fletcher, S. Rangan, V. K. Goyal, and K. RamchandranRobust predictive quantization: analysis and design via convex optimization.IEEE Journal of selected topics in signal processing1(4),pp. 618–632.Cited by:§III.2.
• Harriset al.(2020)C. R. Harris, K. J. Millman, S. J. van der Walt, R. Gommers, P. Virtanen, D. Cournapeau, E. Wieser, J. Taylor, S. Berg, N. J. Smith, R. Kern, M. Picus, S. Hoyer, M. H. van Kerkwijk, M. Brett, A. Haldane, J. F. del Río, M. Wiebe, P. Peterson, P. Gérard-Marchant, K. Sheppard, T. Reddy, W. Weckesser, H. Abbasi, C. Gohlke, and T. E. OliphantArray programming with NumPy.Nature585(7825),pp. 357–362.External Links:Document,LinkCited by:Inverting and Up-channelizing Critically-Sampled Polyphase Filter Banks.
• Harriset al.(2003)F. J. Harris, C. Dick, and M. RiceDigital receivers and transmitters using polyphase filter banks for wireless communications.IEEE transactions on microwave theory and techniques51(4),pp. 1395–1412.Cited by:§I.
• Hui and Neuhoff (2001)D. Hui and D. L. NeuhoffAsymptotic analysis of optimal fixed-rate uniform scalar quantization.IEEE Transactions on Information Theory47(3),pp. 957–977.Cited by:Appendix B,§IV.3.
• Hunter (2007)J. D. HunterMatplotlib: a 2d graphics environment.Computing in science & engineering9(3),pp. 90–95.Cited by:Inverting and Up-channelizing Critically-Sampled Polyphase Filter Banks.
• Leunget al.(2025)C. Leung, S. Andrew, K. W. Masui, C. Brar, T. Cassanelli, S. Chatterjee, V. Kaspi, K. Khairy, A. E. Lanman, M. Lazda,et al.A vlbi software correlator for fast radio transients.The Astronomical Journal170(1),pp. 53.Cited by:§I.
• Lin and Vaidyanathan (1995)Y. Lin and P. VaidyanathanLinear phase cosine modulated maximally decimated filter banks with perfect reconstruction.IEEE Transactions on Signal Processing43(11),pp. 2525–2539.Cited by:§I.
• Malvar (1990)H. S. MalvarModulated qmf filter banks with perfect reconstruction.Electronics letters26(13),pp. 906–907.Cited by:§I.
• Max (1960)J. MaxQuantizing for minimum distortion.IRE Transactions on Information Theory6(1),pp. 7–12.Cited by:§IV.3.
• Morrisonet al.(2020)I. Morrison, J. Bunton, W. van Straten, A. Deller, and A. JamesonPerformance of oversampled polyphase filterbank inversion via fourier transform.Journal of Astronomical Instrumentation9(01),pp. 2050004.Cited by:§I,§IV.3.
• Oppenheim (1999)A. V. OppenheimDiscrete-time signal processing.Pearson Education India.Cited by:§I.
• Price (2021)D. C. PriceSpectrometers and polyphase filterbanks in radio astronomy.InThe WSPC Handbook of Astronomical Instrumentation: Volume 1: Radio Astronomical Instrumentation,pp. 159–179.Cited by:§I.
• Rabiner (1996)L. R. RabinerMultirate digital signal processing.Prentice Hall PTR.Cited by:Appendix B,§II.2.
• Rybicki and Press (1992)G. B. Rybicki and W. H. PressInterpolation, realization, and reconstruction of noisy, irregularly sampled data.Astrophysical Journal, Part 1 (ISSN 0004-637X), vol. 398, no. 1, Oct. 10, 1992, p. 169-176.398,pp. 169–176.Cited by:§IV.1.
• Shaw (2014)J. R. ShawPfb-inverse(Website)External Links:LinkCited by:§III.
• Shewchuk (1994)J. R. ShewchukAn introduction to the conjugate gradient method without the agonizing pain.Carnegie-Mellon University. Department of Computer Science.Cited by:§IV.2.1.
• Smith and Barnwell (1984)M. Smith and T. BarnwellA procedure for designing exact reconstruction filter banks for tree-structured subband coders.InICASSP’84. IEEE International Conference on Acoustics, Speech, and Signal Processing,Vol.9,pp. 421–424.Cited by:§I.
• Tuthillet al.(2015)J. Tuthill, G. Hampson, J. D. Bunton, F. Harris, A. Brown, R. Ferris, and T. BatemanCompensating for oversampling effects in polyphase channelizers: a radio astronomy application.In2015 IEEE Signal Processing and Signal Processing Education Workshop (SP/SPE),pp. 255–260.Cited by:§I.
• Tuthillet al.(2012)J. Tuthill, G. Hampson, J. Bunton, A. Brown, S. Neuhold, T. Bateman, L. de Souza, and J. JosephDevelopment of multi-stage filter banks for askap.In2012 International Conference on Electromagnetics in Advanced Applications,pp. 1067–1070.Cited by:§I.
• Vaidyanathan (2002)P. P. VaidyanathanMultirate digital filters, filter banks, polyphase networks, and applications: a tutorial.Proceedings of the IEEE78(1),pp. 56–93.Cited by:Appendix B,§I,§I.
• Vetterli and Le Gall (2002)M. Vetterli and D. Le GallPerfect reconstruction fir filter banks: some properties and factorizations.IEEE Transactions on Acoustics, Speech, and Signal Processing37(7),pp. 1057–1071.Cited by:§I.
• Virtanenet al.(2020)P. Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. J. van der Walt, M. Brett, J. Wilson, K. J. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C. J. Carey, İ. Polat, Y. Feng, E. W. Moore, J. VanderPlas, D. Laxalde, J. Perktold, R. Cimrman, I. Henriksen, E. A. Quintero, C. R. Harris, A. M. Archibald, A. H. Ribeiro, F. Pedregosa, P. van Mulbregt, and SciPy 1.0 ContributorsSciPy 1.0: Fundamental Algorithms for Scientific Computing in Python.Nature Methods17,pp. 261–272.External Links:DocumentCited by:Inverting and Up-channelizing Critically-Sampled Polyphase Filter Banks.
• Wiener (1949)N. WienerExtrapolation, interpolation, and smoothing of stationary time series: with engineering applications.The MIT press.Cited by:§IV.1.