Imports BasicDSP
Public Class dsp3ex2
    Dim wv As New Signal
    Dim fx As Waveform
    Const WINTIME As Double = 0.03
    Const STPTIME As Double = 0.01
    ' Compute autocorrelation up to len samples (len < wv.Count)
    Function Autocorrelation(ByRef wv As Waveform, ByVal len As Integer) As WaveDouble
        Dim a As New WaveDouble(len + 1, wv.Rate, 0)
        ' calculate normalising factor
        Dim sumsq As Double = 0
        For i As Integer = wv.First To wv.Last
            sumsq += wv(i) * wv(i)
        Next
        ' for each delay, calculate correlation
        For i As Integer = 0 To len
            Dim sum As Double = 0
            For j As Integer = wv.First + i To wv.Last
                sum += wv(j - i) * wv(j)
            Next
            a(i) = sum / sumsq
        Next
        Return a
    End Function
    Public Function FindMaxIndex(ByRef c As WaveDouble, ByVal first As Integer, ByVal last As Integer) As Integer
        Dim max As Double = 0
        Dim idx As Integer = 0
        For i As Integer = first To last
            If c(i) > max Then
                max = c(i)
                idx = i
            End If
        Next
        Return idx
    End Function
    Public Sub ProcessSignal()
        Dim wsize As Integer = WINTIME / wv.Period
        Dim ssize As Integer = STPTIME / wv.Period
        Dim ms2 As Integer = 0.002 / wv.Period
        Dim ms20 As Integer = 0.02 / wv.Period

        fx = New Waveform(0, 1 / STPTIME)
        For pos As Integer = 1 To wv.Count - wsize Step ssize
            Dim win As Waveform = wv.Cut(pos, wsize).Float

            ' calculate the autocorrelation
            Dim a As WaveDouble = Autocorrelation(win, ms20)

            ' find index of maximum magnitude
            Dim idx As Integer = FindMaxIndex(a, ms2, ms20)

            ' add fundamental to end of Fx graph
            If a(idx) > 0.5 Then
                fx.Add(1.0 / (idx * wv.Period))
            Else
                fx.Add(0)
            End If
        Next
    End Sub
    Public Sub DisplayGraphs(ByVal fname As String)
        Dim gr As New Graph(Me.CreateGraphics, zgc, 2, 1)
        gr.PlotSignal(1, wv, "Input from " & fname)
        gr.PlotWaveDouble(2, fx, "Fx Track", "Freq (Hz)", "Time (s)")
        zgc.Refresh()
    End Sub
    Private Sub OpenToolStripMenuItem_Click(ByVal sender As System.Object, ByVal e As System.EventArgs) Handles OpenToolStripMenuItem.Click
        If OpenFileDialog1.ShowDialog() = Windows.Forms.DialogResult.OK Then
            If (wv.LoadWaveFile(OpenFileDialog1.FileName)) Then
                ProcessSignal()
                DisplayGraphs(OpenFileDialog1.FileName)
            End If
        End If
    End Sub

End Class