Фаза сигнала

[phanis]

Здравствуйте уважаемые форумчане. Пытаюсь сделать фазовый пеленгатор. Но пока застрял на теории.
И сам вопрос.
Как определить фазу сигнала. Формула?. Имеется массив из выборок с линейного входа компьютера.

[rustot]

да там чистая тригонометрия.

по простому sin(w*t)*sin(w*t+fi) = cos(fi)/2 - cos(2*w*t+fi)/2, теперь если интегрировать за период, то второй член обращается в ноль, и остается только первый, косинус фазы в чистом виде. поскольку из косинуса вблизи нуля фаза определяется с большими погрешностями, то тем же способом на втором сумматоре можно получить sin(fi)

[phanis]

да там чистая тригонометрия.

по простому sin(w*t)*sin(w*t+fi) = cos(fi)/2 - cos(2*w*t+fi)/2, теперь если интегрировать за период, то второй член обращается в ноль, и остается только первый, косинус фазы в чистом виде. поскольку из косинуса вблизи нуля фаза определяется с большими погрешностями, то тем же способом на втором сумматоре можно получить sin(fi)

Спасибо за быстрый ответ, только я совсем не понял, как это будет работать. Наверное не правильно задал вопрос ..
Вот пример который мне предлагали, может из него станет понятно, то что мне нужно. А нужен угол, желательно в градусах.. Все таки пеленгатор.(будет указывать направление) Здесь пример на 180 точек.

Нужно взять 180 значений синуса и косинуса для углов от 0 до 358 с шагом 2. На самом деле там реально разных значений всего 45, остальное повторы и повторы со знаком минус. И посчитать сумму произведений входного сигнала с АЦП в 180 точках круга на соответствующие синусы и косинусы. Получим 2 числа:

A = sum(U_k*cos(Alfa_k)), Alfa_k = 2*k град., k=0..179
B = sum(U_k*sin(Alfa_k)), Alfa_k = 2*k град., k=0..179

Эти числа - ни что иное, как действительная и мнимая части первого коэффициента преобразования Фурье входного сигнала - соответствующего основной гармонике, т.е. частоте "вращения" антенного переключателя. По ним легко определить азимут - это просто фаза этой гармоники, то есть arg(A+i*B), считается обычно через арктангенс (его тоже задают таблицей на отрезке от -1 до 1):

если А>B и A>-B, то Phi = arctg(B/A);
если А<B и A>-B, то Phi = 90 град. - arctg(A/B);
если А<B и A<-B, то Phi = arctg(B/A) + 180 град.;
если А>B и A<-B, то Phi = 270 град. - arctg(B/A);

Попроще никак??

[rustot]

сдвиг фаз же относительно чего-то, допустим относительно сигнала который вы отправляете. вот этот сигнал обозначим через sin(w*t). одновременно сгенерим сдвинутый по фазе cos(w*t). обратно приняли сигнал sin(w*t+fi).

подаем sin(w*t+fi) и sin(w*t) на умножитель, суммируем полученное произведение в течение периода, получаем N*cos(fi)/2

подаем sin(w*t+fi) и cos(w*t) на умножитель, суммируем полученное произведение в течение периода, получаем N*sin(fi)/2

делим второе на первое и получаем tg(fi), тангенс того сдвига фаз что мы искали. если охота в градусах, то берете арктангенс и умножаете на 180/pi

но тангенс обычно более полезен. допустим тобы нарисовать на бумажке tg(fi)=0.9 не нужно транспортиров, 10 клеток вправо, 9 клетов вверх - получится линия под искомым углом

[phanis]

но тангенс обычно более полезен. допустим тобы нарисовать на бумажке tg(fi)=0.9 не нужно транспортиров, 10 клеток вправо, 9 клетов вверх - получится линия под искомым углом

Да вот теперь не много проясняется. Понимаю, хочется сказать нафига этим занялся, если даже основных тригонометрических понятий нет.
Во общем попробовал, получаются значения от 0 до 89.99999 и от 0 до -89.99999999 , а хотелось бы от 0 до 359.999999999
tg( N*sin(fi)/2 / N*cos(fi)/2) вернусь к клеткам, для меня они воспринимаются проще . Допустим если значения N*sin(fi)/2 N*cos(fi)/2 получатся со знаком минус, или хотя бы одно из них, то куда рисовать по клеточкам??

[rustot]

ну в обратную сторону. cos вправо, -cos влево, sin вверх, -sin вниз

[phanis]

Спасибо разобрался,и сделал именно так ..Все заработало правильно. Кстати там есть еще вариант, когда делитель равен нулю..А на ноль, лучше не делить. Это я тоже учел..
Вот еще вопросик..

Как можно округлять углы например если скачет в районе 350 - 10 градусов среднее должно быть примерно 0 градусов?? простая арифметика не подходит, например (350 + 10 )/2=180 совсем в другую сторону. Грубый пример но тем не менее.

[Shr1mpy]

Как можно округлять углы например если скачет в районе 350 - 10 градусов среднее должно быть примерно 0 градусов?? простая арифметика не подходит, например (350 + 10 )/2=180 совсем в другую сторону. Грубый пример но тем не менее.

Для простой арифметики тебе для начала надо было бы представить 350 как "-(минус) 10"

[phanis]

Все разобрался, спасибо всем кто помогал..

Усреднение сделал не по градусам, а средне арифметически получаемых Сумм синусов и сумм косинусов ,
Я не силен в обозначениях математики. Поэтому приведу результат работы в паскале.. Для режима стерео. Один из каналов.
Усреднение для 100 - получаемых массивов.


Sum(sin(U[0..N]* t))/ Количество сложений Sum(cos(Y[0..N] *t))/ Количество сложений

Код: Выделить всё

function TForm1.GET_PHASE( freq : real; input_buffer : array of smallint; input_buffer_size : integer): integer;
 var
 i: integer;
 SUMM_COS : REAL;
 SUMM_SIN : REAL;
 TEMP_UGOL: REAL;
 Temp_real : real;

 begin

  SUMM_COS :=0;
  SUMM_SIN :=0;
  TEMP_UGOL:=0;

   Temp_real :=(360/(const_freq/(freq+trackbar1.Position)));
  EDit3.Text:= floattostr((freq+trackbar1.Position));


  for i:=0 to  input_buffer_size div 4 -1 do
     begin

     SUMM_COS := SUMM_COS+( input_buffer[i*2] * cos( Gradus(i*Temp_real)) ) ;
     SUMM_SIN := SUMM_SIN+( input_buffer[i*2] * sin( Gradus(i*Temp_real)) ) ;

     end;

    SUMM_SIN:=SUMM_SIN/(input_buffer_size/2);
    SUMM_COS:=SUMM_COS/(input_buffer_size/2);


  GLOBAL_SUMM_COS[GLOBAL_COUNT] := SUMM_COS;
  GLOBAL_SUMM_SIN[GLOBAL_COUNT] := SUMM_SIN;
  GLOBAL_COUNT:=GLOBAL_COUNT+1;
  if (GLOBAL_COUNT>=100) then GLOBAL_COUNT:=0;

      SUMM_COS :=0;
  SUMM_SIN :=0;
    for i:=0 to  99 do
     begin

     SUMM_COS := SUMM_COS+GLOBAL_SUMM_COS[i];
     SUMM_SIN := SUMM_SIN+GLOBAL_SUMM_SIN[i] ;

     end;

      SUMM_SIN:=SUMM_SIN/100;
    SUMM_COS:=SUMM_COS/100;

  if ((SUMM_SIN=0) and (SUMM_COS=0)) then TEMP_UGOL := 0 else
  begin

  if  (SUMM_SIN=0) then
  begin
  if (SUMM_COS>0) then     TEMP_UGOL := 0;
  if (SUMM_COS<0) then     TEMP_UGOL := 180;
   end else

  if  (SUMM_COS=0) then
  begin
  if (SUMM_SIN>0) then     TEMP_UGOL := 90;
  if (SUMM_SIN<0) then     TEMP_UGOL := 270;
   end else
       begin
       TEMP_UGOL:=180/PI*  arctan(SUMM_SIN/SUMM_COS);
       if ((SUMM_COS > 0 )and(SUMM_SIN > 0 )) then TEMP_UGOL := TEMP_UGOL +   0;
       if ((SUMM_COS > 0 )and(SUMM_SIN < 0 )) then TEMP_UGOL := TEMP_UGOL + 360;
       if ((SUMM_COS < 0 )and(SUMM_SIN < 0 )) then TEMP_UGOL := TEMP_UGOL + 180;
       if ((SUMM_COS < 0 )and(SUMM_SIN > 0 )) then TEMP_UGOL := TEMP_UGOL + 180;
       end;

  end;




end;

 RESULT := ROUND(TEMP_UGOL);



 end;

   function TForm1.Gradus(c:real):real;
const
pi=3.1415926535897932384626433832795 ;
begin

 result:=c*pi/180-pi/2;
end;

[rustot]

вы значительно повысите точность, если урежете обрабатываемый буфер до точно кратного половине периода, отбросив лишний хвост. допустим если вы просуммируете данные за _ровно_ 1 период то получится более точное значение, чем если просуммируете за 100 с четвертью периодов

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

[phanis]

Вроде так и сделал. Частота дискретизации 44100 Гц размер буфера 256 (значений), значит частота должна быть кратная 44100/256=172,265625 Гц
К примеру если частота 1722,65625 Гц, то значит в одном буфере я имею 10 периодов синусоиды.

Далее получаю сумму синусов и сумму косинусов всего буфера , делим (сумму синусов и сумму косинусов) на количество периодов синусоид(в одном буфере здесь 10) и получаем среднеарифметическое средн(sin) средне(cos).



Среднеарифметические (sin) и (cos) с каждого буфера складываю необходимое количество раз к примеру 100(буферов) и делю на 100 получается среднеарифметическое значение уже на 1000 периодов нахожу угол через arctg(средее за 1000 периодов синуса / средее за 1000 периодов косинуса) (+/- sin +/-cos)для определения четверти окружности и все. Теперь даже если скачет в районе 350 10 градусов показывает 0 градусов.
В принципе деление это лишняя операция я так понимаю делить нет необходимости.? НЕ на 10 не на 100 просто складывать??
Или может я что то не до понимаю?
И Еще вопрос , из корень квадратный ( A^2+B^2), где А сумма синусов В сумма косинусов что получаем?

[B@R5uk]

Если заранее известно, что сигналом является чистый синус и обработка ведётся на компе, то нет необходимости подстраивать длину буфера под частоту сигнала и наоборот. Достаточно сигнал обработать методом наименьших квадратом взяв в качестве базовых функций синус и косинус частоты, которую, правда, перед этим необходимо определить. Но в любом случае это гораздо удобнее, чем подгонять количество семплов под частоту.

Ещё более продвинутый метод обработки может выдать относительную девиацию частоты (ну или фазы, не принципиально при малых значениях) на промежутке, а так же амплитуды, когда сигнал записан на достаточно большом промежутке времени (в сравнении с периодом дискретизации). Правда тут двумя базовыми функциями не обойдёшься. За то точность определения частоты и фазы (за счёт большой статистики) будет просто заоблачной, единственная погрешность, в ней содержащаяся будет только систематическая -- погрешность частоты кварца звуковой карты.

Я кстати такого рода частотометр хотел замутить для компа (чтобы непрерывно мерял частоту и дисплеил девиацию в реальном времени), просчитал всю математику и метод расчёта, даже немного DirectSound освоил, но так до конца и ниасилил :,-((