Eine Beschleunigung Berechnungen mit numpy Matrizen

Question

Ich habe zwei Matrizen. Beide sind mit Nullen und Einsen gefüllt. Eine davon ist ein großer (3000 x 2000 Elemente), und der andere ist kleiner (20 x 20) Elemente. Ich tue so etwas wie:

newMatrix = (size of bigMatrix), filled with zeros
l = (a constant)

for y in xrange(0, len(bigMatrix[0])):
    for x in xrange(0, len(bigMatrix)):

        for b in xrange(0, len(smallMatrix[0])):
            for a in xrange(0, len(smallMatrix)):

                if (bigMatrix[x, y] == smallMatrix[x + a - l, y + b - l]):
                    newMatrix[x, y] = 1

Was wird schmerzlich langsam. Mache ich etwas falsch? Gibt es eine intelligente Art und Weise, diese Arbeit schneller zu machen?

edit: Grundsätzlich bin ich für jede (x, y) in der großen Matrix, alle Pixel der beiden großen Matrix und die kleine Matrix Überprüfung um (x, y) zu sehen, ob sie 1 sind, wenn sie 1 sind , dann setze ich diesen Wert auf newMatrix. Ich bin eine Art Kollisionserkennung zu tun.

Answer 1

Ich kann es von ein paar Optimierungen denken - Wie Sie 4 verschachtelten Python „für“ Aussagen verwenden, sind Sie so langsam wie Sie sein können.

Ich kann nicht genau herauszufinden, was Sie suchen - aber für eine Sache, wenn Ihre große Matrix „1“ Dichte niedrig ist, können Sie sicher Pythons „any“ Funktion auf bigMtarix die Scheiben verwenden, um schnell zu überprüfen, ob es irgendwelche Satz Elemente gibt es - Sie haben eine mehrfache Geschwindigkeitserhöhung bekommen könnte es:

step = len(smallMatrix[0])
for y in xrange(0, len(bigMatrix[0], step)):
    for x in xrange(0, len(bigMatrix), step):
        if not any(bigMatrix[x: x+step, y: y + step]):
            continue
        (...) 

An diesem Punkt, wenn nach wie vor auf jedes Element interagieren müssen, was Sie tun, ein anderes Paar von Indizes jede Position innerhalb des Schrittes zu gehen -. Aber ich denke, Sie auf die Idee gekommen

Neben der Verwendung Innen numerische Operationen wie diese „jede“ usage, könnte man sicherlich eine gewisse Kontrolle Flow Code fügen Sie die (b, a) Schleife, wenn die erste Anpassungs Pixel gefunden wird, Abrissen. (Wie, das Einfügen einer „Pause“ Anweisung in Ihrer letzten „if“ und ein anderes if..break Paar für die „b“ Schleife.

Ich kann wirklich nicht herausfinden, was Ihre Absicht ist - so kann ich Ihnen mehr specifc Code nicht geben.

Answer 2

Ihr Beispielcode macht keinen Sinn, aber die Beschreibung Ihres Problems klingt wie Sie eine 2D-Faltung eines kleinen BitArray über den großen BitArray zu tun versuchen. Es gibt eine convolve2d Funktion in scipy.signal Paket, das dies tut genau. Just do convolve2d(bigMatrix, smallMatrix) das Ergebnis zu erhalten. Leider ist die scipy Implementierung keinen Sonderfall für boolean-Arrays, so dass die volle Faltung eher langsam ist. Hier ist eine Funktion, die sich die Tatsache zunutze kommt, dass die Arrays enthalten nur Einsen und Nullen:

import numpy as np

def sparse_convolve_of_bools(a, b):
    if a.size < b.size:
        a, b = b, a
    offsets = zip(*np.nonzero(b))
    n = len(offsets)
    dtype = np.byte if n < 128 else np.short if n < 32768 else np.int
    result = np.zeros(np.array(a.shape) + b.shape - (1,1), dtype=dtype)
    for o in offsets:
        result[o[0]:o[0] + a.shape[0], o[1]:o[1] + a.shape[1]] += a
    return result

Auf meinem Rechner läuft es in weniger als 9 Sekunden für eine 3000x2000 von 20x20 Faltung. Die Laufzeit ist abhängig von der Anzahl von Einsen in dem kleineren Array, 20ms für jedes Nicht-Null-Element zu sein.

Answer 3

Wenn Bits verpackt wirklich 8 pro Byte / 32 pro int, und Sie können Ihre smallMatrix zu 20x16, reduzieren
versuchen Sie dann die folgende, hier für eine einzelne Zeile.
(newMatrix[x, y] = 1 wenn jeder Bit des 20x16 um x, y 1 ?? Was suchen Sie wirklich?)

python -m timeit -s '
""" slide 16-bit mask across 32-bit pairs bits[j], bits[j+1] """

import numpy as np

bits = np.zeros( 2000 // 16, np.uint16 )  # 2000 bits
bits[::8] = 1
mask = 32+16
nhit = 16 * [0]

def hit16( bits, mask, nhit ):
    """
        slide 16-bit mask across 32-bit pairs bits[j], bits[j+1]
        bits: long np.array( uint16 )
        mask: 16 bits, int
        out: nhit[j] += 1 where pair & mask != 0
    """
    left = bits[0]
    for b in bits[1:]:
        pair = (left << 16) | b
        if pair:  # np idiom for non-0 words ?
            m = mask
            for j in range(16):
                if pair & m:
                    nhit[j] += 1
                    # hitposition = jb*16 + j
                m <<= 1
        left = b
    # if any(nhit):  print "hit16:", nhit

' \
'
hit16( bits, mask, nhit )
'

# 15 msec per loop, bits[::4] = 1
# 11 msec per loop, bits[::8] = 1
# mac g4 ppc