Итак, когда в начале учебного года более-менее установилось учебное расписание, я связался со своим научным руководителем Сергеем из [Института солнечно-земной физики]( https://ru.wikipedia.org/wiki/%D0%98%D0%BD%D1%81%D1%82%D0%B8%D1%82%D1%83%D1%82_%D1%81%D0%BE%D0%BB%D0%BD%D0%B5%D1%87%D0%BD%D0%BE-%D0%B7%D0%B5%D0%BC%D0%BD%D0%BE%D0%B9_%D1%84%D0%B8%D0%B7%D0%B8%D0%BA%D0%B8 ), в дальнейшем просто ИСЗФ.

Здесь и началось всё веселье

В прошлом году я писал вместе с ним курсовую по обработке радиоизображений с помощью алгоритма [CLEAN]( https://en.wikipedia.org/wiki/CLEAN_(algorithm) ). Данный алгоритм широко используется для устранения шумов с сырых изображений, полученных с радиотелескопов и, в частности, с антенных решёток. Он, с одной стороны, совсем несложный, с другой - нужен для понимания того, как интерпретировать правильно картинки с радиотелескопов. В курсовой я попробовал с нуля запрограммировать алгоритм, поиграться с параметрами и обработать парочку тестовых картинок с [сибирского радиогелиографа]( https://ru.wikipedia.org/wiki/%D0%A1%D0%B8%D0%B1%D0%B8%D1%80%D1%81%D0%BA%D0%B8%D0%B9_%D1%81%D0%BE%D0%BB%D0%BD%D0%B5%D1%87%D0%BD%D1%8B%D0%B9_%D1%80%D0%B0%D0%B4%D0%B8%D0%BE%D1%82%D0%B5%D0%BB%D0%B5%D1%81%D0%BA%D0%BE%D0%BF ) СРГ-48.

Красивые картиночки и код с моей курсовой можно посмотреть на Гитхабе: https://github.com/vit1-irk/clean_lib

### Вернёмся к настоящему

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

Для меня нашлась работёнка в отыскании радиоисточников на Солнце, связанных с явлением гиромагнитного резонанса. Это излучение, порождённое электронами вне атомов, движущимися по замкнутым орбитам. На Солнце оно происходит в на короне и на границе короны с хромосферой, причём сигнал идёт сразу на нескольких гармониках. На частоте 34 ГГц подобные источники практически ни разу не находили, за исключением 2017 года, и целью было проверить, присутствуют ли ещё подобные образования на Солнце в другие года. Проверка не сильно сложная, однако, руки до сих пор ни у кого не дошли.

Использовать надо было данные с [Нобеямской радиообсерватории]( https://ru.wikipedia.org/wiki/%D0%9D%D0%BE%D0%B1%D0%B5%D1%8F%D0%BC%D1%81%D0%BA%D0%B0%D1%8F_%D1%80%D0%B0%D0%B4%D0%B8%D0%BE%D0%BE%D0%B1%D1%81%D0%B5%D1%80%D0%B2%D0%B0%D1%82%D0%BE%D1%80%D0%B8%D1%8F ) в Японии, где были данные аж с 1999 года по 2017. Более поздние тоже имеются, но они порченные, кривые и не подходят для обработки.

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

## Самый первый поиск источников

Начальная фильтрация - грубый поиск образований высокой яркостной температуры (выше 150 000K), при этом на корреляционной кривой сигнал не должен превышать шумовой порог. Изображения делались только раз в день, поэтому многие из источников, вероятно, были пропущены. Ищутся яркие точки, которые долго остаются на диске Солнца, по несколько часов.

https://ii-net.tk/ii/ii-point.php?q=/f/f/alicorn.blog/dGGy1m9rj59gqEfznqEe

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

Отфильтровалось несколько десятков дней, для каждого был построен видеоряд (см. в самый конец, где код).

После этого потребовалось выяснить конфигурацию магнитных полей на фотосфере Солнца в нужных активных областях, магнитограммы надо было доставать со спутников [Solar Dynamics Observatory]( https://en.wikipedia.org/wiki/Solar_Dynamics_Observatory ) и [Hinode]( https://en.wikipedia.org/wiki/Hinode_(satellite) ) .

## Магнитограммы

#### Магнитограмма SDO

https://ii-net.tk/ii/ii-point.php?q=/f/f/alicorn.blog/6lUXgHjXzef2rpEbD6eR

#### Увеличенная

https://ii-net.tk/ii/ii-point.php?q=/f/f/alicorn.blog/4QZiaLpBuCdh9qzQjR20

#### Магнитограмма Hinode

https://ii-net.tk/ii/ii-point.php?q=/f/f/alicorn.blog/RuUGVfkzn97VAPY2dtoB

## Фишки Питона

Пишу код и провожу вычисления в интерактивной среде разработки
Jupyter Notebook (а точнее - в сборке Jupyter Lab). Очень удобная, позволяет избегать ошибок и строить графики с другими результатами прямо в коде.

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

# пример записи в файл
with open("150k-curves.obj", "wb") as dump:
    pickle.dump(curves_filtered, dump)

# пример считывания из файла
with open("magnetic-plots.obj", "rb") as dump:
    magnetic_plots = pickle.load(dump)

## Построение крутых видосиков

А вот код для построения видеоряда в течение дня. Сохранять картинки в png и склеивать в видео оказалось накладно, поэтому отрисовка идёт напрямую в ffmpeg (через matplotlib FFMpegWriter).

Рендеринг первым способом занимает 4 минуты на видео, вторым - 2 минуты. Надеюсь, кому-нибудь пригодится.

    plt.rc('font', size=12)
    plt.rc('figure', titlesize=16)
    
    plt.close()
    fig = plt.figure(figsize=(22, 20))
    gs = fig.add_gridspec(nrows=4, ncols=2, width_ratios=[5, 4], hspace=0.4, wspace=0.1)

    magnplot = fig.add_subplot(gs[0, 0])
    ccplot = fig.add_subplot(gs[1, 0])
    maxplot = fig.add_subplot(gs[2, 0])
    maxplot_17 = fig.add_subplot(gs[3, 0])

    picplot = fig.add_subplot(gs[2, 1])
    picplot_17 = fig.add_subplot(gs[3, 1])

    magnplot.set_title("maximum magnetic field by SDO")
    magnplot.plot(magn["times"], magn["maxvals"], ":y")
    magnplot.plot(magn["times"], magn["maxvals"], "o")

    ccplot.set_title("34 GHz correlation curve")
    ccplot.set_ylim(-std * 0.3, std * 1.1)
    ccplot.plot(cc["times"], cc["data"] - np.mean(cc["data"]))
    ccplot.axhline(std)

    maxplot.set_title("34 GHz, maximum value")
    maxplot.set_ylim(-100, cfg.intensity_threshold * 1.5)
    maxplot.plot(pics_times, maxvals, "-o")
    maxplot.axhline(cfg.intensity_threshold, color="green")

    maxplot_17.set_title("17 GHz, maximum value")
    maxplot_17.plot(pics_times_17, maxvals_17, "-o")

    print("making video")
    path = "/mnt/data/filepath/videos/{0}.mp4".format(pic_start)
    writer = FFMpegWriter(fps=3, extra_args=['-vcodec', 'libx264'])

    with writer.saving(fig, path, dpi=110):
        for a in range(0, len(pics)):
            pic = pics[a]
            pic_17 = pics_17[a]
            maxval = maxvals[a]
            maxval_17 = maxvals_17[a]

            tl1 = ccplot.axvline(pic["times"], color="red")
            tl2 = maxplot.axvline(pic["times"], color="red")
            tl3 = maxplot_17.axvline(pic_17["times"], color="red")
            tl4 = magnplot.axvline(pic["times"], color="red")

            picplot.set_title("34 GHz, maxval = " + str(maxval))
            picplot_17.set_title("17 GHz, maxval = " + str(maxval_17))

            picplot.imshow(pic["data"], interpolation=None, origin='low')
            picplot_17.imshow(pic_17["data"], cmap="winter", interpolation=None, origin='low')

            writer.grab_frame()

            tl1.remove()
            tl2.remove()
            tl3.remove()
            tl4.remove()

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

https://ii-net.tk/ii/ii-point.php?q=/f/f/alicorn.blog/jPqH5A3gnFAsHzLaOuMN

## Что дальше

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

Ссылка в блоге: https://blog.alicorn.tk/posts/the-beginning-sun.html