標簽(空格分隔): 算法學(xué)習(xí)
……(留在日后填寫,原理以及實現(xiàn)細節(jié),這里先直接貼程序吧)
在做衛(wèi)星圖像質(zhì)量提高的時候,相機噪聲通常都是很大的,在項目里面采用了5-3提升方向小波,為了提高程序的運行速度,將原來的9個方向縮減為現(xiàn)在的3個方向,這么做的好處是什么呢?當(dāng)采用9個方向的方向提升小波時,需要對圖像進行預(yù)處理插值,9個方向時需要分別插出1/2像素點和1/4像素點。依次類推,采用5個方向預(yù)測時,則只需要插出1/2像素點。預(yù)測3個方向時,則不需要提前對圖像進行插值。 在寫c語言的時候,我也重新整理了一下3方向提升小波,下面我按照思路將程序一塊塊的貼上來:
方向提升小波的正變換對應(yīng)于小波正變換,會得到圖像的LL,LH,HL,HH四個分量,分別對應(yīng)于圖像的低頻分量,水平方向上的高頻分量,垂直方向上的高頻分量以及45度角上的高頻分量。
cout << "----------------------wavelet denoise start---------------------" << endl; Mat * img_in = new Mat(); Mat * img_out = new Mat(); *img_in = Mat::zeros(rows_ori, cols_ori, CV_64FC1); *img_out = Mat::zeros(rows_ori, cols_ori, CV_64FC1); // double * img_ori_ptr ---> Mat::img_ori int ptr_bias = 0; for (int i = 0; i < rows_ori; ++i) { for (int j = 0; j < cols_ori; ++j) { img_in->at<double>(i, j) = *(image_in + ptr_bias); ptr_bias++; } } int rows = img_in->size().height; int cols = img_in->size().width; clock_t start, end; start = clock(); cout << "------step1: Forward Transform-----" << endl; ////////////////////////////////小波一級分解/////////////////////////////// Mat* LL = new Mat(); LL->create(rows / 2, cols / 2, CV_64FC1); Mat* LH = new Mat(); LH->create(rows / 2, cols / 2, CV_64FC1); Mat* HL = new Mat(); HL->create(rows / 2, cols / 2, CV_64FC1); Mat* HH = new Mat(); HH->create(rows / 2, cols / 2, CV_64FC1); int *dir_mat11 = new int[rows * cols]; int *dir_mat12 = new int[(rows / 2) * cols]; ADL97(img_in, 1, LL, LH, HL, HH, dir_mat11, dir_mat12, blk_size);上述函數(shù)分別完成了將double型的圖像指針轉(zhuǎn)換到mat矩陣中,并提前為LL,LH,HL,HH分配好內(nèi)存空間,其中dirMat變量用來存放預(yù)測的方向矩陣,方向矩陣中的取值為,-1,0,1分別對應(yīng)于三個方向(左上-右下,上-下,右上-左下)。
方向預(yù)測的過程是對圖像塊進行處理,在程序中對512*512的圖像進行方向預(yù)測時默認的塊大小為32*32,小波二級分解時圖像塊的大小變?yōu)?6*16,三級分解變?yōu)?*8。在此處編程時,應(yīng)用到了一個常用的技巧,即就是取出圖像塊,對圖像塊進行處理,再將處理后的圖像塊放置回原圖對應(yīng)的位置。 下面的getBlock和setBlock函數(shù)就是進行上述取塊放塊操作:
void getBlock(Mat* img, double * blk_img, int* blk_position) { /*對圖像進行分塊 輸入: img:輸入的待分塊圖像 win_size : 塊的大小 輸入圖像的長和寬必須為win_size的整數(shù)倍大小 輸出: blk_img:圖像塊 double* */ int p_rows_f = blk_position[0];// 首行 int p_rows_l = blk_position[1];// 末行 int p_cols_f = blk_position[2];// 首列 int p_cols_l = blk_position[3];// 末列 int r_rows = p_rows_l - p_rows_f; int r_cols = p_cols_l - p_cols_f; for (int i = p_rows_f, r_i = 0; i < p_rows_l; i++, r_i++) { for (int j = p_cols_f, r_j = 0; j < p_cols_l; j++, r_j++) { blk_img[r_i*r_cols + r_j] = img->at<double>(i,j); } } }void setBlock(int* dir_mat, int* blk_img, int* blk_position,int img_cols){ /* 輸入: blk_img:圖像塊 double*塊中的數(shù)據(jù)位當(dāng)前塊的預(yù)測方向,取值應(yīng)為(-1,0,1,); blk_position:當(dāng)前塊四個頂點的左邊索引 img_cols:輸入圖像的列數(shù) 輸出: dir_mat:由所有塊拼接而成的方向矩陣,大小與輸入圖像一致; */ //將小圖像塊拼接稱為整個圖像 int p_rows_f = blk_position[0];// 首行 int p_rows_l = blk_position[1];// 末行 int p_cols_f = blk_position[2];// 首列 int p_cols_l = blk_position[3];// 末列 int r_rows = p_rows_l - p_rows_f; int r_cols = p_cols_l - p_cols_f; for (int i = p_rows_f, r_i = 0; i<p_rows_l; i++, r_i++) { for (int j = p_cols_f, r_j = 0; j<p_cols_l; j++, r_j++) { dir_mat[i * img_cols + j] = blk_img[r_i * r_cols + r_j]; /*cout << i << endl;*/ } }}下面貼出方向預(yù)測的主函數(shù):
void getDirMat(Mat* img, int flag, int decom_level, int* dir_mat,int blk_size){ /* 輸入: img:待處理圖像 Mat* flag:當(dāng)前處理標志 flag == 1:列變換 flag == 2:行變換 decom_level:分解級數(shù) 輸出: dir_mat:方向矩陣 int* 大小與輸入圖像一致 */ int win_size = blk_size / (pow(2, (decom_level - 1))); int rows = img->size().height; int cols = img->size().width; int count_rows = rows / win_size; int count_cols = cols / win_size; //當(dāng)前塊的位置 double* blk_img = new double[win_size * win_size]; int blk_position[4] = { 0 }; int temp; int* temp_mat = new int[win_size * win_size]; //將圖像分解成win*win的塊進行計算 for (int i = 0; i < count_rows; i++) { for (int j = 0; j < count_cols; j++) { //當(dāng)前塊的四個角坐標索引 blk_position[0] = i * win_size; blk_position[1] = (i + 1) * win_size; blk_position[2] = j * win_size; blk_position[3] = (j + 1) * win_size; //讀取到當(dāng)前塊 getBlock(img, blk_img, blk_position); //對當(dāng)前塊進行方向預(yù)測 temp = dirEstimation2(blk_img, flag, win_size); for (int k = 0; k < win_size * win_size; k++) { temp_mat[k] = temp; } setBlock(dir_mat, temp_mat, blk_position, cols); } } delete[] blk_img; delete[] temp_mat;}其中子函數(shù)dirEstimation2對輸入的圖像塊用奇數(shù)行預(yù)測偶數(shù)行,偶數(shù)行就包含了高頻信息,分別計算是哪個方向預(yù)測出的高頻矩陣,取高頻矩陣二范數(shù)最小,即能說明該方向的預(yù)測是準確的,塊中所有點都為此方向。下面貼出方向預(yù)測的代碼:
void matTrans(int* mat, int* mat_T, int rows, int cols){ for (int i = 0; i < cols; i++) { for (int j = 0; j < rows; j++) { mat_T[i*rows + j] = mat[j*cols + i]; } } }int minIdx(double* power, int data_size){ double min = power[0]; int idx = 0; for (int i = 1; i < data_size; i++) { if (power[i] < min) { min = power[i]; idx = i; } } return idx;}double calcuPower(double* high_freq, int win_size){ double sum = 0; for (int i = 0; i < win_size / 2 ; i++) { for (int j = 0; j < win_size ; j++) { sum = sum + pow(high_freq[i*win_size + j] ,2); } } return sum;}void downSampling(double* img, int win_size, double* high_freq,int flag){ if (flag == 1)//列變換(行下2采樣) { for (int i = 0; i < win_size / 2; i++) { for (int j = 0; j < win_size; j++) { high_freq[i*win_size + j] = img[(2 * i)*win_size + j]; } } } if (flag == 2)//行變換(列下2采樣) { for (int i = 0; i < win_size; i++) { for (int j = 0; j < win_size / 2; j++) { high_freq[i * win_size / 2 + j] = img[i*win_size + (2*j)]; } } }}int dirEstimation2(double* blk_img, int flag, int win_size){ double power[3] = { 0 }; double* high_freq_mat = new double[win_size * (win_size / 2)]; memset(high_freq_mat, 0, sizeof(double)*win_size*(win_size / 2)); //把圖像塊的內(nèi)容暫存在blk的Mat* 里面 Mat* blk_mat_img = new Mat(); blk_mat_img->create(win_size, win_size, CV_64FC1); //邊界擴展后的圖像塊 Mat* blk_expanded_rows = new Mat(); blk_expanded_rows->create(win_size + 2, win_size , CV_64FC1); Mat* blk_expanded_cols = new Mat(); blk_expanded_cols->create(win_size + 2, win_size + 2, CV_64FC1); //把double*的圖像塊數(shù)據(jù)存放在Mat矩陣里面 for (int i = 0; i < win_size; i++) { for (int j = 0; j < win_size; j++) { blk_mat_img->at<double>(i, j) = blk_img[i*win_size + j]; } } if (flag == 1)//當(dāng)前為列操作 { for (int k = -1; k < 2; k++) { //循環(huán)之前,需要將blk_mat_img的數(shù)據(jù)重新賦值給blk_expanded,并填充邊框 copyMakeBorder(*blk_mat_img, *blk_expanded_rows, 1, 1, 0, 0, BORDER_REFLECT_101); copyMakeBorder(*blk_expanded_rows, *blk_expanded_cols, 0, 0, 1, 1, BORDER_REFLECT); //數(shù)據(jù)準備完畢,進行方向預(yù)測 for (int i = 1; i < win_size + 1; i = i + 2) { for (int j = 1; j < win_size + 1; j++) { blk_expanded_cols->at<double>(i, j) = blk_expanded_cols->at<double>(i, j) - 0.5*( blk_expanded_cols->at<double>(i-1,j+k) + blk_expanded_cols->at<double>(i+1,j-k) ); } } for (int i = 0; i < win_size; i++) { for (int j = 0; j < win_size; j++) { blk_img[i*win_size + j] = blk_expanded_cols->at<double>(i + 1, j + 1); } } downSampling(blk_img, win_size, high_freq_mat, flag); power[k + 1] = calcuPower(high_freq_mat, win_size); } } if (flag == 2) { Mat* blk_img_T = new Mat(); blk_img_T->create(win_size, win_size, CV_64FC1); //圖像直接轉(zhuǎn)置 for (int i = 0; i < win_size; i++) { for (int j = 0; j < win_size; j++) { blk_img_T->at<double>(i, j) = blk_mat_img->at<double>(j, i); } } for (int k = -1; k < 2; k++) { //循環(huán)之前,需要將blk_mat_img的數(shù)據(jù)重新賦值給blk_expanded,并填充邊框 copyMakeBorder(*blk_img_T, *blk_expanded_rows, 1, 1, 0, 0, BORDER_REFLECT_101); copyMakeBorder(*blk_expanded_rows, *blk_expanded_cols, 0, 0, 1, 1, BORDER_REFLECT); //數(shù)據(jù)準備完畢,進行方向預(yù)測 for (int i = 1; i < win_size + 1; i = i + 2) { for (int j = 1; j < win_size + 1; j++) { blk_expanded_cols->at<double>(i, j) = blk_expanded_cols->at<double>(i, j) - 0.5*(blk_expanded_cols->at<double>(i - 1, j + k) + blk_expanded_cols->at<double>(i + 1, j - k)); } } for (int i = 0; i < win_size; i++) { for (int j = 0; j < win_size; j++) { blk_img[i*win_size + j] = blk_expanded_cols->at<double>(i + 1, j + 1); } } downSampling(blk_img, win_size, high_freq_mat, 1); power[k + 1] = calcuPower(high_freq_mat, win_size); } delete blk_img_T; } int dir = 0; dir = minIdx(power, 3) - 1; delete[] high_freq_mat; delete blk_mat_img; delete blk_expanded_rows; delete blk_expanded_cols; return dir;}下面是方向小波核心程序
void singleLifting(Mat* img, int flag, int step, double coe, int* dir_mat){ /*小波變換單步提升,核心步驟 輸入: img:輸入圖像 Mat* flag: 行列處理標志 flag == 1列變換, flag ==2 行變換 step: 小波處理步驟標志 step == 1預(yù)測 step == 2更新 coe:小波提升系數(shù) dir_mat: 方向矩陣,大小與輸入圖像一致 輸出: img: 輸出圖像,在函數(shù)內(nèi)部對圖像進行處理 */ int rows = img->size().height; int cols = img->size().width; if (flag == 1)//列變換 { //預(yù)處理,圖像對稱延拓一圈 Mat *extend_img = new Mat(); extend_img->create(rows + 1, cols + 1, CV_64FC1); copyMakeBorder(*img, *extend_img, 1, 1, 1, 1, BORDER_REFLECT_101); int extend_rows = extend_img->size().height; int extend_cols = extend_img->size().width; //第一步,確認輸入數(shù)據(jù)的正確性 //輸出擴展后的圖像數(shù)據(jù)和進行預(yù)測的方向矩陣 if (step == 1)//預(yù)測 { for (int i = 1; i < extend_rows - 2; i = i + 2)//剔除原圖最后一行,再2:2:rows-1 { for (int j = 1; j < extend_cols - 1; j++)//2:cols-1 { extend_img->at<double>(i, j) = extend_img->at<double>(i, j) + coe*( extend_img->at<double>(i - 1, j + dir_mat[(i - 1)*cols + j - 1]) + extend_img->at<double>(i + 1, j - dir_mat[(i - 1)*cols + j - 1]) ); } } for (int i = 0; i < rows; i++) { for (int j = 0; j < cols ; j++) { img->at<double>(i, j) = extend_img->at<double>(i + 1, j + 1); } } delete extend_img; } if (step == 2)//更新 { for (int i = 2; i < extend_rows - 1; i = i + 2)//剔除原圖第一行,再2:2:rows-1 { for (int j = 1; j < extend_cols - 1; j++)//2:cols-1 { extend_img->at<double>(i, j) = extend_img->at<double>(i, j) + coe*(extend_img->at<double>(i - 1, j + dir_mat[(i - 1)*cols + j - 1]) + extend_img->at<double>(i + 1, j - dir_mat[(i - 1)*cols + j - 1])); } } for (int i = 0; i < rows; i++) { for (int j = 0; j < cols; j++) { img->at<double>(i, j) = extend_img->at<double>(i + 1, j + 1); //outFile << img->at<double>(i, j) << endl; } } delete extend_img; } } if (flag == 2)//行變換 { Mat* img_T = new Mat(); img_T->create(cols, rows, CV_64FC1); int* dir_mat_T = new int[rows*cols]; for (int i = 0; i < cols; i++) { for (int j = 0; j < rows; j++) { img_T->at<double>(i, j) = img->at<double>(j, i); } } ////cvTranspose這一opencv內(nèi)置函數(shù)用于矩陣的轉(zhuǎn)置 matTrans(dir_mat, dir_mat_T, rows, cols); ////預(yù)處理,圖像對稱延拓一圈 Mat *extend_img = new Mat(); extend_img->create(rows + 1, cols + 1, CV_64FC1); copyMakeBorder(*img_T, *extend_img, 1, 1, 1, 1, BORDER_REFLECT_101); int extend_rows = extend_img->size().height; int extend_cols = extend_img->size().width; if (step == 1) { for (int i = 1; i < extend_rows - 2; i = i + 2)//剔除原圖最后一行,再2:2:rows-1 { for (int j = 1; j < extend_cols - 1; j++)//2:cols-1 { extend_img->at<double>(i, j) = extend_img->at<double>(i, j) + coe*(extend_img->at<double>(i - 1, j + dir_mat_T[(i - 1)*rows + j - 1]) + extend_img->at<double>(i + 1, j - dir_mat_T[(i - 1)*rows + j - 1])); //此時的dir_mat_T是cols行,rows列 } } for (int i = 0; i < cols; i++) { for (int j = 0; j < rows; j++) { img_T->at<double>(i, j) = extend_img->at<double>(i + 1, j + 1); } } //將圖像轉(zhuǎn)置回去 //cvTranspose(img_T, img); for (int i = 0; i < rows; i++) { for (int j = 0; j < cols; j++) { img->at<double>(i, j) = img_T->at<double>(j, i); } } delete extend_img; delete img_T; delete[] dir_mat_T; } if (step == 2) { for (int i = 2; i < extend_rows - 1; i = i + 2)//剔除原圖第一行,再2:2:rows-1 { for (int j = 1; j < extend_cols - 1; j++)//2:cols-1 { //extend_img[i*extend_cols + j] = extend_img[i*extend_cols + j] extend_img->at<double>(i, j) = extend_img->at<double>(i, j) + coe*(extend_img->at<double>(i - 1, j + dir_mat_T[(i - 1)*rows + j - 1]) + extend_img->at<double>(i + 1, j - dir_mat_T[(i - 1)*rows + j - 1])); } } for (int i = 0; i < cols; i++) { for (int j = 0; j < rows; j++) { img_T->at<double>(i, j) = extend_img->at<double>(i + 1, j + 1); } } //將圖像轉(zhuǎn)置回去 //cvTranspose(img_T, img); for (int i = 0; i < rows; i++) { for (int j = 0; j < cols; j++) { img->at<double>(i, j) = img_T->at<double>(j, i); } } delete extend_img; delete img_T; delete[] dir_mat_T; } }}新聞熱點
疑難解答