Применение метода наименьших квадратов в лазерной корреляционной анемометрии Текст научной статьи по специальности «Математика»

УЧЕНЫЕ ЗАПИСКИ Ц А Г И

Т о м XI

19 8 0

№ 3

УДК 533.6.071.08:532.57:621.375.8

ПРИМЕНЕНИЕ МЕТОДА НАИМЕНЬШИХ КВАДРАТОВ В ЛАЗЕРНОЙ КОРРЕЛЯЦИОННОЙ АНЕМОМЕТРИИ

Р. Ф. Аврамченко

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

1. Лазерный анемометр с фотонной корреляцией, предназначенный для измерения скорости частиц, описан в ряде работ, например, [1, 2]. Вкратце принцип его работы следующий. Луч лазера с помощью расщепителя делится на два луча, пересекающихся в объеме, где требуется произвести измерение скорости частиц. Рассеянный движущимися частицами свет отчасти поступает на детектор фотонов — фотокатод регистрирующего устройства. Образующийся поток фотоэлектронов поступает на вход коррелятора - счетчик фотоэлектронов. Счетчик фиксирует число фотоэлектронов за следующие подряд друг за другом интервалы длительностью Т. На выходе коррелятора формируется сигнал, описываемый дискретной случайной функцией

где (7 — произвольная нормирующая величина, в работе полагается О — 5; 5 — полное число интервалов длительностью Т, обрабатываемых в корреляторе; — момент включения коррелятора; п (£в + 5 Т) — случайное число электронов, вылетевших из фотокатода в интервале (^в + яГ, + (я 1) Т)\ аналогично определяется и п (£в 4- йТ -Ь /гТ); &Г — временной интервал между выборками длительностью Т; & = 1, . . — 2.

Применение метода наименьших квадратов для оценки допплеровской частоты в сигнале связано, в первую очередь, с определением средних значений функций (1 )—К(ЬТ).

2. Будем рассматривать поток одинаковых с оптической точки зрения частиц, движущихся с одинаковой скоростью. Пусть за время от включения до выключения коррелятора через измерительный объем проходит 1 частиц. При этом условии влиянием концевых эффектов — включением или выключением коррелятора в момент нахождения частицы в измерительном объеме — можно пренебречь. Будем также предполагать, что все частицы пролетают через центр измерительного объема.

5—к—1

(О

Обозначим через ;а случайное время между вылетами последнего фотоэлектрона от а-й частицы и первого фотоэлектрона от а+1-й частицы (а=1, . .АО* Ограничимся в дальнейшем потоками, для которых вероятностью Р(£а^ктТ) можно пренебречь. При этом в интервалах межд'у частицами с вероятностью, близкой к единице, будет содержаться группа более чем из кт интервалов Т с нулевым числом фотоэлектронов и, следовательно, в (1) не будут содержаться произведения числа электронов, порожденных разными частицами. Это обстоятельство позволяет, во-первых, все частицы рассматривать автономно и сумму

5-/е-1 N Ьт Ьт

^ заменить двойной суммой ^ ^ , а во-вторых, сумму ^ заменить сум-

5=0 а= 1 5=0 5=0

со

мои 7 , так как сумма всех добавляемых слагаемых с близкой к единице вероят-

5=0

ностыо равна нулю.

Таким образом, функцию К(кТ) при сделанных предположениях можно представить в виде

N со

^(*7')=т2 2я(в«+*7’)л(8-+*г+*7')* (2)

а=15=0

где 0а—начало интервала длительностью Т, в который попадает первый фотоэлектрон от а-й частицы.

Обозначим через іа момент появления первого фотоэлектрона для а-й частицы. Для каждой частицы можно произвольно выбрать начало отсчета време-

д

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

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

тами вылета первого фотоэлектрона и пролета частицы через центр с максимальной амплитудой освещенности. При этом ва Н-74.

Тогда для отдельной частицы интенсивность (математическое ожидание) потока электронов из фотокатода {х (і) можно представить в виде [3]

__

\ь(І)=АеЄ 1а~ (1 4- СОБ со£), (3)

где Ае — максимальная амплитуда среднего тока; а — длительность импульса освещенности на уровне освещенности, равной е~1/2 от максимальной; со — допплеровская частота тока.

Введем в рассмотрение случайное время о£а (0<;5^а<7’), равное времени между моментом пролета частицы через центр измерительного объема и ближайшим перед ним моментом смены интервалов длительностью Т— моментом переключения счетчика фотоэлектронов. Для рассматриваемого здесь случайного потока частиц величина Ыа распределена равномерно на интервале [0, Т] с плотностью распределения /?йа

Если из іа вычесть 5*а, то для любой частицы, породившей хотя бы один электрон, можно написать неравенства

где /а^0 — целые случайные числа (/а Т — Ыа = 0а).

Вероятность /г(^а/о£а) того, что до момента 0а частица не породила пи одного фотоэлектрона, но в интервале [0а, 0а-|- Т) появился хотя бы один фотоэлектрон, определяется формулой

FVJbt J = e

] — е

' т

В формуле (2) случайные числа /г(0а+и п(Ьа-\- бТ кТ) статистически

независимы. Поэтому среднее их произведения равно произведению их средних значений. При этом, в частности,

(1а+з + 1)Т

й(0а+$Г) = п((/а+б')77Ц,)= | (а (( - Ыа) йг

(/а+^) т

АР Те

Wo.+W-ota.V 2 а* /

1 Ф (to Т) cos со

где Ф (со Т) = sin

Ло7Л 2

I 2 )(й7-

Приведенных предпосылок достаточно для определения среднего значения функции К (kT):

К

СО

N V

Среднее ПО и и bta

N

со оэ 7*

v f /1 ((Za + S) Tlbta)n ((/„ + s-^rk) Tjbta) ^(*«/&*,)

d (o&t)

T

T l(le+J)2T+[(i«+J+*)n*

V V f e *?

NA; T 00

~л, *

/ T /a

a — s—„ (1 +COScof) (it 2«2

.A%1)T ~^(1+cos,at)dt

X\e

— e

X

X jl+Ф (“7’) COS to |/a + s 4- —

T — bt.

X

X ^1 4- ф (©Г) cos O) ^ ^/a -f 5 -1- k + -y- j Г— otaJ d (ota).

Влияние оta на экспоненты в последнем равенстве пренебрежимо мало, и эта величина была опущена в показателях экспонент. Интеграл но bta произведения двух последних сомножителей в скобках дает в результате

1 4- — Ф2 (со Т) cos tok Т 2 v '

Т

плюс совокупность слагаемых, содержащих синусы с аргументом со (/а + $) Т. Сумма по любому из индексов 1а или 5 произведения экспоненциальных функций на знакопеременную функцию синуса существенно меньше суммы одних лишь экспонент.

Поэтому для К (kT) можно написать более простое выражение

К (Л Т)

NA* Т* S

_ Ша +s)*+(la+s+k)*\ Т*

е

Jmd jemi

S=2О la = —оо

со

V V

2 а*

X

(1а +1) Т

- j' ^{t)dt - f ij.{t)dt

— OO —CO

Є — Є

1 4- — Ф2 (w T) cos ink T ].

Сделаем теперь, во-первых, замену переменной == /?га — 5, во-вторых, замену суммирования интегрированием (для экспоненциальных функций при условии а > 2тс/со, которое на практике обычно реализуется, такая замена не дает заметной ошибки) и, наконец, представим разность экспонент в круглых скобках в виде производной.

Используя (3), имеем

U Т (*а-И) Т

— J }x(t)dt

(ma—s)T--------

„ , 0/>2

j ;j.( t)dt ~^e J e (l + cosco/1) dt

— 00 —OO

= Є

(ша-(S-І)) r________(1_

■A _ Г г ‘2a2

(1 + cos wt) dt

Заметим, что в (4) гармонические знакопеременные составляющие в процессе суммирования будут осредняться. Поэтому в последнем выражении функцию cos о)/ можно без существенной потери точности отбросить. Тогда приближенно получим

— А

(ma—s)T _

I* 2а“ ы

\ е dt

(«в-и-1))Г t»

( 2а2

J е

—со

(ma -s)T

dt

d

ds

-A

r-

2 a?

dt

После замены сумм в (4) интегралами и после интегрирования по переменной 5 функцию К{кТ) приведем к виду

К (k Т) =

Т I

ma+(ma+/e)2) Т* 2а?

®a Т

-к f

1 — е

dma х

I 1 4~ Ф2 Т) cos со.ЛТ

(5)

Максимум первой экспоненты в подынтегральном выражении (5) достигается к

при тя=з:: — . Для практически наиболее интересного диапазона значений

-А,_ 1 в *»'

кТ в окрестности этого максимума экспонента £ 03 при сред-

нем числе фотоэлектронов от одной частицы порядка нескольких десятков будет намного меньше единицы и поэтому не окажет существенного влияния на функцию 1( (кТ). После вычисления интеграла в (5) окончательно получим

К (k Т)

ni NT

2 у it aS

{кту

4 а-

1 4 Ф2 (<о Т) cos & к Т).

(6)

где ne = V2ъаАе — среднее число фотоэлектронов от одной частицы.

Величина Ф(соГ)<1. Поэтому функция /<(/г Т) — знакоопределенная колебательная функция с экспоненциально убывающими максимальными и минимальными амплитудами. В дальнейшем наряду с &71 будет применяться и обозначение = & Т, а вместо я2 Л/Т/2 У%. аЗ — К*. Кроме того, возможность получения аналитического решения связана с предположением о том, что Ф(а>Г)^]. Поэтому далее будем полагать, что со 1.

3. В соответствии с методом наименьших квадратов функция

.2 2

'Г(“. Ко) = ^ ( /<(хА) — К0е № (

СОБ

/г = 1

дифференцируется по неизвестным параметрам сигнала со и Ко и производные приравниваются нулю:

_^1

дК

дЧ

дш

гп

= 2КоУ \КЫ)-Кое

/г=1

2с2

\

1 4- ---------- СОБ СОХа

2 /г

/ / 1 \ ХГ^(1 + ТС05“^-

[г ^I ! \ \

- = -2^ уКЫЗ-Шое -4с*-^ + — С08<0ТЙ )

т/г \

У «По>ТЛ 1=0,

4 с1

X

4с’ X

С ОБ сотЛ | = 0

С тсМ

Здесь использована замена а = == —— , где М — число интерференцион-

ных полос в измерительном объеме с уровнем освещенности, не меньшим £'—1/2 от максимальной. Обычно с>1.

Из последних двух уравнений получаем

т

X К

/г = 1

4с2

СОБ сот^

т

/г = 1

2 с2

СОБ сот£

£

/е=1

1г=\

т:? о)2 к

4 с2

СОБ СОТд

к = 1

X *

2 с2

V с3

-п- СОБ С0Т£

1 +

1 \ -у- соя сот/г I

2

"'к

4 с-

СОБ озтЛ

X

— Т/г БШ сотд, = О

(7)

Оценкой допплеровской частоты со в данной реализации сигнала на выходе

Л

коррелятора служит решение уравнения (7). Обозначим эту оценку через о>.

В дальнейшем с целью получения аналитического решения будут использоваться следующие два ограничения:

кт Г> (3 ~ 4) а, Т < тдв (8)

позволяющие с необходимой точностью заменить сумму У] несобственным ин-

2 к

тегралом; в (8) тд =—-—допплеровский период сигнала. При этих условиях

можно убедиться в том, что в (7)

Е

/г=1

9 А

со3 д

2 сл (і 1 ^ \ ш ^ А \ А

----------------- | і + _ С08 \ ^ 8Ш

І 4 “о" СОБ 0>Т£

9 А

СО2

к

2 С2

1 4

А \2 СОБ оз'

т V

/е=1

Поэтому (7) можно представить в более простом виде

б

2 Л

4 с3

/г=1

Г / 2 ^ \ _

/ т*(о / 1 А \ # Д

1 С* А / I 1 4 ^ С ОБ сот к 1 Т Б1П со и к

_ \ со ' \ /

= 0.

(9)

д д д

4. Найдем дисперсию оценки со. Введем обозначения со = со-|-5со и К(^и) =

= К (^) + (т*), с помощью которых из (9) получим уравнение для вариаций

А к =1

О со = — --------------------

4 с2

Х1 ш 1

— I | I — СОБ СОТ* | — 1к БІГ! ют*

4 с-

к =1

С*

----- БІП сот л — т? СОБ сох л

Є» со / 2 * 4 *

Используя при_вычислении знаменателя дроби (6) и (8), найдем, что он равен — т?/2 а3 /С*/4 У 2 Т, так что

А 4/2 71 Л

о со =

т?12аЧ<,

к =1

г6 4 а2 ( О (і + —

А с2 “ / 1 2

СОБ сот/г — тд, БІП со~д.

Обозначим произведение экспоненты на выражение в квадратных скобках через (т^). Тогда

(8ш)2 = £ £ 8/< (**> »/с ы т ы V ы.

к кт пг т

(10)

Кроме того, для краткости записи формулу (2) представим в виде

К оо

/СЫ = -±- £ £

Вй=15А=°

Можно убедиться в том, что

6/с (Т*) В/С (Ті) = /с (ТА) /с (т,) - к Ы к ы,

(11)

129

где

N N оо оо СРЄДНЄЄ ПО /в, Ыа И П

ЯШШ-4т Е XII "'л'л'ЛлЧ-,

“А=1“г=15й=0 *і=о

Выполним в (12) осреднение сначала только по числу фотоэлектронов п:

// СО

Е X

^=1 ^=0

Л/ со

Е Е

а. = 15.= 0

4*

N

“ V, (я? т - — я? «о _ я с я2 п8 - - —

£2 ^ \ -У* ^ ' V Л 7г *&т/г Ч

зк=0

£ к к к к і

+ (п п\ п т — п п2 п т _т ) V ^/г & /г 5ЛТЛ і ^т/г £ і У

(13)

При получении этого выражения ставилась задача разделения произведений случайных величин с субиндексом & от произведений с субиндексом /. С этой целью в общей сумме среднее от слагаемых с произведениями зависимых величин заменялось на произведения их средних значений. Так, например, вместо

п5 пз х п .П8-х. в случае, когда 5/ = $£, в общую сумму подставляется

^ ^ I II

— 9 —

п5 п т /2^ т-., а ниже эта величина вычитается. Подобная же операция выпол-

/г к к к I

нена и для случаев, когда 5/ = 5^ 4- $1 + ^ = 5/ 4- т/ = 4- тд,. После осред-

нения по числу частиц индексы и а/ опущены. Выражения в круглых скобках в (13) относятся к одной и той же частице, так как значения п для разных частиц статистически независимы. Раскроем эти выражения:

-2 -

п1 Пч х Пс т — Пч т П< - ) = ( Пс ' п, „ — п„ п, , п, _

* к 5к',к И' к' х 5 к єк'к 5/г к 'к * кхі

+ ( "** "V* ~ Л [%] Ч** + 15 N ПЧ--ь] =

Аналогично для второго выражения в круглых скобках находим: п _ , /2 _ _ х

Хя„ ; для третьего: /гс /г. , /г. _ ; для четвертого: яс /2 _ _ лс _ _ -4-

_ к ^ і _ ^ ЧГ”ТЛ 3к*1СЧ

+ %Ч^ + 'Ч”у*

Подставив эти выражения в (13), затем (13) в (12), (12) в (11), а (11) в (10), получим

32 Г-Ы

- /г

т оо

среднее ПО / , о/

.*=1 5Л=0

Л~п8п^х,) У2Ы)

кт кт со ,

среднее ПО 1а,

+ Е Е Е л .л

/г=1і = І5л=0^ к к к кі к к к к *1 к к к Ь I

Явные слагаемые в квадратных скобках обозначим по порядку через /ь /о,... ,/7. При их определении будем руководствоваться теми же преобразованиями и допущениями, которые имели место при выводе формулы (б). Для 1Х получим

кт ОО *«'

СО ОО СО Т

к =1 5^= О (/*+*)’

X

Хб

X I 1 + СОБ

2 а2

X (^1 +С05

и + 5 + ^ + ~2~) Т —

^а + 5+ 2 ) Т ) Х

(кТу^ г ^2

’*■ (-5—!г)х

(та-*)Г Г2

X (1 + \ СОБ свт^ — ~к 5?П сот*.

2 (А.

.г

2 а

2(И

*(К)

После повсеместной замены /а + 5 = /72а, зависящей от 5 окажется лишь выражение в фигурных скобках, интеграл от которого по 5 принимается равным 1, поэтому из приведенного выше выражения получим

А = 2 а;

оо

= 2 У71 аА2е | е

к= 0

00 к?Т СО

3 * а1

(/ Л = 0 та=-

Ш Г2 4дз \КкТу<» 1 ^

\ с* “ ^

= 2 /я яЛ2 ^ }/" 7С 3/3"

тп„ Т* * кТ>

а а (та Т)(1(кТ)

1 4- “ГГ СОБ со&Т' ^ — кТ БШ со& Т

У

а (кТ) =

, а3 °Ъ + О (

«2 * а\

3 ]/3

При вычислении этого интеграла использовано условие с2 > 1.

Интегралы /2 и /3 вычисляются таким же образом, причем оказывается, что Т

/2 = /о = О • Однако 7 « тд < я, так что полагаем /2 = /3 = 0.

Далее

ктп кт ОО

/4 = £Х £ п^п*^п*^‘ =л'г

’Л °£7г Л /

& =1г = 1 О

{кТ?+{1Ту

00 оо СО

II I ■

к=01 =0 тп п — —оо

тат*+(та+к)*Т*+(тл+1)*Т*

Та?

X

X е 4 а* (тЛ) х¥ (т 1) X I 1 + ~2~ соб со &Г + ~2 С08 + ~2 008 ш (^ ~ 0 ^ ) X

?Г=„лЗ <*> °° к*-4М+Р

X Л (таТ) (I (&Г) (I (г Г) =

12 а2

Г2

(£7") (/7*) БШ С0^7 81П со/7* X

Уъ к =о г-0

X (1 4- "§■ сое со&Г 4- “^т сое согТ 4- соя со (к — 1) Т^ (I (кТ)с1(1Т) =

— 3 со 00 £2—4£/ + /2

е Г Г Г 12 1кТ)ПТ) — со&<й(к — 1) ТХЛ(кТ)с1(1Т).

ъ «/ •/ 4

/3

/е=#0 /=0

Заменим в последнем выражении высокочастотный множитель cos2со (k — i)T его средним значением, равным 1/2. Кроме того, перейдем к полярным координатам: kT — р cosy, iT = р sin ср. Тогда из этого выражения получим

]/'2tz аА3 °~ ^~2 — Р2sin У cosУ) 8/3

к/ r v1---t siu у <.иэ

J f е 12 а* р3 sin ср cos <р dp dy = 11,4 а5.

е

р=0 <р=0

Аналогичным образом получаем

4 = 4 = —2,13 Л; а?\ In = 8,7 А. аъ.

Подставим /1? . . . , /7 в (14). Тогда относительное среднеквадратичное отклонение оценки допплеровской частоты будет определяться формулой

V

(5со)2 1,1

с У пеЫ

В работе [4] получена аналогичная формула для предельной точности оценки частоты со по реализации сигнала от одной частицы. При равных, но достаточно больших числах фотоэлектронов в сигналах эта формула отличается от (15) заменой коэффициента 1.1 на 1.

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

Решение вопроса о выборе более предпочтительного метода зависит от того, какие потери по точности связаны с применением метода наименьших квадратов. Формула (15) показывает, что потери эти равны 10%, т. е. достаточно малы, чтобы безусловно сделать выбор в пользу этого метода.

ЛИТЕРАТУРА

1. Pike F. R. The application of photon correlation spectroscopy to laser Doppler measurements. J. Phys. D: Appl. Phys., vol. 5, 1972.

2. BlagroveR. J., Foord R., Jakemen E. Determination of diffusion coefficients of haemocyanin at low concentration by intensity fluctuation spectroscopy of scattered laser light. „Nature*, vol. 227, 1970.

3. Г родзовский Г. Л. Анализ точности измерений ЛДИСа. „Труды ЦАГИ“, вып. 1750, 1976.

4. ГродзовскийГ. Л., Овсянников Р. Н., Шумил-к и н В. Г. Совместное оценивание параметров сигнала ЛДИС методом максимального правдоподобия. „Ученые записки ЦАГИ“, т. 10, № 5, 1979.

Рукопись поступала 161VII 1978 г.