#1 楼
好吧,让我们大声笑@ AnonSubmitter85给您一个很好的答案,但让我展示一下在matlab中完成此操作的方式,这也许很容易移植到C:
首先我要在
256
采样的500hz
中创建44100hz
采样,看看如何累积相位,并在第一个循环结束时将相位置于0和2pi之间。 />现在,让我们转到第二个循环,以创建更多的256
样本,并使其连续进行。 > 这是上面代码的图::
#2 楼
看到这么多答案有点奇怪,但没有一个答案能用C给出实际答案,也没有解释如何以及为什么这样做。一般的想法是保持一个逐步增加步长的相位根据频率和采样率计算得出。这样,您将永远不会出现相位不连续性。
执行此操作时,必须非常小心浮点变量中的累积舍入误差,因为尽管相对舍入误差主要保持在不管数字有多大,绝对舍入误差都会随着数字的大小而增加。 (在每位计算机科学家应该了解的有关浮点算法的文章中有深入而详尽的解释。)
如果您只是继续在阶段中增加步长,它将很快达到幅值太大,因此必须加以检查。由于我们处理的是单个波形,因此可以通过环绕或模运算将相位限制在0.0和1.0之间。执行此操作的最简单且可能最有效的方法是为该阶段使用无符号整数,然后让C编译器负责包装。 C令人讨厌,因为未定义许多算术运算,但是以所需的方式定义了无符号整数环绕,因此我们可以利用它。
下面的程序输出波形使用上述技术在标准输出上显示。相位保持为无符号整数,并与基于请求的频率和采样率计算出的另一个整数递增。每当相位超过整数大小设置的限制时,它将自动换为零。
然后将整数相位缩放为正弦函数中所需的索引。可以很容易地将其更改为使用波表或内插查找。
由于可移植性和兼容性,某些部分异常复杂。在固定设置下,可以简化很多操作。同样,在某些地方可能需要添加更好的舍入。
#include <math.h>
#include <stdio.h>
#define M_TWOPI 6.283185307179586476925286766559
/*
* The phase must be an unsigned integer.
* 'maxphase' is for example 65536 if phase_t is 16 bits.
*/
typedef unsigned long phase_t;
double maxphase = (double)((phase_t)0-(phase_t)1)+1.0;
double fs = 44100;
phase_t hz_to_delta( double hz )
{
return maxphase*hz/fs+0.5;
}
float sample_phase( phase_t phase )
{
return sin( phase/maxphase*M_TWOPI );
}
int main( void )
{
long i;
phase_t delta, iphase = 0;
delta = hz_to_delta( 500.0 );
for( i=0; i<fs; ++i )
printf( "%e\n", sample_phase( iphase += delta ) );
delta = hz_to_delta( 1000.0 );
for( i=0; i<fs; ++i )
printf( "%e\n", sample_phase( iphase += delta ) );
return 0;
}
评论
$ \ begingroup $
仅供参考,此帖子最初是关闭的,因为要求进行严格的编码违反SE政策。可以通过在两个循环之间添加一行代码并在第二个循环中修改一行来固定OP的代码。我不会编写代码,而只是在结处设置一个用于相位偏移的变量,然后将其添加到第二个循环的参数中。我只是为了区分“变化”(从标题)和“变化”(从代码示例)而进行区分。
$ \ endgroup $
–雪铁龙·道格(Cedron Dawg)
20-2-12在18:16
$ \ begingroup $
@CedronDawg我知道,尽管我仍然认为“整数阶段”的技巧对于实际的DSP实现而言已经足够基础,至少值得一提。
$ \ endgroup $
–烟斗
20-2-12在18:30
$ \ begingroup $
应该始终考虑它,但并不总是确定性的。您付出的额外除法和类型转换可能不值得,也不是一堆额外的函数调用。附带一提,您的代码将受益于用存储的倒数的乘积替换除数。
$ \ endgroup $
–雪铁龙·道格(Cedron Dawg)
20-2-12在18:37
#3 楼
最简单的方法是注意,根据定义,频率是相位的导数。因此,您可以定义每个样本的频率,然后对其进行积分以获得相位。这样可以保持相位连续。例如,在matlab中,它看起来像这样:评论
$ \ begingroup $
“根据定义,频率是相位的导数” –我不同意这一点。根据定义,频率是每个时间单位中信号重复的次数。根据定义,相位的导数是相位速度。碰巧频率等于ω/2π,尽管可以说它是正弦波特有的。
$ \ endgroup $
–leftaround关于
20-2-13在15:11
$ \ begingroup $
我不能说我曾经听说过相速度一词。也许瞬时频率是我如何使用它的更准确的术语。我想每个领域都有自己的秘密。
$ \ endgroup $
– AnonSubmitter85
20-2-13在15:22
$ \ begingroup $
实际上,我指的是角频率(德语术语是Winkelgeschwindigkeit,点亮。“角速度”)。因此,好的,它在英语中被称为“频率”,但我认为频率和角频率仍应保持分开。如果没有其他问题,则至少是由于2·π因子。
$ \ endgroup $
–leftaround关于
20-2-13在16:03
$ \ begingroup $
@leftaroundabout:有趣的事实:相速度是物理学中的一个术语,用于在空间或介质(例如,传输线中的电信号)中传播的波。相速度=角速度/波数,并且具有距离/时间而不是角度/时间的尺寸。 en.wikipedia.org/wiki/Phase_velocity =波浪的明显“波峰”移动的速度。有一些相关的术语,例如“群速度”,有关动画,请参见Wikipedia。
$ \ endgroup $
– Peter Cordes
20-2-14在3:16
#4 楼
除了增加相位(而不是增加时间并乘以频率,可能会引起锯齿)之外,还请注意,触发函数的输入可能需要进行范围限制(或包装),以防止触发函数降低量程的精度实施。按增量相位递增后,我通常将相位(通过圆形环绕)限制在-2pi和2pi之间,或者将-pi限制为+ pi,等等。
评论
$ \ begingroup $
感谢您的回答。范围缩小精度损失是什么意思?
$ \ endgroup $
– deltafft
20-2-12在11:30
$ \ begingroup $
如果自变量相对于浮点尾数中的位数变大,则相位的小数位(小相位角增量)将越来越量化甚至消失。
$ \ endgroup $
– hotpaw2
20-2-12在15:06
#5 楼
通常,可以将纯稳定的实际正弦信号建模为:参数$ \ theta t + \ phi $通常被称为信号的“相位”,这可能会引起混淆,因为$ \ phi $也使用相同的术语,因此我们将其称为相位函数。对于稳定频率,相位函数是线性的,瞬时频率与恒定频率相同。在规定的情况下,必须确保相位函数的值在结点处匹配,并且信号将是“连续的”(在离散的上下文中实际上是无意义的术语。)
信号,您需要一个不是线性的相位函数。如我的评论所述,最简单的方法是使用二次方。概括系数名称,可以将其表示为:
$$ x(t)= A \ cos(C_0 + C_1 t + C_2 t ^ 2)$$
这导致频率以线性速率变化。其他功能也是可能的。对此进行匹配并进行编码是您的责任。
#6 楼
对于嵌入式设备,高级别的答案是不够的,这是我的方法,一个带有相位计数的查找表。可以使用不同时序循环上的插值来平滑查找之间的值。通过使用固定的查找,可以免费顺畅地适应频率的阶跃变化您只需要波形查找表的1/4即可生成合理的正弦波,但是此示例具有内存不足的全波查找为了简化代码
频率由下一个刻度超时确定
#define F1 500 // in hz
#define F2 1000
#define SAMPLES 255
#define CUTOVER 10*SAMPLES // in phase cycles
#define PHASE_STEP(f) 1.0/(f*SAMPLES)
uint32_t sine[SAMPLES] = {0, ... }; //fixed point sine samples
uint32_t output; // used to drive DAC/DSP output value
void isr()
{
static uint8_t phase = 0;
output = sine[phase];
phase++;
}
void main()
{
// can be a timer or otherwise
for(unsigned int cycle = 0; cycle < CUTOVER; cycle++) {
isr()
sleep(PHASE_STEP(F1));
}
while (1) {
isr()
sleep(PHASE_STEP(F2));
}
}
评论
$ \ begingroup $
当然,查找表可以工作,但是不可避免地在节点上存在一些不连续性,这可能是一个问题。 CORDIC也具有不连续性,因此不适合浮点运算。多项式的逼近可以得到非常好的结果,而出乎意料的很少的计算-OlliW的“快速函数”似乎是对此的第一个很好的参考,但是我已经推导了具有相同基本原理的高阶版本,并且得到了很好的结果。
$ \ endgroup $
–格雷厄姆
20-2-13在14:03
$ \ begingroup $
@Grahan的优点当然,基于查找的示例已得到极大简化,并且在实践中可以完成许多优化和优化的数据结构,但是从根本上说,如果函数可以多项式化(例如,仅正弦输出),那么在定点上这是不昂贵的也一样
$ \ endgroup $
–疯狂
20-2-14在1:22
#7 楼
也许您也会对某些ASCII图形感兴趣。在这里查看我完整的C解决方案,以使用
printf
s在C中生成图。 -the-exponential-function / 53899222#53899222 结果是这样的:
根据需要更改以下内容:
double equation(const double x)
{
return exp(-0.1 * x) * sin(0.5 * x);
}
编辑:
要获得频率不断增长的正弦曲线,可以使用以下方程式作为起点:
double equation(const double x)
{
const double xx = x - LIMIT_X_MIN;
return sin(0.05 * xx * xx);
}
评论
$ \ begingroup $
我很好奇,但我一定做错了,因为该程序在CPU核心最大的情况下对我来说挂起。我从字面上复制粘贴了您的代码,尝试了更少的宽度,高度,点数,相同的内容。我使用Archlinux并使用gcc test.c -lm -o test进行编译,没有错误或消息。在没有X的控制台中,最大化的普通100x40终端窗口是相同的。
$ \ endgroup $
–有关的公民
20-2-12在21:39
$ \ begingroup $
您显示的代码会随时间改变信号的幅度,而不是要求的频率。
$ \ endgroup $
–ambitiose_sed_ineptum
20-2-12在21:59
$ \ begingroup $
@ambitiose_sed_ineptum该代码是作为我以前的回答之一提供的示例。人们需要改变方程式以获得正确的结果。在这种情况下,那将去除指数。
$ \ endgroup $
–约翰
20-2-13在10:48
$ \ begingroup $
@aconcernedcitizen好抓住。我那里有一个错误。请参阅修改后的代码。唯一的变化是您需要x = get_next_x(x)来代替get_next_x(x)。
$ \ endgroup $
–约翰
20-2-13在10:59
$ \ begingroup $
现在可以使用。至于降票,无论谁做,它都可能认为它不能直接解决OP,而不是将其作为补充。只是一个猜测。
$ \ endgroup $
–有关的公民
20-2-13在11:27
#8 楼
自动化的一点点clear all;
close all;
m_phase = 0;
frequency = 500;
fs = 44100
signal = [];
samples = 2^8;
initial = 1;
for k=1:3
frequency = k*frequency
phaseInc = 2*pi*frequency/fs;
if(k == 1)
initial = 1;
else
initial = ((k-1)*samples) + 1;
endif
for i=initial:k*samples
signal(i) = sin(m_phase);
m_phase = m_phase + phaseInc;
end
m_phase = mod(m_phase, 2*pi);
endfor
subplot(211);
plot(signal);
N = 512;
X = fftshift(fft(signal,N));
df = fs / N;
sampleIndex = -N/2:N/2-1;
relative_f = sampleIndex * df
subplot(212);
plot(relative_f, abs(X));
axis([min(relative_f) max(relative_f)])
评论
$ \ begingroup $
但是我不知道我们可以从FFT推论什么?
$ \ endgroup $
– jomegaA
20-2-13在20:19
$ \ begingroup $
我认为有些东西无法显示FFT中的频率
$ \ endgroup $
–ederwander
20-2-13在22:06
评论
@Peter K我认为将其作为编码问题将其关闭可能为时过早。 OP尝试的解决方案不是改变频率,而是改变频率。为此,需要做的就是通过选择正确的相位偏移来确保结点处的信号连续。可以通过小的重叠和滑动平均值来实现平滑度。对于变化的频率,答案更简单。只需在trig参数中使用非线性函数即可。二次方是最简单的二次方。例如。 $ A \ cos(C_0 + C_1 \ theta + C_2 \ theta ^ 2)$。@CedronDawg好!重新打开......
如果要查找各种流行的线性调频(如二次线性调频或三次线性调频),则平滑变化的版本称为“线性调频”。每个都有有趣的有用属性。
这似乎与Math StackExchange上的问题有点相关:术语术语:离散a * sin(t)的连续性。