Пакет анализа данных GLESP для карт реликтового излучения на полной сфере и его реализация в рамках системы обработки FADPS Текст научной статьи по специальности «Математика»

© Специальная астрофизическая обсерватория РАН, 2005

Пакет анализа данных GLESP для карт реликтового излучения на полной сфере и его реализация в рамках системы обработки FADPS

О.В. Верходанов1, А.Г. Дорошкевич2, П.Д. Насельский3’4, Д.И. Новиков5, В.И. Турчанинов6, ИД. Новиков2,4,7, П.Р. Кристенсен4, Л.-И.Чианг4

1 Специальная астрофизическая обсерватория РАН, Нижний Архыз, 369167, Россия

2 Астрокосмический центр ФИАН, Профсоюзная 84/32, Москва

3 Ростовский государственный университет, Зорге 5, Ростов-на-Дону, 344090

4 Институт им, Нильса Бора, Блегдамсвай 17, DK-2100 Копенгаген, Дания

5 Империал колледж, Лондон, Великобритания

6 Институт прикладной математики им, Келдыша РАН, 125047 Москва

7 Обсерватория Копенгагенского университета, DK-2100, Дания

Поступила в редакцию 23.07.2004; принята к печати 10.09.2004-

Разработана новая схема пикселизации GLESP (Gauss-LEgendre Sky Pixelization) для построения карт реликтового излучения всего неба. Схема основана на вычислении интегралов в коэффициентах сферических гармоник методом гауссовых квадратур в нулях полиномов Лежандре и позволяет получать строгое ортогональное разложение карты на сферические гармоники. Разработаны соответствующие программы и проведено сравнение с другими методами. Пакет программ GLESP создан в соответствии с базовыми принципами построения гибкой системы обработки данных FADPS. Описываются основные процедуры и структура пакета.

Ключевые слова: космология, космическое микроволновое излучение, наблюдения, а на-

лиз данных

GLESP PACKAGE FOR FULL SKY CMB MAPS DATA ANALYSIS AND ITS REALIZATION IN THE FADPS DATA PROCESSING SYSTEM, by О. V. Verkho-danov, A. G. Doroshkevich, P. D. Naselsky, D. I. Novikov, V. I. Turchaninov, I. D. Novikov, P. R. Christensen L.-Y.P.Chiang. A new scheme of sky pixelization GLESP (Gauss-LEgendre Sky Pixelization) is developed for CMB maps. The scheme is based on the Gauss-Legendre polynomials zeros and allows one to create strict orthogonal expansion of the map. A corresponding code has been implemented and comparison with other methods has been done. The package has been realized using basic principles of the FADPS data reduction system. Structure and main procedures of the package are described.

Key words: cosmology: cosmic microwave background — cosmology: observations — methods:

data analysis

1. Введение

Процесс анализа данных космического микроволнового излучения (CMB: cosmic microwave background) включает несколько шагов:

1) регистрация данных в формате временных рядов,

2) пикселизация,

3) преобразование карты — сферические гарМОНИКИ,

4) разделение компонент,

5) анализ статистики сигнала,

6) вычисление спектра С (£) - мощности сигнала в зависимости от номера гармоники (мультиполя) ,

7) оценка космологических параметров.

В данной работе мы рассмотрим шаги (2), (3) и (6).

Начиная с космического эксперимента СОВЕ, использующего так называемую квадрилате-

рализованную небесную кубическую проекцию (Quadrilateralized Sky Cube Projection, см. детали у Чана и О’Нила, 1976; О’Нила и Лаубшера, 1976; Грейсена и Калабретты, 1993), проблема пикселизации неба в задачах исследования реликтового излучения (РИ) вызывает большой интерес. Были предложены, как минимум, три метода пикселизации небесной сферы для карт СМВ: пикселизация по граням икосаэдра (Тег-марк, 1996), пикселизация IGLOO (Критенден и Турок, 1998, далее КТ98) и HEALPix1 по методу Горского и др. (1999) (включая последнюю модификацию в 2003). Рассмотрим два важных вопроса, уже упомянутых Тегмарком (1996): а) какой метод является оптимальным для выбора положений центров пикселов, форм и размеров, чтобы получить (наилучшим образом) компактное однородное покрытие неба пикселами с равными площадями, Ь) какой лучший путь вычисления интегралов сверток карт суммированием пикселов.

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

В этой работе мы сосредоточим внимание на задачах обработки данных на сфере и затем определим схему пикселизации. Напомним, что пикселизация данных СМВ на сфере является только частью основной проблемы, которая заключается в определении коэффициентов сферических гармоник при разложении сигнала как для анизотропии, так и для поляризации. Эти коэффициенты, которые мы обозначаем alm, используются в последующих шагах анализа измеряемого сигнала и, в частности, при определени спектра мощности C(£) для анизотропии и поляризации (см. обзор Хиво-на и др., 2002), а в ряде методов — для разделения компонент (Столяров и др., 2002; Насельский и др., 2003а) и анализа статистики фаз (Чианг и др., 2003, Насельский и др., 2003а,Ь, 2004; Коулс и др., 2004).

Здесь мы предлагаем особый метод вычисления коэффициентов aim■ Он основывается на так называемой квадратуре Гаусса и описывается во втором разделе. Для этой схемы пикселизации положение центров пикселов вдоль полярного угла (координата в) соответствует так называемым нулям квадратуры Гаусса-Лежандра и, как будет показано ниже, этот метод существенно увеличивает точность вычислений.

Таким образом, метод вычисления коэффици-

ентов aim диктует выбор метода пикселизации. Мы называем наш метод GLESP (Gauss-Legendre Sky Pixelization: пикселизация неба по Гауссу-Лежандру (Дорошкевич и др., 2003; Верходанов и др., 2003, 2004)). Для реализации подхода GLESP было создано соответствующее программное обеспечение, объединенное в одноименный пакет, который позволяет исследовать данные СМВ, включая определение спектров мощности анизотропии и поляризации, Ci, функционалы Минковского и другие статистики.

Эта работа посвящена описанию основной идеи метода GLESP, оценки точности различных шагов и финальных результатов, описанию пакета GLESP и его тестированию. Мы не обсуждаем в данной работе проблему интегрирования данных из временных рядов внутри пиксела ограниченного размера. Простейшая схема интегрирования по площади пиксела эквивалентна взвешиванию по отношению к центру пиксела. Пакет GLESP может делать это тем же методом, что и HEALPix ИЛИ IGLOO.

Кроме того, здесь мы описываем основные процедуры пакета и их взаимодействие. Пакет реализован с использованием принципов гибкой системы обработки данных (FADPS: flexible data processing system), действующей на радиотелескопе РАТАН-600 (Верходанов и др., 1993; Верходанов, 1997а).

2. Основные идеи и базовые соотношения

Стандартное разложение вариаций температуры на небе AT (в, ф) на сферические гармоники записывается следующим образом:

ОО m=i

AT (в,ф) = ЕЕ almYlm(в, ф)

i=2 m=

Ут(в,ф)

(2£ + 1) {£ — т)\ 4-7Г (£ + то)!

pm(x)e

гтф

(1)

(2)

где рт (х) — присоединенные полиномы Лежандра. Для непрерывной функции АТ (х, ф) коэффициенты разложения а£т выражаются как

/1 р 2п

((х J АТ(х,ф)Уе*т (х,ф)(ф, (3)

где У1т обозначет комплексное сопряжение Ует-Для вычисления интеграла (3) мы будем использовать квадратуру Гаусса — метод, предложенный Гауссом еще в 1814 г. и развитый позднее Кристоф-

х

x

/

Рис. 1: Взвешивающие коэффициенты Гаусса-Лежандра (wj) в зависимости от нулей полиномов Лежандра (xj = cos в^, являющихся центрами колец, используемых в GLESP (для случая N = 31). Положения нулей отмечены вертикальными линиями.

нии (3) представляет собой интеграл по полиному

x

(Пресс и др., 1992):

/ at(x, ф)ут(х, ф)<1х =

j-i

N

j=1

wj AT(xj ,ф)У1*т (xj ,ф) ,

(4)

где — собственная взвешивающая функция квадратуры Гаусса. Здесь взвешивающая функция = ш(х^), а также АТ(х^,ф)У^т(х^,ф) вычисляются в точках х^, которые являются сеткой корней полиномов Лежандра:

рИ (х )=0 , (5)

N

линома. Хорошо известно, что уравнение (5) имеет N нулей в интервале —1 < x < 1. Для метода Гаусса-Лежандра взвешивающие коэффициенты в уравнении (4) равны (Пресс и др., 1992):

2 /

wi = 1^я\.Рм(хз)Г2 , ®

где знак ' обозначает производную. Они могут

xj

процедуры “gauleg” (Пресс и др., 1992, пар. 4.5).

Трапецеидальные пикселы в GLESP ограни-

вф

в

xj = cosej. Таким образом, интервал —1 < x < 1 покрывается кольцами пикселов (число колец N, У1,-т(в, ф) = (—1) У т(в,ф),а1т = (— 1) а^ —т -(8)

Рис. 2: Положение цент,ров колец (х^ = со80^) в зависимости от номера кольца для двух схем пикселизации. НЕАЬРгх показан сполошной линией, ОЬЕЗР — штриховой. Рисунок демонстрирует случай когда N = 31.

см. детали в разделе 3). Угловое разрешение, достигаемое при измерении данных СМВ, определяет верхний предел суммирования в уравнении (1), £ < ^тах- Чтобы избежать ограничений Найкви-ста, мы используем число колец пикселов, равное N > 2£тах. Чтобы сделать пикселы в экваториаль-

ф

но квадратными, число пикселов Nф'ax в этом направлении выбирается Nф'ax « 2N. Число пикселов по другим кольцам, №ф, должно определяться из условия организации размеров пикселов равными размерам экваториальных с такой точностью, с какой это возможно.

На рис. 1 показаны весовые коэффициенты и положение центров пикселов для случая N = 31. На рис. 2 сравниваются особености пикселизаци-онных схем, используемых в НЕЛТРгх и ОЬЕЗР (см. раздел 4). На рис. 3 сравнивается форма пикселов и их распределение на сфере в эллиптической проекции для НЕЛТРгх и ОЬЕЗР.

В (1) коэффициенты а£т являются комплекс-АТ

кете СЬКЗР, построенном исходя из определения

АТ

£тах

АТ (в,ф)=^2агоУго(в,ф) +

£=2

£тах £

+ ЕЕ (а£тУ£т(^1 ф) + а1,-‘тУ£, — ‘т(@,Ф)) 1 (7)

i=2 т=1

ГД0

0

ИИ»»И»«Г«»1»П»Ч ,

'««МИМ*

ГМ'Ш ■ Mil k

ИЯ

mm

■■ ■ і|іі>

111»

P

P

*

І Тії .fill

—.■Htllf ПРДІ f

VPPPMK

риті

Рис. 3: Схематическое представление двух типов пикселизации на сфере: HEALPix (вверху) и GLESP (внизу). Различные цвет,а пикселов 'использованы, для выделения их

Таким образом,

АТ(в,ф) = —= ]Г Re(aefi)Pe (cos в) +

2 m=1

2 (і + m)!

x [Re(a£m) cos(mф) — Im(a£m) sin(m^j

где Pp (cos в) являются хорошо известными присоединенными полиномами Лежандра (см., например, Грандтптейн и Рыжик, 2000). В пакете GLESP мы используем нормализованные присоединенные полиномы Лежандра

fm(x) =

І2£ + 1 (£-m)\ 2 (£ + то)!

Pm(x)

(10)

где x = cos в, а в — полярный угол. Эти полиномы frn(x) могут быть рассчитаны с использованием двух хорошо известных соотношений. Первое из

них дает frn (x) для заданного m и всех £ > m:

fm(x) =.

и2 - 1

£2 -т2^-1

2£+1(£-1)2 - m2 2£ - 3 £2 - to2 Ie-'2 '

Это соотношение стартует с

(И)

(9) fm(x) =

(-l)m /(2m+ 1)Н

у/2 V (2Г77, - 1)!!

_ x2)m/2

(І — x )

/m m

/™+1 = .,л 2///

Второе рекуррентное соотношение дает /£т (х) для заданного £ и всех т < I:

\/(£- т -!)(£ + т + 2)/™+2(.г) +

2х(т + 1) +1

+ ,/1^2 ^ И +

+ у/(£-т)(£ + т + 1)/?(х) =0. (12)

Это соотношение стартует с тех же /£(х) и /0 (х), которые должны находиться с помощью (11).

Как уже обсуждалось (Пресс и др., 1992, пар. 5.5), первое рекуррентное соотношение (11) является формально нестабильным, если число итераций возрастает до бесконечности. К сожалению, нет теоретических рекомендаций, какая максимальная итерация может использоваться в квази-стабильной области. Однако это соотношение может быть использовано, когда мы заинтересованы в так называемом доминантном решении (Пресс и др., 1992, пар. 5.5), которое является приближенно стабильным. Второе соотношение (12) является стабильным для всех 1шт.

3. Свойства СЬЕБР

Следуя предыдущим рассуждениям, мы определяем новую схему пикселизации СЬКВР следующими пунктами I

• В полярном направлении х = сов9 мы задаем х2,] = 1,2,..^ как сеть корней уравнения (5).

• Каждый корень х^ определяет положение кольца с N2 центрами пикселов с ф-коодинатами

фг ■

• Каждый пиксел имеет вес (см. (6)).

Наши вычисления, реализованные для пиксе-лизационной схемы СЬКВР. сделаны при следующих условиях.

ординатных линий 9 и ф. Таким образом, с большой точностью пикселы трапецеидальные.

ф

це их число программа позволяет выбирать произвольно. Число пикселов зависит от £тах, принятого для обработки данных РИ.

•

число колец вдоль оси х = соб(9) должно быть

N > 2£тах + 1-

приблизительно квадратными, число пикселов вдоль азимутальной оси ф берегся равным N'фax = т1(2п/(9и + 0.5), где к = + 1)/2, а (9и =

0.5(9^+! - 9и-1 )■

•

определяется как БРгХе1 = (9их(ф, где (9и и (ф = 2п/Nф'x являются размерами пиксела в экваториальном кольце.

• Число пикселов N2 в кольце при х = х^

вычисляется как №ф = ш1;(27г^/1 — х2/Бр^е1 +0.5);

ring number

Рис. 4: Отношение размер пиксела/площадь пиксела на экваторе в зависимости от номера кольца в GLESP для числа колец N = 300 и N = 5000 соответственно.

• Так как чиело Nj отличается от 2fc, где k — целое, для вычисления быстрого преобразования Фурье вдоль азимутального направления мы используем код FFTW (Фриго и Джонсон, 1997). Этот код позволяет применять не только подход 2”, но также и другие числа оснований, и обеспечивает даже большую скорость вычислений.

По этой схеме размеры пикселов одинаковы в пределах каждого кольца и максимальное различие между пикселами различных колец составляет ~1.5% вблизи полюсов (рис. 4). Увеличение разрешения уменьшает абсолютную ошибку разности площадей, потому что неэквивалентность полярных и экваториальных пикселов пропорциональна N-2

На рис. 5 показана схема пикселизации для карт с высоким разрешением (например, £тах > 500), дающую почти равную ширину по dB для

Рис. 5: Размер пиксела вдоль полярного угла

(£тах = 250;.

Wp(в,ф)=22 Wp(£, т)У1т(в, ф)-

(14)

Соответствующая корреляционная функция (КТ98) для пикселизированного сигнала

(ATpATq) = J2 C(£)Wp(£,m)W*(£,m).

(15)

4.1. Точность оценки функции окна

Дискретность пикселизированной карты определяет свойства сигнала и ограничивает точность, достигаемую в любой пикселизационной схеме. Чтобы оценить эту точность, мы можем использовать разложение (КТ98):

АТтар(9, ф) = ^ БрАТрШр(9,ф) =

большинства колец.

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

мы пикселизации. Схема IGLOO, примененная к дискретам широты GLESP, даст слишком различные размеры пикселов. Пикселы не будут ни равно расположенными по высоте, ни однородными по площади, как того требует IGLOO.

4. Функция окна GLESP пиксела

Для применения схемы GLESP мы должны учесть влияние размера, формы и расположения пиксела на сфере на сигнал в пикселе и его вклад в спектр мощности C (£). Температура в пикселе дается как (Горски и др., 1999; КТ98):

ATp

/дпр

Wp(в, ф)АТ(в,ф№ ,

(13)

где ^р(9, ф) ^ функция окна для р-ого пиксела с площадью АПр. Для функции окна Шр(9,ф) = 1 внутри пиксела и Шр(9,ф) = 0 снаружи (Горски и др., 1999) из уравнений (1) и (13):

АТр = Е агтШр(£, т) ,

где

Wp(£,m)= Wp (в,ф)Ут,(в,фЦП

Еаттр У1т(в,ф),

(16)

map ____

а1т

= J dfiAT^(в^^ (в,ф) =

(17)

= £ Sp ATp Wp (£, m) , p

где Sp — площадь p-ого пиксела. Эти соотношения обобщают уравнение (3), учитывая свойства функции окна. Схема GLESP использует свойства интегрирования Гаусса Лежандра в полярном направлении, в то время как азимутальная пиксели™ зация для каждого кольца подобна схеме IGLOO, и мы получаем (см (4)):

Wp (£, m)

\/27гА-

■ exp

p

Np

Ф

sin | nm/Np

nm/Np

рХр+§.5Дхр

x fim(x)dx,

(18)

где Ахр = (хр+1 — хр-1 )/2 с хр р-ым узлом по

Гауссу-Лежандру, а Nф — число пикселов в азимутальном направлении. Этот интеграл может быть переписан как:

рХр+0.5Дх„

fm(x)dx ~

l + (-l)fc (fc + 1)

k+1

= ■ <»>

где fm (xp) обозначает k-ую производную в x = xp. Тогда для Axp ^ 1 получаем разложение (18):

W^{£,m) = W<?H£,m) ^1 + , (20)

w

p

X

Xp —0.5Дхр

-\/2тг

exp

sin nm

*/к

Np

Ф

ъ/Nl

fmx)

(21)

и W°(£, m) те зависят от Axp. Точность этой оцен-

5Wp(£,m) _ Wp(2)(£,m) - W^;(£,m)

(0)/

Wp(£,m)

(f " )m (Axp )2

Wp(0)(£,m)

24fi

(22)

В соответствии с последней модификацией НЕАЬР1х точность воспроизводства функции окна пиксела порядка 10-3. Чтобы получить такую же точность для Шр(£,т), нам необходимо иметь

Axp < 0.15

f т fi

(f " )l

(23)

Используя приближенную связь между функциями Лежандра и Бесселя для больших £ (Гранд-штейн и Рыжик, 2000) /т гс Jm (£х) получаем:

Ажр < 0Л5хр/\Jm(m + 1) , и для Axp ~ п/N мы имеем из (24)

(24)

SWp(£, то) 2

Wp(£, то) “

£

-От

N

тах —3

10"

(25)

получаем что явля-

Напрнмер, для N = 2£

SWp(£,m)/Wp(£,m) ~ 2.3

ется совершенно разумной точностью для £тах ~ 3000-6000.

5. Структура пакета GLESP

Пакет организован и развивается на двух уровнях. Первый, объединяющий функции Фортрана F77 и С, подпрограммы и процедуры-связки между разными языками, состоит из главных процедур: ‘signal’, которая обращает значения aim в карту на сфере; ‘alm\ которая вычисляет aim по данным карты; lcl2alm\ которая моделирует aim-коэффициенты для заданного набора Cf, и ’almBcF, которая вычисляет Ci для aim. Процедуры тестирования кода, управления параметрами, анализа на гауссовость (в том числе и

aim

и однородность распределения фаз, а также другие, тоже включены в пакет. Работа этих процедур основана на другом блоке подпрограмм, вычисляющих пикселизацию по Гауссу-Лежандру для заданного разрешения, пересчета угла на небе в номер соответствующего пиксела и обратно.

Второй уровень пакета содержит программы, реализующие процедуры первого уровня. В добавление к сказанному, имеются еще процедуры, вычисляющие аналитические функции — шаблоны, генерирующие карты со сферическими функциями, такими, как, например, У20, У21 и У22. Име-

aim

коэффициентов, конвертирования карты GLESP в формат HEALPix и обратно.

На рис. 6 демонстрируется организация пакета GLESP. Круг определяет зону влияния GLESP, основанную на библиотеке пикселизации. Она может включать несколько процедур и действующих программ. Основная программа второго уровня ‘cl_ тар ’, показанная как большой прямоугольник, взаимодействует с процедурами первого уровня. Эти процедуры показаны маленькими прямоугольниками и вызывают внешние библиотеки преобразования Фурье и вычисления полиномов Лежандра. Пакет читает и записывает данные как в виде ASCII-таблиц, так и в виде FITS-формата. Более 15-ти программ пакета GLESP работают в зоне GLESP.

К настоящему моменту в пакете также реализованы процедуры параллельных вычислений. Процедуры динамической визуализации с использованием библиотеки OPEN GL развиваются в Институте астрономии в Кембридже.

5.1. Вазовые принципы организации пакета GLESP

Пакет реализован в соответствии с идеологией системы обработки FADPS (Верходанов и др., 1993; Верходанов, 1997). Он удовлетворяет следующим принципам:

•

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

•

но от пакета.

•

ной строки и управляется внешними параметрами. Кроме того, она может работать в диалоговом режиме и в ряде случаев управляться также файлом настройки.

•

организован стандартным путем и оформляется в виде FITS- или F-форматов (Верходанов и Кононов, 2002) или ASCII-таблиц, доступных для других пакетов.

•

(вводить и выводить данные) с другими процеду-

p

х=х

p

2

Рис. 6: Структура пакета GLESP.

рами FADPS и базой данных CATS.

5.2. Основные операции

В пакете GLESP имеется четрые группы операций с данными:

1. Разложение карты по сферическим гар-

aim

2. Сглаживание карты гауссовой диаграммой направленности (сШтар).

3. Нахождение суммы/разности/среднего нескольких карт (difmap).

4. Скалярное умножение/деление (difmap).

5. Вращение карты (difmap).

6. Преобразование координатной сетки из галактических координат в экваториальные с соответствующим вращением карты (difmap).

7. Обрезание карты по заданным температурным пределам (mapcut).

8. Выделение области на карте с занулени-ем значений внутри/снаружи зоны (mapcut).

9. Вырезание сечений из карты (mapcut).

10. Подготовка стандартных шаблонов карт (mappat).

11. Перевод карты, заданной в ASCII-кодах в FITS (mappat).

12. Чтение таблицы точечных источников и наложение их на карту (mappat).

13. Распечатка значений в пикселах карты в ASCII-формате (mapcut).

14. Поиск минимального/максимального значения в пикселе (йг/тор).

15. Простейший статистический анализ карты (Л'фтар).

16. Поиск коэффициента корреляции двух карт (Л'фтар).

17. Оценка размера пиксела карты (пШ).

18. Построение изображений карт ^2Ау).

• Операции над коэффициентами а£т :

1. Синтез карты по данным а£т (сШтар).

2. Поиск сумм/разностей (difalm).

3. Скалярное умножение/деление (difalm).

4. Векторное умножение/деление (йг/о!т).

5. Добавление фазы во все гармоники (М-}а1т).

6. Вырезание заданной моды из гармоник (йг/о!т).

7. Расчет спектра мощности С (а1т2(11).

8. Расчет фаз гармоник (alm2dl).

9. Выборка гармоник с заданной фазой (alm2dt).

10. Сравнение двух наборов коэффициентов ает (скескаЬт).

11. Вычисление а£т для производных карт по обеим коодинатам (йаЪт).

• Операции со спектром мощности Сц :

1. Вычисление спектра мошности Сц (alm2dt).

2. Моделирование карты для заданного Сц (с12тар).

aim

Ci

• Операции с фазами гармоник ф1т и амплитудами | aim | :

1. Вычисление фаз ф1т (alm2dl).

2. Вычисление амплитуд |aim| (alm2dl).

aim

atealm).

4. Выбор гармоник с данной фазой (alm2dl).

5. Добавление заданной фазы ко всем гармоникам (difalm).

5.3. Основные программы

Сейчас реализованы следующие процедуры, организованы® как отдельные программы, работающие в частотно-пространственной или пиксельной области:

aim

коэффициентам.

aim

стар конвертирует карты из HF.AT.Pix формата в GLESP.

aim

aim

Ci

aim

Ci

dalm вычисляет первую и вторую производные

aim

difalm производит арифметические действия

aim

difmap производит арифметические действия над картами, а также делает координатные преобразования.

fzfig строит цветные изображения карт в формате GIF.

f2map преобразует карту GLESP в формат HF.AT.Pix.

aim

коэффициентов из FITS в F.

aim

коэффициентов из F в FITS.

mapcut обрезает амплитуды и вырезает координаты на карте GLESP, строит одномерные и двумерные сечения для моделирования наблюдений на РАТАН-600 для исследования РИ (Парийский, 2001).

mappat создает стандарные образцы карт, читает ASCII-данные для построения карт, читает положения точечных источников из ASCII-файлов, содержащих вывод базы данных CATS (Верходанов и др., 1997b).

psep производит фазовый анализ и разделение

aim

2003а).

5.4. Формат данных

Данные GLESP представлены в двух форматах,

aim

aim

деке, описывающий номера I и m-мод, соответственно, в том же формате, что и HF.AT.Pix. и ве-

aim

метра описаны по правилам записей с тремя полями как Binary Table для F-формата (Верходанов и Кононов, 2002).

Данные карт описываются тремя полями как Binary Table в F формате, содержащими вектор положений xi = cos в, вектор значений числа пикселов в каждом кольце Nа также массив температур в каждом пикселе, записываемый от Северного полюса.

6. Тестирование и точность кода GLESP

Три теста позволяют нам проверить код. Первый из них заключается в рассмотрении аналитических карт

У2’° = \/т^(3ж2 “1} ’

^2 1 = —\ ---X\Jl — X2 cos ф,

V 8п

/Тб" /-------

Уг -1 = —\ —ху 1 — х2 sin ф ,

12,2 = V^i^1 ~~ ^ cos(2<^’

У2’~2 = ~ ^ sin(2<^ ’

по которым вычисляются величины коэффициен-

aim

величин aim с точностью лучше 10—7.

Второй тест воспроизводит аналитическую карту AT(x, ф) = У1т(x, ф) для заданных aim. Оба эти теста проверяют точность вычисления карт и сферических гармоник независимо.

aim

AT(x, ф)

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

aim

исходные величины.

Точность кода можно оценить введением набора aim = 1 и восстановлением его после обращения. Этот тест показал, что, используя соотношение (11), мы можем воспроизвести введенные aim, с точностью ~ 10-7, ограниченной только точностью представления плавающего числа. Для соотношения (12) точность составляет ~ 10-5.

Рис. 7: Сравнение точности вычислений в пакет,ах НЕАЬРгх версии 1.20 (4 итерации) и в СЬЕЗР (вычисляется без итераций). Количество пикселов приблизительно одинаковое (~ 6 х 107,), а время вычисления пропорционально числу итераций.

На рис. 7 сравнены точность вычисления С с использованием пакетов НЕАЬР1х и С!ЬЕ8Р".

Необходимо отметить, что, в отличие от кода НЕАЬР1х, методика СЬЕвР не требует итерационного подхода при вычислении а^т— коэффициептов, и поэтому значительно быстрее. Наше определение а^то-коэффициентов в точности такое же, как и в НЕАЬР1х, которое используется для оценки спектра мощности:

1

21 + 1

|aio |2 + 2^2 |aim|2

(26)

7. Репикселизация

Чтобы преобразовать карту с распределенным сигналом из одной сетки пикселов в другую (например, из НЕАЫЧх в СЬЕвР), мы должны использовать один из двух методов:

1) рассчитать коэффициенты а^то, после чего восстановить новую карту;

2) использовать процедуры репикселизации для текущего распределения яркости.

Любая процедура репикселизации приведет к потере информации и внесет новые ошибки и невязки. Пакет СЬЕвР дает возможность для репикселизации карт с использованием двух различных методов в области АТ(в,ф): первый заключается в осреднении входящих величин в соответствующем пикселе, второй связан со сплайттовой интерполяцией внутри сетки пикселов.

В первом методе мы рассматриваем входящие пикселы, которые попадают в наш новый пиксел с величинами АТ(в*,ф*) с дальнейшим осреднением с помощью весовой функции. Реализованная весовая функция — простое осреднение с равными весами. Этот метод является широко распространенным при присвоении заданных значений пикселу с соответствующим номером.

1(1+1)С, /2п/Т

2 Провел вычисления и построил рисунок Владислав Столяров, Интститут астрономии, Кембридж

Рис. 8: Спектры мощности, рассчитанные для исходной карты НЕАЬРгх (верхняя кривая) с £тах = 1000, N,^6 = 1024, размером пиксела = 11.80260', и числом пикселов Щог=12 582 912, и для результирующей репикселизированной карты СЬЕЗР (нижняя кривая) с ближайшим возможным размером пиксела = 11.80380', N^=12 581 579. Отклонения спектра мощности на высоких I показывают от,ношение функций окон для НЕАЬРгх, и СгЬЕБР.

Во втором методе репикселизации мы используем интерполяцию сплайнами. Если у нас есть карта АТ(в*, ф*), записанная в узлах, отличных от сетки Гаусса-Лежандра, то мы можем репиксели-зировать ее в нашу новую сетку АТ(в[ ,ф[), используя приблизительно равное число пикселов и стандартную интерполяционную схему, основанную на кубических сплайнах. Этот подход достаточно быстрый, так как сплайн рассчитывается один раз

для одного вектора табулированных данных (например, в одном кольце), а величины интерполируемой функции для любого входного аргумента получаются единственным вызовом отдельной процедуры (см. подпрограммы “зрЫпё’ для расчета вторых производных интерполируемых функций и “splinf для получения интерполированной величины кубическим сплайном в книге Пресса и др. (1992)).

Наша сплайновая интерполяция состоит из трех шагов:

• построение сетки равноотстоящих узлов по оси ф для всех полярных углов;

• изменение сетки по оси х = cos(ff) в соответствии с узлами сетки GLESP;

ф для колец GLESP, проходящих через х-точки.

Рис. 8 демонстрирует отклонение точности спектра мощности в случае репикселизации из карты HEALPix в карту GLESP с одинаковым разрешением. Как видно из рисунка, для диапазона t < tmox/2 репикселизация воспроизводит корректно все свойства спектра мощности. Для t > tmax/2 должны быть проведены дополнительные исследования, чтобы учесть функцию окна пиксела. Эти исследования продолжаются.

8. Резюме

Мы предлагаем новую схему GLESP для пик-селизации неба, основанную на нулях квадратуры Гаусса-Лежандра. Она дает строгое разложение по ортогональным функциям и позволяет по-

aim

лучше 10-7 без применения итерации. Мы реализовали два подхода вычисления полиномов Лежандра, используя L- и M-методы расчетных схем.

Среди основных преимуществ нашей схемы выделим:

• aim

рационного подхода,

ного размера диаграммы, который означает оптимальное число пикселов и размер пиксела.

Соответствующие программы были разработаны на алгоритмических языках Фортран 77 и С для процедур анализа данных СМВ на небесной сфере.

aim

основной целью данного пакета. Они используются при разделении компонент и тестов на не-гауссовость (Чианг и др., 2003, Насельский и

др., 2003а,b). GLESP ориентирован на быстрое и

aim

ния. Используя точно рассчитанные коэффициен-

aim

лизации для заданных центров пикселов: GLESP, HEALPix, IGLOO или икосаэдрную.

Благодарности. Эта работа была частично финансирована Датским фондом научных исследований (Danmark Grundforskningsfond) через его поддержку открытия Центра теоретической астрофизики. Авторы благодарны Владиславу Столярову за тестирование возможности параллельных вычислений и динамическую визуализацию с помощью библиотеки OPEN GL для нулевой версии GLESP. ОВВ благодарит РФФИ за частичную поддержку данной работы грантом No 02-07-90038. Некоторые результаты данной работы были получены с помощью пакета HEALPix (Горски и др., 1999).

Список литературы

Верходанов и др. (Verkhodanov O.V., Erukhimov B.L., Monosov M.L., Chernenkov V.N., Shergin V.S.), 1993, Astrofiz. Issled. (Izv. SAO), 36, 132 Верходанов (Verkhodanov O.V.), 1997, in: “Astronomical Data Analysis Software and Systems VI”, eds.: G.Hunt & H.E.Payne, ASP Conf. Ser., 125, 46 Верходанов и др. (Verkhodanov O.V., Trushkin S. A., An-dernach H., Chernenkov V.N.), 1997, in: “Astronomical Data Analysis Software and Systems VI”, eds.: G.Hunt

& H.E.Payne, ASP Conf. Ser., 125, 322 Верходанов и Кононов (Verkhodanov O.V., Kononov V.K.), 2002, Бюлл. Спец. астрофиз. обсерв., 53, 119 Верходанов и др. (Verkhodanov O.V., Doroshkevich A.G., Naselsky P.D., Turchaninov V.I., Novikov I.D., Cristensen P.R.), 2003, in: Book of Abstracts, Danish Phys. Soc. Ann. Meeting, 2003, Hotel Nyborg Strand, June 12-13, HCOe Tryk, Kobenhavn, AA30P Верходанов О.В., Дорошкевич А.Г., Насельский П.Д., Новиков Д.И., Турчанинов В.И., Новиков И.Д., Кристенсен П.Р., 2004, в сб.: Труды Астроном. Инст. им. Штернберга, T.LXXV, Тез. Всероссийской Астрон. Конф. ВАК-2004 “Горизонты Вселенной”, МГУ, ISSN 0371-6769, 184 Горски и др. (Gorski К.М., Hivon Е., & Wan-delt B.D.), 1999, in: “Evolution of Large-Scale Structure: from Recombination to Garching”

(http://www.eso.org/science/healpix)

Градштейн и Рыжик (Gradsteyn I.S., Ryzhik I.М.), 2000, Tables of Intergrals, Series and Products, Sixth Edition, Academic Press Грейсен и Калабретта (Greisen E.W., Calabretta М.), 1993, Bull. American Astron. Soc., 182, 09.01 Дорошкевич и др. (Doroshkevich A.G., Naselsky P.D., Verkhodanov O.V., Novikov D.I., Turchaninov V.I., Novikov I.D., Christensen P.R., Chiang L.-Y.), 2003, accepted to Internat. J. Modern Phys. D, 14, No. 7 (astro-ph/0305537)

Коулс и др. (Coles P., Dineen P., Earl J., Wright D.), 2004, MNRAS, 350, 989 (astro-ph/0310252) Критенден и Турок (Crittenden R.G., Turok N.G.), 1998, Exactly Azimuthal Pixelizations of the Sky, Report-no: DAMTP-1998-78 (astro-ph/9806374)

(KT98)

Насельский и др. (Naselsky P.D., Verkhodanov O.V., Chiang L.-Y., Novikov I.D.), 2003a, ApJ(submitted) (astro-ph/0310235)

Насельский и др. (Naselsky P.D., Doroshkevich A.G., Verkhodanov O.V.), 2003b, ApJ, 599, L53 (astroph/0310542)

Насельский и др. (Naselsky P.D., Doroshkevich A.G., Verkhodanov O.V.), 2004, MNRAS, 349, 695 (astroph/0310601)

О’Нил и Лаубшер (O’Neil E.M., Laubscher R.E.), 1976, Extended studies of a quadrilateralized spherical cube Earth data base, Computer Sciences Corp., EPRF Technical Report Парийский (Parijskij Yu N.), 2001, Current Topics in As-trofundamental Physics: the Cosmic Microwave Background. Proceedings of the NATO Advanced Study Institute, held December 5-16, 1999, in Erice, Etta

Majorana Centre, Italy. Edited by Norma G. Sanchez. Published by Kluwer Academic Publishers, P. O. Box

17, 3300 AA Dordrecht, ISBN 0-7923-6855-X, 219 Пресс и др. (Press W.H., Teukolsky S.A., Vetterling W.T., Flannery B.P.), 1992, Numerical Recipes in FORTRAN, Second Edition, Cambridge University Press (http://www.nr.com)

Столяров и др. (Stolyarov V., Hobson M.P., Ashdown M.A.J., Lasenby A.N.), 2002, MNRAS, 336, 97 Тегмарк (Tegmark М.), 1996, ApJ, 470, L81 Фриго и Джонсон (Frigo М., Johnson S.G.), 1997, The Fastest Fourier Transform in the West, Technical Report MIT-LCS-TR-728 (http://www.ffiw.org)

Хивон и др. (Hivon, E., Gorski. K.M., Netterfield, C.B.), Crill, B.P., Prunet, S., & Hansen, F., 2002, ApJ, 567,

2

Чан и О’Нил (Chan F.K., O’Neil E.M.), 1976, Feasibility study of a quadrilateralized spherical cube Earth data base, Computer Sciences Corp., EPRF Technical Report

Чианг и др. (Chiang L.-Y., Naselsky P.D., Verkhodanov O.V., Way M.J), 2003, ApJ, 590, L65 (astroph/0303643)