Как шумит радио (цифровая обработка и слуховой приём)

[alex_m]

чем шире "ворота" тем больше шума, СШ ухудшается

Такое ощущение, что этот треск не связан с естественным шумом.
Он и на относительно мощных сигналах присутствует.

как реализован демодулятор?

        public void Demodulate(Complex* iq, float* audio, int length)
        {
            for (var i = 0; i < length; i++)
            {
                // Polar discriminator
                var f = iq[i] * _iqState.Conjugate();

                // Limiting
                var m = (float)Math.Sqrt(f.Real * f.Real + f.Imag * f.Imag);
                if (m > 0.0f)
                {
                    f /= m;
                }

                // Angle estimate
                var a = (float)Math.Atan2(f.Imag, f.Real);

                // Scale
                audio[i] = a * FMGain; // FMGain = 0.00001f;

                _iqState = iq[i];
            }
            if (_mode == FmMode.Narrow)
            {
                for (var i = 0; i < length; i++)
                {
                    audio[i] *= NarrowAFGain;   // NarrowAFGain = 0.5f;
                }
            }
        }

[Genadi Zawidowski]

Записи не радиолюбительского обмена, но в любительском диапазоне (тема про QRM с Питерского форума).
Код демодулятора привожу, просто достаточное соотношение сигнал/шум.

	omega2ftw_k1 = POWF(2, NCOFTWBITS);

...

#define OMEGA2FTWI(angle) ((ncoftwi_t) ((FLOAT_t) (angle) * omega2ftw_k1 / (FLOAT_t) M_TWOPI))	// angle in radians -pi..+pi to signed version of ftw_t

...

// Демодуляция FM
static RAMFUNC ncoftwi_t demodulator_FM(
	FLOAT32P_t vp1,
	const uint_fast8_t rxi				// 0/1: main_RX/sub_RX
	)
{
	// Здесь, имея квадратурные сигналы vp1.IV и vp1.QV, начинаем демодуляции
	//
	// tnx Vladimir Vassilevsky
	// http://www.dsprelated.com/showmessage/71491/2.php
	//
	static ncoftwi_t prev_fi [NTRX] = { 0, };

	if (vp1.IV == 0 && vp1.QV == 0)
		vp1.QV = 1;
	const ncoftwi_t fi = OMEGA2FTWI(ATAN2F(vp1.QV, vp1.IV));	//  returns a value in the range –pi to pi radians, using the signs of both parameters to determine the quadrant of the return value.

	const ncoftwi_t d_fi = (ncoftwi_t) (fi - prev_fi [rxi]);
	prev_fi [rxi] = fi;

	return d_fi;
}

...

		{
			// Демодуляция NBFM
			/*const FLOAT_t gain = */agc_forvard_float(dspmode, vp0f, rxi);	// для отображения S-метра
			saved_delta_fi [rxi] = demodulator_FM(vp0f, rxi);	// погрешность настройки - требуется фильтровать ФНЧ
			//const int fdelta10 = ((int64_t) saved_delta_fi [rxi] * ARMSAIRATE * 10) >> 32;	// Отклнение частоты в 0.1 герц единицах
			// значение для прослушивания
			// 0.707 == M_SQRT1_2
			const FLOAT_t sample = saved_delta_fi [rxi] * (FLOAT_t) M_SQRT1_2;
			BEGIN_STAMP2();
			r = sample * rxoutdenom;	// масштабирование к разрядности аудио-кодека 1_31 -> 1_15
			END_STAMP2();
			if (glob_nb_off == 0)
			{
				// "шумодав"
				testholdmax3(iir_nfmnbbpf(sample));
				const int nbopen = getholdmax3() < nbfence;
				r = (nbopen != 0) ? r : (r / 16);
			}
		}

[EU1SW]

// Angle estimate var a = (float)Math.Atan2(f. Imag, f.Real); // Scale audio[i] = a * FMGain; // FMGain = 0.00001f;

не вижу дифференциатора
вот гляньте на заведомо рабочий код, вычисляется разница фазы между соседними семплами, скачки фазы между семплами больше чем на PI должны быть устранены, реально это работает очень хорошо, проверено
вместо argof используете ессно своё

float phase, dphase;



    for (int i=0; i<input_size; i++)



    {



        phase=argof(input,i);



        dphase=phase-last_phase;



        if(dphase<-PI) dphase+=2*PI;



        if(dphase>PI) dphase-=2*PI;



        output[i]=dphase/PI;



        last_phase=phase;



    }



    return last_phase;

[Genadi Zawidowski]

не вижу дифференциатора

Точно!

[alex_m]

вот гляньте на заведомо рабочий код

Что должно возвращать argof(input,i) ?

попробовал подставить значение Q, не работает, вот результат:

[EU1SW]

Фазу, ессно... У вас же в явном виде в коде аналогично все идет // англ естимейт
Или просто цельнотянутое, не вникая? )))

[alex_m]

Фазу, ессно... У вас же в явном виде в коде аналогично все идет

Чтото не понял. Код, который вы привели - это демодулятор или только часть от него?
Если только часть, то фазу от чего на вход подавать? Брать из значений входного iq потока?

И что такое output[i]? Это выходной аудиопоток или массив значений Q?
Если массив значений Q, то что дальше делать с обработанным таким образом IQ потоком?

[Genadi Zawidowski]

atan2(I, q) дает текущую фазу. вычитаем из предыдущей - получаем скорость изменения фазы (т.е. частоту). Нормируем. Изящный способ решения проблем с необходимостью оставаться в пределах 0..2 пи предложен по ссылке, которая оставлена в моем исходнике. Хотя можно и явно проверять на выход и корректировать на 2пи - пока разбираетесь с алгоритмом.

А было вот так:

const s32 MAX_S32 = 0x7FFFFFFF;
const f32 PI = 3.141592653;
const f32 radians_to_int = MAX_S32/PI;


{
...

static s32 prev_fi;

f32 fi = atan2(Q,I);

s32 fi_int = (s32)(fi*radians_to_int);

s32 d_fi = fi_int - prev_fi_int;

prev_fi_int = fi_int;

So on, so forth...

Vladimir Vassilevsky

DSP and Mixed Signal Design Consultant

http://www.abvolt.com

Как результат в d_fi значение соответствующее частоте. MAX_S32 с плюсом или минусом соответствуют отклонению от нулевой састоты на 1/2 частоты дискретизации.
Для тестирования можете настраиваться на немодулированную несущую и наблюдать в d_fi значение.

[alex_m]

Теперь понятно. Переделал демодулятор, надеюсь что теперь правильно понял:

        public void Demodulate(Complex* iq, float* audio, int length)
        {
            var gain = FMGain;
            if (_mode == FmMode.Narrow)
                gain *= NarrowAFGain;
            for (var i = 0; i < length; i++)
            {
                // Polar discriminator
                var f = iq[i] * _iqState.Conjugate();

                // Limiting
                var m = (float)Math.Sqrt(f.Real * f.Real + f.Imag * f.Imag);
                if (m > 0.0f)
                {
                    f /= m;
                }

                // Angle estimate
                var phase = (float)Math.Atan2(f.Imag, f.Real);

                var dphase = phase - last_phase;
                if (dphase < -Math.PI) dphase += TWO_PI;
                if (dphase > Math.PI) dphase -= TWO_PI;
                dphase /= (float)Math.PI;
                last_phase = phase;

                // Scale
                audio[i] = dphase * gain;
                _iqState = iq[i];
            }
        }

Такой вариант работает, но не могу сказать что стало лучше.
Приаттачил примеры со старым вариантом и с новым.

[Genadi Zawidowski]

Выкинуть работу с _iqState, протестировать на немодулированном сигнале (получаем постоянную составляющую, пропорциональную расстройке). Выкинуть Limiting.