#include <math.h>
#include <Windows.h>
#include <iostream>

using namespace std;

const int N = 3 , MAX_SIZE=200;
double q[N][N], r[N][N], m[N][N],h[N],tmp[N][N],tmpr[N][N],tmpq[N][N], sum,x[MAX_SIZE],y[MAX_SIZE],f[N],b[N],sol[N];
int n;

double Function(double x) {
	return x;
}

void dupeR() {
    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            tmpr[i][j] = r[i][j];
    }
void dupeQ(){
    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            tmpq[i][j] = q[i][j];
}

double norm(double vec[],int size,int from) {
    double s = 0;
    for (int i = from; i < size; i++) {
        s += vec[i] * vec[i];
    }
    return sqrt(s);
}

double basis(double x, int deg) {
    return pow(x, deg);
}

 double resF(double x){
    double sum = 0;
    for (int i = 0 ; i < N ; i++) sum+= basis(x,i)*sol[i];
    return sum;
}

void graph() {
    float p;
    HDC hDC = GetDC(GetConsoleWindow());
    HPEN Pen = CreatePen(PS_SOLID, 2, RGB(255, 255, 255));
    SelectObject(hDC, Pen);
    MoveToEx(hDC, 0, 200, NULL);
    LineTo(hDC, 500, 200);
    MoveToEx(hDC, 300, 0, NULL);
    LineTo(hDC, 300, 470);
    for (p = -8.0f; p <= 6.0f; p += 0.001f) // O(300,200) - center
    {
        MoveToEx(hDC, 70 * p + 300, -70 * resF(p) + 200, NULL);//10 - scale
        LineTo(hDC, 70 * p + 300, -70 * resF(p) + 200);
    }
}

void graphTrue() {
    float p;
    HDC hDC = GetDC(GetConsoleWindow());
    HPEN Pen = CreatePen(PS_SOLID, 2, RGB(0, 255, 0));
    SelectObject(hDC, Pen);
    for (p = -8.0f; p <= 6.0f; p += 0.001f) // O(300,200) - center
    {
        MoveToEx(hDC, 70 * p + 300, -70 * Function(p) + 200, NULL);//10 - scale
        LineTo(hDC, 70 * p + 300, -70 * Function(p) + 200);
    }
}

void QR() {

    for (int t = 0; t < N-1; t++) {
        for (int i = 0; i < N; i++)
            for (int j = 0; j < N; j++)
                if (i == j) {
                    tmp[i][j] = 1;
                }
                else {
                    tmp[i][j] = 0;
                }
        for (int i = 0; i < N; i++) if (i != 0) sum += m[i][1] * m[i][1];
        if (sum != 0) {
                for (int i = 0; i < N; i++) if (i < t) h[i] = 1; else h[i] = r[i][t];
            h[t] = h[t] + norm(h, N,t);
            double normh = norm(h, N,t);
            for (int i = t; i < N; i++) h[i] = h[i] / normh;
            for (int i = t; i < N; i++)
                for (int j = t; j < N; j++)
                {
                    tmp[i][j] += -2 * h[i] * h[j];
                }

            cout << endl;
            cout << "R:";
            for (int i = 0; i < N; i++) {
                cout << endl;
                for (int j = 0; j < N; j++) {
                    cout << tmp[i][j];
                    cout << "  ";
                }
            }

            dupeR();
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++) {
                    r[i][j] = 0;
                    for (int k = 0; k < N; k++) {
                        r[i][j] += tmp[i][k] * tmpr[k][j];
                    }
                }

            cout << "Q:";
            cout << endl;

            for (int i = 0; i < N; i++) {
                cout << endl;
                for (int j = 0; j < N; j++) {
                    cout << r[i][j];
                    cout << "  ";
                }
            }


            dupeQ();
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++) {
                    q[i][j] = 0;
                    for (int k = 0; k < N; k++) {
                        q[i][j] += tmp[i][k] * tmpq[k][j];
                    }
                }



        }

    }




    /*
    Транспонирование Q и вычисление A = QR
    dupeQ();
    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            q[i][j] = tmpq[j][i];


    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++) {
            tmp[i][j] = 0;
            for (int k = 0; k < N; k++) {
                tmp[i][j] += q[i][k] * r[k][j];
            }
        }

    cout << "QR:";
    cout << endl;

    for (int i = 0; i < N; i++) {
        cout << endl;
        for (int j = 0; j < N; j++) {
            cout << tmp[i][j];
            cout << "  ";
        }
    }
    */

}


void birds() {
	float p;
	HDC hDC = GetDC(GetConsoleWindow());
	HPEN Pen = CreatePen(PS_SOLID, 2, RGB(0, 200, 200));
	SelectObject(hDC, Pen);

		for (int i=0 ; i <= n ; i++)
	{
		MoveToEx(hDC, 70 * x[i] + 300, -70 * y[i] + 200, NULL);//10 - scale
		LineTo(hDC, 70 * x[i] +304, -70 * y[i] + 200);
		MoveToEx(hDC, 70 * x[i] + 302, -72 * y[i] + 200, NULL);//10 - scale
		LineTo(hDC, 70 * x[i] + 302, -68 * y[i] + 200);
	}
 }
int main() {

    //cin >> n; // Количество точек < MAX_SIZE
    //cin >> a;
    //cin >> b;
    /*for (int i = 0; i < n; i++) {
        cin >> x[i];
        cin >> y[i];
    }*/
    n = 20;
    for (int i = 0; i < n; i++) {
        x[i] = 0.25*i;
        y[i] = Function(x[i]) + pow(-1.0,i)*0.01*(rand()%50) ;
    }

    for (int i = 0; i < N; i++) {
        sum = 0;
        for (int j = 0; j < n; j++) sum += y[j] * basis(x[j], i);
        f[i] = sum;
    }

    for (int i = 0; i < N;i++)
        for (int j = 0; j < N; j++)
        {
            sum = 0;
                for (int k = 0; k < n; k++) sum += basis(x[k], i)*basis(x[k], j);
                m[i][j] = sum;
        }

    for (int i= 0; i<N ; i++)
        for (int j = 0; j < N; j++)
        if (i == j) {
            q[i][j] = 1;
            r[i][j] = m[i][j];
        }
        else {
            q[i][j] = 0;
            r[i][j] = m[i][j];
        }
    QR();

    for (int i = 0; i < N; i++)
    {
        sum = 0;
        for (int k = 0; k < N; k++) sum += q[i][k]*f[k];
        b[i] = sum;
    }
    cout << endl;
    cout << "coef:";

    for (int i = N-1; i >=0; i--)
    {
        sum = b[i];
        for (int k = i+1; k < N; k++) sum += -sol[k]*r[i][k] ;
        sol[i] = sum/r[i][i];

    }

    for (int i = 0; i < N; i++) {
        cout << sol[i];
        cout << ' ';
    }

    graph();
   graphTrue();
	birds();
    cin >> m[0][0];
}