2.1.22.4.4 ocmath_savitsky_golay
Description
Savitzky-Golay smoothing filter.
Syntax
int ocmath_savitsky_golay( const double * pY, double * pYs, uint nSize, int nLeft, int nRight = -1, int nPolyDeg = 2, int nDervOrder = 0, int nPadding = EDGEPAD_EXTRAPOLATE )
Parameters
- pY
- [input] pointer to Y vector data
- pYs
- [output] pointer to smoothed data, or derivatives, depending on nDervOrder
- nSize
- [input] vector size of both pY and pYs
- nLeft
- [input] number of data points to the left to be used in filter convolution,
- nRight
- [input] number of data points to the right to be used in filter convolution, default (-1) will assume nLeft.
- Total number of points used in the polynomial fits are (nLeft + nRight + 1) and it must be odd, namely nLeft + nRight must be even.
- nPolyDeg
- [input] The polynomial order. Higher order will preserve sharper features. nPolyDeg must be less then (nLeft + nRight + 1).
- nDervOrder
- [input] order of derivative desired (0 = smoothing). To generate a fourth derivative, a minimum quartic (order 4) smoothing must be used.
- For a third derivative, a minimum cubic (order 3) smoothing is needed. Similarly, a second derivative requires a minimum quadratic (order 2) smoothing.
- nPadding
- [input] data are padded with nLeft and nRight points on both ends. Possible values are
- EDGEPAD_NONE no padding
- EDGEPAD_REFLECT pad reflect, end points are repeated such that on the left, [-1] = [0], [-2] = [1], [-3] = [2] and etc
- EDGEPAD_REPEAT pad with [0] values the left and with [nSize-1] on the right
- EDGEPAD_EXTRAPOLATE linear extrapolation using nLeft points on the left and nRight points on the right
- EDGEPAD_PERIODIC pad periodic, [-1] = [nSize -1], [-2] = [nSize - 2]
Return
OE_NOERROR for success
Examples
EX1
//Assume in the current graph, curve's XY data is in the first data plot. This piece
//of code get the XY data of the curve from the first data plot and Savitzky-Golay
//smoothing filter on it. The result is output in a new worksheet and the new curve
//will plot in the original data plot with color red.
void ocmath_savitsky_golay_ex1()
{
GraphLayer gl = Project.ActiveLayer();
if (!gl)
{
out_str("Active layer is not a graph.");
return;
}
//get XY data from the first dataplot
DataPlot dp = gl.DataPlots(0);
DataRange dr;
vector vx, vy;
if(dp.GetDataRange(dr))
{
DWORD dwPlotID;
if(dr.GetData(DRR_GET_DEPENDENT | DRR_NO_FACTORS, 0, &dwPlotID, NULL, &vy, &vx) < 0)
{
printf("get data failed GetData");
return;
}
}
vector vSmooth;
vSmooth.SetSize(vy.GetSize());
ocmath_savitsky_golay(vy, vSmooth, vy.GetSize(), 3, 3, 2, 1); // Left=Right=7, quadratic, 1st order derivative
//new a worksheet to output the result
Worksheet wks;
wks.Create("Smooth");
wks.SetSize(-1, 2);
wks.SetColDesignations("XY");
DataRange drOut;
drOut.Add("X", wks, 0, 0, -1, 0);
drOut.Add("Y", wks, 0, 1, -1, 1);
drOut.SetData(&vSmooth, &vx);
//plot the curve after smoothing filter
int nPlot = gl.AddPlot(drOut, IDM_PLOT_LINE);
dp = gl.DataPlots(nPlot);
dp.SetColor(1);
}
Remark
Savitzky-Golay smoothing filter.
This time-domain smoothing is based on least squares polynomial fitting across a moving window within the data.
This smoothing method preserves higher moments in the data, thus reducing the distortion of essential features of the data like peak heights and line widths in a spectrum. It can be used to just do data smoothing or to calculate derivatives.
Please note that this routine requires uniform spacing in the time domain for the data.
This function does not handle missing values so work with a copy of data with missing values removed if needed.
cf. Numerical Recipes (Press et al. 1992, Sec.14.8)
See Also
ocmath_adjave_smooth
Header to Include
origin.h
Reference
|