7

積分画像を計算するために、Repa を使用して累積合計関数を実装しようとしています。私の現在の実装は次のようになります。

cumsum :: (Elt a, Num a) => Array DIM2 a -> Array DIM2 a
cumsum array = traverse array id cumsum' 
    where
        elementSlice inner outer = slice array (inner :. (0 :: Int))
        cumsum' f (inner :. outer) = Repa.sumAll $ elementSlice inner outer

問題は elementSlice 関数にあります。matlab または numpy では、これは array[inner,0:outer] として指定できます。だから私が探しているのは、次のようなものです：

slice array (inner :. (Range 0 outer))

ただ、現在 Repa ではスライスを範囲指定することはできないようです。Haskell での効率的な並列ステンシル畳み込みで説明されているように、パーティションの使用を検討しましたが、反復ごとに変更される場合、これはかなり重量級のアプローチのように思えます。また、スライスをマスキングすることも検討しました (バイナリ ベクトルを乗算する) - しかし、行列内のすべての点にマスク ベクトルを割り当てるため、大きな行列ではパフォーマンスが非常に悪いように思えました...

私の質問 - Repa に範囲でのスライスのサポートを追加する計画があるかどうか誰か知っていますか? または、おそらく別のアプローチで、この問題にすでに対処できるパフォーマンスの高い方法はありますか?

4

2 に答える 2

5

サブレンジの抽出は、fromFunction で表現するのに十分簡単なインデックス スペース操作ですが、API にもっと適切なラッパーを追加する必要があります。

let arr = fromList (Z :. (5 :: Int)) [1, 2, 3, 4, 5 :: Int] 
in  fromFunction (Z :. 3) (\(Z :. ix) -> arr ! (Z :. ix + 1))

> [2,3,4]

結果の要素は、提供されたインデックスをオフセットし、それをソースから検索することによって取得されます。この手法は、より高いランクの配列に自然に拡張されます。

並列フォールドとスキャンの実装に関しては、そのためのプリミティブをライブラリに追加することでこれを行います。マップに関して並列リダクションを定義することはできませんが、遅延配列の全体的なアプローチを引き続き使用できます。それは合理的に直交する拡張です。

于 2011-05-31T06:42:04.763 に答える
3

実際、主な問題は、Repa にスキャン プリミティブがないことだと思います。(ただし、非常によく似たライブラリAccelerateが行います。) スキャンには、プレフィックス スキャンとサフィックス スキャンの 2 つのバリエーションがあります。与えられた 1D 配列

[a_1, ..., a_n]

プレフィックススキャンが返されます

[0, a_0, a_0 + a_1, ..., a_0 + ... + a_{n-1} ]

サフィックススキャンが生成する間

[a_0, a_0 + a_1, ..., a_0 + a_1 + ... + a_n ]

cumsumこれが、累積合計（ ）関数であなたが目指していることだと思います。

接頭辞と接尾辞のスキャンは、自然に多次元配列に一般化され、ツリー削減に基づいて効率的に実装されます。このトピックに関する比較的古い論文は、"Scan Primitives for Vector Computers"です。また、Conal Elliott は最近、Haskell での効率的な並列スキャンの導出に関するいくつかの ブログ 記事を書いています。

インテグラル イメージ (2D アレイ上) は、水平方向と垂直方向の 2 つのスキャンを実行することで計算できます。スキャン プリミティブがないため、非常に非効率的に実装しました。

horizScan :: Array DIM2 Int -> Array DIM2 Int
horizScan arr = foldl addIt arr [0 .. n - 1]
  where 
    addIt :: Array DIM2 Int -> Int -> Array DIM2 Int
    addIt accum i = accum +^ vs
       where 
         vs = toAdd (i+1) n (slice arr (Z:.All:.i))
    (Z:.m:.n) = arrayExtent arr

--
-- Given an @i@ and a length @len@ and a 1D array @arr@ 
-- (of length n) produces a 2D array of dimensions n X len.
-- The columns are filled with zeroes up to row @i@.
-- Subsequently they are filled with columns equal to the 
-- values in @arr.
--
toAdd :: Int -> Int -> Array DIM1 Int -> Array DIM2 Int
toAdd i len arr = traverse arr (\sh -> sh:.(len::Int)) 
               (\_ (Z:.n:.m) -> if m >= i then arr ! (Z:.n) else 0)

積分画像を計算する関数は、次のように定義できます。

vertScan :: Array DIM2 Int -> Array DIM2 Int
vertScan = transpose . horizScan . transpose

integralImage = horizScan . vertScan

Scan が Accelerate に実装されていることを考えると、それを Repa に追加するのはそれほど難しいことではありません。既存の Repa プリミティブを使用した効率的な実装が可能かどうかはわかりません。

于 2011-05-30T02:02:42.983 に答える