haskell - 使用 repa 计算图像直方图

Question

我正在从磁盘加载 RGB 图像 JuicyPixels-repa 。不幸的是，图像的数组表示是 Array F DIM3 Word8其中内部维度是 RGB 像素。这与现有的repa有点不兼容。 RGB 图像为 Array U DIM2 (Word8, Word8, Word8) 的图像处理算法.

我想计算图像的 RGB 直方图，我正在搜索带有签名的函数:

type Hist = Array U DIM1 Int
histogram:: Array F DIM3  Word8 -> (Hist, Hist, Hist)

如何折叠 3d 数组以获得每个颜色 channel 的 1d 数组？

编辑:

主要问题不是我无法从 DIM3 进行转换至DIM2对于每个 channel (通过切片很容易完成)。问题是我必须迭代源图像 DIM2或DIM3并且必须累积到 DIM1不同的数组Shape (Z:.256) 和范围。所以我不能使用repa的foldS因为它将维度减少了一倍，但程度相同。

我还尝试过traverse但它会迭代目标图像的范围，提供从源图像获取像素的函数，这会导致代码效率非常低，为每个颜色值计算相同的像素。

一个好方法是简单折叠 Vector使用直方图类型作为累加器，但不幸的是我没有 U (未装箱)或V基于(向量)的数组，我可以从中有效地获得 Vector 。我有一个Array F (外部指针)。

Answer 1

好的，我找到了几分钟。下面，我介绍了四种解决方案，并使最差的解决方案(中间两个，涉及 O(n) 数据转换)对您来说非常容易。

让我们承认愚蠢的解决方案

从显而易见的事情开始是合理的。您可以使用 Data.List.foldl 遍历行和列，从初始零数组构建直方图(以下是未经测试/部分代码):

foldl (\(histR, histG, histB) (row,col) ->
            let r = arr ! (Z:.row:.col:.0)
                g = arr ! (Z:.row:.col:.1)
                b = arr ! (Z:.row:.col:.2)
            in (incElem r histR, incElem g histG, incElem b histB)
         (zero,zero,zero)
         [ (row,col) | row <- [0..nrRow-1], col <- [0..nrCol-1] ]
 ...
 where (Z:.nrRow:.nrCol:._) = extent arr

我不确定这有多高效，但怀疑它会进行太多的边界检查。切换到 unsafeIndex 应该是合理的，假设延迟数组 hist* 表现良好，因为您选择实现 incElem。

您可以构建您想要的数组

使用traverse，您实际上可以将 JP-Repa 样式数组转换为带有元素元组的 DIM2 数组:

main = do
    let arr = R.fromFunction (Z:.a:.b:.c) (\(Z:.i:.j:.k) -> i+j-k)
        a =4 :: Int
        b = 4 :: Int
        c= 4 :: Int
        new = R.traverse arr
                       (\(Z:.r:.c:._) -> Z:.r:.c) -- the extent
                       (\l idx -> (l (idx:.0)
                                  ,l (idx:.1)
                                  ,l (idx :. 2)))
    print (R.computeS new :: R.Array R.U DIM2 (Int,Int,Int))

您能否指出您提到的使用这种格式的代码主体？修补 JP-Repa 以包含此类功能会很简单。

您可以构建您提到的未装箱矢量

您提到一个简单的解决方案是折叠未装箱的向量，但遗憾的是 JP-repa 不提供未装箱的数组。幸运的是，转换很简单:

toUnboxed :: Img a -> VU.Vector Word8
toUnboxed = R.toUnboxed . R.computeUnboxedS . R.delay . imgData

我们可以修补 Repa

这实际上只是一个问题，因为 Repa 没有我认为的正常遍历功能。 Repa 的遍历更多的是一种数组构造，它恰好为另一个数组提供索引功能。我们希望以以下形式进行遍历:

newTraverse :: Array r sh e -> a -> (a -> sh -> e -> a) -> a

但粗略地说，这实际上只是一个畸形的折叠。因此，让我们重命名它并重新排序参数:

foldAllIdxS :: (sh -> a - > e -> a) -> a -> Array r sh e -> a

这与(预先存在的)foldAllS 操作形成鲜明对比:

foldAllS :: (a -> a -> a) -> a -> Array r sh a -> a

请注意我们的新折叠有两个关键特征。结果类型不需要与元素类型匹配，因此我们可以从直方图元组开始。其次，我们的折叠版本传递索引，它允许您选择要更新元组中的哪个直方图(如果有)。

您可以懒惰地使用最新的JuicyPixels-Repa

要获取您喜欢的 Repa 数组格式，或者获取未装箱的向量，您可以使用新上传的 JuicyPixels-Repa-0.6。

someImg <- readImage path :: IO (Either String (Img RGBA))
let img = either (error "Blah") id someImg
    uvec = toUnboxed img
    tupleArr = collapseColorChannel img

现在您可以按照您最初的期望折叠向量或直接使用元组数组。

我还尝试了一个丑陋的尝试来充实第一个非常天真的解决方案:

histograms :: Img a -> (Histogram, Histogram, Histogram, Histogram)
histograms (Img arr) =
    let (Z:.nrRow:.nrCol:._) = R.extent arr
        zero = R.fromFunction (Z:.256) (\_ -> 0 :: Word8)
        incElem idx x = RU.unsafeTraverse x id (\l i -> l i + if i==(Z:.fromIntegral idx) then 1 else 0)
    in Prelude.foldl (\(hR, hG, hB, hA) (row,col) ->
             let r = R.unsafeIndex arr (Z:.row:.col:.0)
                 g = R.unsafeIndex arr (Z:.row:.col:.1)
                 b = R.unsafeIndex arr (Z:.row:.col:.2)
                 a = R.unsafeIndex arr (Z:.row:.col:.3)
             in (incElem r hR, incElem g hG, incElem b hB, incElem a hA))
          (zero,zero,zero,zero)
          [ (row,col) | row <- [0..nrRow-1], col <- [0..nrCol-1] ]

我对这段代码的性能过于警惕(每个索引 3 次遍历......我一定很累)，无法将其放入 JP-Repa，但如果您发现它运行良好，请告诉我。