SignalUtils.Integrate Method

void Integrate(TVec Src, TVec Dst, ref TIntegrateState State, Double Dt)

Integrate signal.

#	Name	Type	Description
1	Src	TVec	source TVec
2	Dst	TVec	source TVec
3	State	TIntegrateState (ref)
4	Dt	Double	scalar

Result: stored in self (calling object)

Remarks:

Use Simpson's formula to integrate the Src and place the result in Dst. Src can be real or complex. dT defines the sampling period (dT = 1/FS). The State variable holds the initial conditions.

Dst[j] = 1/6sum_(i=0)^j(Src[i-2] + 4 Src[i-1] + Src[i]) dT, j = 0,1,...,n-1

The integrate routine does not preserve linear phase. An analytical solution for a linear phase integrator does not exists. In general there are two types of applications for integration: 1. The signal has a more or less fixed mean value. (DC offset). An example of such a signal is the signal comming from the accelerometer measuring vibrations. The numerical integration will work well on signals whose mean value is exactly zero (otherwise it will rise or fall in infinity). Because this is often not the case, the signal must be passed through a DC filter first. The DC filter can be IIR or FIR type. An alternative to numerical integration and the Integrate routine is a linear phase integration filter designed with the Dew.Signal.Units.OptimalFir.RemezImpulse routine. Linear phase integrator also removes the DC offset. 2. The signal does not have a mean value. An example of such a signal is the signal comming from the accelerometer measuring the acceleration and deceleration of a driving car. In this case the Integrator routine can be used directly to obtain speed and/or distance from the acceleration data.

Examples

using Dew.Math;
using Dew.Math.Units;
using Dew.Signal;
using Dew.Signal.Units;
using Dew.Math.Tee;
using Dew.Math.Editors;
using Dew.Signal.Tee;

private void button1_Click(object sender, EventArgs e)
{
    int n,i;
    double FS = 1; //sampling frequency
    TIntegrateState State2 = new TIntegrateState();
    TDiffState State = new TDiffState();

    Vector h = MtxExpr.Ramp(30, TMtxFloatPrecision.mvDouble, 0, 1);
    Vector b = new Vector(h.Length);
    Vector c = new Vector(h.Length);

//single block

    SignalUtils.Integrate(h, b, ref State2, 1.0/FS);
    SignalUtils.Differentiate(b, c, ref State, 1.0 / FS);
    MtxVecTee.DrawIt(c,"Processed in one block",false);

//streaming

    //reset initial conditions
    State2 = new TIntegrateState();
    State = new TDiffState();

    n = h.Length / 10; //integer division (!)
    for (i = 0; i < 10; i++)
    {
        h.SetSubRange(i*n,n);
        b.SetSubRange(i*n,n);
        c.SetSubRange(i*n,n);
        SignalUtils.Integrate(h,b, ref State2, 1.0/FS);

//            Should be:  b = [0 , 1 , 3, 6, 10, 15, 21,... ]
//            But becomes: b = [0, 0.1666, 1.1666, 3.166, 6.166, 10.166, 15.166, 21.166,... ]
//            because of Simpson

    SignalUtils.Differentiate(b,c, ref State,1.0/FS);

//            Should be: c = [0,1,2,3,4,5,6....]
//            But becomes: c = [0, 0.08333, 0.5833, 1.5, 2.5, 3.5, 4.5....]

    }
    c.SetFullRange();
    MtxVecTee.DrawIt(c,"Processed per partes, but same result",false);
}