Внезапный скачок частоты синусоидальной функции

Я пишу источник сигнала для цифровой обработки сигналов и обнаружил странное поведение. Вот код:

    float complex *output = malloc(sizeof(float complex) * 48000);
    FILE *f = fopen("test_signal.cf32", "wb");
    int64_t freq = -3355;
    float phase_increment = 2 * (float) M_PI * (float) freq / 48000;
    float phase = 0.0F;
    for (int i = 0; i < 150 * 5; i++) {
        for (int j = 0; j < 9600; j++) {
            output[j] = cosf(phase) + sinf(phase) * I;
            phase += phase_increment;
        }
        fwrite(output, sizeof(float complex), 9600, f);
    }
    fclose(f);

Это создаст файл со сложным сигналом, смещенным на -3355 Гц от центра. Поэтому, когда я запускаю БПФ над этим файлом, я ожидаю сигнал на частоте -3355 Гц. С некоторым незначительным частотным распределением вокруг этого числа из-за квантования. Но на самом деле я получаю следующее:

Внезапный скачок частоты синусоидальной функции

Есть довольно значительный скачок частоты (~2000 Гц) в районе 50 секунд.

Кто-нибудь знает, почему?

Мои эксперименты показали:

в output[j] = cosf(phase) + sinf(phase) * I; является ли эта большая I другой переменной, чем маленькая i, и если да, то что это за переменная I?

Sedenion 17.05.2022 19:25

@Sedenion, это воображаемый юнит, импортированный из complex.h

tstanisl 17.05.2022 19:29

Тема генерации бесконечного синуса обсуждалась на stackoverflow.com/q/69729326/4989451. Это не ответ на ваш вопрос, но вы можете найти его полезным

tstanisl 17.05.2022 19:32

MS справочная страница для sin() говорит: «Если Икс больше или равно 263 или меньше или равно -263, происходит потеря значимости в результате».

Weather Vane 17.05.2022 19:35
Стоит ли изучать PHP в 2023-2024 годах?
Стоит ли изучать PHP в 2023-2024 годах?
Привет всем, сегодня я хочу высказать свои соображения по поводу вопроса, который я уже много раз получал в своем сообществе: "Стоит ли изучать PHP в...
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
В JavaScript одним из самых запутанных понятий является поведение ключевого слова "this" в стрелочной и обычной функциях.
Приемы CSS-макетирования - floats и Flexbox
Приемы CSS-макетирования - floats и Flexbox
Здравствуйте, друзья-студенты! Готовы совершенствовать свои навыки веб-дизайна? Сегодня в нашем путешествии мы рассмотрим приемы CSS-верстки - в...
Тестирование функциональных ngrx-эффектов в Angular 16 с помощью Jest
В системе управления состояниями ngrx, совместимой с Angular 16, появились функциональные эффекты. Это здорово и делает код определенно легче для...
Концепция локализации и ее применение в приложениях React ⚡️
Концепция локализации и ее применение в приложениях React ⚡️
Локализация - это процесс адаптации приложения к различным языкам и культурным требованиям. Это позволяет пользователям получить опыт, соответствующий...
Пользовательский скаляр GraphQL
Пользовательский скаляр GraphQL
Листовые узлы системы типов GraphQL называются скалярами. Достигнув скалярного типа, невозможно спуститься дальше по иерархии типов. Скалярный тип...
3
4
48
2
Перейти к ответу Данный вопрос помечен как решенный

Ответы 2

Ответ принят как подходящий

Как только phase становится достаточно большим, приращения к нему квантуются путем округления с плавающей запятой, может быть, достаточно, чтобы округление до ближайшего с тай-брейком округления до четной мантиссы не уравновешивалось, а вместо этого означало, что вы последовательно продвигаете фазу меньше. См. статья в Википедии о плавающей запятой одинарной точности IEEE

Аналогия с основанием 10 со значениями, ограниченными 4 значащими цифрами, будет 101.2 + 0.111 только до 101.3, а не 101.311. (Компьютеры используют двоичные числа с плавающей запятой, поэтому на самом деле такие вещи, как 101.125, точно представимы, а 101.1 — нет.)

Я подозреваю, что это то, что происходит. Вы не переполнитесь до +Inf, но, учитывая ограниченную относительную точность float, добавление небольшого числа к огромному числу в конечном итоге вообще не переместит его: ближайшим к phase + phase_increment представимым числом с плавающей запятой будет исходное phase.

Быстрый способ проверить эту гипотезу — использовать double phase (без внесения каких-либо других изменений) и посмотрите, сдвинется ли при этом точка, в которой вы получаете шаги частоты. (Чтобы быть намного лучше, возможно, за пределами того, что вы генерируете).

(О, только что заметил, что вы уже сделали это, и это действительно помогло. Если вы имеете в виду, что это полностью устранило проблему, то, вероятно, это так.)

Вы все еще можете сохранить одинарную точность sinf и cosf, хотя их изменение было бы другой задачей: как указывает @Weather Vane, некоторые реализации sin теряют значительную точность при уменьшении диапазона. Но вы ожидаете, что это будет более шумный сигнал, а не постоянное изменение частоты (которое потребует масштабного коэффициента для фазы).

Вы также можете позволить ему работать дольше с помощью float и посмотреть, получите ли вы еще один шаг вниз по частоте или изменение до постоянного тока (постоянная фаза, частота = 0).

Или с помощью float разверните вручную, чтобы у вас было два счетчика, смещенных на phase_increment, и вы увеличили каждый на 2*phase_increment. Таким образом, относительная величина phase по сравнению с тем, что добавляется, не становится такой экстремальной. Однако они по-прежнему будут иметь ту же абсолютную величину, большую по сравнению с исходной phase_increment, поэтому некоторое округление все же будет.


Я не знаю, есть ли лучшие стратегии для создания сложных синусоидальных волн, я просто пытаюсь объяснить поведение, которое вы видите.

Для производительности этого метода вы, вероятно, захотите sincosf вычислять оба значения одновременно, если ваш компилятор не оптимизирует ваши вызовы для этого. (И, надеюсь, автоматизировать векторизацию, вызвав SIMD sincosf.)

Понятно. Благодарю вас! Кстати: -O2 создаст sincosf для этого кода.

Andrey Rodionov 17.05.2022 20:01

То, что вы наблюдаете, - это проблемы с точностью с плавающей запятой. Точность поплавков ухудшается с величиной.

К счастью, для фаз или подобных угловых значений есть простой обходной путь: оберните значение вокруг 2π, чтобы величина никогда не была слишком большой.

Добавьте следующий код после phase += phase_increment; и посмотрите, поможет ли это.

if ( phase < (float)( 2.0 * M_PI ) )
    continue;
phase -= (float)( 2.0 * M_PI );

Да, редукция работает хорошо: ``` if (phase < -(float)( 2.0 * M_PI )) { Phase += (float)( 2.0 * M_PI ); } if ( фаза > (с плавающей запятой) ( 2.0 * M_PI )) { фаза -= (с плавающей запятой) ( 2.0 * M_PI ); }```

Andrey Rodionov 17.05.2022 20:25

Другие вопросы по теме