numpy的2D陣列（座標）需要分配到3D陣列，一些特定倉

Question

我有以下代碼來計算與8分箱和計數的預索貝爾X和Y濾波圖像HOG圖像：

for y in xrange(0, 480): 
    for x in xrange(0, 640): 
     base_angle = np.arctan2(sobel_Y[y,x], sobel_X[y,x]) * 180/np.pi 
     if base_angle < 0: base_angle += 360 
     angle = int(round(base_angle/45)) 
     if angle == 8: angle = 0 
     hog[y,x,angle] += np.sqrt(sobel_X[y,x]**2 + sobel_Y[y,x]**2) 

我試圖修改它以避免循環：

base_angle = np.arctan2(sobel_Y, sobel_X) * 180/np.pi 
base_angle[base_angle < 0] += 360 
angle =(base_angle/45).round().astype(np.uint8) 
angle[angle == bins] = 0 
hog[:,:,angle] += np.sqrt(sobel_X**2 + sobel_Y**2) 

但是，最後一個表達式不能正確計數。我基本需要的是根據角度數組中的索引，在豬陣列的每個（y，x）點處將幅度（np.sqrt ...表達式）添加到豬羣數組中。任何解決方案

Answer 1

使用

magnitude = np.sqrt(sobel_X**2 + sobel_Y**2) 
Y, X = np.ogrid[0:angle.shape[0], 0:angle.shape[1]] 
hog[Y, X, angle] += magnitude 

更新hog。

---

import numpy as np 

def using_for_loop(hog, sobel_Y, sobel_X): 
    for y in xrange(0, sobel_Y.shape[0]): 
     for x in xrange(0, sobel_X.shape[1]): 
      base_angle = np.arctan2(sobel_Y[y, x], sobel_X[y, x]) * 180/np.pi 
      if base_angle < 0: 
       base_angle += 360 
      angle = int(round(base_angle/45)) 
      if angle == 8: 
       angle = 0 
      hog[y, x, angle] += np.sqrt(sobel_X[y, x] ** 2 + 
             sobel_Y[y, x] ** 2) 
    return hog 

def using_indexing(hog, sobel_Y, sobel_X): 
    base_angle = np.arctan2(sobel_Y, sobel_X) * 180/np.pi 
    base_angle[base_angle < 0] += 360 
    angle = (base_angle/45).round().astype(np.uint8) 
    angle[angle == bins] = 0 
    magnitude = np.sqrt(sobel_X ** 2 + sobel_Y ** 2) 
    Y, X = np.ogrid[0:angle.shape[0], 0:angle.shape[1]] 
    hog[Y, X, angle] += magnitude 
    return hog 

bins = 8 
sobel_Y, sobel_X = np.meshgrid([1, 2, 3], [4, 5, 6, 7]) 
# hog = np.zeros(sobel_X.shape + (bins,)) 
hog = np.random.random(sobel_X.shape + (bins,)) 
answer = using_for_loop(hog, sobel_Y, sobel_X) 
result = using_indexing(hog, sobel_Y, sobel_X) 
assert np.allclose(answer, result) 

---

注意，如果

In [62]: angle.shape 
Out[62]: (4, 3) 

然後

In [74]: hog[:,:,angle].shape 
Out[74]: (4, 3, 4, 3) 

這不是正確的形狀。相反，如果你定義

In [75]: Y, X = np.ogrid[0:angle.shape[0], 0:angle.shape[1]] 

然後hog[Y, X, angle]具有相同的形狀magnitude：

In [76]: hog[Y, X, angle].shape 
Out[76]: (4, 3) 

In [77]: magnitude = np.sqrt(sobel_X ** 2 + sobel_Y ** 2) 

In [78]: magnitude.shape 
Out[78]: (4, 3) 

當然，這不是一個證明hog[Y, X, angle]是正確的表達，但它是一個相當於物理學家的dimensionality check，這表明我們至少可以在正確的軌道上。

---

隨着NumPy fancy indexing，當Y，X，和angle全部具有相同的形狀（或廣播到相同的形狀），

hog[Y, X, angle] 

也將是相同的形狀。在hog[Y, X, angle]的(i,j)個元素將等於

hog[Y[i,j], X[i,j], angle[i,j]] 

這就是爲什麼hog[Y, X, angle] += magnitude作品。