view release on metacpan or  search on metacpan

lib/Math/BSpline/Basis.pm  view on Meta::CPAN

    return $class->$orig($munged_args);
};



has 'degree' => (
    is       => 'ro',
    required => 1,
);



has 'knot_vector' => (
    is      => 'lazy',
    builder => sub {
        my ($self) = @_;
        my $p      = $self->degree;

        return [
            (map { 0 } (0..$p)),
            (map { 1 } (0..$p)),
        ]
    }
);



# I use the same variable names as in the NURBS book, although some
# of them are very generic. The use of $p, $U, $P, and $n is
# consistent throughout the relevant chapters of the book.
sub find_knot_span {
    my ($self, $u) = @_;
    my $p          = $self->degree;
    my $U          = $self->knot_vector;
    my $n          = (@$U - 1) - $p - 1;

    # We expect $u in [$U->[$p], $U->[$n+1]]. We only support
    # values outside this range for rounding errors, do not assume
    # that the result makes sense otherwise.
    return $n if ($u >= $U->[$n+1]);
    return $p if ($u <= $U->[$p]);

    # binary search
    my $low  = $p;
    my $high = $n + 1;
    my $mid  = int(($low + $high) / 2);
    while ($u < $U->[$mid] or $u >= $U->[$mid+1]) {
        if ($u < $U->[$mid]) { $high = $mid }
        else                 { $low  = $mid }
        $mid = int(($low + $high) / 2);
    }

    return $mid;
}



# The variable names are inspired by the theory as laid out in the
# NURBS book. We want to calculate N_{i,p}, that inspires $N and
# $p. U is the knot vector, left and right are inspired by the
# terms in the formulas used in the theoretical derivation.
sub evaluate_basis_functions {
    my ($self, $i, $u) = @_;
    my $p              = $self->degree;
    my $U              = $self->knot_vector;
    my $n              = (@$U - 1) - $p - 1;

    if ($u < $U->[$p] or $u > $U->[$n+1]) {
        return [map { 0 } (0..$p)];
    }

    my $N     = [1];
    my $left  = [];
    my $right = [];
    for (my $j=1;$j<=$p;$j++) {
        $left->[$j]  = $u - $U->[$i+1-$j];
        $right->[$j] = $U->[$i+$j] - $u;
        my $saved = 0;
        for (my $r=0;$r<$j;$r++) {
            my $temp = $N->[$r] / ($right->[$r+1] + $left->[$j-$r]);
            $N->[$r] = $saved + $right->[$r+1] * $temp;
            $saved   = $left->[$j-$r] * $temp;
        }
        $N->[$j] = $saved;
    }

    return $N;
}



sub evaluate_basis_derivatives {
    my ($self, $i, $u, $d) = @_;
    my $p                  = $self->degree;
    my $U                  = $self->knot_vector;
    my $n                  = (@$U - 1) - $p - 1;
    my $result             = [];

    $d = min($d, $p);

    if ($u < $U->[$p] or $u > $U->[$n+1]) {
        for (my $k=0;$k<=$d;$k++) {
            push(@$result, [map { 0 } (0..$p)]);
        }
        return $result;
    }

    my $ndu   = [[1]];
    my $left  = [];
    my $right = [];
    for (my $j=1;$j<=$p;$j++) {
        $left->[$j]  = $u - $U->[$i+1-$j];
        $right->[$j] = $U->[$i+$j] - $u;
        my $saved = 0;
        for (my $r=0;$r<$j;$r++) {
            $ndu->[$j]->[$r] = $right->[$r+1] + $left->[$j-$r];
            my $temp = $ndu->[$r]->[$j-1] / $ndu->[$j]->[$r];
            $ndu->[$r]->[$j] = $saved + $right->[$r+1] * $temp;
            $saved           = $left->[$j-$r] * $temp;
        }
        $ndu->[$j]->[$j] = $saved;