- PVSM.RU - https://www.pvsm.ru -
В предыдущей статье [1] я уделил внимание исключительно алгоритму и «физике» метода. В этой статье будет описано немного про получение результатов и их представлении.
Целью были результаты, потому ничего странного не было в том, что бы выбрать готовый продукт который уже умеет АРСС. Выбор был сделан в пользу oldCronos [4]. Причинами, сподвигнувшими меня на такой выбор были:
Для всех пар из множества была посчитана функция максимального правдоподобия.
Все расчёты проводились на маленьком «кластере» из 16 ПК на базе AMD Athlon 64 X2 4400+ c 2 Гб оперативной памяти. Впрочем, это и всё что понадобилось.
Для 120 дней (а это в общей сложности 480 точек для чистых данных и 120 для прореженных) и для 9 существенных точках на гринвичском меридиане имеем 800 моделей. Расчёты заняли существенное время — 1 неделя и 3 дня.
Как же всё это было сделано, для начала подключили все нужные header-файлы
#include <boost/mpi.hpp>
namespace mpi = boost::mpi;
#include <boost/serialization/access.hpp>
#include <setjmp.h> // jmp_buf, longjmp, setjmp
#include <gsl/gsl_errno.h> // GSL: gsl_set_error_handler
// oldCronos headers
#include "ARMAModel.hpp"
#include "TimeSeries.hpp"
boost.serialization нужен для передачи через boost.mpi собственноручно написанного простенького класса, который будет хранить порядки модели и значение.
class task {
private:
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive &ar, const unsigned int version) {
ar &p;
ar &q;
ar &llh;
}
public:
int p;
int q;
double llh;
task() {}
explicit task(int p, int q, double llh)
: p(p), q(q), llh(llh)
{}
};
Следующим важным моментом стало задание своего обработчика ошибок для GSL. Это необходимо чтобы GSL не делал abort() при расходимости методов оптимизации — в частности, используемого ЛМ метода.
jmp_buf pos;
int p;
int q;
bool error;
void my_handler(const char *reason, const char *file, int line, int gsl_errno)
{
std::cout << "t" << p << " " << q << " # " << file << ":" << line << ' ' << reason << std::endl;
error = true;
longjmp(pos, 1);
}
auto old_handler = gsl_set_error_handler(&my_handler); // setting own gsl_error handler
…
setjmp(pos);
…
gsl_set_error_handler(old_handler); // returning default gsl error handler
Создать контейнер с набором не составляет проблем, потому сразу отправляем задания на удалённые машины и ждём ответ. И снова отправляем.
Полный код можно увидеть здесь [5].
Концом заданий является специально сформированный объект класса task, где p или q равны -1. В результате из файлы примерно такого содержимого [6], получаем такое [7].
С помощью элементарного кода на LaTeX и используя пакет pgfplots можно получить красивые картинки (1 [8] и 2 [9]):
documentclass{standalone}
usepackage{pgfplots}
pgfplotsset{compat=newest}
pgfplotsset{every axis/.append style={
font=tiny,
tiny,
tick style={ultra thin}}}
begin{document}
begin{tikzpicture}
begin{axis}[height=65mm
,width=65mm,
scatter/@pre marker code/.code={
pgfmathparse{100-33.3*ln(pgfplotspointmetatransformed+1)}
letopacity=pgfmathresult
defmarkopts{mark=*,color=red!opacity,fill=red!opacity}
expandafterscopeexpandafter[markopts]}
,scatter/@post marker code/.code={endscope}
,scatter src=explicit
,xlabel={$p$}
,ylabel={$q$}]
addplot[only marks,scatter] table[x index=0, y index=1,meta index=2] {<data>};
end{axis}
end{tikzpicture}
end{document}
С помощью ImageMagick, или, к примеру, этого [10] сервиса были получены картинки в формате png.
Рис. 1: Максимальное правдоподобие для 45°S, 0°.
Значение для пары сначала линейно отображается на шкалу от 0 до 1000, а затем элементарно преобразовывается формулой .
Рис. 2: Максимальное правдоподобие для 0°, 0°. Прореженные данные.
А с помощью простенького приложения [11] на Qt можно посмотреть IONEX формат и сохранить изображение на память.
А с помощью LaTeX это всё можно приукрасить:
documentclass{standalone}
usepackage{pgfplots}
begin{document}
begin{tikzpicture}
begin{axis}[height=45mm
,width=45mm
,scale only axis
,enlargelimits=false
,axis on top
,point meta min=0
,point meta max=6.5e100
,xtick={-120,0,120}
,ytick={-60,0,60}
,xlabel={$lambda,,{}^circ$}
,ylabel={$varphi,,{}^circ$}
,xlabel shift={-1mm}
,ylabel shift={-5mm}]
addplot graphics[xmin=-180,xmax=180,ymin=-87.5,ymax=87.5] <image>};
end{axis}
end{tikzpicture}
end{document}
И в результате получим pdf с картинкой, осями и подписями.
Автор: m0nhawk
Сайт-источник PVSM.RU: https://www.pvsm.ru
Путь до страницы источника: https://www.pvsm.ru/algoritmy/6670
Ссылки в тексте:
[1] статье: http://habrahabr.ru/post/142434/
[2] Программирование как оно было: #programming
[3] Представление результатов: #representation
[4] oldCronos: http://www.stat.cmu.edu/~abrock/oldcronos/
[5] здесь: https://gist.github.com/2569840
[6] содержимого: https://gist.github.com/2569781
[7] такое: https://gist.github.com/2569785
[8] 1: #first_picture
[9] 2: #second_picture
[10] этого: http://www.convertfiles.com/
[11] приложения: http://m0nhawk.github.com/ionex-library/
Нажмите здесь для печати.