在给定时间/样本量下,频率f1和f2之间指数变化的正弦波

Eth*_*ith 5 c python algorithm math waveform

我正在尝试实现生成正弦波的Python方法,它以指数方式在两个频率之间斜升.使用以下Python代码在[this question]中解决了线性变化:

from math import pi, sin

def sweep(f_start, f_end, interval, n_steps):    
    for i in range(n_steps):
        delta = i / float(n_steps)
        t = interval * delta
        phase = 2 * pi * t * (f_start + (f_end - f_start) * delta / 2)
        print t, phase * 180 / pi, 3 * sin(phase)

sweep(1, 10, 5, 1000)
Run Code Online (Sandbox Code Playgroud)

如何将这种线性累积相位/ delta方法改变为指数频率扫描并且对人耳是平滑的.

Bas*_*els 5

这类问题的诀窍在于理解频率调制相位调制之间的关系,这两者密切相关.具有恒定频率f和幅度的正弦A可以描述为(公式,而不是python代码):

x(t) = A sin(2 * pi * f * t)
Run Code Online (Sandbox Code Playgroud)

但是写一个不同的方法是首先定义一个阶段phi作为时间的函数:

phi(t) = 2 * pi * f * t
x(t) = A sin(phi(t))
Run Code Online (Sandbox Code Playgroud)

这里要注意的是频率f是相位的导数,除以2*pi : f = d/dt(phi(t)) / (2*pi).

对于频率随时间变化的信号,您可以类似地定义瞬时频率 f_inst:

x(t) = A sin(phi(t))
f_inst = d/dt(phi(t)) / (2*pi)
Run Code Online (Sandbox Code Playgroud)

您想要做的是与此相反,您有一个给定的瞬时频率(您的对数扫描),您需要将其转换为一个阶段.由于推导的反面是积分,你可以像这样计算适当的相位(仍然是公式):

phi(t) = 2 * pi * Integral_0_to_t {f_inst(t) dt}
x(t) = A sin(phi(t))
Run Code Online (Sandbox Code Playgroud)

你在这里做的是某种信号的相位调制(零频率),以获得所需的瞬时频率.在numpy中这很容易做到:

from pylab import *
n = 1000 # number of points
f1, f2 = 10, 30 # frequency sweep range in Hertz

t = linspace(0,1,1000)
dt = t[1] - t[0] # needed for integration

# define desired logarithmic frequency sweep
f_inst = logspace(log10(f1), log10(f2), n)
phi = 2 * pi * cumsum(f_inst) * dt # integrate to get phase

# make plot
plot(t, sin(phi))
xlabel('Time (s)')
ylim([-1.2, 1.2])
grid()
show()
Run Code Online (Sandbox Code Playgroud)

结果图片:

对数频率扫描

但是(正如戴夫提到的那种欺骗中所指出的那样),你可能不希望进行对数扫描,而是指数扫描.你的耳朵具有频率的对数感知,因此平滑/线性音阶(想想钢琴上的键)因此呈指数间隔.这可以通过简单地重新定义瞬时频率来实现f_inst(t) = f1 * exp(k * t),其中k选择了这样的频率f_inst(t2) = f2.

如果要同时使用幅度调制,只需更改AA(t)与公式相关的时间即可.


and*_*oke 5

Bas的答案很棒,但实际上并没有给出解析解,所以这部分是......

据我所知,你想要的sin(Aexp(Bt))地方AB常量.我会假设时间开始0并继续C(如果它在其他时间开始,从两者中减去).

然后,巴斯说,我认为,如果我们有sin(g(t))频率f是这样的2 * pi * f = dg / dt.我们希望f0时间0fC时间都是如此C.

如果你通过数学,这很容易(它确实是 - 去年的学校水平),你会得到:

B = 1/C * log(fC/f0)
A = 2 * pi * f0 / B
Run Code Online (Sandbox Code Playgroud)

这里是一些使用1000个样本在5秒内从1到10Hz的代码:

from math import pi, sin, log, exp

def sweep(f_start, f_end, interval, n_steps):
    b = log(f_end/f_start) / interval
    a = 2 * pi * f_start / b
    for i in range(n_steps):
        delta = i / float(n_steps)
        t = interval * delta
        g_t = a * exp(b * t)
        print t, 3 * sin(g_t)

sweep(1, 10, 5, 1000)
Run Code Online (Sandbox Code Playgroud)

这使:

一个漂亮的情节

(并且您可以添加常量 - sin(g_t + k)以获得您想要的起始阶段).

更新

为了表明你所看到的问题是一个抽样的假象,这里有一个做过采样的版本(如果你把它设置为一个参数):

from math import pi, sin, log, exp

def sweep(f_start, f_end, interval, n_steps, n_oversample=1):
    b = log(f_end/f_start) / interval
    a = 2 * pi * f_start / b
    for i in range(n_steps):
        for oversample in range(n_oversample):
            fractional_step = oversample / float(n_oversample)
            delta = (i + fractional_step) / float(n_steps)
            t = interval * delta
            g_t = a * exp(b * t)
            print t, 3 * sin(g_t)

sweep(16000.0, 16500.0, 256.0/48000.0, 256)      # looks strange
sweep(16000.0, 16500.0, 256.0/48000.0, 256, 4)   # looks fine with better resolution
Run Code Online (Sandbox Code Playgroud)

如果您检查代码,您将看到n_oversample4的所有设置(第二次调用)都会为时间步长添加更高的分辨率.特别是,代码何时oversample = 0(即fractional_step = 0)之前相同,所以第二个图包括第一个图中的点,加上额外的"填写"缺失数据并使一切看起来不那么令人惊讶.

这是原始和开始附近的过采样曲线的特写,详细显示了正在发生的事情:

在此输入图像描述

最后,这种事情是完全正常的,并不表示任何类型的错误.当从数字波形生成模拟信号时,您将获得"正确"的结果(假设硬件正常工作). 如果不清楚,这个优秀的视频将解释事情.