feat: GuitarAnalyzer (spectrogram power in freq range)

rhythm.py

@@ -284,3 +284,71 @@ class BassAnalyzer:
     def _viterbi_ampl(self, Scp2, path):
         max_amplitudes = np.array([np.abs(Scp2[i, path[i]]) for i in range(Scp2.shape[0])])
         return max_amplitudes
+
+
+class GuitarAnalyzer:
+    """
+    Rhythm analysis from songs. 
+    Provides a beat amplitude signal from the audio signal.
+
+    Performs short-time Fourier Transform on the signal, 
+    then returns the f1..f2 band power for each window.
+
+    For low-frequency instruments like percussion, 
+    use BassAnalyzer instead.
+    """
+    W = 1024  #: window size (must be even, so that right padding W/2-1 works)
+    shift_sec = 0.008  #: window shift in sec ('delta_tau') between subsequent windows
+    target_band_f1 = 800.0   #: lower bound of target freq band in Hz
+    target_band_f2 = 4300.0  #: upper bound of target freq band in Hz
+
+    def __init__(self, fs, sig):
+        """
+        :param fs: sampling rate
+        :param sig: audio signal normalized to [-1,1]
+        """
+        self.f = np.pad(sig, (self.W//2, self.W//2-1))  #: signal padded (W-FFT to determine scalogram parameters)
+        self.fs = fs
+
+        self.D = int(self.shift_sec * fs)  #: spectrogram step
+        self.L = self.f.shape[0]
+        self.M = (self.L-self.W) // self.D + 1  #: number of time steps
+
+    def spectrogram_power_amplitudes(self):
+        """
+        Compute spectrogram power from a target frequency range.
+        NOTE: downsampled from the original 'fs'.
+        :returns: (fsd, sig): sampling rate, amplitude signal
+        """
+        Spf = self._spectrogram()
+        ampl = self._spectrogram_power(Spf)
+        return ampl
+
+    def _spectrogram_power(self, Spf):
+        # Spf
+        # fs, W
+        fs, W = self.fs, self.W
+        k1, k2 = int(self.target_band_f1/fs*W), int(self.target_band_f2/fs*W)  # freq band range in W-FFT
+        #
+        # spectrum power in f1..f2 bands
+        #
+        #hp_slice = highpass(np.sum(np.abs(Spf_slice[:, k1:k2]), axis=1), fps=fs/Dp, cf=2.0, tw=0.2)
+        hp_slice = np.sum(np.abs(Spf[:, k1:k2]), axis=1)-np.mean(np.sum(np.abs(Spf[:, k1:k2]), axis=1))
+        return hp_slice
+
+    def _spectrogram(self):
+        """W-FFT (STFTs) to determine scalogram parameters"""
+        # *f
+        # M, W, D
+        f = self.f
+        M, W, D = self.M, self.W, self.D
+
+        #
+        # compute spectrogram: 'Spf' (M x W)  <-  from 'f'
+        #
+        iwss = np.linspace(W//2, W//2 + (M-1)*D, M, dtype=int)  # 'D'-spaced start time indices of windows on 's'
+        Spf = np.zeros((M, W), dtype=np.complex128)
+        for i, iw in zip(range(M), iwss):
+            iws, iwe = iw-W//2, iw+W//2
+            Spf[i,:] = fft(f[iws:iwe])
+        return Spf