## 3 Data Processing
### Time-Frequency Filtering
In the next step, a time-frequency bandpass filter \(WD^{*}(t,f)\) is constructed around one of the components of interest. The asterisk \(\star\) refers to the filtered component: _lin_, \(HH2\) or \(HH3\). The filter is constructed as:  
\[\begin{array}{l}WD^{*}(t,f)=WD_{1}^{*}(t,f).WD_{2}(t)\\ \text{with :}\\ WD_{1}^{*}(t,f)\\ =\left\{ \begin{array}{cc}1&\left|f^{*}-f_{lin}(t)\right|<\frac{FT(t)}{2}\\ 0&\left|f^{*}-f_{lin}(t)\right|>\frac{FT(t)}{2}+\frac{\pi}{4}BW(t)\\ \frac{1}{2}+\frac{1}{2}\text{cos}\left(\frac{4\left(\left|f^{*}-f_{lin}(t) \right|-\frac{FT(t)}{2}\right)}{BW(t)}\right)&\text{elsewhere}\\ WD_{2}(t)=\left\{\begin{array}{cc}0&t(0.06~{}t_{\text{meas}}~{}\text{OR}~{}t).09~{}t_{\text{meas}}\\ 1&\text{elsewhere}\\ \end{array}\right.\\ \text{and}\\ f_{lin}(t)=f_{start}+k(t)~{}(f_{end}-f_{start})\\ FT(t)=FT_{1}+k(t)~{}(FT_{2}-FT_{1})\\ BW(t)=BW_{1}+k(t)~{}(BW_{2}-BW_{1})\\ f^{*}=f\text{ for linear component extraction}\\ f^{*}=f/2\text{ for }\text{ see end HH component extraction}\\ f^{*}=f/3\text{ for third HH component extraction}\\ k(t)=\frac{t-0.05~{}t_{\text{meas}}}{0.90~{}t_{\text{meas}}}\end{array}  \tag{2}\]  
where \(f_{start}\) and \(f_{end}\) are the start and end frequency of the chirp, respectively. \(t_{meas}\) is the total measurement time, \(k(t)\) represents the chirp advance and \(f_{lin}(t)\) is the instantaneous chirp excitation frequency at time \(t\). The filter is Tukey shaped with time-dependent flat top length \(FT(t)\) and edge bandwidth \(BW(t)\). The term \(WD_{2}(t)\) limits the filter in time domain such that the response is only retained when the excitation signal is active. This is advised because the abrupt start and end of the burst chirp excitation signal, at 2 ms and 38 ms respectively, results in a broadband response of the sample (see Figure 4a). As an example, Figure 4b shows the filter for extracting the second harmonic component \(WD^{HH2}\) (with \(f^{*}=2f\)). The excitation signal properties are: \(f_{start}\) = 10 kHz, \(f_{end}\) = 60 kHz and \(t_{meas}\) = 40 ms and the filter settings are: flat top length \(FT_{1}\) = 3 kHz (at start), \(FT_{2}\) = 7 kHz (at end) and edge band width \(BW_{1}\) = 3 kHz (at start), \(BW_{2}\) = 7 kHz (at end). These filter settings are chosen such that the window fully envelops only the \(HH2\) component (see Figure 4c). As a final step, inverse STFT is performed to transform the bandpass filtered signal from the time-frequency domain back to the time domain:  
\[ \begin{array}{l}V_{Z}^{*}(x,y,t(m))=\sum\limits_{k=1}^{K}\overline{V}_{Z}( x,y,t(l),f(k))~{}~{}WD^{*}(t(l),f(k))~{}e^{2\pi i\frac{\pi k}{K}}\\ \to V_{Z}^{*}(x,y,t(n))=\frac{H}{E_{sw}}\sum\limits_{l=1}^{L}V_{Z}^{*}(x,y,t(n- lH))~{}w(n-lH)~{}v(n-lH)\end{array}  \tag{3}\]  
Note that the time divisions \(V_{Z}^{*}\) are shifted in time (\(m=n-lH\)) and added in order to obtain the final filtered signal \(V_{Z}^{*}(x,y,t(n))\). The function \(v\) is a window which helps to fulfill the correct overlap-and-add (OLA) conditions: \(\sum\limits_{l=1}^{L}w(n-lH)~{}v(n-lH)\) must be constant over time index \(n\) and \(E_{uv}=\sum\limits_{m=1}^{M}w(m)\ v(m)\). Here, a Hanning window is used and the OLA requirement was verified using the \"_OLAExam_\" software [37]. It is often desired to process the obtained velocity component \(V_{Z}^{*}(x,y,t)\) in the frequency domain instead of the time domain. Therefore, FFT is performed:  
\[V_{Z}^{*}(x,y,f(k))=\frac{1}{N}\sum\limits_{n=1}^{N}V_{Z}^{*}(x,y,t(n))\ e^{-2\pi i \frac{nk}{N}} \tag{4}\]