Помощь - Поиск - Пользователи - Календарь
Полная версия этой страницы: Усреднение периодической величины
Форум разработчиков электроники ELECTRONIX.ru > Цифровая обработка сигналов - ЦОС (DSP) > Алгоритмы ЦОС (DSP)
Brains
Доброго дня, форумчане!

У меня есть периодическая величина, значение которой нужно усреднить средним арифметическим (X1+X2+..+Хn)/n - в простейшем случае (Х1+Х2)/2.
Величина - это угол сдвига фаз, меняющийся от 0 до 359 градусов. Измерение ведется путем захвата таймером интервала времени между фронтами прямоугольных импульсов. При условии синфазности сигналов (должно измеряться значение 0 градусов), появляется неразрешимая пока мною проблема: из-за джиттера, являющегося следствием шума на входе, начинается гонка фронтов сигналов, при этом я получаю захваченные таймером интервалы, соответствующие углу в градусах, к примеру в такой последовательности: 1, 359, 1, 359.... Если попытаться усреднить эту последовательность, то получится значение 180 градусов, вместо ожидаемого нуля.
Вариант перехода к формату -180->0->+180 уже рассматривался, и он имеет такую же проблему при усреднении, возникающую в окрестностях 180 градусов (-179,+179 после усреднения дает 0 вместо 180).
Алгоритм нужен для целочисленного МК, желательно простой арифметическо-логический.
Какие есть варианты решения этой проблемы?
GetSmart
Вообще-то угол - это величина из двумерного пространства. Там числа должны быть двумерными/комплексными. Но математика там будет посложнее целочисленной и сложности эти нужны не всегда.

Проще можно так. Если конечно фаза не скачет по всему кругу. Перед прибавлением очередного значения к сумме, новое значение корректируется на -360 или 0 и выбирается ближайшее. Потом, результат можно так же скорректировать на -360, если он превысит 360 (360*N). Или на +360, если он будет отрицательным, а нужен положительный.
Brains
Цитата
Перед прибавлением очередного значения к сумме, новое значение корректируется на -360 или 0 и выбирается ближайшее.


GetSmart, извините, не уловил смысл того что нужно сделать. Можете пояснить на примере: что нужно сделать с цифрами, если сначала приходит значение 1, затем 359, чтобы получить результат 0?
MKS
Можно от угла перейти к декартовым координатам (x,y) и усреднять их а потом обратно к углу вернуться (это реализуемо в целочисленой арифметике). Или детектировать скачок фазы при разрыве и компенсировать его.
GetSmart
Цитата(Brains @ Nov 25 2011, 01:10) *
GetSmart, извините, не уловил смысл того что нужно сделать. Можете пояснить на примере: что нужно сделать с цифрами, если сначала приходит значение 1, затем 359, чтобы получить результат 0?

359 можно представить как -1, если прибавить к исходному значению -360. Вообще, можно прибавлять к углу хоть -360, хоть +360 и это будет тот же угол.

Так что, если первое число было 1, а второе 359, то смотрим разницу между ними = 358. А если к 359 прибавить -360, то разница между числами будет 2 (1 - (359-360)). А 2 меньше чем 358. Значит его и используем в среднем арифметическом. Соответственно, 1+(-1) = 0. 0 это и есть среднее значение. Разницу только нужно сравнивать по модулю (без знака).
Brains
Цитата
Можно от угла перейти к декартовым координатам (x,y) и усреднять их а потом обратно к углу вернуться (это реализуемо в целочисленой арифметике).


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


GetSmart, спасибо, я попробую этот метод. Если переложить в псевдокод получится:
Код

if(|Xn - Xn-1| < |Xn - 360 - Xn-1|)
      Acc += Xn - 360;
else
      Acc += Xn;
Xn-1 = Xn;

Понятно стало, что кроме аккумулятора и счетчика накопления потребуется еще переменная предыдущего отчета Xn-1, необходимая для предложенного вами кастинга. Но еще не понял как тогда ее инициализировать.
GetSmart
Цитата(Brains @ Nov 25 2011, 02:09) *
Понятно стало, что кроме аккумулятора и счетчика накопления потребуется еще переменная предыдущего отчета Xn-1, необходимая для предложенного вами кастинга. Но еще не понял как тогда ее инициализировать.

Во-первых, "кастинг" должен быть именно с аккумулятором. Во-вторых, сравнивать нужно модуль, то есть abs(acc/n - X), либо abs(acc-X*n). Обычно проще и быстрее умножать, чем делить.
В-третьих, нужно контролировать диапазон аккумулятора, чтобы он был диапазоне 0..360*n. Иначе корректировать на -360*n, либо на +360*n.
Ну и арифметика со знаком ессно.

А если после каждого добавление аккума алгоритм считает среднее значение от акка, то так будет даже проще. Защита от переполнения не понадобится.
Brains
Цитата
Во-первых, "кастинг" должен быть именно с аккумулятором. Во-вторых, сравнивать нужно модуль, то есть abs(acc/n - X), либо abs(acc-X*n). Обычно проще и быстрее умножать, чем делить.

Про это понял.

Цитата
В-третьих, нужно контролировать диапазон аккумулятора, чтобы он был диапазоне 0..360*n. Иначе корректировать на -360*n, либо на +360*n.


А эту фразу нет. Вы имеете ввиду, что после каждого прибавления (и увеличения n на единицу) в аккумуляторе должно быть число не больше чем 360*n ?



Цитата
А если после каждого добавление аккума алгоритм считает среднее значение от акка, то так будет даже проще. Защита от переполнения не понадобится.


Этого, к сожалению, достичь нельзя, нет ресурсов для деления на любое число, а потому рассчитываю n как степень двойки, подбираю разрядность аккума для худшего случая, и после накопления 2^n отсчетов делю аккум сдвигами.
GetSmart
Цитата(Brains @ Nov 25 2011, 02:28) *
А эту фразу нет. Вы имеете ввиду, что после каждого прибавления (и увеличения n на единицу) в аккумуляторе должно быть число не больше чем 360*n ?

Да. Причём акк может стать и отрицательным. И тогда его нужно увеличить на 360*n. Либо допускать отрицательные значения до -180*n. Тогда диапазон будет от -180 до +359(.9999). А каждое новое значение угла должно быть в диапазоне 0..359(.9999 если fixed point). И когда уже в акке накопится 2^k значений, и если акк отрицательный, то добавить к нему 360*(2^k) чтобы сдвигать исключительно положительные числа.
fontp
QUOTE (Brains @ Nov 24 2011, 21:47) *
Вариант перехода к формату -180->0->+180 уже рассматривался, и он имеет такую же проблему при усреднении, возникающую в окрестностях 180 градусов (-179,+179 после усреднения дает 0 вместо 180).
Алгоритм нужен для целочисленного МК, желательно простой арифметическо-логический.
Какие есть варианты решения этой проблемы?


Можно разрыв арктангенса сделать на противоположной стороне круга по отношению к текущим данным, т.е.
например к первому отсчету
fi0 + 180 % 360
а все другие значения приводить уже к этому новому интервалу (который определяется точкой разрыва) и просто суммировать.

Если есть такие сильные помехи, что выбивают угол больше чем на 180 градусов, то от них всё равно нельзя избавиться любым способом
Alexey Lukin
Есть способ лучше: построить гистограмму углов, сгладить её круговой свёрткой и найти максимум. Он даст "основной угол". А затем все значения углов из несглаженной гистограммы суммируются относительно "основного угла".

А с аккумулятором будет нестабильно: последовательность типа {0, 0, 90, 180, 270, 0, 0, 90, 180, 270, ...} будет ваш аккумулятор водить кругами... Впрочем, это, наверное, редкий на практике случай.
GetSmart
Цитата(fontp @ Nov 25 2011, 12:54) *
Можно разрыв арктангенса сделать на противоположной стороне круга по отношению к текущим данным, т.е.
например к первому отсчету
fi0 + 180 % 360
а все другие значения приводить уже к этому новому интервалу (который определяется точкой разрыва) и просто суммировать.

Воруют ©
Перельман и китайцы sm.gif

PS. Это там деление, или мне показалось? sm.gif
fontp
QUOTE (GetSmart @ Nov 25 2011, 12:53) *
Воруют ©
Перельман и китайцы sm.gif

PS. Это там деление, или мне показалось? sm.gif


Это там остаток по модулю 360.
Сначада интервал был 0-360. Вычисляем то что там написано Х. В этой точке делаем на круге разрыв.
Т.е. все последующие значения фазы приводятся к интервалу [X-360, Х] и тупо суммируются.
Здесь, есть конечно неустойчивость к забитости первой точки помехой - вдруг она ляжет криво по отношению к среднему. Собственно разрыв на круге нужно сделать бы в точке среднего, но оно нам как раз неизвестно заранее

Поэтому действительно как-то так
QUOTE (Alexey Lukin @ Nov 25 2011, 11:17) *
Есть способ лучше: построить гистограмму углов, сгладить её круговой свёрткой и найти максимум. Он даст "основной угол". А затем все значения углов из несглаженной гистограммы суммируются относительно "основного угла".

GetSmart
Цитата(fontp @ Nov 25 2011, 18:16) *
Это там остаток по модулю 360.

То есть остаток по модулю имеем без деления? Как-то с полуслова Вы меня не понимаете.
ТС-у делить в напряг. Он же написал.


Цитата(fontp @ Nov 25 2011, 18:16) *
Сначада интервал был 0-360. Вычисляем то что там написано Х. В этой точке делаем на круге разрыв.
Т.е. все последующие значения фазы приводятся к интервалу [X-360, Х] и тупо суммируются.
Здесь, есть конечно неустойчивость к забитости первой точки помехой - вдруг она ляжет криво по отношению к среднему. Собственно разрыв на круге нужно сделать бы в точке среднего, но оно нам как раз неизвестно заранее

Я так и не понял, чем это отличается от предложенного мной?

Цитата
Поэтому действительно как-то так

А свёртка это вообще жесть. Причём свёртка будет каждого элемента в гистограмме. Хотя результат вобщем достоверный. Но сложение векторов в декартовых координатах будет явно проще.
fontp
QUOTE (GetSmart @ Nov 25 2011, 16:49) *
Я так и не понял, чем это отличается от предложенного мной?


Где вы предлагали смещать интервал? Не вижу. Вы же рассматривали стандартный интервал [0, 360] или аналогичный стандартный [-180, 180],а я предложил убрать разрыв от данных подальше. Проблема же в точке разрыва.

QUOTE (GetSmart @ Nov 25 2011, 16:49) *
То есть остаток по модулю имеем без деления? Как-то с полуслова Вы меня не понимаете.
ТС-у делить в напряг. Он же написал.


Да нет там деления, это просто запись, там 2 if
GetSmart
Цитата(fontp @ Nov 25 2011, 19:26) *
Где вы предлагали смещать интервал? Не вижу. Вы же рассматривали стандартный интервал [0, 360] или аналогичный,
а я предложил убрать разрыв от данных подальше. Проблема же в точке разрыва.

в самом первом моём посте.
Цитата(GetSmart @ Nov 25 2011, 00:30) *
Проще можно так. Если конечно фаза не скачет по всему кругу. Перед прибавлением очередного значения к сумме, новое значение корректируется на -360 или 0 и выбирается ближайшее. Потом, результат можно так же скорректировать на -360, если он превысит 360 (360*N). Или на +360, если он будет отрицательным, а нужен положительный.

Это и есть взять ближайшее в окресности +-180 от текущего.
Сразу оптимизированно без делений и модулей. А потом контролировать текущий аккум как раз чтобы не приходилось делать лишние деления. Как и просил ТС
Цитата
желательно простой арифметическо-логический.
fontp
QUOTE (GetSmart @ Nov 25 2011, 16:49) *
А свёртка это вообще жесть. Причём свёртка будет каждого элемента в гистограмме. Хотя результат вобщем достоверный. Но сложение векторов в декартовых координатах будет явно проще.


А почему жесть? Если напрягают вычисления можно взять очень грубую гистограмму. Можно и свертку не делать. На этом этапе нужно всего лишь определить среднее грубо, чтобы знать где окружность разрезать

QUOTE (GetSmart @ Nov 25 2011, 17:30) *
в самом первом моём посте.
Это и есть взять ближайшее в окресности +-180 от текущего.


Если так - то извиняйте. Если это так, то значит не понял Вашу довольно мутную фразу

QUOTE (GetSmart @ Nov 25 2011, 17:30) *
Сразу оптимизированно без делений и модулей. А потом контролировать текущий аккум как раз чтобы не приходилось делать лишние деления. Как и просил ТС


Уже объяснял, что нет никаких там делений. Взять по модулю 360 - это 2 if (в контексте + 180 даже 1). И где бы окружность не разрезать, приведение к любому интервалу это 2 if + пара арифметических операций
GetSmart
Цитата(fontp @ Nov 25 2011, 19:35) *
Если это так, то значит не понял Вашу довольно мутную фразу

Бывает. Написав сразу оптимизированный алгоритм, забыл раскрыть тему в камментах. На английском sm.gif


Цитата(fontp @ Nov 25 2011, 19:35) *
А почему жесть? Если напрягают вычисления можно взять очень грубую гистограмму. Можно и свертку не делать. На этом этапе нужно всего лишь определить среднее грубо, чтобы знать где окружность разрезать

Это непонятно. Определить нужно минимум сглаженной функции на круговой гистограмме. И если делать свёртку, то кол-во только умножений будет (кол-во всего углов)*(кол-во элементов гистограммы). Причём реальный минимум может оказаться там, где элемент на гистограмме пустой. А как это всё оптимизировать и угрубить без ущерба - не очень ясно. Возьмём например 5 элементов 0,90,180,270 и 45.
Alexey Lukin
Сложность свёртки с окном Ханна — 7 умножений на отсчёт, независимо от размера фильтра. Число отсчётов равно числу бинов в гистограмме, а не числу углов. Свёртка делается один раз, для окончательно построенной гистограммы. Если максимум на пустом элементе — ничего страшного, на точности алгоритма это не скажется, в т.ч. для ваших пяти элементов.

Цитата(fontp @ Nov 25 2011, 18:35) *
Если напрягают вычисления можно взять очень грубую гистограмму. Можно и свертку не делать.

Можно и не делать, но лучше сделать — хотя бы с прямоугольным окном, покрывающим где-нибудь 1/3 всей гистограммы. Сложность — 0 умножений, 2 сложения на бин гистограммы.
AndrewN
QUOTE (Brains @ Nov 24 2011, 22:47) *
У меня есть периодическая величина, значение которой нужно усреднить средним арифметическим (X1+X2+..+Хn)/n - в простейшем случае (Х1+Х2)/2.

На мой взгляд, проблема находится в математической постановке задачи. Формулировка "есть периодическая величина" или неполна (точнее "есть периодическая случайная величина") или вообще величина неслучайная, т.е. функция.

Усреднять можно всё что угодно, например синус - тоже 2-пи периодическая функция, но она не имеет разрыва в +/-пи/2. Пила (рэмп) имеет разрыв в 2 пи. Но усреднение само по себе устраняет информацию о разрывах - как любая интегральная оценка.

В первом варианте (случайная величина) среднее просто должно быть поставлено в зависимость от _желаемой_ (известной до усреднения) точки разрыва, а значение случайной величины в окрестности точки разрыва нужно модифицировать так, чтобы устранить разрыв. Естественно, что разрыв нужно ещё как-то разумно выделить - например, по величине "производной" (конечно, это не производная), т.е. по величине конечной разности.

Например, используются N значений для усреднения, тогда нужно вычислить (N-1)+(N-2)+...+ 1 = N*(N-1)/2 разностей и проверить каждую на превышение заранее заданного порога, и если он превышен, то из большей величины вычесть период. Причём нужно учитывать, что для остальных пар значение разностей может измениться!

Плохое свойство периодических (или ограниченных) случайных величин в том, что дисперсия ограничена интервалом и не всегда можно наверняка выделить разрыв. Например, в генераторах равномерно распределённых псевдослучайных чисел (от 0 до M) никто же не пытается искать разрывы, хотя их принцип действия (алгоритм) очень от разрывов зависит, например в конгруэнтных генераторах.

Если дисперсия много меньше интервала, то разрыв можно легко опознать и устранить. Если наоборот, дисперсия сравнима или "больше" (она, конечно, не может стать больше, это некая натяжка, когда величина меняется сразу на несколько периодов) то сложно - невозможно! надёжно выделить разрыв.
Alexey Lukin
Цитата(AndrewN @ Nov 25 2011, 21:19) *
Формулировка "есть периодическая величина" или неполна (точнее "есть периодическая случайная величина") или вообще величина неслучайная, т.е. функция.

Под "периодической величиной" автор темы имеет в виду не периодический по времени сигнал, а закольцованную область значений сигнала. При этом отсчёты сигнала по времени могут быть случайны и их порядок не важен.
AndrewN
QUOTE (Alexey Lukin @ Nov 25 2011, 21:31) *
Под "периодической величиной" автор темы имеет в виду не периодический по времени сигнал, а закольцованную область значений сигнала. При этом отсчёты сигнала по времени могут быть случайны и их порядок не важен.

Числа появляются в процессе измерения. Они вполне себе могут быть или коррелированы или просто функция времени.

Неважность порядка я учел в количестве пар из N измерений, если учитывать порядок, то достаточно (N-1) пары.

И, кроме всего прочего, понятие разрыва имеет смысл всё-таки для функции, но не для "истинно" случайной величины...
AndrewN
QUOTE (Brains @ Nov 24 2011, 21:47) *
Алгоритм нужен для целочисленного МК, желательно простой арифметическо-логический.

Очень много хороших идей было высказано, что в итоге подвело меня к выводу, что самое вычислительно легкое (считая, что точность данных измерения - 1 градус, соответственно, таблица синусов в 90 входов) решение - усреднять тригонометрические функции. Они непрерывны, в отличие от угла. Осталось доказать, что среднее косинуса равно косинусу среднего угла.
GetSmart
Цитата(Brains @ Nov 25 2011, 02:28) *
Про это понял.

Я кое-что не учёл.
Кроме "кастинга" со смещениями -360 и 0 нужно ещё +360. То есть из трёх углов берётся ближайший к текущему.

Это если первым (или текущим) будет угол >180, а новый появится допустим 0. В этом случае брать нужно угол 0+360.
И диапазон аккума лучше наверное держать в положительных числах.

Upd.
Лучше даже так.
Держать диапазон акка в -180..+179. А диапазон новых углов 0..359. Тогда "кастинг" будет только с двумя значениями - реальным и сдвинутым на -360.
Это значит, если первым был угол 200, то акк сразу скорректирует его в -160.
Alexey Lukin
Цитата(AndrewN @ Nov 26 2011, 00:09) *
Осталось доказать, что среднее косинуса равно косинусу среднего угла.

Удачи! biggrin.gif
Brains
Спасибо всем участникам обсуждения темы!

В обсуждении появлялись фразы о том, что нужно закладываться на случаи разбега фазы по всем квадрантам, или погрешности измерения в 1 градус.
Для того, чтобы прояснить детали задачи, я нарисовал две картинки, объясняющие процесс:

Нажмите для просмотра прикрепленного файла Нажмите для просмотра прикрепленного файла

Джиттер в 1 градус был указан мною только для примера. На самом деле достигается погрешность 0,1 градуса. Но дело не в этом: какой бы мизерной ни была бы погрешность, всегда в случайном процессе (шум) возникнут условия "гонки".

Особая благодарность GetSmart! Ваш метод буду пробовать "в железе".
qxov
Можно перемножить два сигнала и считать скважность результата.
Для просмотра полной версии этой страницы, пожалуйста, пройдите по ссылке.
Invision Power Board © 2001-2025 Invision Power Services, Inc.