2022-09-19

WARNING: This article may be obsolete
This post was published in 2022-09-19. Obviously, expired content is less useful to users if it has already pasted its expiration date.

This article is categorized as "Garbage" . It should NEVER be appeared in your search engine's results.

多普勒效应

介绍

🔗 [多普勒效应 - 维基百科，自由的百科全书] https://zh.wikipedia.org/wiki/多普勒效应

example:

假设我们正在街边静止不动，有一辆救护车从远处驶过来，然后又走了，现在要用R/python生成我们听到的救护车声音（多普勒）：

R语言版本

R版本：

library(tuneR)
# only available on macOS
# setWavPlayer("afplay")

Fs = 8000
t = seq(0, 10 - 1 / Fs, by = 1 / Fs)
cos_theta = (50 * t - 250) / ((100 ^ 2 + (50 * t - 250) ^ 2) ^ (1 / 2))
loudness_series = 1 / ((100 ^ 2 + (50 * t - 250) ^ 2) ^ (1 / 2))
# pitch shift function
k = 0.05
pitch_shift_func = 1 / (1 + k * cos_theta)
# plot(pitch_shift_func)

xt_total = rep(0, 10 * Fs)
for (n in c(1, 2, 4, 8)) {
  f1 = 440 * n
  f2 = 440 * 4 / 3 * n
  
  delta = 1 / Fs
  f_delta_series_1 = f1 * pitch_shift_func[1:(Fs / 2)]
  xt = (1 / n) * cos(2 * pi * cumsum(delta * f_delta_series_1))
  
  f_delta_series_2 = f2 * pitch_shift_func[(Fs / 2 + 1):(Fs)]
  xt = (1 / n) * c(xt, cos(2 * pi * cumsum(delta * f_delta_series_2)))
  
  for (x in 1:9) {
    f_delta_series_1 = f1 * pitch_shift_func[(Fs * x + 1):(Fs * x + Fs / 2)]
    xt = c(xt, cos(2 * pi * cumsum(delta * f_delta_series_1)))
    
    f_delta_series_2 = f2 * pitch_shift_func[(Fs * x + Fs / 2 + 1):(Fs * x + Fs)]
    xt = c(xt, cos(2 * pi * cumsum(delta * f_delta_series_2)))
  }
  xt_total = xt_total + xt
}

xt_total = xt_total * loudness_series

bits = 16
wave = floor(2 ^ (bits - 1) * xt_total)
wave = as.integer((wave - min(wave)) / (max(wave) - min(wave)) * (32767 +
                                                                    32768) - 32768)
# play(Wave(wave, samp.rate = Fs, bit = bits))
writeWave(Wave(wave, samp.rate = Fs, bit = bits),'./doppler_demo.wav')

Python版本（chatgpt）

python版本（chatgpt帮我把上面的R code翻译成了python版本）

import numpy as np
from scipy.io import wavfile

Fs = 8000
t = np.arange(0, 10, 1 / Fs)
cos_theta = (50 * t - 250) / ((100 ** 2 + (50 * t - 250) ** 2) ** (1 / 2))
loudness_series = 1 / ((100 ** 2 + (50 * t - 250) ** 2) ** (1 / 2))

# pitch shift function
k = 0.05
pitch_shift_func = 1 / (1 + k * cos_theta)

xt_total = np.zeros(10 * Fs)
for n in [1, 2, 4, 8]:
    f1 = 440 * n
    f2 = 440 * 4 / 3 * n

    delta = 1 / Fs
    f_delta_series_1 = f1 * pitch_shift_func[:int(Fs / 2)]
    xt = (1 / n) * np.cos(2 * np.pi * np.cumsum(delta * f_delta_series_1))

    f_delta_series_2 = f2 * pitch_shift_func[int(Fs / 2):Fs]
    xt = np.concatenate((xt, (1 / n) * np.cos(2 * np.pi * np.cumsum(delta * f_delta_series_2))))

    for x in range(1, 10):
        f_delta_series_1 = f1 * pitch_shift_func[(Fs * x + 1):(Fs * x + int(Fs / 2) + 1)]
        xt = np.concatenate((xt, np.cos(2 * np.pi * np.cumsum(delta * f_delta_series_1))))

        f_delta_series_2 = f2 * pitch_shift_func[(Fs * x + int(Fs / 2) + 1):(Fs * x + Fs + 1)]
        xt = np.concatenate((xt, np.cos(2 * np.pi * np.cumsum(delta * f_delta_series_2))))

    xt_total[:len(xt)] += xt

xt_total *= loudness_series

bits = 16
wave = np.floor(2 ** (bits - 1) * xt_total)
wave = np.int16((wave - np.min(wave)) / (np.max(wave) - np.min(wave)) * (32767 + 32768) - 32768)

wavfile.write('./doppler_demo_python.wav', Fs, wave)

音频文件

下面的2个音频文件几乎一样：

使用上面的R语言代码生成：

使用上面的python代码生成：