「残声」概念验证

zansei.py

#!/usr/bin/env python
#coding:utf-8
# Author:  Sg4Dylan --<sg4dylan#gmail.com>
# Created: 06/18/2021

import numpy as np
from scipy import signal
from scipy.io.wavfile import write
import librosa
import resampy
from tqdm import tqdm

# SSB Modulation
def freq_shift_mono(x, f_shift, d_sr):
    N_orig = len(x)
    N_padded = 2**int(np.ceil(np.log2(np.abs(N_orig))))
    S_hilbert = signal.hilbert(np.hstack((x, np.zeros(N_padded-N_orig, x.dtype))))
    S_factor = np.exp(2j*np.pi*f_shift*d_sr*np.arange(0, N_padded))
    return (S_hilbert*S_factor)[:N_orig].real

def freq_shift_multi(x, f_shift, d_sr):
    return np.asarray([freq_shift_mono(x[i], f_shift, d_sr) for i in range(len(x))])

def zansei_impl(x, sr, m=8, decay=1.25):
    # Pre-highpass
    b,a = signal.butter(11,3000/(sr/2),'high')
    d_src = signal.filtfilt(b,a,x)
    # !!!
    d_sr = 1/sr
    f_dn = freq_shift_mono if len(x.shape) == 1 else freq_shift_multi
    d_res = np.zeros_like(x)
    for i in tqdm(range(m)):
        d_res += f_dn(d_src, sr*(i+1)/(m*2), d_sr) * np.exp(-(i+1)*decay)
    # Post-highpass
    b,a = signal.butter(11,16000/(sr/2),'high')
    d_res = signal.filtfilt(b,a,d_res)
    # Dynamic-protect
    adp_power = np.mean(np.abs(d_res))
    src_power  = np.mean(np.abs(x))
    adj_factor = src_power/(adp_power+src_power)
    return (x + d_res)*adj_factor


path_of_audio = 'test.wav'
path_of_result = 'mod.wav'

audio_y, audio_sr = librosa.load(path_of_audio,mono=False,sr=None)
audio_y = resampy.resample(audio_y, audio_sr, audio_sr*2, filter='kaiser_fast')
audio_sr *= 2

shifted = zansei_impl(audio_y, audio_sr)
write(path_of_result, audio_sr, np.asarray(shifted).astype(np.float32).T)