CASCADED ALL-PASS FILTERS WITH RANDOMIZED CENTER FREQUENCIES AND PHASE POLARITY FOR ACOUSTIC AND SPEECH MEASUREMENT AND
DATA AUGMENTATION

Let us start from its building block. The following equation defines the transfer function of a simple discrete-time all-pass filter, characterized by the center frequency $f_{k}$ and the bandwidth $b_{k}$[27].

	$\displaystyle H_{k}(z)$	$\displaystyle=\frac{z^{-1}-z_{k}^{\ast}}{1-z_{k}^{\ast}z^{-1}},\ \ \ \ \left(\mbox{Note:}\ \ H_{k}^{\ast}(z)H_{k}(z)=1\right)$		(1)
	$\displaystyle z_{k}$	$\displaystyle=\exp\left(-\frac{\pi b_{k}}{f_{s}}+j\frac{2\pi f_{k}}{f_{s}}\right),$		(2)

where $f_{s}$ represents the sampling frequency and $j=\sqrt{-1}$. A term $a^{\ast}$ represents the complex conjugate of $a$. By cascading these building blocks also yields an all-pass filter. Its inverse z-transform provides an impulse response. We call the impulse response unit-CAPRICEP afterward.

The transfer function of the cascaded all-pass filters is the product of constituent transfer function $H_{k}(z)$.

$\displaystyle H_{C}(z)$

$\displaystyle=\prod_{k=1}^{K}H_{k}(z),$

(3)

where $K$ represents the total number of the cascaded filters. A systematic procedure determines each center frequency $f_{k}$ using two set of random numbers $r_{1}[k]$ and $r_{2}[k]$. We introduce random numbers because each pole corresponds to a damped sinusoid. Intuitively speaking, the sum of many sinusoids with the same amplitude with the random phase yields items behave like samples from normal distribution[28]. We hope frequency randomization yields noise-like behavior of the impulse response of $H_{C}(z)$.

The following equation is an implementation of the systematic procedure of filter assignment.

$\displaystyle f_{k}$

$\displaystyle=r_{1}[k]\,F_{d}\sum_{k=1}^{K}r_{2}[k],$

(4)

where $F_{d}$ represents the average interval between aligned center frequencies (i.e. $f_{k+1}>f_{k}$). The random number $r_{1}[k]$ has the value $1$ and $-1$ with equal probability. Note that $r_{1}[k]=-1$ makes $H_{k}(z)$ anti-causal and makes its impulse response time-reversed. This randomization is for design flexibility described in the next section. Without introducing $r_{1}[k]$, the response is causal and does not provide flexible design.

The random number $r_{2}[k]$ obeys a probability density distribution function (pdf) $g(x)$ defined on $[0,1]$. We do not randomize $f_{k}$ directly. Instead, we randomized the frequency interval between neighboring all-pass filters. References[19, 20] reported simple random allocation does not yield well-behaved velvet noise²²2 They named the noise “velvet noise” because the noise sounds smoother than Gaussian white noise. Simple random allocation destroyed the smoothness[19, 20]. . Intuitively speaking, the original velvet noise and FVN use the triangular pdf for interval distribution. They also use the additional constraint, fixed-width partitioning of the time axes. This fixed-width partitioning forced us to introduce non-linear frequency axis warping for shaping temporal distribution[5].

We introduce the Beta distribution for $r_{2}[k]$ because it is a simple adjustable distribution having two design parameters. The following is the definition of the Beta distribution.

$\displaystyle b(x)$

$\displaystyle=\frac{x^{\alpha-1}(1-x)^{\beta-1}\Gamma(\alpha+\beta)}{\Gamma(\alpha)\Gamma(\beta)},$

(5)

where $\Gamma(\alpha)$ is Gamma function. When $\alpha\equiv\beta$, $b(x)$ seamlessly deforms from concave to convex via uniform distribution with $\alpha$.

For $b_{k}$ (Eq. 2), we make it proportional to $F_{d}$ using a coefficient $c_{\mathrm{mag}}$. In short, we set $b_{k}=c_{\mathrm{mag}}F_{d}$ (at least for the first time). In summary, $H_{C}(z)$ has three design parameters, $F_{d},c_{\mathrm{mag}}$, and $\alpha$. This additional parameter $\alpha$ that FVN does not have enables us to design the shape of the waveform distribution. Let $h_{c}[n]$ represent the impulse response of $H_{C}(z)$. It is unit-CAPRICEP.

Shaping waveform distribution is an optimization problem. Let us define the target shape on the discrete-time axis $W[n]$. We design the sample variance $V(h_{c}[n])$ by minimizing a distance measure between $V(h_{c}[n])$ and $W[n]$. We use Wasserstein distance[29] as the cost function to minimize because it is an energy allocation problem on the discrete-time axis, in the other words.

We selected the rectangular shape as the first target because it is suitable for making less-fluctuating test signals using periodic allocation (i.e. OLA: overlap and add) of unit-CAPRICEPs. We conducted several coarse multivariate optimizations before final optimization. It gave us the setting $c_{\mathrm{mag}}=2^{\frac{1}{4}}$ and $\alpha=8$. Then, we optimized the duration of the rectangle $T_{\mathrm{ERD}}$ using a grid search.

Figure 1 shows an example of shape design. The sample variance calculation used 5000 unit-CAPRICEPs. The optimized $T_{\mathrm{ERD}}$ is 1.736 times of the nominal duration of $1/F_{d}$. Note that the decay outside of the rectangle is exponential. The total number of cascaded filters is 1102 for $f_{s}=44100$ Hz, $F_{d}=40$ Hz.

Left plot of Fig. 2 shows the distribution of the maximum cross correlation (absolute value) of $400\times 399$ pairs of different unit-CAPRICEPs. The thick line shows probability $p(x_{c}|\ \ |x_{c}|\leq\theta_{c})$ where $x_{c}$ represents the maximum cross correlation and $\theta_{c}$ represents the threshold correlation value shown in the horizontal axis. The thin line shows $p(x_{c}|\ \ |x_{c}|>\theta_{c})$. The median of the distribution is 0.0905 and close to zero.

Another useful target shape is the raised cosine because OLA with 50% overlap provides a constant level. Cascading a short (0.5 ms) raised cosine-shaped unit-CAPRICEP and the rectangular one yielded the shape shown in the right plot of Fig. 2. The outstanding spike at the center shown in Fig. 1 disappears. This composite shape is appropriate for acoustic measurement applications.

CAPRICEP fulfilled the design goal as follows: a) Using IIR all-pass filter as the building block made the decaying exponential that enables well-behaving localization (section 2.1). b) Cascading many filters provides many degrees of freedom (section 2.2 and 2.3). c) Many degrees of freedom made cross-correlation between different unit-CAPRICEPs close to zero (section 2.3).

Note that making $T_{\mathrm{ERD}}$ frequency-dependent is straightforward for CAPRICEP, while complicated and indirect for FVN[5]. There is much additional flexibility in design using $\beta$ and $b_{k}$, which are fix-valued in the current implementation.

A set of four periodic sequences made from unit-CAPRICEPs enables the simultaneous measurement. Define a weight matrix $\mathrm{B}_{4}$, which consists of four orthogonal binary rows shown below.

	$\displaystyle\mathrm{B}_{4}$	$\displaystyle=\left[\begin{array}[]{rrrrrrrr}1&1&1&1&1&1&1&1\cr 1&-1&1&-1&1&-1&1&-1\cr 1&1&-1&-1&1&1&-1&-1\cr 1&1&1&1&-1&-1&-1&-1\end{array}\right]$		(10)
		$\displaystyle=[\mathbf{b}^{(1)}\mathbf{b}^{(2)}\mathbf{b}^{(3)}\mathbf{b}^{(4)}]^{T}.$		(11)

The following equation defines the $m$-th CAPRICEP sequence $s_{4}^{(m)}[n],$ $(m=1,\ldots,4)$. This procedure is a periodic OLA.

$\displaystyle s_{4}^{(m)}[n]$

$\displaystyle=\sum_{k=0}^{K-1}h_{c}^{(m)}[n-kn_{\mathrm{o}}]\,b^{(m)}_{\mathrm{mod}(k,8)+1}\ ,$

(12)

where $b_{j}^{(m)}$ represents the $j$-th element of $\mathbf{b}^{(m)}$. $K$ represents the number of repetitions. Note that 8 is the number of columns of $\mathrm{B}_{4}$.

Make a test $x_{\mathrm{testR}}[n]$ by adding the first three CAPRICEP sequences. Convolution of the test signals and the time-reversed unit-CAPRICEP yields compressed signal $q^{(m)}_{\mathrm{R}}[n]$. They consist of recovered pulses corresponding to repetitions. Since a set of unit-CAPRICEPs are not strictly orthogonal, correlations remain.

$\displaystyle q^{(m)}_{\mathrm{R}}[n]$

$\displaystyle=h_{\mathrm{fvn}}^{(m)}[-n]\ast x_{\mathrm{testR}}[n],$

(13)

where “$\ast$” represents the convolution operation.

The inner product of $\mathbf{b}^{(i)}$ and $\mathbf{b}^{(j)}$ of (11) is zero when $i\neq j$. Periodic time shift and averaging using the binary weights in $\mathbf{b}^{(1)},\ldots,\mathbf{b}^{(4)}$ remove these cross-correlations. The following equation yields the orthogonalized signal $r_{\mathrm{itrR}}^{(m)}[n]$.

$\displaystyle r_{\mathrm{itrR}}^{(m)}[n]$

$\displaystyle=\frac{1}{8}\sum_{k=0}^{8-1}q^{(m)}_{\mathrm{R}}[n-kn_{\mathrm{o}}]\,b^{(m)}_{\widetilde{k+1}}\ ,$

(14)

where $8$ is the length of the weight vectors $\mathbf{b}^{(1)},\ldots\mathbf{b}^{(4)}$. The notation $\widetilde{k+1}$ represents cyclic indexing between 1 and 8.

Synchronous averaging given below yields the unit impulse response. It also reduces the effect of random fluctuation of the signal.

$\displaystyle r^{(m)}_{\mathrm{R}}[n]$

$\displaystyle=\frac{1}{\#(\Omega)}\sum_{n_{\mathrm{ini}}+8kn_{\mathrm{o}}\in\Omega}r_{\mathrm{itrR}}^{(m)}[n+n_{\mathrm{ini}}+8kn_{\mathrm{o}}],$

(15)

where $0\leq n<n_{\mathrm{o}}$, and the symbol $\Omega$ represents the region where cross-correlation is effectively vanished. The function $\#(\Omega)$ yields the number of pulses separated by $8n_{\mathrm{o}}$ and located inside the region. Then, further averaging the first three responses $r^{(m)}_{\mathrm{R}}[n],(m=1,2,3)$ provides the final averaged response $r_{\mathrm{R}}[n]$.

$\displaystyle r_{\mathrm{R}}[n]$

$\displaystyle=\frac{1}{3}\left(r^{(1)}_{\mathrm{R}}[n]+r^{(2)}_{\mathrm{R}}[n]+r^{(3)}_{\mathrm{R}}[n]\right).$

(16)

Figure 3 shows examples to illustrates these procedures. This example uses CAPRICEP designed to have the equivalent rectangular duration $T_{\mathrm{ERD}}=0.2$ s. We truncated the length of CAPRICEP 0.8 s for implementation using FFT-based fast convolution. IIR-based implementation requires a two-pass procedure consisting of time axis reversal, which is not relevant for real-time applications. The top three plots of Fig. 3 show the recovered signals $q_{\mathrm{R}}^{(m)}[n],$ ($m=1,2,3)$. The fourth plot shows $q_{\mathrm{R}}^{(1)}[n]$ using dB scale. The background noise around 45 dB is the result of the cross-correlation between unit-CAPRICEPs. The bottom plot shows the orthogonalized signal $r^{(1)}_{\mathrm{R}}[n]$ using the dB scale. Note that the orthogonalization procedure removes components due to cross-correlation between unit-CAPRICEPs.

The orthogonalized signal using the fourth unit-CAPRICEP $r_{\mathrm{itrR}}^{(4)}[n]$ does not consist of components in the test signal $x_{\mathrm{testR}}[n]$. Therefore, $r_{\mathrm{itrR}}^{(4)}[n]$ consists of the background noise, the effects of temporal variation in the response, and randomly generated non-linear time-variant response.

The consecutive eight repetition segments provide the eight possible combinations of three unit-CAPRICEPs’ polarities. The averaged response $r_{\mathrm{R}}[n]$ is the average of all these combinations, deviations from this average consists of the non-linear time-invariant responses (and the random component mentioned above). Detailed derivation in [6] provides the procedure of the simultaneous measurement and simulation tests quantitatively validated the procedure using FVN. This derivation holds for CAPRICEP because FVN and CAPRICEP are members of TSP. We conducted the same simulation by replacing FVNs with CAPRICEPs. The results using CAPRICEP also validated the procedure.