Nice, but trying to play with this I found that it requires PDL::IO::Image which depends on Alien::FreeImage which won't build on modern g++ (it requires -std=c++14 and redefines powf64) and it has security issues (see issues).

Other than that, thanks! Some things to look back at when I encounter similar tasks.

I suspected result on Linux be wilder than "wrong_setup" because terminals are almost certain to use palette, not pure RGB. "Security"? Life corrects as always.

##################################
########### Deal with image
##################################

use PDL::IO::FlexRaw;

my $image;
{
    my $im = Imager-> new( file => $image_file )
        -> to_rgb8
        -> convert( preset => 'noalpha' )
        -> scale( xscalefactor => 1, 
                  yscalefactor => 1 / $distortion );

    my $w_r   = $cols / $im-> getwidth;
    my $h_r   = $rows / $im-> getheight;
    my $ratio = $w_r < $h_r ? $w_r : $h_r;
    $im = $im-> scale( scalefactor => $ratio ) if $ratio < 1;

    $im-> write( data => \my $buf, type => 'raw' );
    my $dims = [ 3, $im-> getwidth, $im-> getheight ];
    open my $fh, '<:raw', \$buf;
    $image = readflex( $fh, [{ Type => 'byte', Dims => $dims }])
}
[download]

TIMTOWTDI. One reason to use PDL::IO::Image was to play with its rescale filters (which amounted to nothing) to average glyph bitmap to single pixel. In the end it's arithmetic mean with Imager::scale(qtype => 'mixing',... and itself not necessary neither, could be unpack '%32C*', ...

(Another reason for PDL::IO::Image was to use the best what happened to image IO and basic manipulation, for PDL if not Perl. Very sad if it's abandoned. PDL::IO::GD and PDL::IO::Pic don't compare)

Computations above take noticeable time, as opposed to simple algorithm in OP. If final section is formatted as:

##################################
########### Ready to go
##################################

use Time::HiRes 'time';

my $t = time;

$image   = t_lab-> apply(( !t_srgb )-> apply( $image ));
$palette = t_lab-> apply(( !t_srgb )-> apply( $palette ));

my $indexed = delta_e_2000( 
    $palette, $image-> dummy( 1, $plt_size )
)-> minimum_ind;

print time - $t, "\n";

print @pixels[ @$_ ], "\n" for @{ $indexed-> unpdl };
[download]

then time is 4.56991696357727 for, ultimately, a single Perl statement; using e.g. DejaVuSansMono and console/font size approximately as used to produce picture above. Not cheap. For 10 brushes (i.e. glyphs from this font) and all combinations of FG/BG colors (16*16), the palette size upper limit is 2560. It was reduced to 'only' 1854, excluding trivial duplicates. Distorted image size is 142 by 89.

For each of its 12638 pixels, calculate 1854 distances with

https://en.wikipedia.org/wiki/Color_difference#CIEDE2000

and pick the best palette index, i.e. ~23e6 times. Maybe we are lucky we have PDL, and this horrible formula was already implemented in C.

Unfortunately, it's immediately obvious it'll take same long time if image is a flat single color rectangle (or simple logo, etc.), while it must be 12638 times faster!

For our test image, distorted and scaled to 142x89, only half (6265 out of 12638) RGB triplets are unique (Imager::getcolorcount). It follows, run time should be twice as low.

I didn't find anything w.r.t. 'memoization' and PDL, maybe not tried hard enough? Duplication in arrays and expensive computations are inevitable. It only requires to reduce fat data to unique slim set, keeping track how to fatten result to original shape. Not 'memoization' per definition, but in spirit.

##################################
########### Ready to go
##################################

use Time::HiRes 'time';

my $t = time;

my ( undef, $w, $h ) = $image-> dims;

my $fat_rgb     = $image-> clump( 1, 2 );
my $perm_2_sort = $fat_rgb-> qsortveci;
my $sorted      = $fat_rgb-> dice_axis( 1, $perm_2_sort );
my $enumerated  = $sorted-> enumvecg;
my $selector    = $enumerated-> uniqind;
my $slim_rgb    = $sorted-> dice_axis( 1, $selector );

$slim_rgb = t_lab-> apply(( !t_srgb )-> apply( $slim_rgb ));
$palette  = t_lab-> apply(( !t_srgb )-> apply( $palette ));

my $slim_indexed = delta_e_2000( 
    $palette, $slim_rgb-> dummy( 1, $plt_size )
)-> minimum_ind;

my $indexed = $slim_indexed
    -> index( $enumerated )             # fatten
    -> index( $perm_2_sort-> qsorti )   # permute to original
    -> reshape( $w, $h );               # reshape :)

print time - $t, "\n";

print @pixels[ @$_ ], "\n" for @{ $indexed-> unpdl };
[download]

And time now is exactly as predicted: 2.25233697891235. (PDL has uniq and uniqind pair. But only uniqvec i.e. no uniqvecind. Hence workaround with vector enumeration)

###############

In fact, RGB (or Lab) points in palette are spaced so very sparsely. If further speed is required, some precision can be sacrificed, and even slimmer dataset fed to expensive formula. With

$fat_rgb-> sever-> inplace-> and2( 0b1111_1110 );
[download]

i.e. RGB in image reduced to 7 bits (if GD always does it (to alpha), why can't we?), time becomes 1.3685450553894. With 2 bits zeroed, it's 0.765566825866699; and result is practically the same if observer was not warned. With reduction to 5 bits (again 2x faster) the picture is beginning to feel slightly uncomfortable compared to 'lossless' result.

###############

W.r.t. need to sever above, I'm confused if clump creates dataflow. But what follows looks like plain bugs:

pdl> $x = zeroes 3,2,2; $y = $x-> clump( 1,2 )
pdl> $x .= 3
pdl> p $x, $y, $y-> uniqvec

[
 [
  [3 3 3]
  [3 3 3]
 ]
 [
  [3 3 3]
  [3 3 3]
 ]
]

[
 [3 3 3]
 [3 3 3]
 [3 3 3]
 [3 3 3]
]

[
 [3 3 3]
]

pdl> $x = zeroes 3,2,2; $y = $x-> clump( 1,2 )
pdl> $y .= 3
pdl> p $x, $y, $y-> uniqvec

[
 [
  [3 3 3]
  [3 3 3]
 ]
 [
  [3 3 3]
  [3 3 3]
 ]
]

[
 [3 3 3]
 [3 3 3]
 [3 3 3]
 [3 3 3]
]

[
 [0 0 0]
]
[download]

So far, only uniqvec looks weird. But:

pdl> $x = zeroes 3,2,2; $y = $x-> clump( 1,2 )
pdl> p $y

[
 [0 0 0]
 [0 0 0]
 [0 0 0]
 [0 0 0]
]

pdl> $x .= 3
pdl> p $x, $y, $y-> uniqvec

[
 [
  [3 3 3]
  [3 3 3]
 ]
 [
  [3 3 3]
  [3 3 3]
 ]
]

[
 [0 0 0]
 [0 0 0]
 [0 0 0]
 [0 0 0]
]

[
 [0 0 0]
]
[download]

What?! But I only printed the $y, compared to 1st demo.

Using FlexRaw to read the PV of a scalar as an ndarray is very clever, nice work! However, a more idiomatic way might be to replace the last 4 lines with (untested):

$image = zeroes byte, $im->getwidth, $im->getheight;
$im->write( data => $image->get_dataref, type => 'raw' );
$image->upd_data;
[download]

Another way to read images into PDL with lots of capability is PDL::OpenCV::Imgcodecs.