2010-05-19 22 views
9

Ich möchte den folgenden MATLAB-Code vektorisieren. Ich denke, es muss einfach sein, aber ich finde es trotzdem verwirrend.Vektorisieren von Summen verschiedener Diagonalen in einer Matrix

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 

Der Code berechnet eine Tabelle von Partituren, S, für einen dynamischen Programmier-Algorithmus, von einer anderen Auswertungstabelle, C.
Die Diagonalsummierung dient zur Erzeugung von Werten für einzelne Teile der Daten, die zur Erzeugung von C für alle möglichen Stücke (der Größe r) verwendet werden.

Vielen Dank im Voraus für alle Antworten! Sorry, wenn dies eine offensichtliche ...

Hinweis
Der eingebaute in CONV2 als convnfft erwies sich als schneller sein sollte, weil mein Auge (r) ist recht klein (5 < = r < = 20) . convnfft.m gibt an, dass r> 20 sein sollte, damit sich ein Vorteil manifestieren kann.

Antwort

10

Wenn ich richtig verstehe, versuchen Sie, die diagonale Summe jedes Unterarrays von C zu berechnen, wo Sie die letzte Zeile und Spalte von C entfernt haben (wenn Sie die Zeile/Spalte nicht entfernen sollten, müssen Sie eine Schleife erstellen) zu m-r + 1, und Sie müssen das gesamte Array C an die Funktion in meiner Lösung weiterleiten).

Sie können wie so über eine convolution, um diesen Vorgang tun:

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

Wenn C und r groß sind, können Sie einen Blick auf CONVNFFT aus dem Matlab Exchange haben Datei Berechnungen zu beschleunigen.

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

Grundsätzlich B enthält die Elemente aller Gleitschuhe:

+0

Super, danke. Ich werde morgen auf CONVNFFT schauen. Auf der Teilmenge der Daten, die ich zum Testen verwende (etwa 500 mal kleiner als echte Daten), erreicht die eingebaute conv2-Funktion eine 69,652-fache Reduktion der Anzahl der Aufrufe und 34,56-fache Verringerung der Ausführungszeit im Vergleich zu den Schleifen (23,5 vs 0,68 s). –

2

Ich würde denken, dass Sie möglicherweise C in eine 3D-Matrix neu anordnen müssen, bevor Sie es entlang einer der Dimensionen summieren. Ich werde in Kürze eine Antwort schreiben.

EDIT

Ich schaffe es nicht einen Weg zu finden, es sauber vektorisieren, aber ich habe die Funktion accumarray finden, die eine Hilfe sein könnte. Ich werde es genauer betrachten, wenn ich zu Hause bin.

EDIT # 2

Gefunden eine einfachere Lösung durch lineare Indizierung, aber dieser speicherintensiv sein könnte.

Bei C (1,1) sind die Indizes, die wir summieren wollen, 1+ [0, m + 1, 2 * m + 2, 3 * m + 3, 4 * m + 4, ...] oder (0: r-1) + (0: m: (r-1) * m)

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

erstellen S_offset, ein (mr) durch (nr) mit R-Matrix ist, so dass S_offset (: ,:, 1) = 0, S_offset (:,:, 2) = m + 1, S_offset (:,:, 3) = 2 * m + 2, und so weiter.

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

S_base schaffen, eine Matrix aus Basis-Array Adressen, von denen der Offset berechnet.

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

Schließlich verwenden S_base+S_offset die Werte von C. Adresse

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

Sie können natürlich bsxfun und andere Methoden verwenden, um sie effizienter zu gestalten; hier entschied ich mich für Klarheit. Ich muss dies noch vergleichen, um zu sehen, wie es mit der Doppelschleife-Methode vergleicht; Ich muss zuerst nach Hause zum Abendessen gehen!

+0

Hm, mit der Diagonale von einem 2D-Index an diesem Index die dritte Dimension zu werden? –

+1

@JS: Schöne Idee. 'im2col' könnte dir hier ein paar Zeilen Code sparen. – Jonas

+0

@Jonas: Ich habe eine Lösung hinzugefügt, die genau das tut! – Amro

3

Basierend auf der Idee der JS und als Jonas in den Kommentaren erwähnt, kann dies in zwei Zeilen IM2COL mit einigen Array-Manipulation unter Verwendung von durchgeführt werden der Größe r-by-r über die Matrix C. Dann nehmen wir die Elemente auf der Diagonale von jedem dieser Blöcke B(1:r+1:end,:), berechnen ihre Summe und formen das Ergebnis auf die erwartete Größe um.


Vergleicht man dies auf die Faltung-basierte Lösung von Jonas, dies führt keine Matrix-Multiplikation, nur Indizierung ...

+0

Wenn Sie die Erinnerung leisten können, ist im2col in der Regel ein guter Weg zu gehen. – Jonas

+0

im2col gibt Werte statt Indizes zurück, also sollte die zweite Zeile 'reshape (Summe (B (1: r + 1: Ende, :))), Größe (C) -r + 1);'. –

+0

@reve_etrange: du hast absolut recht, danke, dass du darauf hingewiesen hast – Amro

1

Ist das, was Sie suchen? Diese Funktion fügt die Diagonalen hinzu und fügt sie in einen Vektor ein, ähnlich wie die Funktion 'sum' alle Spalten einer Matrix aufsummiert und sie in einen Vektor setzt.

function [diagSum] = diagSumCalc(squareMatrix, LLUR0_ULLR1) 
% 
% Input: squareMatrix: A square matrix. 
%  LLUR0_ULLR1: LowerLeft to UpperRight addition = 0  
%      UpperLeft to LowerRight addition = 1 
% 
% Output: diagSum: A vector of the sum of the diagnols of the matrix. 
% 
% Example: 
% 
% >> squareMatrix = [1 2 3; 
%     4 5 6; 
%     7 8 9]; 
% 
% >> diagSum = diagSumCalc(squareMatrix, 0); 
% 
% diagSum = 
% 
%  1 6 15 14 9 
% 
% >> diagSum = diagSumCalc(squareMatrix, 1); 
% 
% diagSum = 
% 
%  7 12 15 8 3 
% 
% Written by M. Phillips 
% Oct. 16th, 2013 
% MIT Open Source Copywrite 
% Contact [email protected] fmi. 
% 

if (nargin < 2) 
    disp('Error on input. Needs two inputs.'); 
    return; 
end 

if (LLUR0_ULLR1 ~= 0 && LLUR0_ULLR1~= 1) 
    disp('Error on input. Only accepts 0 or 1 as input for second condition.'); 
    return; 
end 

[M, N] = size(squareMatrix); 

if (M ~= N) 
    disp('Error on input. Only accepts a square matrix as input.'); 
    return; 
end 

diagSum = zeros(1, M+N-1); 

if LLUR0_ULLR1 == 1 
    squareMatrix = rot90(squareMatrix, -1); 
end 

for i = 1:length(diagSum) 
    if i <= M 
     countUp = 1; 
     countDown = i; 
     while countDown ~= 0 
      diagSum(i) = squareMatrix(countUp, countDown) + diagSum(i); 
      countUp = countUp+1; 
      countDown = countDown-1; 
     end 
    end 
    if i > M 
     countUp = i-M+1; 
     countDown = M; 
     while countUp ~= M+1 
      diagSum(i) = squareMatrix(countUp, countDown) + diagSum(i); 
      countUp = countUp+1; 
      countDown = countDown-1; 
     end 
    end 
end 

Prost

Verwandte Themen