Ошибка в алгоритм получения интерполяции с помощью Лагранжа?

Question

[BlinCT]

Всем привет
С такими расчетами я столкнулся впервые, сам его не писал полностью, но вот странность что при определенных значениях в расчетах я получаю Nan.
Мне нужно получить расчитаные точки для интерполяции, и нашел в инете примерно как он должен выглядеть.
Ниже привожу с какими точками я тестирую и сам аргаритм.

18;24
21;21
23;18
24;16
23;11
20;8

И вот сам алгоритм

QList<QPointF> tmpResultPoints;
            QList<QList<QPointF> > tmpResultPoints2;

            auto currentX = list.m_List.begin()->m_X;
            auto endX = list.m_List.last().m_X;

            while(currentX <= endX)
            {
                qreal workX = currentX;
                qreal calcY = 0.0;

                for(int j = 0; j < list.m_List.count(); ++j)
                {
                    qreal Hi = 0.0;
                    bool bFirst = true;

                    for(int i = 0; i < list.m_List.count(); ++i)
                    {
                        if(i != j)
                        {
                            qreal calc = (workX - list.m_List.at(i).m_X) / (list.m_List.at(j).m_X - list.m_List.at(i).m_X);

                            if(bFirst)
                            {
                                Hi = calc;
                                bFirst = false;
                            }
                            else
                            {
                                Hi *= calc;
                                bFirst = false;
                            }
                        }
                    }

                    qreal Px = list.m_List[j].m_Y * Hi;
                    calcY += Px;
                }

                tmpResultPoints.append(QPointF(workX, calcY));
                currentX += m_SizeStep;
            }

Помогите пожалуйста кто понимает в таких расчетах. Не понимаю где ошибка.
Если надо я могу даже тестовую рабочую прогу скинуть.

Comments

[Wataru]

Кстати, по коду, его можно сильно упростить, не надо вам bFirst вообще. Вы просто инициализируйте Hi = 1.0 и в цикле всегда домножайте на calc. Тогда и calc можно выкинуть вообще и просто сразу домножать. Или вообще сразу инициализируйте как list.m_List[j].m_Y.

Answer 1

[Wataru]

У вас там во входных данных 2 точки с одинаковым x (23) и разным y (18, 11). В результате получается деление на 0.

Нельзя инетрполировать такие данные полиномом. Ибо это функция от x - для каждого x одно значение y.

Ошибка не в программе, а в некорректных входных данных.

Можно интерполировать параметрически, если хотите. Заведите параметр t и ищите две функции x(t), y(t) - скармливайте этой программе 2 набора данных 1,x1;2,x2;3,x3...n;xn и 1,y1;2,y2;...n,yn.

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

Comments

[BlinCT]

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

[Wataru]

Что именно? вы хотите вариант с параметрической функцией для x и y?

Просто весь ваш код оборачиваете в функцию. Передаете ей list 2 раза, только первый раз m_X равен номеру точки, m_Y равен m_X точки. Второй раз m_X тот же номер, а m_Y равен m_Y точки.

[BlinCT]

Wataru, если я правильно вас понимаю то лист с набором точек X и Y надо разбить на 2 листа с индексом и точками соответственно?
Получается на выходе из каждой функции я получу 2 листа с расчетами, потмо просто их обьеденить?
Хотя погодите, как пример первый лист должен выглядеть так?
0, 18
1, 21
2, 23
3, 24
4, 23
5, 20

[Wataru]

BlinCT, Да, первый лист выглядит так. Второй:
0;24
1;21
2;18
3;16
4;11
5;8

Результаты нельзя просто объединить. Вам над взять y от первого листа и потом y от второго. Это и будет список x,y координат точек. При этом они могут идти вертикально, справа-налево или вообще закручиватья в спираль. Что совсем не полином, но по крайней мере, эта гладкая кривая пройдет через все ваши изначальные точки.

[BlinCT]

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

qreal Model::translateY(qreal y, qreal coeff)
    {
        if(y > 0)
        {
            return -y * coeff;
        }
        else if(y < 0)
        {
            return (y * coeff) * (-1);
        }

        return 0.0;
    }

    qreal Model::translateX(qreal x, qreal coeff)
    {
        return x * coeff;
    }

Правильный ли это расчет? У меня просто по Y значения ооочень большие получаются.
Второй аргумент это высота или ширина области деленая на максимальное значение на шкале.
В моем тестовом случаи на вертикальнйо шкале это 25 а на горизонтальной это 300.

[Wataru]

Стоп. Вам надо, чтобы линия проходила через вот эти вот все точки, так? Тогда никакого пересчета делать нельзя.

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

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

[BlinCT]

Wataru, а вы не могли бы написать как кубическим бы выглядил алгаритм расчетов?
Я в теории понимаю что если точек больше 3 то считается 2 отрезка потом точка смещается на 1 и снова 2 отрезка. Но вот написать это очень не просто. Я этот вариант нарыл на форуме и то до конца его не понимаю.
заранее спасибо.

[Wataru]

BlinCT, Вот, например, статья. Там с кодом. https://habr.com/ru/articles/749288/

Вам точно также придется 2 раза прогнать для данных {i, xi} и {i, yi}, потому что этот метод все еще не работает с вертикально расположенными точками.

И вообще, гуглите: интерполяция кубический сплайн на плоскости.

[BlinCT]

Wataru, спасибо за этот пример. Сейчас изучаю и попробую переписать на кубтческий сплайн. Но все таки вопрос еще про зависимость точек от размера плоскости. То что я выше написал это правильно, так надо расчитывать?

[Wataru]

BlinCT, Я не понимаю, что у вас за размер. У вас есть ключевые точки, кривая через них пройдет и никак иначе.

[BlinCT]

Wataru, предположим вертикальная шкала от 0 до 25, горизонтальная от 0 до 300.
Есть набор точек, расчеты в конечном итоге будут рисоватся по пикселям. А теперь взяли и растянули график по горизонтале, и уже имеем не соответствие плоскости где мы рисуем и точек которые надо нарисовать.
То есть изменилась ширина.
Я об этмо расчете.

[Wataru]

BlinCT, Ну тут у вас все точки растягиваются вместе с графиком. По идее просто все числа домножаются на коээфициент растяжения по оси OX и по оси OY. Если вы график увеличили c 25 x 300 до 100 x 600, то вы все x координаты умножаете на 4, а y - на 2.

[BlinCT]

Wataru, нет, у меня ничего автоматически не растягивается. Мне самому приходится расчитать точки, и потом на полотне сказать где эта точка находится. То есть если мне надо по шкале точку со значением 26 нарисовать, а при этом по этой ширине полотно у меня в 725 пикселей, то мне надо расчитать где эта точка должна находится. Вопрос был именно об этом. Как правильно этот коэфициент расчитать. И нужно ли это делать до алгаритма интерполяции или это лучше делать после?

[Wataru]

BlinCT, Ну вот вы рассчитали точки, когда полотно было 25x300. Когда его размер 100x600 вам надо все координаты домножить на 100/25 по x и 600/300 по y.