Crazy_Challenger

导航

 

SG平滑算法是由Savizkg和Golag提出来的。基于最小二乘原理的多项式平滑算法,也称卷积平滑。为啥叫多项式平滑呢?且看下去。

下面使用五点平滑算法来说明平滑过程

原理很简单如图:

把光谱一段区间的等波长间隔的5个点记为X集合,多项式平滑就是利用在波长点为Xm-2,Xm-1,Xm,Xm+1,Xm+2的数据的多项式拟合值来取代Xm,,然后依次移动,直到把光谱遍历完。

代码 order必须要比滑动窗口要小,且滑动窗口不能为偶数。
void SavitskyGolaySmoothing(float *arr,int window_size, int order, int rows, int cols){
if(window_size % 2 == 0){
throw std::logic_error("only odd window size allowed");
}
if(order >= window_size){
throw std::logic_error("Order must < window_size");
}

cv::Mat A = cv::Mat::zeros(window_size, order, CV_32FC1);
cv::Mat A_T, A_INV, B;
cv::Mat kernel;
cv::Mat result = cv::Mat::zeros(window_size, 1, CV_32FC1);

int step = int((window_size - 1)/2);
for(int i = 0; i < window_size; i++){
    for(int j = 0 ; j < order; j++){
        float x = pow(-step + i, j);
        A.at<float>(i,j) = x;
    }
}
A_T = A.t();
A_INV = (A_T * A).inv();
B = A * A_INV * A_T;

B.row(step).copyTo(kernel);

float *wrap_data = new float[step*2 + cols];
for(int row =0; row < rows; row++){
    //Extend start data, size is step
    for(int n = 0; n < step; n++)
    {
        wrap_data[n] = arr[row * cols];
    }
    //Copy input data
    for(int col =0; col < cols; col++){
        wrap_data[col + step] = arr[row * cols + col];
    }
    //Extend end data, size is step
    for(int n = 0; n < step; n++){
        wrap_data[cols  + step + n] = arr[row * cols + cols -1];
    }
    for(int m = step; m < step + cols; m++){

        for(int n = -step, j = 0; n <=step; n++, j++){
            result.at<float>(0, j) = wrap_data[m + n];
        }

        arr[row * cols + m - step] = cv::Mat(kernel * result).at<float>(0 ,0);
    }
}
delete []wrap_data;

}

posted on 2019-09-04 20:36  Crazy_Challenger  阅读(3720)  评论(0)    收藏  举报