Deconvolution/2D+: Difference between revisions

=={{header|Perl 6}}==

Translation of Tcl

<lang perl6># Deconvolution of N dimensional matricies.

sub deconv_ND ( @g, @f ) {

my @gsize = size_of @g;

my @fsize = size_of @f;

my @hsize = @gsize >>-<< @fsize >>+>> 1;

my @toSolve;

for loopcoords(@gsize) -> $coords {

@toSolve.push( [ row( @g, @f, @gsize, $coords, @fsize, @hsize ) ]);

}

my @solved = rref( @toSolve );

# Uncomment if you want to see the rref system of equations.

# pretty_print( @solved );

my @h;

my $index = 0;

for loopcoords(@hsize) -> @coords {

insert( @h, @coords, @solved[$index][*-1] );

$index++;

}

return @h;

# Inserts a value in the correct spot in an N dimensional array.

sub insert ( $array is rw, @coords is copy, $value ) {

my $level = @coords.shift;

if +@coords {

insert( $array[$level], @coords, $value );

} else {

$array[$level] = $value;

}

}

}

# Returns a list containing the number of elements in

# each level of an N dimensional array.

sub size_of ( $m is copy ) {

my @size;

while $m ~~ Array { @size.push(+$m); $m = $m[0]; }

return @size;

}

# Construct a row (equation) for each value in @g to be sent

# to the simultaneous equation solver.

# @Xsize = Dimensions of @X, # of elems per level.

# @Xcoords = Path to each element of @X given as a series of indicies.

sub row ( @g, @f, @gsize, @gcoords, @fsize, $hsize ) {

my @row;

for loopcoords( $hsize ) -> @hcoords {

my @fcoords;

for ^@hcoords -> $index {

my $window = @gcoords[$index] - @hcoords[$index];

@fcoords.push($window) and next if 0 <= $window < @fsize[$index];

last;

}

@row.push( +@fcoords == +@hcoords ?? fetch( @f, |@fcoords ) !! 0 );

}

@row.push( fetch( @g, |@gcoords ) );

return @row;

# Returns the value found in @array with the

# coordinates given in the list of @indicies.

sub fetch (@array, *@indicies) {

my $index = @indicies.shift;

return @array[*-1] ~~ Array

?? fetch( @array[$index], @indicies )

!! @array[$index];

}

}

# Constructs an array of arrays of coordinates to each

# element in an N dimensional array.

sub loopcoords ( @hsize ) {

my @hcoords;

for ^([*] @hsize) -> $index {

my @coords;

my $j = $index;

for @hsize -> $dim{

@coords.push( $j % $dim );

$j div= $dim;

}

@hcoords.push( [@coords] );

}

return @hcoords;

}

# Reduced Row Echelon Form simultaneous equation solver.

# Can handle over-specified systems of equations.

# (n unknowns in n + m equations)

sub rref ($m is rw) {

return unless $m;

my ($lead, $rows, $cols) = 0, +$m, +$m[0];

# Trim off over specified rows if they exist.

# Not strictly necessary, but can save a lot of

# redundant calculations.

if $rows >= $cols {

$m = trim_system($m);

$rows = +$m;

}

for ^$rows -> $r {

$lead < $cols or return $m;

my $i = $r;

until $m[$i][$lead] {

++$i == $rows or next;

$i = $r;

++$lead == $cols and return $m;

}

$m[$i, $r] = $m[$r, $i] if $r != $i;

my $lv = $m[$r][$lead];

$m[$r] >>/=>> $lv;

for ^$rows -> $n {

next if $n == $r;

$m[$n] >>-=>> $m[$r] >>*>> $m[$n][$lead];

}

++$lead;

}

return $m;

# Reduce a system of equations to n equations with n unknowns.

# Looks for an equation with a true value for each position.

# If it can't find one, assumes that it has already taken one

# and pushes in the first equation it sees. This assumtion

# will alway be successful except in some cases where an

# under-specified system has been supplied, in which case,

# it would not have been able to reduce the system anyway.

sub trim_system ($m is rw) {

my ($vars, @t) = +$m[0]-1, ();

for ^$vars -> $lead {

for ^$m -> $row {

@t.push( $m.splice( $row, 1 ) ) and last if $m[$row][$lead];

}

}

while (+@t < $vars) and +$m { @t.push( $m.splice( 0, 1 ) ) };

return @t;

}

}</lang>

Use with a pretty printer as follows:

<lang perl6># Pretty printer for N dimensional arrays. Assumes that

# if the FIRST element in any particular level is an array,

# then ALL the elements at that level are arrays.

sub pretty_print ( @array, $indent = 0 ) {

my $tab = 2;

if @array[0] ~~ Array {

say ' ' x $indent,"[";

pretty_print( $_, $indent + $tab ) for @array;

say ' ' x $indent, "]{$indent??','!!''}";

} else {

say ' ' x $indent, "[{say_it(@array)} ]{$indent??','!!''}";

}

sub say_it ( @array ) {

return join ",", @array>>.fmt("%4s");

}

}

my @f = (

[

[ -9, 5, -8 ], [ 3, 5, 1 ],

],

[

[ -1, -7, 2 ], [ -5, -6, 6 ],

],

[

[ 8, 5, 8 ], [ -2, -6, -4 ],

]

);

my @g = (

[

[ 54, 42, 53, -42, 85, -72 ],

[ 45,-170, 94, -36, 48, 73 ],

[ -39, 65,-112, -16, -78, -72 ],

[ 6, -11, -6, 62, 49, 8 ],

],

[

[ -57, 49, -23, 52,-135, 66 ],

[ -23, 127, -58, -5,-118, 64 ],

[ 87, -16, 121, 23, -41, -12 ],

[ -19, 29, 35,-148, -11, 45 ],

],

[

[ -55,-147,-146, -31, 55, 60 ],

[ -88, -45, -28, 46, -26,-144 ],

[ -12,-107, -34, 150, 249, 66 ],

[ 11, -15, -34, 27, -78, -50 ],

],

[

[ 56, 67, 108, 4, 2, -48 ],

[ 58, 67, 89, 32, 32, -8 ],

[ -42, -31,-103, -30, -23, -8 ],

[ 6, 4, -26, -10, 26, 12 ],

]

);

=begin skip_output

say "@g =";

pretty_print( @g );

say '-' x 79;

say "@f =";

pretty_print( @f );

say '-' x 79;

=end skip_output

say "# {+size_of(@f)}D array:";

say "@h =";

pretty_print( deconv_ND( @g, @f ) );</lang>

Output:

<pre># 3D array:

@h =

[

[

[ -6, -8, -5, 9 ],

[ -7, 9, -6, -8 ],

[ 2, -7, 9, 8 ],

],

[

[ 7, 4, 4, -6 ],

[ 9, 9, 4, -4 ],

[ -3, 7, -2, -3 ],

],

]</pre>