#include <string>
#include <fstream>
#include <iostream>
#include <vector>
#include <math.h>
#include <random>
#include <algorithm>

using namespace std;

//size of each batch of secants
size_t SAMPLE_SIZE=100;
//N
size_t SAMPLE_NUM=100000;

bool test(vector<double> beta, double crit_value){
    //returns true if Nullhypothesis is rejected
    double statistics=0;
    for(auto k=0;k<SAMPLE_SIZE;++k){
        statistics+=log(abs(sin(beta[k])));
    }
    if(statistics>=crit_value){
        return true;
    }else{
        return false;
    }
}

bool saveArray( const double* pdata1,const double* pdata2, size_t length, const std::string& file_path )
{
    //function for writing data into .txt document
    std::ofstream os(file_path.c_str(), std::ios::binary | std::ios::out);
    if ( !os.is_open() )
        return false;
    for(auto k=0;k<length;++k){
        os<< *(pdata1+k)<<"     "<<*(pdata2+k) <<endl;
    }
    os.close();
    return true;
}

int main()
{
    //random numbers generator
    uniform_real_distribution<double> unif(0,1);
    default_random_engine re;

    //sample for finding empirical value of T
    vector < vector< double > > sample(SAMPLE_NUM,vector<double>  (SAMPLE_SIZE,(2)));
    vector <double> T(SAMPLE_NUM);
    for(auto k=0;k<SAMPLE_NUM;++k){
        for(int l=0;l<SAMPLE_SIZE;++l){
            //create random angle and compute empirical value of T
            sample[k][l]=unif(re)*M_PI;
            T[k]+=log(abs((sin(sample[k][l]/2))));
        }

    }

    //sorting list of values of T
    sort (T.begin(), T.end());

// find empirically ciritical values for test statistics
    //belonging to a significance of 10% resp. 1%
    double T_10crit=T[9*SAMPLE_NUM/10];
    double T_1crit=T[99*SAMPLE_NUM/100];


    //independently drawn sample respecting the nullhypothesis
    vector < vector< double > > Hyp_null(SAMPLE_NUM,vector< double > (SAMPLE_SIZE));
    for(auto k=0;k<SAMPLE_NUM;++k){
        //cout<<k<<": "<<endl;
        for(auto l=0;l<SAMPLE_SIZE;++l){
            Hyp_null[k][l] =abs(unif(re)*M_PI)/2;
        }
    }

    //independently drawn sample respecting the alternative hypothesis
    vector < vector< double > > Hyp_alt(SAMPLE_NUM,vector< double > (SAMPLE_SIZE));
    for(auto k=0;k<SAMPLE_NUM;++k){
        for(auto l=0;l<SAMPLE_SIZE;++l){
            Hyp_alt[k][l] =acos(unif(re));
        }

    }

//variables storing rejection rate
    int results_null1=0;
    int results_alt1=0;
    int results_null10=0;
    int results_alt10=0;
    for(auto k=0;k<SAMPLE_NUM;++k){
        //finding rejection rate of the two tests (significance level 1%, 10% respectively)
        results_null1+=test(Hyp_null[k],T_1crit);
        results_alt1+=test(Hyp_alt[k],T_1crit);
        results_null10+=test(Hyp_null[k],T_10crit);
        results_alt10+=test(Hyp_alt[k],T_10crit);
    }

    //print out rejection rates
                                                         //empirical significance(should be close to 1%*N
    cout<<"test results:"<<endl<<"1% significance: null hypothesis holds: "<<results_null1;
    //N*Power (Macht) of the test, should be as high as possible (less than or equal to N though)
    cout<<"   , alternative hypothesis holds: "<<results_alt1<<endl;
                                                             //empirical significance(should be close to 10%*N
    cout<<"test results:"<<endl<<"10% significance: null hypothesis holds: "<<results_null10;
    //N*Power (Macht) of the test, should be as high as possible (less than or equal to N though)
    cout<<"   , alternative hypothesis holds: "<<results_alt10<<endl;

//rearrange stuff for printing and visualizing
    vector<double> data_null(SAMPLE_SIZE*SAMPLE_NUM);
    vector<double> data_alt(SAMPLE_SIZE*SAMPLE_NUM);
    for(auto u=0;u<SAMPLE_NUM;++u){
        for(auto l=0;l<SAMPLE_SIZE;++l){
            data_null[u*SAMPLE_SIZE+l]=Hyp_null[u][l];
            data_alt[u*SAMPLE_SIZE+l]=Hyp_alt[u][l];
        }

    }
    saveArray(data_null.data(),data_alt.data(),SAMPLE_SIZE*SAMPLE_NUM,"both_hyp.txt");

    vector<double> statistics_null(SAMPLE_NUM);
    vector<double> statistics_alt(SAMPLE_NUM);
    for(auto l=0;l<SAMPLE_NUM;++l){
        for(auto k=0;k<SAMPLE_SIZE;++k){
            statistics_null[l]+=-log(sin(Hyp_null[l][k]));
            statistics_alt[l]+=-log(sin(Hyp_alt[l][k]));
        }
    }
    saveArray(statistics_null.data(),statistics_alt.data(),SAMPLE_NUM,"statistics.txt");

    return 0;
}