Усрднение значений фазы Опции

[amusin]

Есть измерения фазы в диапазоне 0...359 градусов.
Допустим, показания болтаются около отметки 0 градусов:
половина значений в первом квадранте, половина - в четвертом.
Среднее арифметическое даст примерно 180 градусов, что неверно.

Пока видится одно решение:
сохраняем массив измерений,
доворачиваем измерения настолько, чтобы сектор "болтанки" не лежал на угле 0 градусов.
Решение не нравится тем, что надо копить данные в ОЗУ.

Есть более красивые варианты?

Сообщение отредактировал amusin - Aug 20 2010, 06:38

[Oldring]

Есть более красивые варианты?

Есть множество вариантов осреднения, когда хранить в памяти не нужно. Ну например:

1. mean( unwrap( phi ) )

2. arg( mean( exp( 1j * phi ) ) )

[Fat Robot]

усрядняете cos(phi) и sin(phi). в конце - 4-х квадрантный atan отношения.

[amusin]

Спасибо большое! Алгоритм работы unwrap подойдет
http://www.mathworks.com/access/helpdesk/h...ref/unwrap.html

усрядняете cos(phi) и sin(phi). в конце - 4-х квадрантный atan отношения.

Спасибо, хороший вариант для ряда задач, но не в моем случае - задача на микроконтроллере.

[_4afc_]

Есть измерения фазы в диапазоне 0...359 градусов.
Допустим, показания болтаются около отметки 0 градусов:
половина значений в первом квадранте, половина - в четвертом.
Среднее арифметическое даст примерно 180 градусов, что неверно.

Пока видится одно решение:
сохраняем массив измерений,
доворачиваем измерения настолько, чтобы сектор "болтанки" не лежал на угле 0 градусов.
Решение не нравится тем, что надо копить данные в ОЗУ.

Есть более красивые варианты?

Вариант 1:

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

Если бы мне нужна была точность в 1 градус - я бы завёл массив из 360 счётчиков и увеличивал их необходимое время.
Далее хоть максимумы берите, хоть распределения и дисперсии смотрите.
-----------------------------------------------------------------------------------------
Вариант 2:

2.0 заводим 3 переменных S,N,F

2.1 изменить представление диапазона с 0...359, на -179...180
2.2 за 0 принять результат предыдущего измерения фазы (F)
2.3 преобразовать измеренное значение в представление 2.1+2.2

2.4 Цикл S=S+измеренное; N++.

2.5 F=F+S/N;

Но лучше ещё добавить скользящее окно.
Тогда к этим трём переменным добавится массив длинной в окно (кольцевой буфер), а N превратится в указатель текущего отсчёта в кольцевом буфере.

Соответственно S=S+измеренное-B[N];B[N]=измеренное;(N++)%длинна окна. F=F+S/длинна окна;

Для борьбы с накоплением ошибки - я бы добавил стекание типа F=F-F/1024;

PS: Никогда не усреднял фазу, но вроде должно сработать

[Oldring]

Вариант 2:

2.0 заводим 3 переменных S,N,F

2.1 изменить представление диапазона с 0...359, на -179...180
2.2 за 0 принять результат предыдущего измерения фазы (F)
2.3 преобразовать измеренное значение в представление 2.1+2.2

2.4 Цикл S=S+измеренное; N++.

2.5 F=F+S/N;

Но лучше ещё добавить скользящее окно.
Тогда к этим трём переменным добавится массив длинной в окно (кольцевой буфер), а N превратится в указатель текущего отсчёта в кольцевом буфере.

Соответственно S=S+измеренное-B[N];B[N]=измеренное;(N++)%длинна окна. F=F+S/длинна окна;

Для борьбы с накоплением ошибки - я бы добавил стекание типа F=F-F/1024;

PS: Никогда не усреднял фазу, но вроде должно сработать

Мне кажется, что если я после двух прочтений так и не понял, что именно вы хотите сделать - это всё слишком сложно, и поэтому, малореализуемо.

[_4afc_]

Для борьбы с накоплением ошибки - я бы добавил стекание типа F=F-F/1024;

Вернее S=S-S/1024; конечно же

Мне кажется, что если я после двух прочтений так и не понял, что именно вы хотите сделать - это всё слишком сложно, и поэтому, малореализуемо.

Не вижу сложностей в алгоритме, возможно вас пугает детальное описание - имелось ввиду следующее:

Я предлагаю работать не со значениями фазы, а с разницей между измеренным значением фазы (А) и значением фазы вычисленное в предыдущий момент (F) и выводимое например на индикатор.
Т.е. B=f(A,F). Функция f - не просто B=A-F. Она ещё приводит значение B в диапазон -179...+180.
Числа B усредняем скользящим окном. Накопление ошибок в скользящем окне устраняем стеканием.F - выводим на экран.

[fontp]

Вернее S=S-S/1024; конечно же



Не вижу сложностей в алгоритме, возможно вас пугает детальное описание - имелось ввиду следующее:

Я предлагаю работать не со значениями фазы, а с разницей между измеренным значением фазы (А) и значением фазы вычисленное в предыдущий момент (F) и выводимое например на индикатор.
Т.е. B=f(A,F). Функция f - не просто B=A-F. Она ещё приводит значение B в диапазон -179...+180.
Числа B усредняем скользящим окном. Накопление ошибок в скользящем окне устраняем стеканием.F - выводим на экран.

Это последовательное восстановление каждого следующего значения фазы по предыдущим по непрерывности
"Допустим что болтаются возле 0" - допущение произвольное. Обычно нет механизма, который может гарантировать
Если есть - то конечно, достаточно привести к интервалу (-pi pi)
В самом общем случае, если фаза не в окрестности нуля а где-угодно на окружности

if (fabs(fi0-fi) > PI) {// фазовый разрыв
if (fabs(fi0-fi-2*PI) <= M_PI) {
fi += 2*PI;
}else {
fi -= 2*PI;
}
fi0 = Усреднение(fi) // например интегратором с утечкой

Это очень простой последовательный алгоритм, однако по отношению к шуму хуже чем накоплять экспоненту,
а потом уже брать угол/ Последний алгоритм Олдринга - параллельный и по-существу оптимальный по отношению к шуму (второй Олдринга, первый Олдринга в общем случае не работает , только когда фаза болтается где-то возле 0, а не где-то возле пи) "

[Oldring]

первый Олдринга в общем случае не работает , только когда фаза болтается где-то возле 0, а не где-то возле пи) "

Чушь. Нужно только чтобы СКО было мало по сравнению с pi.
Что подразумевается, исходя из примененного автором темы термина "болтается около".

Если же СКО отсчетов фазы велико, что второй метод тоже не даст точного "осреднения фазы". Оба метода дают совпадающие результате лишь в пределе малого шума. Что там на самом деле нужно, "осреднение фазы" или "комплексное осреднение", зависит от решаемой задачи, но обычно для малого шума разница не принципиальна.

[fontp]

Чушь. Нужно только чтобы СКО было мало по сравнению с pi.
Что подразумевается, исходя из примененного автором темы термина "болтается около".

допустим,
mean( unwrap( phi ) )=0, если фаза болтается в окрестности pi и около

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

[Oldring]

Не вижу сложностей в алгоритме, возможно вас пугает детальное описание - имелось ввиду следующее:

Я предлагаю работать не со значениями фазы, а с разницей между измеренным значением фазы (А) и значением фазы вычисленное в предыдущий момент (F) и выводимое например на индикатор.
Т.е. B=f(A,F). Функция f - не просто B=A-F. Она ещё приводит значение B в диапазон -179...+180.
Числа B усредняем скользящим окном. Накопление ошибок в скользящем окне устраняем стеканием.F - выводим на экран.

Например, тут нужен некоторый механизм запуска, формирующий начальное значение F.
"Постоянно обновляемый индикатор" - это тоже домыслы. Человек вроде хотел простое средне по блоку отсчетов.
Да и описано сложно для понимания словами, а не строгим алгоритмом. Куда идет "В" вы так и не описали. Вы с ним собрались зачем-то "работать", но результат нигде не используется в вашем неформальном описании, нужно догадываться, что именно вы хотели сказать.

допустим,
mean( unwrap( phi ) )=0, если фаза болтается в окрестности pi и около

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

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

А именно чушь вы написали и еще раз. Достаточно запустить в Матлабе команду

Код

mean(unwrap(pi+1e-2*randn(1000,1)))

чтобы в этом убедиться каждому самостоятельно и независимо от моего субъективного мнения.

И еще раз повторю. Знать, какое среднее фазы, априорно совершенно не нужно. Нужно априорно знать, что разброс не чрезмерно велик.

[fontp]

CODE

mean(unwrap(pi+1e-2*randn(1000,1)))

чтобы в этом убедиться каждому самостоятельно и независимо от моего субъективного мнения.

И еще раз повторю. Знать, какое среднее фазы, априорно совершенно не нужно. Нужно априорно знать, что разброс не чрезмерно велик.

Посмотрел реализацию unwrap. Здесь был неправ, этот именно тот последовательный алгоритм, что я привёл выше. И так он и будет работать
Я думал, что это просто приведение к стандартному диапазону.
4afc на самом деле предложил примерно то же самое

[thermit]

amusin:
Есть измерения фазы в диапазоне 0...359 градусов.
Допустим, показания болтаются около отметки 0 градусов:
половина значений в первом квадранте, половина - в четвертом.
Среднее арифметическое даст примерно 180 градусов, что неверно.

Пока видится одно решение:
сохраняем массив измерений,
доворачиваем измерения настолько, чтобы сектор "болтанки" не лежал на угле 0 градусов.
Решение не нравится тем, что надо копить данные в ОЗУ.

Есть более красивые варианты?

Фазу привести в диапазон -180 ... 180
Если микроконтроллер c алу в n разрядов - удобнее положить например -(2^(n-1)) ... 2^(n-1) - 1 (-pi ... pi)

Усреднение fmean(n) = fmean(n-1) + ( f(n) - fmean(n-1) )*k

умножение на k можно заменить арифметическим правым сдвигом.

[Oldring]

Посмотрел реализацию unwrap. Здесь был неправ, этот именно тот последовательный алгоритм, что я привёл выше. И так он и будет работать

Про "распараллеливание" речь вообще не шла. У человека последовательный микроконтроллер, отсчеты он получает последовательно, и даже не хочет их все хранить для дальнейшей обработки. Человек искал именно последовательный алгоритм. Если же кому-то нужно именно хорошо распараллеливаемый алгоритм - unwrap можно привести к таковому, вычислив первую конечную разность и переставив индексы в двойной сумме, в которую он входит.

Усреднение fmean(n) = fmean(n-1) + ( f(n) - fmean(n-1) )*k

Это не усреднение, а фильтрация IIR фильтром первого порядка.

[thermit]

Oldring:
Это не усреднение, а фильтрация IIR фильтром первого порядка.

Усреднение - тоже фильтрация однородным фиром. Этот однородный фир проще заменить ииром без потери качества.
Хотя, есть любители преодолевать трудности...

[Oldring]

Усреднение - тоже фильтрация однородным фиром. Этот однородный фир проще заменить ииром без потери качества.
Хотя, есть любители преодолевать трудности...

С некоторой небольшой потерей качества.
И проще как раз не заменять фильтром, а просто суммировать, если не требуется постоянное плавное обновление отсчетов. Всё равно нужна сравнимая разрядность суммы.

[_4afc_]

Ребят. Я в научных терминах не силён. Но предложенный мной алгоритм - усредняет.
Код примерно следующий:

Код

#define L 16

signed short W[L]; //must be cleared
signed short F=0;  //0..359
signed int S=0;    //-179L...180L
signed int N=0;    //0...L-1
signed short GetPhaseDiff(signed short A)
{
static signed short Old;
signed short B;
B=A-Old;
if (B>180)  B=B-360;
if (B<-179) B=B+360;
Old=A;
return B;
}

void PhaseCalc(signed short A)
{
signed short B;
B=GetPhaseDiff(A);  //-179...180
//-----------------window filter of B
S=S-W[N]+B;W[N]=B;
N=N+1; if (N>=L) N=0;
S=S-S/65536;               //stekanie
if (S>(180*L)) S=180*L;   //perepolnenie
if (S<(-179*L)) S=-179*L; //perepolnenie
//-----------------F update
F=F+S/L;
if (F<0) F=F+360;
if (F>359) F=F-360;
}

На рисунке - красным идеальная фаза, синим - снятая датчиком, зелёным - отфильтрованная окном с длиной 16.

Эскизы прикрепленных изображений

 

Прикрепленные файлы

 PhaseFilter.zip ( 742 байт ) Кол-во скачиваний: 18

 

[Oldring]

Код примерно следующий:

Просто кошмар!
Во-первых, на многих микроконтроллерах int 16-битный, поэтому алгоритм будет как-то работать, только если число точек в окне не превышает порядка 150. Кстати, такого буфера окна автор вопроса хотел избежать.

Во-вторых, F - это кумулятивная сумма от среднего, но не само среднее. Зачем она вам вобще нужна?

В-третьих, любому программисту известно, что "стекание" будет работать только при превышении модулем S величины 65535.

В-четвертых, "переполнение" S может наблюдаться только в том случае, если каким-то образом поработает немного "стекание", а иначе оно возникнуть не может.

[_4afc_]

Во-первых, на многих микроконтроллерах int 16-битный, поэтому алгоритм будет как-то работать, только если число точек в окне не превышает порядка 150. Кстати, такого буфера окна автор вопроса хотел избежать.

Прикол в том, что результат GetPhaseDiff() уже можно фильтровать любым фильтром.
Собственно GetPhaseDiff() и есть ответ на вопрос заданный вначале темы.
В данном случае B отлично отфильтрован 16 точечным скользящим окном.
Переменная S может быть любой разрядности, было бы желание.

Во-вторых, F - это кумулятивная сумма от среднего, но не само среднее. Зачем она вам вобще нужна?

Я таких умных словей как "кумулятивная сумма от среднего" не знаю.
Если окно не нужно - его можно не использовать. Тупо складывайте значения B и делите на их количество - узнаете среднее изменение фазы за всё время наблюдения.

В-третьих, любому программисту известно, что "стекание" будет работать только при превышении модулем S величины 65535.

Обычно S у меня с фиксированной точкой как и В. Например 16.16.

В-четвертых, "переполнение" S может наблюдаться только в том случае, если каким-то образом поработает немного "стекание", а иначе оно возникнуть не может.

В реальной жизни чего только не бывает, при реализации алгоритмов ЦОС на конкретном железе.
Приведённый код - показывает идею, и как использовать результат GetPhaseDiff().

[Oldring]

Обычно S у меня с фиксированной точкой как и В. Например 16.16.

В коде написано int. Арифметика - целочисленная. Где вы там вообразили запятую - значения не имеет.

чего только не бывает, при реализации алгоритмов ЦОС на конкретном железе.

На работающем - не бывает.
На неработающем - лучше не работать, всё равно хорошо не получится.

[Stanislav]

Есть более красивые варианты?

Осредните сперва компоненты сигнального вектора по отдельности, а потом уж находите его фазу.
Способ "осреднения" можно выбрать, исходя из статистических характеристик сигнала и требуемой СКО погрешности оценки фазы.
Осреднение значения фазового угла после его вычисления, в общем случае, будет давать оценку, далекую от оптимальной в смысле СКО на заданной выборке.