riquadro di delimitazione dell'array numpy

Question

Supponiamo di disporre di un array numpy 2D con alcuni valori casuali e zeri circostanti.

Esempio "rettangolo inclinato":

import numpy as np 
from skimage import transform 

img1 = np.zeros((100,100)) 
img1[25:75,25:75] = 1. 
img2 = transform.rotate(img1, 45) 

Ora vorrei trovare il più piccolo rettangolo di delimitazione per tutti i dati diversi da zero. Per esempio:

a = np.where(img2 != 0) 
bbox = img2[np.min(a[0]):np.max(a[0])+1, np.min(a[1]):np.max(a[1])+1] 

Quale sarebbe il modo più veloce per raggiungere questo risultato? Sono sicuro che c'è un modo migliore poiché la funzione np.where richiede un bel po 'di tempo se io sono, per esempio, utilizzando set di dati 1000x1000.

Edit: dovrebbe funzionare anche in 3D ...

Answer 1

È possibile o meno dimezzare il tempo di esecuzione utilizzando np.any per ridurre le righe e le colonne che contengono valori diversi da zero ai vettori 1D, piuttosto che trovare gli indici di tutti i valori diversi da zero utilizzando np.where:

def bbox1(img): 
    a = np.where(img != 0) 
    bbox = np.min(a[0]), np.max(a[0]), np.min(a[1]), np.max(a[1]) 
    return bbox 

def bbox2(img): 
    rows = np.any(img, axis=1) 
    cols = np.any(img, axis=0) 
    rmin, rmax = np.where(rows)[0][[0, -1]] 
    cmin, cmax = np.where(cols)[0][[0, -1]] 

    return rmin, rmax, cmin, cmax 

alcuni parametri di riferimento:

%timeit bbox1(img2) 
10000 loops, best of 3: 63.5 µs per loop 

%timeit bbox2(img2) 
10000 loops, best of 3: 37.1 µs per loop 

estendere questo approccio al caso 3D solo comporta l'esecuzione della riduzione lungo ciascuna coppia di assi:

È facile generalizzare questo N dimensioni utilizzando itertools.combinations per scorrere ciascuna combinazione unica di assi per eseguire la riduzione su:

import itertools 

def bbox2_ND(img): 
    N = img.ndim 
    out = [] 
    for ax in itertools.combinations(range(N), N - 1): 
     nonzero = np.any(img, axis=ax) 
     out.extend(np.where(nonzero)[0][[0, -1]]) 
    return tuple(out) 

---

Se si conoscono le coordinate degli angoli del riquadro di delimitazione originale, l'angolo di rotazione e il centro di rotazione, è possibile ottenere le coordinate degli angoli del rettangolo di selezione trasformati direttamente calcolando la corrispondente affine transformation matrix e punteggiano con l'ingresso coordinate:

def bbox_rotate(bbox_in, angle, centre): 

    rmin, rmax, cmin, cmax = bbox_in 

    # bounding box corners in homogeneous coordinates 
    xyz_in = np.array(([[cmin, cmin, cmax, cmax], 
         [rmin, rmax, rmin, rmax], 
         [ 1, 1, 1, 1]])) 

    # translate centre to origin 
    cr, cc = centre 
    cent2ori = np.eye(3) 
    cent2ori[:2, 2] = -cr, -cc 

    # rotate about the origin 
    theta = np.deg2rad(angle) 
    rmat = np.eye(3) 
    rmat[:2, :2] = np.array([[ np.cos(theta),-np.sin(theta)], 
          [ np.sin(theta), np.cos(theta)]]) 

    # translate from origin back to centre 
    ori2cent = np.eye(3) 
    ori2cent[:2, 2] = cr, cc 

    # combine transformations (rightmost matrix is applied first) 
    xyz_out = ori2cent.dot(rmat).dot(cent2ori).dot(xyz_in) 

    r, c = xyz_out[:2] 

    rmin = int(r.min()) 
    rmax = int(r.max()) 
    cmin = int(c.min()) 
    cmax = int(c.max()) 

    return rmin, rmax, cmin, cmax 

Questo funziona a essere leggermente più veloce dell'utilizzo np.any per le piccole esempio matrice:

%timeit bbox_rotate([25, 75, 25, 75], 45, (50, 50)) 
10000 loops, best of 3: 33 µs per loop 

Tuttavia, dato che la velocità di questo metodo è indipendente dalla dimensione della matrice di input, può essere molto più veloce per gli array più grandi.

L'estensione dell'approccio alla trasformazione in 3D è leggermente più complicato, in quanto la rotazione ora ha tre componenti diversi (uno sull'asse x, uno sull'asse y e uno sull'asse z), ma l'elemento base metodo è lo stesso:

def bbox_rotate_3d(bbox_in, angle_x, angle_y, angle_z, centre): 

    rmin, rmax, cmin, cmax, zmin, zmax = bbox_in 

    # bounding box corners in homogeneous coordinates 
    xyzu_in = np.array(([[cmin, cmin, cmin, cmin, cmax, cmax, cmax, cmax], 
         [rmin, rmin, rmax, rmax, rmin, rmin, rmax, rmax], 
         [zmin, zmax, zmin, zmax, zmin, zmax, zmin, zmax], 
         [ 1, 1, 1, 1, 1, 1, 1, 1]])) 

    # translate centre to origin 
    cr, cc, cz = centre 
    cent2ori = np.eye(4) 
    cent2ori[:3, 3] = -cr, -cc -cz 

    # rotation about the x-axis 
    theta = np.deg2rad(angle_x) 
    rmat_x = np.eye(4) 
    rmat_x[1:3, 1:3] = np.array([[ np.cos(theta),-np.sin(theta)], 
           [ np.sin(theta), np.cos(theta)]]) 

    # rotation about the y-axis 
    theta = np.deg2rad(angle_y) 
    rmat_y = np.eye(4) 
    rmat_y[[0, 0, 2, 2], [0, 2, 0, 2]] = (
     np.cos(theta), np.sin(theta), -np.sin(theta), np.cos(theta)) 

    # rotation about the z-axis 
    theta = np.deg2rad(angle_z) 
    rmat_z = np.eye(4) 
    rmat_z[:2, :2] = np.array([[ np.cos(theta),-np.sin(theta)], 
           [ np.sin(theta), np.cos(theta)]]) 

    # translate from origin back to centre 
    ori2cent = np.eye(4) 
    ori2cent[:3, 3] = cr, cc, cz 

    # combine transformations (rightmost matrix is applied first) 
    tform = ori2cent.dot(rmat_z).dot(rmat_y).dot(rmat_x).dot(cent2ori) 
    xyzu_out = tform.dot(xyzu_in) 

    r, c, z = xyzu_out[:3] 

    rmin = int(r.min()) 
    rmax = int(r.max()) 
    cmin = int(c.min()) 
    cmax = int(c.max()) 
    zmin = int(z.min()) 
    zmax = int(z.max()) 

    return rmin, rmax, cmin, cmax, zmin, zmax 

sono essenzialmente solo modificato la funzione sopra utilizzando le espressioni matrice di rotazione da here - non ho avuto tempo di scrivere un test-case ancora, in modo da utilizzare con cautela.

Answer 2

Ecco un algoritmo per calcolare il riquadro per N matrici tridimensionali,

def get_bounding_box(x): 
    """ Calculates the bounding box of a ndarray""" 
    mask = x == 0 
    bbox = [] 
    all_axis = np.arange(x.ndim) 
    for kdim in all_axis: 
     nk_dim = np.delete(all_axis, kdim) 
     mask_i = mask.all(axis=tuple(nk_dim)) 
     dmask_i = np.diff(mask_i) 
     idx_i = np.nonzero(dmask_i)[0] 
     if len(idx_i) != 2: 
      raise ValueError('Algorithm failed, {} does not have 2 elements!'.format(idx_i)) 
     bbox.append(slice(idx_i[0]+1, idx_i[1]+1)) 
    return bbox 

che può essere utilizzato con 2D, 3D, ecc array come segue,

In [1]: print((img2!=0).astype(int)) 
    ...: bbox = get_bounding_box(img2) 
    ...: print((img2[bbox]!=0).astype(int)) 
    ...: 
[[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 1 1 1 1 1 1 1 1 0 0 0 0 0 0] 
[0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0] 
[0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0] 
[0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0] 
[0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0] 
[0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0] 
[0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0] 
[0 0 0 0 0 0 1 1 1 1 1 1 1 1 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] 
[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]] 
[[0 0 0 0 0 0 1 1 0 0 0 0 0 0] 
[0 0 0 0 0 1 1 1 1 0 0 0 0 0] 
[0 0 0 0 1 1 1 1 1 1 0 0 0 0] 
[0 0 0 1 1 1 1 1 1 1 1 0 0 0] 
[0 0 1 1 1 1 1 1 1 1 1 1 0 0] 
[0 1 1 1 1 1 1 1 1 1 1 1 1 0] 
[1 1 1 1 1 1 1 1 1 1 1 1 1 1] 
[1 1 1 1 1 1 1 1 1 1 1 1 1 1] 
[0 1 1 1 1 1 1 1 1 1 1 1 1 0] 
[0 0 1 1 1 1 1 1 1 1 1 1 0 0] 
[0 0 0 1 1 1 1 1 1 1 1 0 0 0] 
[0 0 0 0 1 1 1 1 1 1 0 0 0 0] 
[0 0 0 0 0 1 1 1 1 0 0 0 0 0] 
[0 0 0 0 0 0 1 1 0 0 0 0 0 0]] 

Sebbene sostituendo il np.diff e le chiamate np.nonzero da uno np.where potrebbero essere migliori.

Answer 3

Sono stato in grado di spremere un po 'più di prestazioni sostituendo np.where con np.argmax e lavorando su una maschera booleana.

 
def bbox(img): 
    img = (img > 0) 
    rows = np.any(img, axis=1) 
    cols = np.any(img, axis=0) 
    rmin, rmax = np.argmax(rows), img.shape[0] - 1 - np.argmax(np.flipud(rows)) 
    cmin, cmax = np.argmax(cols), img.shape[1] - 1 - np.argmax(np.flipud(cols)) 
    return rmin, rmax, cmin, cmax

Per me era più veloce di 10μs rispetto alla soluzione bbox2 sopra lo stesso benchmark. Ci dovrebbe anche essere un modo per usare semplicemente il risultato di argmax per trovare le righe e le colonne diverse da zero, evitando la ricerca extra fatta usando np.any, ma questo potrebbe richiedere qualche indicizzazione difficile che non ero in grado di lavorare in modo efficiente con codice vettoriale semplice