#ifndef _WIN32
#include "../include/lminmaxgd.h"
#else
#include "..\include\lminmaxgd.h"

#endif
#include <time.h>
using namespace std;
#include <iostream>
#include <time.h>
#include "..\src\gd.cpp"
#include "..\src\poly.cpp"
#include "..\src\serie.cpp"
#include "..\src\smatrix.cpp"
#include "..\src\tools.cpp"

//smatrix calc_prefilter(int n, TEG_system S[]); //function for pre filter calculation (optimal communication topology)
//smatrix calc_feedback(int n, smatrix &Com, smatrix &A);
//smatrix calc_equal_trans(smatrix &H1, smatrix &H2);


int main()
{
    srand (time(NULL));

    //Random system 1
    //------------------------------------------
    poly temp_poly_1(rand()%20,rand()%20);
    gd r_temp_1(rand()%20,rand()%20);

    for (int i = 0;i<rand()%5;i++)
        temp_poly_1.add(gd(rand()%20,rand()%20));

    //Random system 2
    //------------------------------------------
    poly temp_poly_2(rand()%20,rand()%20);
    gd r_temp_2(rand()%20,rand()%20);

    for (int i = 0;i<rand()%5;i++)
        temp_poly_2.add(gd(rand()%20,rand()%20));

    serie system_1;
    serie system_2;
    serie grow_rate;

    system_1.init(epsilon, temp_poly_1, r_temp_1);
    system_2.init(epsilon, temp_poly_2, r_temp_2);


    cout<<"Single System:"<<endl;
    cout<<"--------------------------------"<<endl;
    cout<<"System 1:"<<system_1<<endl;
    cout<<"System 2:"<<system_2<<endl<<endl;


    serie sat_system_1, sat_system_2;

    sat_system_1 = frac(system_1,system_1);
    sat_system_1 = frac(sat_system_1,system_1);
    sat_system_1 = prcaus(sat_system_1);
    system_1 = otimes(sat_system_1,system_1);

    sat_system_2 = frac(system_2,system_2);
    sat_system_2 = frac(sat_system_2,system_2);
    sat_system_2 = prcaus(sat_system_2);
    system_2 = otimes(sat_system_2,system_2);

    cout<<"Single System saturated (h = h+):"<<endl;
    cout<<"--------------------------------"<<endl;
    cout<<"System 1:"<<system_1<<endl;
    cout<<"System 2:"<<system_2<<endl<<endl;


    if ( ((float)r_temp_1.getd())/r_temp_1.getg() > ((float)r_temp_2.getd())/r_temp_2.getg()){
        grow_rate.init(epsilon,e,r_temp_1);
    }else{
        grow_rate.init(epsilon,e,r_temp_2);
    }

    float slope, max_slope;
    int index_1 =0, index_2 = 0;
    gd *p1, *p2;
    p1 = system_1.getq().getdata();
    p2 = system_2.getq().getdata();
    for (int i =0;i<system_1.getq().getn();i++)
    {
        slope =((float)p1->getd())/p1->getg();
        if (slope> max_slope) {
            max_slope = slope;
            index_1 = i;
        }
        p1++;
    }

    for (int i =0;i<system_2.getq().getn();i++)
    {
        slope =((float)p2->getd())/p2->getg();
        if (slope> max_slope) {
            max_slope = slope;
            index_2 = i;
        }
        p2++;
    }

    cout<<"Elements for Filter calculation:"<<endl;
    cout<<"--------------------------------"<<endl;
    cout<<"Index: "<<index_1<<"  Element: "<<system_1.getq().getpol(index_1)<<endl;
    cout<<"Index: "<<index_2<<"  Element: "<<system_2.getq().getpol(index_2)<<endl<<endl;

    serie filter_1, filter_2;

    filter_1 = frac(grow_rate,system_1.getq().getpol(index_1));
    filter_2 = frac(grow_rate,system_2.getq().getpol(index_2));


    filter_1  = prcaus(filter_1);
    filter_2  = prcaus(filter_2);

    cout<<"Filters F:"<<endl;
    cout<<"----------"<<endl;
    cout<<"F1:"<< filter_1 << endl;
    cout<<"F2:"<< filter_2 << endl<<endl;

    getchar();

    smatrix H(2,2), A_0(2,2), H_loop;

    serie h_test, h_loop;

    poly h_poly(13,13);
    h_poly.add(gd(20,15));
    gd h_r(35,23);

    h_test.init(epsilon,h_poly,h_r);

    h_loop = h_test;
    h_loop = frac(h_test,h_test);
    h_loop = frac(h_loop,h_test);
    h_loop = prcaus(h_loop);
    h_loop = otimes(h_loop,h_test);

   // H(0,0) = otimes(system_1,filter_1);
   // H(1,1) = otimes(system_2,filter_2);

   // H(0,0) = otimes(system_1,filter_1.getq());
   // H(1,1) = otimes(system_2,filter_2.getq());

   H(0,0) = h_loop;
   H(1,1) = h_loop;


    cout<<"Unconected System:"<<endl;
    cout<<"------------------"<<endl;
    cout<<H<<endl<<endl;



    A_0(0,1) = e;
    A_0(1,0) = e;
   A_0(0,0) = e;
    A_0(1,1) = e;


    H_loop = otimes(H,A_0);
    H_loop = star(H_loop);

    H_loop = otimes(H_loop,H);
    cout<<"Connected System H_loop:"<<endl;
    cout<<"------------------------"<<endl;
    cout<<H_loop<<endl;

    cout<<"Connected System H_loop Periodic Behavier:"<<endl;
    cout<<"------------------------------------------"<<endl;
    cout<<"H_loop[0,0]=("<<H_loop(0,0).getq()<<")("<<H_loop(0,0).getr()<<")*"<<endl;
    cout<<"H_loop[0,1]=("<<H_loop(0,1).getq()<<")("<<H_loop(0,1).getr()<<")*"<<endl;
    cout<<"H_loop[1,0]=("<<H_loop(1,0).getq()<<")("<<H_loop(1,0).getr()<<")*"<<endl;
    cout<<"H_loop[1,1]=("<<H_loop(1,1).getq()<<")("<<H_loop(1,1).getr()<<")*"<<endl;


    return 0;
}




