#!/usr/bin/env perl

use strict;
use warnings;

use Statistics::Test::RandomWalk;
use  Math::Random::MT;

my $NUM_TRIALS ||=100;
my $NUM_BINS ||= 10;

my ($quant, $got, $expected);
my $MTgen = Math::Random::MT->new();

sub     perls_rand { rand() }
sub     mt_rand { $MTgen->rand }
my %results = ();
my $best_results = [ undef, undef ];
foreach( 
        (["RandomWalk test : Perl's rand()", \&perls_rand, $NUM_TRIALS, $NUM_BINS],
         ["RandomWalk test : Math::Random::MT's rand()", \&mt_rand, $NUM_TRIALS, $NUM_BINS],
        )
){
        my $aresult = do_a_test(@$_);
        $results{$aresult->[3]} = $aresult;
        print $aresult->[0]."\n\n";
        if( ! defined($best_results->[0]) or $best_results->[0] > $aresult->[1] ){
                $best_results->[0] = $aresult->[1]; # mean perc diff
                $best_results->[1] = $_->[0];
        }
}
print "\n\n$0 : done, here is a summary:\n";
foreach(values %results){
        print "$0 : ".$_->[3]." with mean %diff from expected to be ".sprintf("%.3f",$_->[1])." %\n";
}
print "\n\n$0 : best results for ".$best_results->[1]." with mean %diff from expected to be ".sprintf("%.3f",$best_results->[0])." %\n";
print "$0 : end.\n";
exit(0);

# Does a trial given
# a test label (string)
# a sub ref to a rand()
# the number of iterations
# the number of bins
# see http://search.cpan.org/dist/Statistics-Test-RandomWalk/lib/Statistics/Test/RandomWalk.pm
# returns an arrayref with a report, mean percentage difference and change
# and iters and nbins (see the return statement below)
sub     do_a_test {
        my ($test_label, $rand_sub, $niters, $nbins) = @_;
        my $tester_MT = Statistics::Test::RandomWalk->new();
        defined($tester_MT) or die "tester_MT new() failed.";
        $tester_MT->set_data( [map {$rand_sub->()} 1..$niters] );
        my ($quant, $got, $expected) = $tester_MT->test($nbins);
        my ($subreport, $mean_diff, $mean_change) = @{my_report_make($quant, $got, $expected)};
        my $report = "$test_label\n                  after $niters iterations and using $nbins bins:\n"
                . $subreport
                . "======================= end ========================"
        ;

        return [
                $report,
                $mean_diff,
                $mean_change,
                $test_label, $niters, $nbins
        ]
}
# returns an arrayref of
#   a string report which can be printed
#   the mean percentage difference ( 200 * abs(A-B)/(A+B) %)
#   the mean percentage change ( 100 - 100*A/B % )
sub     my_report_make {
        my ($quant, $got, $expected) = @_;
        my $N = scalar(@$quant);
        my ($i, $percent_diff, $percent_change);
        my $sum_diff = 0;
        my $sum_change = 0;
        my $ret = "Quantile |         Got  |     Expected  | % diff   | % change\n==============================================================\n";
        for($i=0;$i<$N;$i++){
                $percent_diff = percentage_difference($got->[$i], $expected->[$i]);
                $percent_change = percentage_change($got->[$i], $expected->[$i]);
                $ret .= sprintf("%-11.3f%12.3f %14.3f%10.3f%12.3f\n", $quant->[$i], $got->[$i], $expected->[$i], $percent_diff, $percent_change);
                $sum_diff += $percent_diff;
                $sum_change += $percent_change;
        }
        my $mean_diff = $sum_diff/$N;
        my $mean_change = $sum_change/$N;
        $ret .= sprintf("\ntotal percentage difference between got and expected: %.7f %%\n>>>> mean percentage difference: %.7f %% <<<<\n", $sum_diff, $mean_diff);
        $ret .= sprintf("\ntotal percentage change between got and expected: %.7f %%\n>>>>  mean percentage change: %.7f %% <<<<\n", $sum_change, $mean_change);
        return [$ret, $mean_diff, $mean_change]
}
sub     percentage_difference {
        my ($a, $b) = @_;
        return 200 * abs($a-$b) / ($a+$b)
}
sub     percentage_change {
        my ($a, $b) = @_;
        return 100 - 100 * $a / $b
}