用repa计算图像直方图

Question

我正在从磁盘加载RGB图像JuicyPixels-repa.不幸的是,图像的阵列表示是Array F DIM3 Word8内部维度是RGB像素的位置.这与repaRGB图像的现有图像处理算法有点不相容Array U DIM2 (Word8, Word8, Word8).

我想计算图像的RGB直方图,我正在使用签名搜索函数:

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

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

编辑:

主要的问题不在于我无法从转换DIM3到DIM2每通道(容易切片完成).问题是我必须迭代源图像,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没有我认为的正常traverse功能.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,但是如果你发现它运行良好,请告诉我.