I am using Non-negative Matrix Factorization (NMF) on an audio file of a working motor to decompose it into five NMF components. I need to explore the quefrency/cepstral domain of each component, but I am unsure how to proceed. Since the NMF components are computed after applying the Short-Time Fourier Transform (STFT) to the original signal and are represented as a 2D matrix, I am uncertain about how to transform them into the quefrency domain. I attempted to generate a cepstrogram for one of my NMF components using the following code: # STFT Parameters n_fft = 2048 hop_length = 256 win_length = n_fft window = 'hann' # Compute STFT S = librosa.stft(y, n_fft=n_fft, hop_length=hop_length, win_length=win_length, window=window) X, X_phase = librosa.magphase(S) # Normalization X_norm = librosa.util.normalize(X) # NMF Calculation n_components = 5 model = NMF(n_components=n_components, init='random', beta_loss='kullback-leibler', solver='mu', max_iter=3000, tol=1e-5, l1_ratio=0.8) W = model.fit_transform(X_norm) H = model.components_ # Compute the cepstrogram from log magnitude spectrum c = np.outer(W[:, 4], H[4]) #c is already a absolute value # Compute log magnitude log_magnitude = np.log(c + 1e-10) # Add small value to avoid log(0) # Compute cepstrum using IFFT cepstrum = np.fft.ifft(log_magnitude, axis=0).real # Compute real cepstrum across time frames # Number of unique cepstral coefficients NUP = cepstrum.shape[0] // 2 # Keep only the first half # Extract the meaningful part of the cepstrum C = cepstrum[:NUP, :] # Compute quefrency and time vectors q = np.arange(NUP) / sr # Quefrency (s) t = np.arange(C.shape[1]) * (hop_length / sr) # Time (s) # ===== APPLY CONDITIONING ===== q_mask = q > 0.0005 # Ignore quefrencies below 0.5 ms C = C[q_mask, :] # Filter cepstrum coefficients q = q[q_mask] * 1000 # Convert quefrency from seconds to milliseconds plt.figure(figsize=(10, 6)) plt.imshow(C, aspect='auto', origin='lower', extent=[t[0], t[-1], q[0], q[-1]]) plt.xlabel('Time (s)') plt.ylabel('Quefrency (ms)') plt.title('Cepstrogram from STFT (Magnitude Spectrum)') plt.colorbar(label='Amplitude') plt.show() I need to determine whether this approach is correct. If it is, how can I extract quefrency domain features such as Linear Frequency Cepstral Coefficients (LFCC), Constant Q Cepstral Coefficients (CQCC), Teager Energy Cepstral Coefficients (TECC), and Gabor Filterbank Cepstral Coefficients (GTFB)? Continue reading...
