用repa计算图像直方图。

时间:2022-04-02 18:01:39

I'm loading an RGB image from disk with JuicyPixels-repa. Unfortunately the Array representation of the image is Array F DIM3 Word8 where the inner dimension is the RGB pixels. That's a bit incompatibe with existing repa imageprocessing algorithms where an RGB image is Array U DIM2 (Word8, Word8, Word8).

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

I want to calculate the RGB histograms of the image, I'm searching a function with the Signature:

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

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

how can I fold my 3d array to get a 1d array for each colorchannel?

我如何折叠我的3d数组以得到每个colorchannel的1d阵列?

Edit:

编辑:

The main problem is not that I'm not able to convert from DIM3 to DIM2 for each channel (easy done with slicing). The problem is that I have to iterate the source image DIM2 or DIM3 and have to accumulate to an DIM1 array of a different Shape (Z:.256) and extent. So I can't use repa's foldS as it reduces the dimension by one, but with the same extent.

主要的问题不是我不能将每个通道的DIM3转换为DIM2(简单的切片)。问题是,我必须迭代源图像DIM2或DIM3,并且必须累积到一个不同形状的DIM1数组(Z:.256)和范围。所以我不能使用repa的折叠,因为它减少了一个维度,但同样的程度。

I also experimented with traverse but it iterates over the extent of the destination image, providing a function to get pixels from the source image, that would lead to very inefficient code, counting the same pixels for each colorvalue.

我还尝试了遍历,但它遍历目标图像的范围,提供一个从源图像获取像素的函数,这会导致非常低效的代码,计算每个颜色值的相同像素。

A good way would be a simple folding over a Vector with the histogram type as accumulator, but unfortunately I have no U (unboxed) or V (vector) based array, from which I can efficiently get a Vector. I have an Array F (foreign pointer).

一个好的方法是用直方图类型作为累加器来简单地折叠一个矢量,但不幸的是,我没有U (unboxed)或V(矢量)基数组,从中我可以有效地得到一个向量。我有一个数组F(外指针)。

1 个解决方案

#1


7  

Ok, I found a few minutes. Below, I cover four solutions and have made the worst solutions (the middle two, involving O(n) data conversion) really easy for you.

好的,我找到了几分钟。下面,我介绍了四个解决方案,并做出了最糟糕的解决方案(中间两个,涉及O(n)数据转换)对您来说非常简单。

Lets Acknowledge the Dumb Solution

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

It's reasonable to start with the obvious. You could use Data.List.foldl to traverse the rows and columns, building up your histograms from initial zero arrays (untested/partial code follows):

从最明显的开始是合理的。您可以使用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

I'm not sure how efficient this will be, but suspect that it will do too much bounds checking. Switching to unsafeIndex should do reasonably, assuming the delayed arrays, hist*, do well due to however you'd pick to implement incElem.

我不确定这会有多有效,但怀疑它会做太多的边界检查。切换到unsafeIndex应该是合理的,假设延迟的数组,hist*,由于您选择执行incElem而做得很好。

You Can Build the Array You Want

您可以构建您想要的数组。

Using traverse you can actually convert JP-Repa style arrays into DIM2 arrays with tuples for elements:

使用遍历,您实际上可以将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))

Could you point me to the body of code you talked about that uses this format? It would be simple to patch JP-Repa to include a function of this type.

你能给我讲讲你说过的使用这种格式的代码吗?将JP-Repa添加到该类型的函数中是很简单的。

You can build the Unboxed Vector You Mentioned

您可以构建您提到的Unboxed向量。

You mentioned an easy solution is to fold over unboxed vectors, but lamented that JP-repa doesn't provide an unboxed array. Luckily, conversion is simple:

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

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

We Could Patch Repa

我们可以修补Repa

This is really only a problem because Repa doesn't have what I consider a normal traverse function. Repa's traverse is more of an array construction that happens to provide an indexing function into another array. We want traverse in the form:

这只是一个问题,因为Repa没有我认为的正常的遍历函数。Repa的遍历更多的是一个数组构造,它恰好为另一个数组提供了索引函数。我们想要的形式是:

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

but of coarse this is actually just a malformed fold. So lets rename it and reorder the arguments:

但粗糙的这只是一个畸形的褶皱。因此,让我们重命名它并重新排序参数:

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

which contrasts nicely with the (preexisting) foldAllS operation:

这与(现有的)foldAllS操作形成了鲜明的对比:

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

Notice how our new fold has two critical characteristics. The result type is not required to match the element type, so we could start with a tuple of Histograms. Second, our version of fold passes the index, which allows you to select which histogram in the tuple to update (if any).

注意我们的新折叠有两个关键特征。结果类型不需要与元素类型匹配,所以我们可以从一个直方图元组开始。其次,我们的折叠版本通过索引,它允许您在tuple中选择要更新的直方图(如果有的话)。

You can lazily use the latest JuicyPixels-Repa

您可以延迟使用最新的JuicyPixels-Repa。

To acquire your preferred Repa array format, or to acquire an unboxed vector, you can just use the newly uploaded JuicyPixels-Repa-0.6.

为了获得您喜欢的Repa数组格式,或者获得一个unboxed矢量,您可以使用刚刚上传的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

Now you can fold over the vector or use the tuple array directly, as you originally desired.

现在,您可以像您最初希望的那样,直接折叠这个向量或使用tuple数组。

I also took an ugly stab at fleshing out the first, horribly naive, solution:

我也做了一个丑陋的尝试,把第一个,可怕的,天真的,解决方案:

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] ]

I'm too wary of the performance of this code (3 traversals per index... I must be tired) to throw it into JP-Repa, but if you find it works well then let me know.

我对这个代码的性能太过谨慎了(每个索引有3个遍历……)我一定是累了)把它扔进了JP-Repa,但是如果你发现它很好用,那就告诉我。

#1


7  

Ok, I found a few minutes. Below, I cover four solutions and have made the worst solutions (the middle two, involving O(n) data conversion) really easy for you.

好的,我找到了几分钟。下面,我介绍了四个解决方案,并做出了最糟糕的解决方案(中间两个,涉及O(n)数据转换)对您来说非常简单。

Lets Acknowledge the Dumb Solution

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

It's reasonable to start with the obvious. You could use Data.List.foldl to traverse the rows and columns, building up your histograms from initial zero arrays (untested/partial code follows):

从最明显的开始是合理的。您可以使用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

I'm not sure how efficient this will be, but suspect that it will do too much bounds checking. Switching to unsafeIndex should do reasonably, assuming the delayed arrays, hist*, do well due to however you'd pick to implement incElem.

我不确定这会有多有效,但怀疑它会做太多的边界检查。切换到unsafeIndex应该是合理的,假设延迟的数组,hist*,由于您选择执行incElem而做得很好。

You Can Build the Array You Want

您可以构建您想要的数组。

Using traverse you can actually convert JP-Repa style arrays into DIM2 arrays with tuples for elements:

使用遍历,您实际上可以将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))

Could you point me to the body of code you talked about that uses this format? It would be simple to patch JP-Repa to include a function of this type.

你能给我讲讲你说过的使用这种格式的代码吗?将JP-Repa添加到该类型的函数中是很简单的。

You can build the Unboxed Vector You Mentioned

您可以构建您提到的Unboxed向量。

You mentioned an easy solution is to fold over unboxed vectors, but lamented that JP-repa doesn't provide an unboxed array. Luckily, conversion is simple:

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

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

We Could Patch Repa

我们可以修补Repa

This is really only a problem because Repa doesn't have what I consider a normal traverse function. Repa's traverse is more of an array construction that happens to provide an indexing function into another array. We want traverse in the form:

这只是一个问题,因为Repa没有我认为的正常的遍历函数。Repa的遍历更多的是一个数组构造,它恰好为另一个数组提供了索引函数。我们想要的形式是:

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

but of coarse this is actually just a malformed fold. So lets rename it and reorder the arguments:

但粗糙的这只是一个畸形的褶皱。因此,让我们重命名它并重新排序参数:

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

which contrasts nicely with the (preexisting) foldAllS operation:

这与(现有的)foldAllS操作形成了鲜明的对比:

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

Notice how our new fold has two critical characteristics. The result type is not required to match the element type, so we could start with a tuple of Histograms. Second, our version of fold passes the index, which allows you to select which histogram in the tuple to update (if any).

注意我们的新折叠有两个关键特征。结果类型不需要与元素类型匹配,所以我们可以从一个直方图元组开始。其次,我们的折叠版本通过索引,它允许您在tuple中选择要更新的直方图(如果有的话)。

You can lazily use the latest JuicyPixels-Repa

您可以延迟使用最新的JuicyPixels-Repa。

To acquire your preferred Repa array format, or to acquire an unboxed vector, you can just use the newly uploaded JuicyPixels-Repa-0.6.

为了获得您喜欢的Repa数组格式,或者获得一个unboxed矢量,您可以使用刚刚上传的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

Now you can fold over the vector or use the tuple array directly, as you originally desired.

现在,您可以像您最初希望的那样,直接折叠这个向量或使用tuple数组。

I also took an ugly stab at fleshing out the first, horribly naive, solution:

我也做了一个丑陋的尝试,把第一个,可怕的,天真的,解决方案:

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] ]

I'm too wary of the performance of this code (3 traversals per index... I must be tired) to throw it into JP-Repa, but if you find it works well then let me know.

我对这个代码的性能太过谨慎了(每个索引有3个遍历……)我一定是累了)把它扔进了JP-Repa,但是如果你发现它很好用,那就告诉我。