雙邊濾波(Bilateral filter)是一種非線性的濾波方法,是結合圖像的空間鄰近度和像素值類似度的一種折衷處理,同時考慮空域信息和灰度類似性,達到保邊去噪的目的。具備簡單、非迭代、局部的特色。算法
雙邊濾波器的好處是能夠作邊緣保存(edge preserving),通常過去用的維納濾波或者高斯濾波去降噪,都會較明顯地模糊邊緣,對於高頻細節的保護效果並不明顯。雙邊濾波器顧名思義比高斯濾波多了一個高斯方差sigma-d,它是基於空間分佈的高斯濾波函數,因此在邊緣附近,離的較遠的像素不會太多影響到邊緣上的像素值,這樣就保證了邊緣附近像素值的保存。可是因爲保存了過多的高頻信息,對於彩色圖像裏的高頻噪聲,雙邊濾波器不可以乾淨的濾掉,只可以對於低頻信息進行較好的濾波。數組
濾波算法中,目標點上的像素值一般是由其所在位置上的周圍的一個小局部鄰居像素的值所決定。在2D高斯濾波中的具體實現就是對周圍的必定範圍內的像素值分別賦以不一樣的高斯權重值,並在加權平均後獲得當前點的最終結果。而這裏的高斯權重因子是利用兩個像素之間的空間距離(在圖像中爲2D)關係來生成。經過高斯分佈的曲線能夠發現,離目標像素越近的點對最終結果的貢獻越大,反之則越小。其公式化的描述通常以下所述:函數
其中的c即爲基於空間距離的高斯權重,而用來對結果進行單位化。spa
高斯濾波在低通濾波算法中有不錯的表現,可是其卻有另一個問題,那就是隻考慮了像素間的空間位置上的關係,所以濾波的結果會丟失邊緣的信息。這裏的邊緣主要是指圖像中主要的不一樣顏色區域(好比藍色的天空,黑色的頭髮等),而Bilateral就是在Gaussian blur中加入了另外的一個權重分部來解決這一問題。Bilateral濾波中對於邊緣的保持經過下述表達式來實現:3d
其中的s爲基於像素間類似程度的高斯權重,一樣用來對結果進行單位化。對二者進行結合便可以獲得基於空間距離、類似程度綜合考量的Bilateral濾波:指針
上式中的單位化分部綜合了兩種高斯權重於一塊兒而獲得,其中的c與s計算能夠詳細描述以下:code
且有htm
且有blog
上述給出的表達式均是在空間上的無限積分,而在像素化的圖像中固然沒法這麼作,並且也不必如此作,於是在使用前須要對其進行離散化。並且也不須要對於每一個局部像素從整張圖像上進行加權操做,距離超過必定程度的像素實際上對當前的目標像素影響很小,能夠忽略的。限定局部子區域後的離散化公就能夠簡化爲以下形式:get
上述理論公式就構成了Bilateral濾波實現的基礎。爲了直觀地瞭解高斯濾波與雙邊濾波的區別,咱們能夠從下列圖示中看出依據。假設目標源圖像爲下述左右區域分明的帶有噪聲的圖像(由程序自動生成),藍色框的中心即爲目標像素所在的位置,那麼當前像素處所對應的高斯權重與雙邊權重因子3D可視化後的形狀如後邊兩圖所示:
左圖爲原始的噪聲圖像;中間爲高斯採樣的權重;右圖爲Bilateral採樣的權重。從圖中能夠看出Bilateral加入了類似程度分部之後能夠將源圖像左側那些跟當前像素差值過大的點給濾去,這樣就很好地保持了邊緣。爲了更加形象地觀察二者間的區別,使用Matlab將該圖在兩種不一樣方式下的高度圖3D繪製出來,以下:
上述三圖從左到右依次爲:雙邊濾波,原始圖像,高斯濾波。從高度圖中能夠明顯看出Bilateral和Gaussian兩種方法的區別,前者較好地保持了邊緣處的梯度,而在高斯濾波中,因爲其在邊緣處的變化是線性的,於是就使用連累的梯度呈現出漸變的狀態,而這表如今圖像中的話就是邊界的丟失(圖像的示例可見於後述)。
有了上述理論之後實現Bilateral Filter就比較簡單了,其實它也與普通的Gaussian Blur沒有太大的區別。這裏主要包括3部分的操做: 基於空間距離的權重因子生成;基於類似度的權重因子的生成;最終filter顏色的計算。
這就是一般的Gaussian Blur中使用的計算高斯權重的方法,其主要經過兩個pixel之間的距離並使用以下公式計算而來:
與基於距離的高斯權重計算相似,只不過此處再也不根據兩個pixel之間的空間距離,而是根據其類似程度(或者兩個pixel的值之間的距離)。
其中的表示兩個像素值之間的距離,能夠直接使用其灰度值之間的差值或者RGB向量之間的歐氏距離。
其中的類似度分部的權重s主要根據兩個Pixel之間的顏色差值計算面來。對於灰度圖而言,這個差值的範圍是能夠預知的,即[-255, 255],於是爲了提升計算的效率咱們能夠將該部分權重因子預計算生成並存表,在使用時快速查詢便可。使用上述實現的算法對幾張帶有噪聲的圖像進行濾波後的結果以下所示:
1 #include <time.h> 2 #include <stdio.h> 3 #include <stdlib.h> 4 #include <math.h> 5 #include "bf.h" 6 7 8 //-------------------------------------------------- 9 /** 10 雙邊濾波器圖像去噪方法 11 \param pDest 去噪處理後圖像buffer指針 12 \param pSrc 原始圖像buffer指針 13 \paran nImgWid 圖像寬 14 \param nImgHei 圖像高 15 \param nSearchR 搜索區域半徑 16 \param nSigma 噪聲方差 17 \param nCoff DCT變換系數個數 18 19 return void 20 */ 21 //-------------------------------------------------- 22 void BilateralFilterRGB(float *pDest, BYTE *pSrc,int nImgWid, int nImgHei, int nSearchR, int nSigma, int nCoff) 23 { 24 25 BYTE *pSrcR = new BYTE[nImgHei * nImgWid]; // 新建圖像緩衝區R 26 BYTE *pSrcG = new BYTE[nImgHei * nImgWid]; // 新建圖像緩衝區G 27 BYTE *pSrcB = new BYTE[nImgHei * nImgWid]; // 新建圖像緩衝區B 28 29 30 for(int i = 0;i < nImgHei; i++) 31 { 32 for (int j = 0; j < nImgWid; j++) 33 { 34 pSrcR[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 0]; 35 pSrcG[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 1]; 36 pSrcB[i * nImgWid + j] = pSrc[i * nImgWid * 3 + j * 3 + 2]; 37 } 38 } 39 40 41 float* pDestR = new float[nImgHei * nImgWid]; // 去噪後的圖像數據R通道 42 float* pDestG = new float[nImgHei * nImgWid]; // 去噪後的圖像數據G通道 43 float* pDestB = new float[nImgHei * nImgWid]; // 去噪後的圖像數據B通道 44 memset(pDestR,255,nImgHei * nImgWid * sizeof(float)); // 傳遞變量初始化 45 memset(pDestG,255,nImgHei * nImgWid * sizeof(float)); 46 memset(pDestB,255,nImgHei * nImgWid * sizeof(float)); 47 48 49 float *pII = new float[nImgHei * nImgWid ]; // 積分圖 50 float *pWeights = new float[nImgHei * nImgWid ]; // 每個窗的權重之和,用於歸一化 51 52 53 float mDctc[MAX_RANGE]; // DCT變換系數 54 float mGker[MAX_RANGE]; // 高斯核函數的值 55 56 for (int i = 0; i< MAX_RANGE; ++i) 57 { 58 mGker[i] = exp(-0.5 * pow(i / nSigma, 2.0)); // 首先計算0-255的灰度值在高斯變換下的取值,而後存在mGker數組中 59 } 60 61 CvMat G = cvMat(MIN_RANGE, MAX_RANGE, CV_32F, mGker); // 轉換成opencv中的向量 62 CvMat D = cvMat(MIN_RANGE, MAX_RANGE, CV_32F, mDctc); // 把mDctc係數數組轉換成opencv中的向量形式 63 64 cvDCT(&G,&D,CV_DXT_ROWS); // 而後將高斯核進行離散的dct變換 65 mDctc[0] /= sqrt(2.0); // 離散餘弦變換的第一個係數是1/1.414; 66 67 68 bf(pSrcR, pDestR, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR); 69 bf(pSrcG, pDestG, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR); 70 bf(pSrcB, pDestB, mDctc, pII, pWeights, nImgHei, nImgWid, nCoff, nSearchR); 71 72 73 for(int i=0; i < nImgHei;i++) 74 { 75 for (int j=0; j < nImgWid; j++) 76 { 77 pDest[i * nImgWid * 3 + j * 3 + 0] = pDestR[i * nImgWid + j]; 78 pDest[i * nImgWid * 3 + j * 3 + 1] = pDestG[i * nImgWid + j]; 79 pDest[i * nImgWid * 3 + j * 3 + 2] = pDestB[i * nImgWid + j]; 80 } 81 } 82 83 84 free(pSrcR); // 釋放臨時緩衝區 85 free(pSrcG); 86 free(pSrcB); 87 free(pDestR); 88 free(pDestG); 89 free(pDestB); 90 91 }
1 //-------------------------------------------------- 2 /** 3 雙邊濾波器圖像去噪方法 4 \param pDest 去噪處理後圖像buffer指針 5 \param pSrc 原始圖像buffer指針 6 \paran nImgWid 圖像寬 7 \param nImgHei 圖像高 8 \param pII 積分圖緩衝區 9 \param pWeights 歸一化權重係數 10 \param nCoff DCT變換系數 11 \param nPatchWid 圖像塊寬度 12 13 return void 14 */ 15 //-------------------------------------------------- 16 void bf (BYTE *pSrc, float *pDest, float *pDctc, float *pII, float *pWeights, int nImgHei, int nImgWid, int nCoff, int nPatchWid) 17 { 18 if (pDest == NULL || pSrc ==NULL) 19 { 20 return; 21 } 22 23 float *cR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float)); 24 float *sR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float)); 25 float *dcR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float)); 26 float *dsR = (float*) malloc((nCoff - 1) * MAX_RANGE * sizeof(float)); 27 28 29 // 餘弦變換查找表 30 int fx = 0; 31 int ck = 0; 32 int r2 = pow(2.0*nPatchWid + 1.0, 2.0); // 這個是圖像塊的像素點的總個數 33 34 float c0 = pDctc[0]; // 這是DCT的第一個係數 35 float c0r2 = pDctc[0] * r2; // 這個用於初始化歸一化權重之和 36 37 38 for (ck = 1; ck < nCoff; ck++) // 注意到這裏面的係數並不包含C0,因此後面纔有把C0額外的加上 39 { // 一個矩陣,用於存放各類係數和數值,便於查找 40 int ckr = (ck - 1) * MAX_RANGE; // 數組是從0開始的,因此這裏減1 41 42 for (fx = 0; fx<MAX_RANGE; fx++) // fx就是圖像的灰度值 43 { 44 float tmpfx = PI * float(fx) * float(ck) / MAX_RANGE; // ck其實就至關於那個餘弦函數裏面的t,fx至關於u,PI/MAX至關於前面的那個係數,這都是餘弦變換相關的東西 45 46 cR[ckr + fx] = cos(tmpfx); // 存儲餘弦變換,這個用在空間域上 47 sR[ckr + fx] = sin(tmpfx); // 存儲正弦變換 48 dcR[ckr + fx] = pDctc[ck] * cos(tmpfx); // 存儲餘弦變換和係數的乘積,這個則用在強度範圍上 49 dsR[ckr + fx] = pDctc[ck] * sin(tmpfx); // 存儲正弦變換和係數的乘積 50 } 51 } 52 53 54 float *pw = pWeights; // 新建一個歸一化權重的中間變量進行數據的初始化 55 float *pwe = pWeights + nImgWid * nImgHei; // 限定最大位置 56 57 while (pw < pwe) 58 { 59 *pw++ = c0r2; // 賦初值,讓它們都等於第一個DCT係數乘以圖像塊中總的像素點的個數 60 } 61 62 for (int ck = 1; ck < nCoff; ck++) 63 { 64 int ckr = (ck-1)*MAX_RANGE; // 數組是從0開始的,因此這裏減1 65 66 add_lut(pWeights, pSrc, pII,cR + ckr, dcR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid); // add cos term to pWeights 67 add_lut(pWeights, pSrc, pII,sR + ckr, dsR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid); // add cos term to pWeights 68 69 } 70 71 ImgRectSumO(pDest, pSrc, pII, c0, nImgHei, nImgWid, nPatchWid, nPatchWid); //加上C0的變換以後,再初始化濾波後的函數 72 73 for (int ck = 1; ck < nCoff; ck++) 74 { 75 int ckr = (ck-1)*MAX_RANGE; 76 77 add_f_lut(pDest, pSrc, pII, cR + ckr, dcR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid); // add cosine term to dataf 78 add_f_lut(pDest, pSrc, pII, sR + ckr, dsR + ckr, nImgHei, nImgWid, nPatchWid, nPatchWid); // add sine term to dataf 79 80 } 81 82 83 float *pd = pDest + nPatchWid * nImgWid; 84 pw = pWeights + nPatchWid * nImgWid; 85 86 float *pdie = pd + nImgWid - nPatchWid; 87 float *pdend = pDest + (nImgHei - nPatchWid) * nImgWid; 88 89 while (pd < pdend) 90 { 91 pd += nPatchWid; //把邊都給去掉了 92 pw += nPatchWid; //這個也是把邊去掉了 93 94 while (pd < pdie) 95 { 96 *pd++ /= ( (*pw++)*255); // 之因此要除以255,是爲了把它化到[0,1]之間,便於顯示 97 } 98 99 pd += nPatchWid; // 把邊都給去掉了 100 pw += nPatchWid; // 這個也是把邊去掉了 101 pdie += nImgWid; // 過渡到下一行 102 } 103 104 free(cR); 105 free(sR); 106 free(dcR); 107 free(dsR); // 釋放緩衝區 108 }
//-------------------------------------------------- /** 餘弦函數首係數積分圖 \param pDest 去噪處理後圖像buffer指針 \param pSrc 原始圖像buffer指針 \paran nImgWid 圖像寬 \param nImgHei 圖像高 \param pII 積分圖緩衝區 \param c0 DCT變換第一個係數 \param nPatchWid 圖像塊寬度 \param nPatchHei 圖像塊高度 return void */ //-------------------------------------------------- void ImgRectSumO (float* pDest, const BYTE* pSrc, float* pII, float c0,const int nImgHei, const int nImgWid, const int nPatchWid, const int nPatchHei) { if (pDest == NULL || pSrc ==NULL) { return; } int ri1 = nPatchWid + 1; // 中心像素點 int rj1 = nPatchHei + 1; int ri21 = 2 * nPatchWid + 1; // 圖像塊的寬度 int rj21 = 2 * nPatchHei + 1; const BYTE *pSrcImg = pSrc; // 原始圖像數據 const BYTE *pSrcImgEnd = pSrc + nImgHei * nImgWid; // 最大的像素點位置 const BYTE *pSrcImgWid = pSrcImg + nImgWid; // 第一行的最大的像素點位置,下面循環用的的變量 float *pIITmp = pII; // 積分圖數據 float *ppII = pIITmp; // 積分圖中間變量 *pIITmp++ = c0 * float( *pSrcImg++ ); // 第一個係數乘以圖像數據 while (pSrcImg < pSrcImgWid) // 遍歷第一行進行求積分圖,其實這個是圖像數據的積分圖 { *pIITmp++ = (*ppII++) + c0 * float(*pSrcImg++); } while (pSrcImg < pSrcImgEnd) // 遍歷全部的像素點的變量 { pSrcImgWid += nImgWid; // 進行第二行的循環 ppII = pIITmp; // 積分圖中間變量 *pIITmp++ = c0 * float(*pSrcImg++); // 第二行第一個像素點的處理 while (pSrcImg < pSrcImgWid) { (*pIITmp++) = (*ppII++) + c0 * float(*pSrcImg++); // 求第二行的積分圖 } float *pIIWid = pIITmp; // 當前位置像素點 pIITmp = pIITmp - nImgWid; // 上一行的起始點位置 float *pIIup = pIITmp - nImgWid; // 上上一行的起始點位置 while (pIITmp < pIIWid) // 遍歷整個行的每個數據 { (*pIITmp++) = (*pIITmp) + (*pIIup++); // 行與行之間的積分圖 } } float *pres = pDest + ri1 * nImgWid; // 最小的行位置 float *pend = pDest + (nImgHei - nPatchWid) * nImgWid; // 最大的行位置 float *pii1 = NULL; float *pii2 = NULL; float *pii3 = NULL; float *pii4 = NULL; // 求積分圖的四個點 pii1 = pII + ri21 * nImgWid + rj21; // 右下角 pii2 = pII + ri21 * nImgWid; // 左下角 pii3 = pII + rj21; // 右上角 pii4 = pII; // 左上角 while (pres < pend) { float *pe = pres + nImgWid - nPatchHei; pres += rj1; while (pres < pe) // 這個只是求圖像數據的積分圖 { (*pres++) = (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++); } pres += nPatchHei; pii1 += rj21; pii2 += rj21; pii3 += rj21; pii4 += rj21; } }
//-------------------------------------------------- /** 餘弦積分圖加函數 \param pNomalWts 像素點的歸一化權重矩陣 \param pSrc 原始圖像buffer指針 \param pII 積分圖 \paran nImgWid 圖像寬 \param nImgHei 圖像高 \param pII 積分圖緩衝區 \param pCosMtx 餘弦函數矩陣 \param pCosDctcMtx 餘弦函數與DCT變換系數乘積矩陣 \param nPatchWid 圖像塊寬度 \param nPatchHei 圖像塊高度 return void */ //-------------------------------------------------- void add_lut (float* pNomalWts, const BYTE* pSrc, float* pII, float* pCosMtx, float* pCosDctcMtx, const int nImgHei, const int nImgWid, const int nPatchWid, const int nPatchHei) { int ri1 = nPatchWid + 1; // kernel相關 int rj1 = nPatchHei + 1; int ri21 = 2 * nPatchWid + 1; // 這是kernel的寬度和長度 int rj21 = 2 * nPatchHei + 1; const BYTE *pi = pSrc; // 這個是圖像數據 float *pii = pII; // 這個是中間的緩衝變量,積分圖的緩衝區 const BYTE *piw = pi + nImgWid; // 也是賦初值的時候使用的中間變量 float *pii_p = pii; // 直線積分圖指針的指針 *pii++ = pCosMtx[*pi++]; //圖像數據值所對應的查找表中的餘弦或者正弦變換 while (pi < piw) { *pii++ = (*pii_p++) + pCosMtx[*pi++]; // 第一行的積分圖 } const BYTE *piend = pSrc + nImgHei * nImgWid; // 限定循環位置 while (pi < piend) // 遍歷整個圖像數據 { piw += nImgWid; // 第二行 pii_p = pii; // 指向積分圖指針的指針 *pii++ = pCosMtx[*pi++]; // 獲取第一個圖像數據值所對應的查找表中的餘弦或者正弦變換 while (pi < piw) { (*pii++) = (*pii_p++) + pCosMtx[*pi++]; // 求這一行的積分圖 } float *piiw = pii; pii = pii - nImgWid; float *pii_p1 = pii - nImgWid; while (pii < piiw) { (*pii++) = (*pii) + (*pii_p1++); // 行與行之間求積分圖 } } float *pres = pNomalWts + ri1 * nImgWid; // 定位要處理點的像素的那一行 float *pend = pNomalWts + (nImgHei - nPatchWid) * nImgWid; // 限定了邊界 float *pii1 = NULL; float *pii2 = NULL; float *pii3 = NULL; float *pii4 = NULL; pii1 = pII + ri21 * nImgWid + rj21; // 得到積分圖中那個最大的數據右下角 pii2 = pII + ri21 * nImgWid; // 積分圖左下角 pii3 = pII + rj21; // 積分圖右上角 pii4 = pII; // 積分圖左上角 pi = pSrc + ri1 * nImgWid; // 定位要處理的像素的位置 while (pres < pend) // 設定高度作循環 { float *pe = pres + nImgWid - nPatchHei; // 限定了寬度的範圍 pres = pres + rj1; // 定位了要處理的像素點的歸一化係數矩陣的位置 pi = pi + rj1; // 定位了要處理的像素點 while (pres < pe) // 設定寬度作循環 { (*pres++) = (*pres) + pCosDctcMtx[*pi++] * ( (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++) ); // 獲得的就是歸一化係數矩陣之和 } pres += nPatchHei; pi += nPatchHei; pii1 += rj21; pii2 += rj21; pii3 += rj21; pii4 += rj21; // 略過不處理的邊界區域 } }
//-------------------------------------------------- /** 積分圖 \param pDest 像素點的歸一化權重矩陣 \param pSrc 原始圖像buffer指針 \param pII 積分圖 \paran nImgWid 圖像寬 \param nImgHei 圖像高 \param pII 積分圖緩衝區 \param pCosMtx 餘弦函數矩陣 \param pCosDctcMtx 餘弦函數與DCT變換系數乘積矩陣 \param nPatchWid 圖像塊寬度 \param nPatchHei 圖像塊高度 return void */ //-------------------------------------------------- void add_f_lut (float* pDest,const BYTE* pSrc,float* pII,float* pCosMtx,float* pCosDctcMtx,const int nImgHei, const int nImgWid,const int nPatchWid, const int nPatchHei) { int ri1 = nPatchWid + 1; // kernel相關,目標是指向中心像素點 int rj1 = nPatchHei + 1; // 目標是指向中心像素點 int ri21 = 2 * nPatchWid + 1; // 這是kernel的寬度和長度 int rj21 = 2 * nPatchHei + 1; // 其實就是指向搜索區域的右下角,也就是最大位置的那個像素點,就是邊界了 const BYTE *pi = pSrc; // 這個是圖像數據,圖像的灰度值 float *pii = pII; // 這是積分圖 const BYTE *piw = pi + nImgWid; // 指向第一行的末尾位置,用於循環控制 float *pii_p = pii; // 指向積分圖指針的指針 *pii++ = float(*pi) * pCosMtx[*pi]; // 灰度值乘以餘弦係數 ++pi; while (pi < piw) // 第一行遍歷 { *pii++ = (*pii_p++) + float(*pi) * pCosMtx[*pi]; ++pi; } const BYTE *piend = pSrc + nImgHei * nImgWid; while (pi < piend) { piw += nImgWid; pii_p = pii; *pii++ = float(*pi) * pCosMtx[*pi]; ++pi; while (pi < piw) { (*pii++) = (*pii_p++) + float(*pi) * pCosMtx[*pi]; ++pi; } float *piiw = pii; pii = pii - nImgWid; // 上一行起點 float *pii_p1 = pii - nImgWid; // 上上一行的起始點 while (pii < piiw) // 循環一行 { (*pii++) += (*pii_p1++); //其實就是在列的位置上加上上一行 } } float *pres = pDest + ri1*nImgWid; float *pend = pDest + (nImgHei - nPatchWid) * nImgWid; float *pii1 = NULL; float *pii2 = NULL; float *pii3 = NULL; float *pii4 = NULL; pii1 = pII + ri21 * nImgWid + rj21; // 積分圖搜索區域最大位置的元素,堪稱右下角的元素 pii2 = pII + ri21 * nImgWid; // 能夠當作搜索區域中左下角的像素位置 pii3 = pII + rj21; // 搜索區域右上角像素位置 pii4 = pII; // 搜索區域左上角像素位置 pi = pSrc + ri1 * nImgWid; // 定位要處理的像素的位置的那一行 while (pres < pend) { float *pe = pres + nImgWid - nPatchHei; //限定了寬度的範圍 pres = pres + rj1; //定位了要處理的像素點的歸一化係數矩陣的位置 pi = pi + rj1; //定位了要處理的像素點 while (pres < pe) //遍歷整個行 { (*pres++) = (*pres) + pCosDctcMtx[*pi++] * ( (*pii1++) - (*pii2++) - (*pii3++) + (*pii4++) ); //這個實際上是計算整個像素塊 } pres += nPatchHei; pi += nPatchHei; pii1 += rj21; pii2 += rj21; pii3 += rj21; pii4 += rj21; } }