Module 5. Savitzky-Golay smoothing - (2)특정 구간만 다항식 fitting 하기

Savitzky-Golay smoothing은 한 점의 전/후 일정 구간의 데이터를 이용하여 다항식으로 fitting하고, fitting식으로부터 다시 근사값을 구하여 해당 점의 값을 대신하게 됩니다.

따라서, 전체 데이터 구간으로 다항식 fitting하는 것이 아니라, 데이터의 일부 구간에 대해서만 fitting하는 함수로 변경해야합니다.

또한, fitting 결과로 반환되는 다항식의 계수를 이용하여 임의의 X값이 주어질 때 근사값을 구하는 함수가 필요하며, Savitzky-Golay smoothing의 특징인 해당 점에서 n차 미분계수를 구하는 함수를 작성할 수도 있습니다.

이번 글에서는 지난번 작성했던 PolynomialFit 함수를 일부 수정하여, 일부 구간의 데이터만으로 다항식 fitting하는 함수로 수정하도록 하겠습니다.

'---------------------------------------------------
Function PolynomialFit(iX(), iY(), iOrder As Long, Optional iSN As Long = -1, Optional iFN As Long = -1) As Variant()

  '전체 구간 데이터를 이용하는 다항식 Fitting 함수를 약간만 수정하도록 하되, 기존의 매크로를 수정하지 않고 사용할 수 있도록 일부 코드만 변경합니다.

  'X, Y 데이터 및 다항식의 차수가 입력하고, 만약 일부 구간에 대해서만 fitting하겠다면, 옵션으로 시작 및 끝 데이터 구간을 입력하도록 합니다.

  Dim i As Long, j As Long, k As Long, tMatX(), tMat(), tX(), tY(), tA()
  Dim tVal As Double, tAve As Double, n As Long, tN As Long, tSN As Long, tFN As Long

  '만약 시작 및 끝점이 음수로 입력이 된다면, 데이터 전체에 대해 fitting하고, 구간이 양수로 입력이 된다면, 시작과 끝의 데이터 구간을 지정하도록 합니다. 이 부분을 옵션값으로 입력하고, 기본값을 음수로 지정해두면, 기존에 작성한 다른 코드는 변경하지 않고 그대로 사용할 수 있습니다.

  If iSN < 0 Then tSN = LBound(iX, 1) Else tSN = iSN
  If iFN < 0 Then tFN = UBound(iX, 1) Else tFN = iFN

  '데이터 시작과 끝을 변경할 때, 영향을 받는 부분은 iX 및 iY 데이터로부터 계산을 수행하는 부분이 됩니다. 따라서, 아래의 일부 코드만 변경하도록 합니다. 색상은 빨간색으로 표시해두었습니다.

  n = tFN - tSN + 1
  tN = iOrder + 1
  tAve = 0
  For i = tSN To tFN
    tAve = tAve + iX(i, 1)
  Next
  tAve = tAve / n
  
  ReDim tX(1 To n, 1 To 1)
  For i = tSN To tFN
    tX(i - tSN + 1, 1) = iX(i, 1) - tAve
  Next
  
  ReDim tMatX(1 To tN, 1 To n)
  For j = 1 To n
    tMatX(1, j) = 1
    For i = 2 To tN
      tMatX(i, j) = tMatX(i - 1, j) * tX(j, 1)
    Next
  Next
  
  ReDim tMat(1 To tN, 1 To tN)
  For i = 1 To tN
    For j = 1 To tN
      tVal = 0
      For k = 1 To n
        tVal = tVal + tMatX(i, k) * tMatX(j, k)
      Next
      tMat(i, j) = tVal
    Next
  Next
  
  ReDim tY(1 To tN, 1 To 1)
  For i = 1 To tN
    tVal = 0
    For j = tSN To tFN
      tVal = tVal + tMatX(i, j - tSN + 1) * iY(j, 1)
    Next
    tY(i, 1) = tVal
  Next
  
  If Solve_GaussJordan(tMat, tY, tA) Then
    ReDim tMat(1 To tN, 1 To tN)
    For i = 1 To tN
      tMat(1, i) = 1
      tMat(i, i) = 1
    Next
    For i = 2 To tN
      For j = i + 1 To tN
        tMat(i, j) = tMat(i - 1, j - 1) + tMat(i, j - 1)
      Next
    Next
    tVal = 1
    For j = 2 To tN
      tVal = tVal * (-tAve)
      For i = 1 To tN - j + 1
        tMat(i, i + j - 1) = tMat(i, i + j - 1) * tVal
      Next
    Next
    ReDim tY(1 To tN, 1 To 1)
    For i = 1 To tN
      For j = 1 To tN
        tY(i, 1) = tY(i, 1) + tMat(i, j) * tA(j, 1)
      Next
    Next
  Else
    ReDim tY(1 To tN, 1 To 1)
  End If
  
  PolynomialFit = tY
  Erase tMatX, tMat, tX, tY, tA
End Function
'---------------------------------------------------

나머지 코드 부분은 기존에 작성된 PolynomailFit 함수와 동일하기 때문에 별다른 설명은 생략하도록 하겠습니다.