4

これは、numpyでストライドを使用する最初の試みであり、さまざまなフィルターでの単純な反復と比較して速度が向上しましたが、それでもかなり遅いです(完全に冗長または非効率的なものが少なくとも1つまたは2つあるように感じます) 。

だから私の質問は:これを実行するより良い方法や、コードを大幅に高速化するための微調整はありますか?

アルゴリズムは、ピクセルごとに9つの異なるフィルターのローカル評価を実行し、標準偏差が最も小さいフィルターを選択します(画像分析で説明されているように、Nagau and Matsuyma(1980)「複雑な領域写真の構造分析」を実装しようとしました)。本)。その結果、画像は滑らかになり、エッジがシャープになります(私に言わせればかなりクールです!)

import numpy as np
from scipy import ndimage
from numpy.lib import stride_tricks

def get_rotating_kernels():

    kernels = list()

    protokernel = np.arange(9).reshape(3,  3)

    for k in xrange(9):

        ax1, ax2 = np.where(protokernel==k)
        kernel = np.zeros((5,5), dtype=bool)
        kernel[ax1: ax1+3, ax2: ax2+3] = 1
        kernels.append(kernel)

    return kernels


def get_rotation_smooth(im, **kwargs):

    kernels = np.array([k.ravel() for k in get_rotating_kernels()],
                dtype=bool)

    def rotation_matrix(section):

        multi_s = stride_tricks.as_strided(section, shape=(9,25),
            strides=(0, section.itemsize))

        rot_filters = multi_s[kernels].reshape(9,9)

        return rot_filters[rot_filters.std(1).argmin(),:].mean()

    return ndimage.filters.generic_filter(im, rotation_matrix, size=5, **kwargs)

from scipy import lena
im = lena()
im2 = get_rotation_smooth(im)

(コメントget_rotating_kernelですが、とにかくほとんど時間が費やされていないため、実際には最適化されていません)

私のネットブックでは、126秒かかりましたが、レナは結局のところ非常に小さな画像です。

編集:

かなりの数の平方根を節約するために変更rot_filters.std(1)するという提案があり、5秒のオーダーで何かを削りました。rot_filters.var(1)

4

2 に答える 2

1

複雑なピクセルごと+近傍操作の場合、パフォーマンスを向上させるためにcythonの使用を検討できます。これにより、後でCコードに変換されるPython構文に近いforループとしてコードを効率的に記述できます。

インスピレーションを得るために、たとえばscikit-imageコードを見てください。

https://github.com/scikit-image/scikit-image/blob/master/skimage/filter/_denoise.pyx

于 2012-10-19T08:54:58.047 に答える
1

Python+を使用して大幅に最適化するのは難しいと思いますscipy。ただし、ブールインデックスを使用するのではなく、を使用as_stridedして直接生成することで、わずかな改善を行うことができました。rot_filtersこれは、非常に単純なn次元windows関数に基づいています。(2D畳み込み関数がに存在することに気付く前に、この問題scipyを解決するために作成しました。)次のコードは、私のマシンで適度な10%の高速化を提供します。それがどのように機能するかの説明については、以下を参照してください。

import numpy as np
from scipy import ndimage
from numpy.lib import stride_tricks

# pass in `as_strided` as a default arg to save a global lookup
def rotation_matrix2(section, _as_strided=stride_tricks.as_strided):
    section = section.reshape(5, 5)  # sqrt(section.size), sqrt(section.size)
    windows_shape = (3, 3, 3, 3)     # 5 - 3 + 1, 5 - 3 + 1, 3, 3
    windows_strides = section.strides + section.strides
    windows = _as_strided(section, windows_shape, windows_strides)
    rot_filters = windows.reshape(9, 9)
    return rot_filters[rot_filters.std(1).argmin(),:].mean()

def get_rotation_smooth(im, _rm=rotation_matrix2, **kwargs):
    return ndimage.filters.generic_filter(im, _rm, size=5, **kwargs)

if __name__ == '__main__':
    import matplotlib.pyplot as plt
    from scipy.misc import lena
    im = lena()
    im2 = get_rotation_smooth(im)
    #plt.gray()      # Uncomment these lines for
    #plt.imshow(im2) # demo purposes.
    #plt.show()

上記の関数rotation_matrix2は、次の2つの関数と同等です(windowsより一般化されているため、実際には元の関数よりも少し遅くなります)。これは、元のコードが行うこととまったく同じです。9つの3x3ウィンドウを5x5配列に作成し、処理のためにそれらを9x9配列に再形成します。

def windows(a, w, _as_strided=stride_tricks.as_strided):
    windows_shape = tuple(sa - sw + 1 for sa, sw in zip(a.shape, w))
    windows_shape += w
    windows_strides = a.strides + a.strides
    return _as_strided(a, windows_shape, windows_strides)

def rotation_matrix1(section, _windows=windows):
    rot_filters = windows(section.reshape(5, 5), (3, 3)).reshape(9, 9)
    return rot_filters[rot_filters.std(1).argmin(),:].mean()

windowsウィンドウの次元数が同じである限り、任意の次元の配列で機能します。仕組みの内訳は次のとおりです。

    windows_shape = tuple(sa - sw + 1 for sa, sw in zip(a.shape, w))

windows配列は、nd配列のnd配列と考えることができます。外側の配列の形状は、大きい方の配列内のウィンドウの自由度によって決まります。すべての次元で、ウィンドウが取ることができる位置の数は、大きい方の配列の長さからウィンドウの長さを引いたものに1を加えたものに等しくなります。この場合、5x5配列への3x3ウィンドウがあるため、外側の2次元配列は3x3配列になります。

    windows_shape += w

内側の配列の形状は、ウィンドウ自体の形状と同じです。この場合も、これは3x3アレイです。

さあ、歩みを進めましょう。外側のnd配列と内側のnd配列のストライドを定義する必要があります。しかし、それらは同じであることが判明しました!結局のところ、ウィンドウは、個々のインデックスが配列内を移動するのとまったく同じ方法で、より大きな配列内を移動します。

    windows_strides = a.strides + a.strides

これで、ウィンドウを作成するために必要なすべての情報が得られました。

    return _as_strided(a, windows_shape, windows_strides)
于 2012-10-19T16:06:43.080 に答える