package WT;

import java.io.IOException;
import java.util.Random;

import WT.tools.Utilities;

public class Wiener2 {

/**
 *  Bestimmt den Wert der Schauder-Funktionen S_n(t)
 *  @param t  betrachtete Stelle
 *  @param n  Ordnung der Schauder-Funktion
 *
 *  @return Funtkionswert der Schauderfunktion n-ter Ordnung an Stelle t
 */
      private static double schauderFunction(double t, int n){
	if (n==1) { return t;}

	/* Bestimme die Darstellung n=2^m+k mit k=1,...,2^m */
	int m = 0, k;
	double c = Math.pow(2,m);
	while (2*c < n){
	    m++;
	    c = Math.pow(2,m);
	}
	k = (int)(n - c);

        /*  Funktionswert wie in Vorlesung bzw. Musterlösung beschrieben */
	if( t < (k-1)/c || t >= k/c ){ 
	    return 0; }
	
	if( t < ((2*k)-1)/(2.0*c) ){
	    return ( Math.sqrt(c)*(t - (k-1)/c) ); }
	else{
	    return (-Math.sqrt(c)*(t - k/c) );     }
    }

/**
 *  Verwendung des ersten Algorithmus zur Approximation eines Wiener Prozesses X
 *  @param m  2^m = Ordnung der Approximation
 *  @param seed Startwert für Zufallszahlengenerator
 *  @return Realisierungen von X an den Stützstellen
 */
    private static double[] algorithm1(int m, long seed){
	double n_d = Math.pow(2,m);
	int n = (int)n_d, k=0;
	double[] result = new double[n+1];
	/* Start in (0,0) */
	result[0] = 0;

	/* erzeuge n=2^m N(0,1)-Zufallszahlen */
	Random rand = new Random(seed);
	double[] y = new double[n+1];
	for(k = 0; k <= n; k++){
	    y[k] = rand.nextGaussian();
	}

	/* Berechnung der Funktionswerte an den Stützstellen */
	double sum = 0;
	for(int i = 1; i <= n; i++){
	    sum = 0;
	    for(k = 0; k <= n; k++){
		sum = sum + y[k]*schauderFunction(i/n_d, k+1);
	    }
	    result[i] = sum;
	}
	return result;
    }
 
/**
 *  Verwendung des zweiten Algorithmus zur Approximation eines Wiener Prozesses X
 *  @param m  2^m = Ordnung der Approximation
 *  @param seed  Startwert für Zufallszahlengenerator
 *  @return Realisierungen von X an den Stützstellen
 */
    private static double[] algorithm2(int m, long seed){
	double n_d = Math.pow(2,m);
	int n = (int)n_d;
	double[] result = new double[n+1];
	Random rand = new Random(seed);

	/* Start in (0,0) */
	result[0] = 0;

	/* Berechnung der Funktionswerte an den Stützstellen 
	   unter Verwendung von U({-1,1})-Zufallszahlen */
	double sum = 0;
	double z = 0;
	for(int i = 1; i <= n; i++){
	    z = rand.nextDouble();
	    if( z < 0.5){ sum = sum - 1.0/Math.sqrt(n_d);  }
	    else{ 	  sum = sum + 1.0/Math.sqrt(n_d);  }
	    result[i] = sum;
	}
	return result;
    }


/**
 *   Hauptprogramm:
 *   Simulation eines Wiener Prozesses mit n = 2^6, 2^7 bzw. 2^8 Stützstellen.
 *   Approximation der Ordnung n = 2^6, 2^7 bzw. 2^8.
 */
    public static void main(String[] args) {
	int sampleSize = 100, m=8;
	double[] sample = new double[sampleSize];
	double[] values = new double[(int)Math.pow(2,m)+1];
	Random rand_seed = new Random(System.currentTimeMillis());

	/* Algorithmus 1 */
         System.out.println("Algorithmus 1");
	for(int i=0; i < sampleSize; i++){
	    /* Unabhängige Simulationen werden z.B. über Initialisierung mit 
	       verschiedenen random seeds erzielt */
	    values = algorithm1(m,rand_seed.nextLong());
	    sample[i] = Utilities.getMax(values);
	}
	/* Bestimme Mittelwert und Stichprobenvarianz */
	double mean = Utilities.getMean(sample), var = Utilities.getVar(sample);
	System.out.println("Erwartungswert: " +Utilities.out(mean)+ "\n");
	System.out.println("Standardabweichung: "+Utilities.out(var)+"\n");

	/* Teste auf signifikante Abweichungen vom Erwartungswert */
	double test_statistic_mu = (mean - Math.sqrt(2.0/Math.PI))/Math.sqrt(var/sampleSize);
	if(Math.abs(test_statistic_mu) > 1.96){
	    System.out.println("Die Nullhypothese mu = 0.797 wird abgelehnt, "+
                                "Teststatistik: T="+Utilities.out(test_statistic_mu)+".\n");
	}
	else{
	    System.out.println("Die Nullhypothese mu = 0.797 wird nicht abgelehnt, "+
                                "Teststatistik: T="+Utilities.out(test_statistic_mu)+".\n");
	}

	/* Teste auf signifikante Abweichungen von der theoretischen Varianz */
	double mu4 = Utilities.getFourthMoment(sample);
	double sigmaq = Math.sqrt(1.0-2.0/Math.PI);
	double test_statistic_sigmaq = 0.0;
	if( mu4 > sigmaq*sigmaq){
	    test_statistic_sigmaq = (var - sigmaq)/Math.sqrt( (mu4 - sigmaq*sigmaq)/sampleSize );
	}
	else{
	    test_statistic_sigmaq = 5.0;
	}
	if(Math.abs(test_statistic_sigmaq) > 1.96){
	    System.out.println("Die Nullhypothese sigma² = "+Utilities.out(sigmaq)+" wird abgelehnt, "+ 
                                "Teststatistik: T="+Utilities.out(test_statistic_sigmaq)+".\n");
	}
	else{
	    System.out.println("Die Nullhypothese sigma² = "+Utilities.out(sigmaq)+" wird nicht abgelehnt," +
                                "Teststatistik: T="+Utilities.out(test_statistic_sigmaq)+".\n");
	}


	/* Algorithmus 2 */
        System.out.println("Algorithmus 2");
	for(int i=0; i < sampleSize; i++){
	    values = algorithm2(m,rand_seed.nextLong());
	    sample[i] = Utilities.getMax(values);
	}
	/* Bestimme Mittelwert und Stichprobenvarianz */
	mean = Utilities.getMean(sample); var = Utilities.getVar(sample);
	System.out.println("Erwartungswert: "+Utilities.out(mean)+"\n");
	System.out.println("Standardabweichung: "+Utilities.out(var)+"\n");

	/* Teste auf signifikante Abweichungen vom Erwartungswert */
	if(Math.abs(test_statistic_mu) > 1.96){
	    System.out.println("Die Nullhypothese mu = 0.797 wird abgelehnt, "+ 
                               "Teststatistik: T="+Utilities.out(test_statistic_mu)+".\n");
	}
	else{
	    System.out.println("Die Nullhypothese mu = 0.797 wird nicht abgelehnt, "+ 
                                "Teststatistik: T="+Utilities.out(test_statistic_mu)+".\n");
	}

	/* Teste auf signifikante Abweichungen von der theoretischen Varianz */
	mu4 = Utilities.getFourthMoment(sample);
	if( mu4 > sigmaq*sigmaq){
	    test_statistic_sigmaq =(var - sigmaq)/Math.sqrt( (mu4 - sigmaq*sigmaq)/sampleSize );
	}
	else{
	    test_statistic_sigmaq = 5.0;
	}
	if(Math.abs(test_statistic_sigmaq) > 1.96){
	    System.out.println("Die Nullhypothese sigma² = "+Utilities.out(sigmaq)+" wird abgelehnt,"+ 
                                "Teststatistik: T="+Utilities.out(test_statistic_sigmaq)+".\n");
	}
	else{
	    System.out.println("Die Nullhypothese sigma² = "+Utilities.out(sigmaq)+" wird nicht abgelehnt, "+
                                "Teststatistik: T="+Utilities.out(test_statistic_sigmaq)+".\n");
	}
    }
}



//Ergebnisse für den E-Wert und die Varianz des Maximums
//    Algorithmus 1
//    Erwartungswert: 0,780
//    Standardabweichung: 0,417
//
//    Die Nullhypothese mu = 0.797 wird nicht abgelehnt, Teststatistik: T=-0,280.
//    Die Nullhypothese sigma² = 0,603 wird abgelehnt, Teststatistik: T=-3,123.
//
//    Algorithmus 2
//    Erwartungswert: 0,744
//    Standardabweichung: 0,341
//
//    Die Nullhypothese mu = 0.797 wird nicht abgelehnt, Teststatistik: T=-0,280.
//    Die Nullhypothese sigma² = 0,603 wird abgelehnt,Teststatistik: T=-7,407.