File Coverage

blib/lib/DSP/LinPred_XS.pm

Criterion	Covered	Total	%
statement	62	96	64.5
branch	7	16	43.7
condition	8	9	88.8
subroutine	7	10	70.0
pod	5	7	71.4
total	89	138	64.4

package DSP::LinPred_XS;
use 5.008005;
use Mouse;
our $VERSION = "0.02";
use XSLoader;
XSLoader::load(__PACKAGE__, $VERSION);




has 'mu' => (
    is => 'rw',
    isa => 'Num',
    default => 0.001
    );
has 'mu_mode' => (
    is => 'rw',
    isa => 'Int',
    default => 0
    );
has 'h_length' => (
    is => 'rw',
    isa => 'Int',
    default => 100
    );
has 'h' => (
    is => 'rw',
    isa => 'ArrayRef[Num]',
    default => sub{[(0) x 100]}
    );
has 'x_stack' => (
    is => 'rw',
    isa => 'ArrayRef[Num]',
    default => sub{[(0) x 100]}
    );
has 'x_count' => (
    is => 'rw',
    isa => 'Int',
    default => 0
    );
has 'current_error' => (
    is => 'rw',
    isa => 'Num',
    default => 0
    );
has 'dc' => (
    is => 'rw',
    isa => 'Num',
    default => 0
    );
has 'dc_init' => (
    is => 'rw',
    isa => 'Num',
    default => 0
    );
has 'stddev' => (
    is => 'rw',
    isa => 'Num',
    default => 0
    );
has 'stddev_init' => (
    is => 'rw',
    isa => 'Num',
    default => 1
    );
has 'est_mode' => (
    is => 'rw',
    isa => 'Num',
    default => 1
    );


# filter specification
# mu : step size
# h_length : filter size
sub set_filter{
    my $self = shift;
    my $conf = shift;
    if(defined($conf->{filter_length})){
	$self->h_length($conf->{filter_length});
	$self->h([(0) x $conf->{filter_length}]);
	if(defined($conf->{dc_init})){
	    $self->x_stack([($conf->{dc_init}) x $conf->{filter_length}]);
	}else{
	    $self->x_stack([(0) x $conf->{filter_length}]);
	}
    }
    if(defined($conf->{dc_init})){
	$self->dc($conf->{dc_init});
	$self->dc_init($conf->{dc_init});
    }
    if(defined($conf->{est_mode})){
	$self->est_mode($conf->{est_mode});
    }
    if(defined($conf->{stddev_init})){
        $self->stddev($conf->{stddev_init});
        $self->stddev_init($conf->{stddev_init});
    }
}

# reset filter state
sub reset_state{
    my $self = shift;
    my $h_length = $self->h_length;
    $self->h([(0) x $h_length]);
    $self->x_stack([($self->dc_init) x $h_length]);
    $self->current_error(0);
    $self->dc($self->dc_init);
    $self->x_count(0);
    $self->stddev($self->stddev_init);
}

# prediction only
# predict_num : number of output predicted values
# this method returns list reference of predicted values
sub predict{
    my $self = shift;
    my $predict_num = shift;
    my $h = $self->h;
    my $x_stack = $self->x_stack;
    my $estimated;
    for(1 .. $predict_num){
	my $x_est = 0;
	for( my $k = 0; $k <= $#{$h} and $k <= $self->x_count; $k++){
	    $x_est += $h->[$k] * ($x_stack->[$k] - $self->dc);
	}
	$x_est += $self->dc;
	unshift(@$x_stack,$x_est);
	push(@$estimated,$x_est);
	pop(@$x_stack);
    }
    return($estimated);
}

# update only
# x should be array reference

sub update{
    my $self = shift;
    my $x = shift;
    my $h_length = $self->h_length;
    my $h = $self->h;
    my $x_stack = $self->x_stack;

    for ( my $kx=0; $kx <= $#{$x}; $kx++){
	unshift(@$x_stack,$x->[$kx]);
	pop(@$x_stack);
	$self->x_count($self->x_count + 1);
	if($self->est_mode == 1){
	    $self->dc_stddev_update;
	}
	my $x_est = 0;
	for( my $k = 0; $k <= $#{$h} and $k <= $self->x_count;$k++){
	    $x_est += $h->[$k] * ($x_stack->[$k] - $self->dc);
	}
	my $error = $x->[$kx] - ($x_est + $self->dc);
	$self->current_error($error);
	my $h_new = $h;
	my $tmp_coef = 1;
	if($self->est_mode == 1){
	    $tmp_coef = $self->mu * $error / (1 + $self->stddev);
	}else{
	    $tmp_coef = $self->mu * $error;
	}
	if($self->mu_mode == 1){
	    $tmp_coef = 10 * $self->mu / (1 + $self->h_length);
	}

	for(my $k = 0;$k <= $#{$h} and $k <= $self->x_count; $k++){
	    $h_new->[$k] =
		$h->[$k]
		+ $tmp_coef * ($x_stack->[$k] - $self->dc);
	}
	$self->h($h_new);
    }
}

## DC component calculation and update
# using x_stack

sub dc_stddev_update{
    my $self = shift;
    my $x_stack = $self->x_stack;
    my ($sum,$mean,$variance,$stddev) = &get_stat($x_stack);
    $self->dc($mean);
    $self->stddev($stddev);
}


## calculation of mean value of filter
sub filter_dc{
    my $self = shift;
    my $h = $self->h;
    my $mean = 0;
    my $num = $#$h + 1;
    for(0 .. $#$h){
        $mean += $h->[$_];
    }
    return($mean / $num);
}

## calculation of stddev of filter
sub filter_stddev{
    my $self = shift;
    my $h = $self->h;
    my $variance = 0;
    my $num = $#$h + 1;
    for(0 .. $#$h){
        $variance += ($h->[$_])**2;
    }
    return(sqrt($variance / $num));
}



1;
__END__