PerlでReverbを実装する

といっても、後部残響音の方。
シュレーダーによる生成アルゴリズムを実装してみた。

やることは、4つのコムフィルタを並列で処理して、
その出力を2回オールパスフィルタに入力するだけ。

use v5.14;
use strict;
use warnings;
use constant X_MAX => 1024;
use constant Y_MAX =>  256;
use constant Y_MIN => -256;

use constant VALUE_MAX =>  256;
use constant VALUE_MIN => -256;
use constant X_STEP => 2;

use Imager;
use List::Util qw/sum/;
use Math::Trig qw/pi/;

my $margin = 10;
my $width = $margin + (X_MAX + 1) + $margin;
my $height = $margin + (Y_MAX - Y_MIN + 1) + $margin;
my ( $x0, $y0 ) = ( $margin, $margin + (Y_MAX + 1) );


my $f1 = comb_filters(
    [ 0.5, 31 ],
    [ 0.5, 37 ],
    [ 0.5, 41 ],
    [ 0.5, 43 ]
);

my $f2 = all_pass_filter( 0.7, 17 );
my $f3 = all_pass_filter( 0.7, 5 );

my $f = sub {
    my $in = shift;
    return $f3->( $f2->($f1->($in)) );
};

=pod
my $f = comb_filters(
    [ 0.5, 31 ],
    [ 0.5, 37 ],
    [ 0.5, 41 ],
    [ 0.5, 43 ]
);
=cut

my $log = {
    in => [],
    out => []
};

my $in = VALUE_MAX;
foreach ( 0..(X_MAX / X_STEP) ) {
    push @{$log->{in}}, $in;
    push @{$log->{out}}, $f->( $in );
    $in = 0;
}
write_graph( 'aaa.png', {
    red     => $log->{in},
    blue    => $log->{out}
} );

sub write_graph {
    my ( $dst_file, $waves ) = @_;

    my $img = Imager->new(
        xsize => $width, ysize => $height, channels => 4 );
    $img->box( filled => 1, color => 'white' );
    draw_graduation( $img, Imager::Color->new(192, 192, 192) );

    foreach my $color ( keys %{$waves} ) {
        my $x = 0;
        my @tmp = ();
        foreach ( @{$waves->{$color}} ) {
            push @tmp, [ $x, $_ ] if 0 != $_; # <<< 0は無視!
            $x += X_STEP;
        }

        draw_points_with_bar( $img, \@tmp, $color );
    }

    $img->write( file => $dst_file ) or die $img->errstr;
}

sub calc_points {
    my ( $fa, $fb ) = @_;

    my @src = ();
    for (my $i=0; $i<=X_MAX; $i+=8) {
        my $t = $i / X_MAX;
        my $y = ($fa * $t) + ($fb * $t ** 5.0);
        push @src, [
            $i,
            int( ($y * 128) )
        ];
    }

    return \@src;
}

sub draw_graduation {
    my ( $img, $color ) = @_;

    {
        my $gray = Imager::Color->new( 192, 192, 192 );

        my $x = 128;#($w / 4);
        while ( $x <= X_MAX ) {
            $img->line( color => $gray,
                x1 => $x0 + $x, y1 => $y0 + Y_MIN,
                x2 => $x0 + $x, y2 => $y0 + Y_MAX );
            $x += 128;#($w / 4);
        }

        my $y = 128;#($h / 4);
        while ( $y <= Y_MAX ) {
            $img->line( color => $gray,
                x1 => $x0 + 0,     y1 => $y0 - $y,
                x2 => $x0 + X_MAX, y2 => $y0 - $y );
            $img->line( color => $gray,
                x1 => $x0 + 0,     y1 => $y0 + $y,
                x2 => $x0 + X_MAX, y2 => $y0 + $y );
            $y += 128;#($h / 4);
        }
    }

    {
        $img->line( color => 'black',
            x1 => $x0, y1 => $y0 + Y_MIN,
            x2 => $x0, y2 => $y0 + Y_MAX );

        $img->line( color => 'black',
            x1 => $x0 + 0,     y1 => $y0,
            x2 => $x0 + X_MAX, y2 => $y0 );
    }
}

sub plot_points {
    my ( $img, $data, $color ) = @_;
    my $n = 1;

    foreach my $pt ( @{$data} ) {
        my ( $x, $y ) = ( $pt->[0], $pt->[1] );
        $y /= (VALUE_MAX / Y_MAX);
        $y = ( $y < .0 ) ? int($y - .5) : int($y + .5);
        #printf( "%6.3f, %6.3f\n", $x, $y );
        $img->box( xmin => $x0 + $x - $n, ymin => $y0 - $y - $n,
                   xmax => $x0 + $x + $n, ymax => $y0 - $y + $n,
                   color => $color, filled => 0 );
    }
}

sub draw_points {
    my ( $img, $data, $color ) = @_;

    foreach my $pt ( @{$data} ) {
        my ( $x, $y ) = ( $pt->[0], $pt->[1] );
        $y /= (VALUE_MAX / Y_MAX);
        $y = ( $y < .0 ) ? int($y - .5) : int($y + .5);
        #printf( "%6.3f, %6.3f\n", $x, $y );
        $img->setpixel( x => $x0 + $x, y => $y0 - $y, color => $color );
        #$img->box( xmin => $x0 + $x - $n, ymin => $y0 - $y - $n,
        #           xmax => $x0 + $x + $n, ymax => $y0 - $y + $n,
        #           color => $color, filled => 1 );
    }
}

sub draw_polyline {
    my ( $img, $data, $color ) = @_;

    my @points = map {
        my ( $x, $y ) = ( $_->[0], $_->[1] );
        $y /= (VALUE_MAX / Y_MAX);
        $y = ( $y < .0 ) ? int($y - .5) : int($y + .5);
        [ $x0 + $x, $y0 - $y ];
    } @{$data};

    $img->polyline( points => \@points, color => $color );
}

# 棒グラフっぽく
sub draw_points_with_bar {
    my ( $img, $data, $color, $opacity ) = @_;

    my $img_dst = $img; 
    if ( defined($opacity) and $opacity < 1.0 ) {
        $img_dst = Imager->new(
            xsize => $img->getwidth(), ysize => $img->getheight(), channels => 4 );
    }

    foreach my $pt ( @{$data} ) {
        my ( $x, $y ) = ( $pt->[0], $pt->[1] );
        $y /= (VALUE_MAX / Y_MAX);
        $y = ( $y < .0 ) ? int($y - .5) : int($y + .5);

        $img_dst->line(
            x1 => ($x0 + $x), y1 => $y0,
            x2 => ($x0 + $x), y2 => ($y0 - $y),
            color => $color );
    }

    if ( $img != $img_dst ) {
        $img->compose(
            src => $img_dst, opacity => $opacity, combine => 'add' );
    }

    plot_points( $img, $data, $color );
}

sub all_pass_filter {
    my ( $aa, $n, $fill ) = @_;
    $n = 1 if not defined($n);
    $fill = 0 if not defined($fill);

    my @za = map { $fill; } 1..$n;
    my @zb = map { $fill; } 1..$n;

    return sub {
        my $in = shift;
        unshift @za, $in;
        my $out = (pop(@zb) * $aa) + pop(@za) - ($aa * $in);
        unshift @zb, $out;
        return $out;
    };
}

sub comb_filter {
    my ( $aa, $n, $fill ) = @_;
    $n = 1 if not defined($n);
    $fill = 0 if not defined($fill);

    my @za = map { $fill; } 1..$n;

    return sub {
        my $in = shift;
        my $out = pop @za;
        unshift @za, ($in + ($aa * $out));
        return $out;
    };
}

sub comb_filters {
    my @args = @_;

    my @filters = ();
    foreach ( @args ) {
        push @filters, comb_filter( @{$_} );
    }

    return sub {
        my $in = shift;
        return sum( map {
            $_->( $in );
        } @filters );
    };
}

コムフィルタだけ通すと、
インパルス応答はディレイが4つみたいな感じになる。

20161001-1

それをオールパスフィルタで拡散する。
つまり、乱反射を再現してるんだと思われる。

20161001-2

という訳で、本日のネタ元はこちら。

あと、このページの方が分かりやすいかも。
(このページ以外も必見!)

リバーブ2/4 リバーブレーター(残響装置)
http://www.ari-web.com/service/soft/reverb-2.htm

次のアプリに、Reverbが搭載される訳ではないので悪しからず。

おしまい。

Leave a Comment