10

次の MATLAB コードをベクトル化したいと考えています。私はそれが単純でなければならないと思いますが、それでも混乱しています。

r = some constant less than m or n
[m,n] = size(C);
S = zeros(m-r,n-r);
for i=1:m-r+1
    for j=1:n-r+1
        S(i,j) = sum(diag(C(i:i+r-1,j:j+r-1)));
    end
end

このコードは、別のスコア テーブルCから動的計画法アルゴリズムのスコア テーブルSを計算します。対角合計は、 C を生成するために使用されるデータの個々の部分のスコアを生成することです。

ご回答ありがとうございます。これが明らかな場合は申し訳ありません...

私の eye(r) は非常に小さい ( 5 <= r <= 20 ) ため、組み込みのconv2
は convnfft よりも高速であることが判明しました。convnfft.m は、メリットが現れるには r > 20 である必要があると述べています。

4

4 に答える 4

11

私が正しく理解している場合、C の最後の行と列を削除した C のすべてのサブ配列の対角和を計算しようとしています (行/列を削除しない場合は、m-r にループする必要があります)。 +1、および配列 C 全体を以下のソリューションの関数に渡す必要があります)。

次のように、畳み込みを介してこの操作を実行できます。

S = conv2(C(1:end-1,1:end-1),eye(r),'valid');

C と r が大きい場合は、計算を高速化するために、Matlab File Exchange からCONVNFFTを参照することをお勧めします。

于 2010-05-19T11:58:06.940 に答える
3

JSの考え方に基づいており、Jonasがコメントで指摘したように、これはIM2COLいくつかの配列操作を使用して2 行で実行できます。

B = im2col(C, [r r], 'sliding');
S = reshape( sum(B(1:r+1:end,:)), size(C)-r+1 );

基本的Bに、サイズ r 行 r 列のすべてのスライディング ブロックの要素が含まれますC。次に、これらの各ブロックの対角線上の要素を取得し、B(1:r+1:end,:)それらの合計を計算して、結果を予想されるサイズに再形成します。

これをJonasによる畳み込みベースのソリューションと比較すると、これは行列の乗算を実行せず、インデックス付けのみを実行します...

于 2010-05-21T20:56:52.653 に答える
2

次元の 1 つに沿って合計する前に、C を 3D 行列に再配置する必要があると思います。私はすぐに答えを投稿します。

編集

私はそれをきれいにベクトル化する方法を見つけることができaccumarrayませんでしたが、いくつかの助けになるかもしれない関数を見つけました。家に帰ったら詳しく見てみます。

編集＃2

線形インデックスを使用してより簡単な解決策を見つけましたが、これはメモリを大量に消費する可能性があります。

C(1,1) で合計したいインデックスは 1+[0, m+1, 2*m+2, 3*m+3, 4*m+4, ... ] または (0 :r-1)+(0:m:(r-1)*m)

sum_ind = (0:r-1)+(0:m:(r-1)*m);

S_offset S_offset(:,:,1) = 0、S_offset(:,:,2) = m+1、S_offset(:,:,3) = 2 となる (mr) by (nr) by r の行列を作成します。 *m+2 など。

S_offset = permute(repmat( sum_ind, [m-r, 1, n-r] ), [1, 3, 2]);

create S_base、オフセットが計算されるベース配列アドレスの行列。

S_base = reshape(1:m*n,[m n]);
S_base = repmat(S_base(1:m-r,1:n-r), [1, 1, r]);

最後に、S_base+S_offsetC の値に対処するために使用します。

S = sum(C(S_base+S_offset), 3);

もちろん、bsxfunその他の方法を使用してより効率的にすることもできます。ここでは、わかりやすくするためにレイアウトすることにしました。ただし、これをベンチマークして、ダブルループ方式と比較する方法を確認する必要があります。まず家に帰って夕食を食べなきゃ！

于 2010-05-19T09:08:05.717 に答える