vr++ correctly identifies this is a histogram type of problem. This note is just about an alternative approach I stumbled on very recently.

The new PDL::IO::STL, in the imminent PDL 2.079, uses PDL::VectorValued if installed, to deduplicate its set of colour-vectors, then point the triangles' vertex-indices at the correct entry in the vertex-vector. It uses vsearchvec of PDL::VectorValued::Utils for this. Summarised from https://github.com/PDLPorters/pdl/blob/32a99b563358f7ad3ee26fba22fb24b5c02380e8/IO/STL/STL.pm#L150-L158:

use PDL::VectorValued::Utils;
# $pdl is ncoords=3,nvertices=3,ntriangles
my $uniqv = $pdl->uniqvec; # ncoords=3,nunique
my $face_indices = vsearchvec($pdl, $uniqv); # nvertices=3,ntriangles
[download]

From this, if you histogram the vertex-indices to get a count of each, you can then turn this into a pair of ndarrays ($coords, $count), and do

$cube->indexND($coords) .= $count;
[download]

This will work because then you are only writing into each cell once. As the "CAVEAT" in range of PDL::Slices notes:

It's quite possible to select multiple ranges that intersect. In that case, modifying the ranges doesn't have a guaranteed result in the original PDL -- the result is an arbitrary choice among the valid values. For some things that's OK; but for others it's not. In particular, this doesn't work:

pdl> $photon_list = new PDL::RandVar->sample(500)->reshape(2,250)*10
pdl> $histogram = zeroes(10,10)
pdl> $histogram->range($photon_list,1)++;  #not what you wanted
[download]

The reason is that if two photons land in the same bin, then that bin doesn't get incremented twice. (That may get fixed in a later version...)

Another alternative approach would be to try PDL::NDBin.