#!/usr/bin/perl -w
#use strict;

########## Variables declaration ###################################################################

my $dice_seq = '3151166666624';

my @obs_list = split(//,$dice_seq);

my %P_emission = (
	'F' => { 1 => 0, 2 => 0, 3 => 0, 4 => 0, 5 => 0, 6 => 0 },
	'L' => { 1 => -0.51, 2 => -0.51, 3 => -0.51, 4 => -0.51 , 5 => -0.51 , 6 => 1.10 } );
my %P_transition = ( 'FF' => 0.64, 'LL' => 0.59, 'FL' => -2.30, 'LF' => -1.61 );

my %path;


########################à VITERBI  algorithm #######################################################


# Initialization:
my %v = ( 'F' => [0], 'L' => [-0.51] );

# Recursion:
my $idx = 1;

foreach my $tmp_obs (@obs_list) {
  #if($obs_list[$idx]){my $tmp_obs=$obs_list[$idx];}else{my $tmp_obs="END";}
  print "i =  ".($idx)." obs_i = $obs_list[$idx] i-1 = ".($idx-1)." obs_i-1 = $obs_list[$idx-1]\n";
  my $max_state = '';
  my $max_state_v = 0;
	foreach my $state2 qw(F L) {
		my ($max_v, $max_prev_state, $max_path) = (-1,'X','XX');

		foreach my $state1 qw(F L) {
			my $tmp_trans = $state1.$state2;


			####### FILL IN ######### VITERBI ITERATION STEP ###### 
		#	my $tmp_v = 


			print sprintf("[%d]%s -> %s: %.2f * %.2f * %.2f = %.2f\n", 
											$idx, $state1, $state2, $P_emission{$state2}->{$tmp_obs}, $P_transition{$tmp_trans},
											$v{$state1}->[$idx-1], $tmp_v);											

			if( $tmp_v > $max_v ) {
				$max_v = $tmp_v; 
        $max_prev_state = $state1;
				$max_path = $state1.$state2;
			}
		}
		print sprintf("MAX PATH to %s: %s, %s%d = %.2f\n", $state2, $max_path, 
										$state2,$idx+1,$max_v);
    $path{$idx}->{$state2} = $max_prev_state;

		#push(@{$path{$state2}}, $max_path);
		$v{$state2}->[$idx] = $max_v;

	}
	print "\n";
	$idx += 1;
}

## Determine the last state
my $last_state = 'X';
if( $v{'F'}->[$idx-1] > $v{'L'}->[$idx-1] ) {
  $last_state = 'F';
} else {
  $last_state = 'L';
}

print "Final state: $last_state (v= $v{$last_state}->[$idx-1])\n";

## Traceback
my @state_seq = ();
for(my $i=$idx-1;$i>0;$i--) {
  my $prev_state = $path{$i}->{$last_state};
  push(@state_seq,$prev_state);
  $last_state = $prev_state;
}

## Print most probable path
print "\n[Output]\n";
print join("",@obs_list),"\n";
print join("",reverse @state_seq),"\n";


##### OUTPUT ADDITIONAL INFO ###########################################

##### VITERBI MAX PREVSTATE k PROBS : example of hash of arrays #########

print "\n\nVITERBI matrix:\n\n";
my $jseq = join("   \t", @obs_list);
print "\tB\t$jseq\n";
my $arrayref = $v{"F"};
my @Vf = @{ $arrayref };
my $Fvit=join("\t", @Vf);
print "F\t$Fvit\n";
$arrayref = $v{"L"};
my @Vl = @{ $arrayref };
my $Lvit=join("\t", @Vl);
print "L\t$Lvit\n\n";

##### PATH POINTERS MATRIX : example of hash of hashes ###############
print "\n\nTRACEBACK pointers:\n\n";
    my $rHoH = \%path;
    my $curpos=0;
    for my $k1 ( sort {$a <=> $b} (keys %$rHoH) ) {
        print "Step: $k1  \tEmission: $obs_list[$curpos]\n";
        for my $k2 ( keys %{$rHoH->{ $k1 }} ) {
            print "\tmax path to $k2 :  $rHoH->{ $k1 }{ $k2 }\t";
			if($k2 ne $rHoH->{ $k1 }{ $k2 }){print "*\n";}else{print "\n";}
        }
		$curpos++;
		print "\n";
    }