Looking at the screw-generating code reminded me of a coding challenge that defeated me 10 years ago, when I was first using PDL. In gimp-perl (the plugin for GIMP allowing Perl scripts as filters etc) there are various scripts contributed by people, nearly all before I ever came along. One of them is a fun one called "Intensity Landscape", which treats the (e.g.) red values of an image as intensity (i.e. height), then "looks at it from the side and above", and makes an image of what it sees (extracted and cut down from https://gitlab.gnome.org/GNOME/gimp-perl/-/blob/master/examples/iland):

use strict;
use warnings;
use PDL;

# user params
my $floor = 100;
my $component = 0;
my $delta = 6;
my $elevation = 2;
my $camerapos = -1;
$| = 1;

my $alpha = 0;
my $srcdata = rpic($ARGV[0]);
my $destdata = $srcdata->zeroes; # "filled" with black
my (undef, $width, $height) = $srcdata->dims;
my $relord = $width / 255;
$delta = 1 if $delta < 1;

for (my $y = 0; $y < $height; $y++) {
  my $row = $srcdata->slice('', '', "($y)");
  my $drow = $destdata->slice('', '', "($y)");
  my $red = $row->slice("($component)");
  my $dred = $drow->slice("(0)");
  my $dgreen = $drow->slice("(1)");
  my $dblue = $drow->slice("(2)");
  $drow->slice(3) .= 255 if $alpha;
  for (my $x = 0; $x < $width; $x++) {
print "." unless $x%10;
    my $r = at($red, $x);
    next if $r <= $floor;
    my $remain = $r;
    my $currentx = $width - $r * $relord + ($x / $elevation);
    #Apply elevation following the x offset in original picture
    while ($remain > 0 && $currentx < $width) {
      if ($remain > 150) {
        set ($dred, $currentx, 0);
        set ($dgreen, $currentx, $remain);
        set ($dblue, $currentx, $remain);
      }
      if ($remain < 150 && $remain > 50) {
        set ($dred, $currentx, 0);
        set ($dgreen, $currentx, $remain + 55);
        set ($dblue, $currentx, 0);
      }
      if ($remain < 50) {
        set ($dred, $currentx, 0);
        set ($dgreen, $currentx, 0);
        set ($dblue, $currentx, $remain + 200);
      }
      $remain -= $delta;
      $currentx++;
    }
  }
print "\n";#  Gimp::Progress->update($y / $height);
}

$destdata->wpic("OUT$ARGV[0]");
[download]

The script worked and made interesting pictures, but it was terribly slow (15 secs+ for a 300x300 image). This is largely because it's written like C, and Perl-loops over all X and Y coordinates, reading each pixel value, etc. I knew (sort of) that PDL could be used to do better, but for a long time I didn't really have any clue how.

Encouraged by the screw thingy, and realising the similarity of problem might mean the approach for one could be applied to the other, I dug back into it. This is the more PDL-idiomatic version, which now runs in <4sec, and is actually shorter:

use strict;
use warnings;
use PDL;

# user params
my $floor = 100;
my $component = 0;
my $delta = 6;
my $elevation = 2;
my $camerapos = -1;
$| = 1;

my $alpha = 0;
my $srcdata = rpic($ARGV[0]);
my $destdata = $srcdata->zeroes; # "filled" with black
$destdata->slice(3) .= 255 if $alpha;
my $destmv = $destdata->mv(0,-1); # x y rgb
my (undef, $width, $height) = $srcdata->dims;
my $relord = $width / 255;
$delta = 1 if $delta < 1;

my $quant = ($srcdata->slice("($component)")->max / $delta)->floor->sc
+lr;
return if $quant <= 0;
for my $x (0..$width-1) {
  my $col = $srcdata->slice("($component),($x)");
  my $exceed_floor = ($col > $floor);
  my $r = $col->where($exceed_floor); # nvals
  my $destx = ($width - $r * $relord + ($x / $elevation))->long; # nva
+ls
  #Apply elevation following the x offset in original picture
  my $remain_s = zeroes(long, 3, $quant, $r->dim(0)); # xyr quant nval
+s
  my $yslice = $remain_s->slice("(1)") .= $exceed_floor->which->dummy(
+0); # quant nvals
  my $xslice = $remain_s->slice("(0)") .= $yslice->xvals + $destx->dum
+my(0); # quant nvals
  my $rslice = $remain_s->slice("(2)") .= $yslice->xlinvals(0,-1) * $q
+uant*$delta + $r->dummy(0); # quant nvals
  $rslice->whereND($xslice >= $width) .= -1;
  my $gt150_ind = whichND($rslice > 150);
  my $btwn_ind = whichND(($rslice <= 150) & ($rslice >= 50));
  my $lt50_ind = whichND(($rslice < 50) & ($rslice > 0));
  $destmv->slice(',,1:2')->indexND(cat(map $_->indexND($gt150_ind), $x
+slice, $yslice)->mv(-1,0)) .= $rslice->indexND($gt150_ind) if $gt150_
+ind->nelem;
  $destmv->slice(',,1')->indexND(cat(map $_->indexND($btwn_ind), $xsli
+ce, $yslice)->mv(-1,0)) .= $rslice->indexND($btwn_ind) + 55 if $btwn_
+ind->nelem;
  $destmv->slice(',,2')->indexND(cat(map $_->indexND($lt50_ind), $xsli
+ce, $yslice)->mv(-1,0)) .= $rslice->indexND($lt50_ind) + 200 if $lt50
+_ind->nelem;
#  Gimp::Progress->update($x / $height);
}

$destdata->wpic("OUT$ARGV[0]");
[download]