解析opencv中Box Filter的實作并提出進一步加速的方案（源碼共享）。

Box Filter，最經典的一種領域操作，在無數的場合中都有着廣泛的應用，作為一個很基礎的函數，其性能的好壞也直接影響着其他相關函數的性能。是以其優化的檔次和程度是非常重要的，本文對opencv的Box Filter進行了解讀，并給出了對其進一步優化的方案和代碼，在此，也邀請各位高手進一步優化。

　　說明：本文所有算法的涉及到的優化均指在PC上進行的，對于其他構架是否合适未知，請自行試驗。

Box Filter，最經典的一種領域操作，在無數的場合中都有着廣泛的應用，作為一個很基礎的函數，其性能的好壞也直接影響着其他相關函數的性能，最典型莫如現在很好的EPF濾波器：GuideFilter。是以其優化的檔次和程度是非常重要的，網絡上有很多相關的代碼和部落格對該算法進行講解和優化，提出了不少O(1)算法，但所謂的0(1)算法也有優劣之分，0(1)隻是表示執行時間和某個參數無關，但本身的耗時還是有差別。比較流行的就是累積直方圖法，其實這個是非常不可取的，因為稍大的圖就會導緻溢出，這裡我們解析下 opencv的相關代碼，看看他是如何進行優化的。

首先找到opencv的代碼的位置，其在\opencv\sources\modules\imgproc\src\smooth.cpp中。

Box Filter 是一種行列可分離的濾波，是以，opencv也是這樣處理的，先進行行方向的濾波，得到中間結果，然後再對中間結果進行列方向的處理，得到最終的結果。

opencv 行方向處理的相關代碼如下：

template<typename T, typename ST>
struct RowSum :
        public BaseRowFilter
{
    RowSum( int _ksize, int _anchor ) :
        BaseRowFilter()
    {
        ksize = _ksize;
        anchor = _anchor;
    }

    virtual void operator()(const uchar* src, uchar* dst, int width, int cn)
    {
        const T* S = (const T*)src;
        ST* D = (ST*)dst;
        int i = 0, k, ksz_cn = ksize*cn;

        width = (width - 1)*cn;
        for( k = 0; k < cn; k++, S++, D++ )
        {
            ST s = 0;
            for( i = 0; i < ksz_cn; i += cn )
                s += S[i];
            D[0] = s;
            for( i = 0; i < width; i += cn )
            {
                s += S[i + ksz_cn] - S[i];
                D[i+cn] = s;
            }
        }
    }
};

　　這個代碼考慮了多個通道以及多種資料類型的情況，為了分析友善我們重寫下單通道時的核心代碼。

for(Z = 0, Value = 0; Z < Size; Z++)    
    Value += RowData[Z];
LinePD[0] = Value;

for(X = 1; X < Width; X ++)
{
    Value += RowData[X + Size - 1] - RowData[X - 1];    
    LinePD[X] = Value;               
}

　　上述代碼中RowData是指對某一行像素進行擴充後的像素值，其寬度為Width + Radius + Radius, Radius是值濾波的半徑，而代碼中Size = 2 * Radius + 1，LinePD即所謂的中間結果，考慮資料類型能表達的範圍，必須使用 int類型。

對于每行第一個點很簡單，直接用for計算出行方向的指定半徑内的累加值。而之後所有點，用前一個累計值加上新加入的點，然後去除掉移出的哪一個點，得到新的累加值。這個方法在很多文獻中都有提及，并沒有什麼新鮮之處，這裡不多說。

按照正常的思維，在列方向的處理應該和行方向完全相同，隻是處理的方向的不同、處理的資料源的不同以及最後需要對計算結果多一個歸一化的過程而已。事實也是如此，有很多算法都是這樣做的，但是由于CPU構架的一些原因（主要是cache miss的增加），同樣的算法沿列方向處理總是會比沿行方向慢一個檔次，解決方式有很多，例如先對中間結果進行轉置，然後再按照行方向的規則進行處理，處理完後在将資料轉置回去。轉置過程是有非常高效的處理方式的，借助于SSE以及Cache優化，能實作比原始兩重for循環快4倍以上的效果。還有一種方式就正如opencv中列處理過程所示，這正是下面将要重點描述的。

在opencv的代碼中，并沒有直接沿列方向處理，而是繼續沿着行的方向一行一行處理，我先貼出我自己翻譯的有關純C的代碼進行解說：

for (Y = 0; Y < Size - 1; Y++)            //    注意沒有最後一項哦                      
    {
        int *LinePS = (int *)(Sum->Data + ColPos[Y] * Sum->WidthStep);
        for(X = 0; X < Width; X++)    ColSum[X] += LinePS[X];
    }

    for (Y = 0; Y < Height; Y++)
    {
        unsigned char* LinePD    = Dest->Data + Y * Dest->WidthStep;    
        int *AddPos              = (int*)(Sum->Data + ColPos[Y + Size - 1] * Sum->WidthStep);
        int *SubPos              = (int*)(Sum->Data + ColPos[Y] * Sum->WidthStep);

        for(X = 0; X < Width; X++)
        {
            Value = ColSum[X] + AddPos[X];
            LinePD[X] = Value * Scale;                    
            ColSum[X] = Value - SubPos[X];
        }
    }

上述代碼中定義了一個ColSum用于儲存每行某個位置處在列方向上指定半徑内各中間元素的累加值，對于第一行，做特殊處理，其他行的每個元素的處理方式和BaseRowFilter裡的處理方式很像，隻是在書寫上有所差別，特别注意的是對第一行的累加時，最後一個元素并沒有計算在内，這個處理技巧在下面的X循環裡有展現，請大家仔細體味下。

上述代碼這樣做的好處是，能有效的減少CPU的cache miss，但是總的計算量是沒有改變的，是以能有效的提高速度。

針對PC，在opencv内部，其在列方向上采用了SSE進行進一步的優化，我們貼出其處理uchar類型的代碼：

解析opencv中Box Filter的實作并提出進一步加速的方案（源碼共享）。

View Code

　　代碼比較多，我稍微精簡下（注意，精簡後的是不可運作的，隻是更清晰的表達了思路）。

virtual void operator()(const uchar** src, uchar* dst, int dststep, int count, int width)
{
    int i;
    int* SUM;
    bool haveScale = scale != 1;
    double _scale = scale;
    if( width != (int)sum.size() )
    {
        sum.resize(width);
        sumCount = 0;
    }

    SUM = &sum[0];
    if( sumCount == 0 )
    {
        memset((void*)SUM, 0, width*sizeof(int));
        for( ; sumCount < ksize - 1; sumCount++, src++ )
        {
            const int* Sp = (const int*)src[0];
            i = 0;
           
            for( ; i <= width-4; i+=4 )
            {
                __m128i _sum = _mm_loadu_si128((const __m128i*)(SUM+i));
                __m128i _sp = _mm_loadu_si128((const __m128i*)(Sp+i));
                _mm_storeu_si128((__m128i*)(SUM+i),_mm_add_epi32(_sum, _sp));
            }
            for( ; i < width; i++ )
                SUM[i] += Sp[i];
        }
    }
    else
    {
        src += ksize-1;
    }

    for( ; count--; src++ )
    {
        const int* Sp = (const int*)src[0];
        const int* Sm = (const int*)src[1-ksize];
        uchar* D = (uchar*)dst;

        i = 0;
        const __m128 scale4 = _mm_set1_ps((float)_scale);
        for( ; i <= width-8; i+=8 )
        {
            __m128i _sm  = _mm_loadu_si128((const __m128i*)(Sm+i));
            __m128i _sm1  = _mm_loadu_si128((const __m128i*)(Sm+i+4));

            __m128i _s0  = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(SUM+i)),
                                         _mm_loadu_si128((const __m128i*)(Sp+i)));
            __m128i _s01  = _mm_add_epi32(_mm_loadu_si128((const __m128i*)(SUM+i+4)),
                                          _mm_loadu_si128((const __m128i*)(Sp+i+4)));

            __m128i _s0T = _mm_cvtps_epi32(_mm_mul_ps(scale4, _mm_cvtepi32_ps(_s0)));
            __m128i _s0T1 = _mm_cvtps_epi32(_mm_mul_ps(scale4, _mm_cvtepi32_ps(_s01)));

            _s0T = _mm_packs_epi32(_s0T, _s0T1);

            _mm_storel_epi64((__m128i*)(D+i), _mm_packus_epi16(_s0T, _s0T));

            _mm_storeu_si128((__m128i*)(SUM+i), _mm_sub_epi32(_s0,_sm));
            _mm_storeu_si128((__m128i*)(SUM+i+4),_mm_sub_epi32(_s01,_sm1));
        }
        for( ; i < width; i++ )
        {
            int s0 = SUM[i] + Sp[i];
            D[i] = saturate_cast<uchar>(s0*_scale);
            SUM[i] = s0 - Sm[i];
        }

        dst += dststep;
    }
}

在行方向上，ColSum每個元素的更新互相之間是沒有任何關系的，是以，借助于SSE可以實作一次性的四個元素的更新，并且上述代碼還将第一行的特殊計算也用SSE實作了，雖然這部分計算量是非常小的。

在具體的SSE細節上，對于uchar類型，雖然中間結果是用int類型表達的，但是由于SSE沒有整形的除法指令（是否有？），是以上面借用了浮點數的乘法SSE指令實作，當然就多了_mm_cvtepi32_ps 以及 _mm_cvtps_epi32這樣的類型轉換的函數。如果有整形除法，那就能更好了。

SSE的實作，無非也就是用_mm_loadu_si128加載資料，用_mm_add_epi32， _mm_mul_ps之類的函數進行基本的運算，用_mm_storeu_si128來儲存資料，并沒有什麼特别複雜的部分，注意到考慮到資料的普遍性，不一定都是16位元組對齊的，是以在加載和儲存是需要使用u方面的函數，其實作在的_mm_loadu_si128和_mm_load_si128在處理速度上已經沒有特别明顯的差別了。

注意到在每個SSE代碼後面，總還有部分C代碼，這是因為我們實際資料寬度并不一定是4的整數倍，未被SSE處理的部分必須使用C實作，這其實對讀者來說是非常好的事情，因為我們能從這部分C代碼中搞明白上述的SSE代碼是幹什麼的。

　　以上就是opencv的Box Filter實作的關鍵代碼，如果讀者去看具體細節，opencv還有針對手機上的Neon優化，其實這個和SSE的意思基本是一樣的。

　　那是否還有改進的空間呢，從我的認知中，在對垂直方向的處理上，應該沒有了，但是水準方向呢， SSE有沒有用武之地，經過我的思考我認為還是有的。

　　在行方向的計算中，這個for循環是主要的計算。

for(X = 1; X < Width; X ++)
{
    Value += RowData[X + Size - 1] - RowData[X - 1];    
    LinePD[X] = Value;               
}

本人認為雖然這裡的操作是前後依賴的，全局無法并行化，但是觀察這一行代碼：

Value += RowData[X + Size - 1] - RowData[X - 1];

其中RowData[X + Size - 1] - RowData[X - 1]; 并不是前後依賴的，是可以并行的，是以如果提前快速的計算出所有的內插補點，那麼提速的空間還比較大，即需要提前計算出下面的函數：

for(X = 0; X < (Width - 1); X++)
     Diff[X] = AddPos[X] - SubPos[X];

　　這裡用SSE實作則非常簡單了：

unsigned char *AddPos = RowData + Size * Channel;
        unsigned char *SubPos = RowData;
        X = 0;                    //    注意這個指派在下面的循環外部，這可以避免當Width<8時第二個for循環循環變量未初始化            
        __m128i Zero = _mm_setzero_si128();
        for (; X <= (Width - 1) * Channel - 8; X += 8)
        {
            __m128i Add = _mm_unpacklo_epi8(_mm_loadl_epi64((__m128i const *)(AddPos + X)), Zero);        
            __m128i Sub = _mm_unpacklo_epi8(_mm_loadl_epi64((__m128i const *)(SubPos + X)), Zero);        
            _mm_store_si128((__m128i *)(Diff + X + 0), _mm_sub_epi32(_mm_unpacklo_epi16(Add, Zero), _mm_unpacklo_epi16(Sub, Zero)));        //    由于采用了_aligned_malloc函數配置設定記憶體，可是使用_mm_store_si128
            _mm_store_si128((__m128i *)(Diff + X + 4), _mm_sub_epi32(_mm_unpackhi_epi16(Add, Zero), _mm_unpackhi_epi16(Sub, Zero)));
        }
        for(; X < (Width - 1) * Channel; X++)
            Diff[X] = AddPos[X] - SubPos[X];

　　和列方向的SSE處理代碼不同的是，這裡加載的是uchar資料，是以加載的函數就有所不同，處理的方式也有差別，上面幾個SSE函數各位查查MSDN就能了解其意思，還是很有味道的。

　　經過這樣的處理，經過測試發現，速度能夠比opencv的版本在提高30%，也是額外的驚喜。

再有的優化可能提速有限，比如把剩下的一些for循環内部分成四路等等。

在我的I5桌上型電腦中，使用上述算法對1024*1024的三通道彩色圖像進行半徑為5的測試，進行100次的計算純算法部分的耗時為800ms，如果是純C版本大概為1800ms；當半徑為200時，SSE版本約為950ms， C版本約1900ms；當半徑為400時，SSE版本約為1000ms， C版本約2100ms；可見，半徑增大，耗時稍有增加，這主要是由于算法中有部分初始化的計算和半徑有關，但是這些都是不重要的。

在不使用多線程（雖然本算法非常适合多線程計算），不使用GPU，隻使用單線程\CPU進行計算的情況下，個人覺得目前我這個代碼是很難有質的超越的，從某個方面講，代碼中的用于計算的時間占用的時間比從記憶體等待資料的時間可能還要短，而似乎也沒有更好的算法能進一步減少記憶體資料的通路量。

本人在這裡邀請各位高手對目前我優化的這個代碼進行進一步的優化，希望高人不要謙虛。

　　運作界面：

解析opencv中Box Filter的實作并提出進一步加速的方案（源碼共享）。

　　本文的代碼是針對常用的圖像資料進行的優化處理，在很多場合下，需要對int或者float類型進行處理，比如GuideFilter，如果讀者了解了本文解讀的代碼的原理，更改代碼以便适應他們則是一件非常簡單的事情，如果您不會，那麼也請你不要給我留言，或者我可以有償提供，因為從本質上講我喜歡提供漁，而不是魚，漁具已經送給你了，你卻找我要魚，那隻能..............。

本文源代碼下載下傳位址（VS2010編寫）：Box Filter

****************************作者： laviewpbt 時間： 2015.12.17 聯系QQ: 33184777 轉載請保留本行資訊**********************

解析opencv中Box Filter的實作并提出進一步加速的方案（源碼共享）。

繼續閱讀

基于FAST的TSN交換（3）TSN循環隊列轉發（CQF）原理基于FAST的TSN交換（3）TSN循環隊列轉發（CQF）原理

Fast DFS概述以及單機版環境搭建

【特征比對】Fast原理及源碼解析

特征提取算法與圖像描述符[ORB]ORB特征算法

高性能可程式設計實驗平台OpenBox-S28

第18篇：fastjson反序列化漏洞區分版本号的方法總結

哥唱rap額grape塵土落central.附近的nearby買了mail.東方吃east.fast刺那亭fascina

FAST 3.0流水線擴充模型推出，支援TSN交換！FAST 3.0流水線擴充模型推出，支援TSN交換

2years，toofast，toomany#聯考#畢業季

你的眼睛背叛你的心。鄭中基。你的眼睛背叛你的心。詞曲：古倩敏。别裝作仍然溫柔，别裝作一切平如。我被逼曾捱過了多少個年，不

BFCNetworkBreakfast倫敦時裝周時隔三年再次見到了CarolineRush、CarolineIssa老朋

深度學習目标檢測之 R-CNN 系列：從 R-CNN 和 Fast R-CNN 到 Faster R-CNN深度學習目标檢測之 R-CNN 系列：從 R-CNN 和 Fast R-CNN 到 Faster R-CNN 總覽

給AI套上限制的枷鎖。不到一個月chatGPT就從萬人熱捧變成世界各地圍堵的麻煩制造者。先是馬斯克、圖靈獎得主等1200

上交所FAST行情接口對接前言FAST接口登入STEP解析FAST解析完整代碼如下：

⬆️國外怎麼玩瓦羅蘭特國服是海外玩家們想知道的問題之一，海外玩家想要遊玩《無畏契約》也就是國服版本的瓦羅蘭特，需要使用w

Java 中 Fail-fast 與 Fail-safe 機制的深入了解與應用