Signal generation

Introduction

If you want to test your system, you always need an input signal. In the tools.signalgenerator you can find some basic signals like sine wave, cosine wave or square wave etc.

Sine wave

The sine wave or sinusoid is a mathematical function that describes a smooth repetitive oscillation. $$y(t) = A \cdot \sin(\omega t + \phi) + D$$
Where: 

• A, the amplitude, is the peak deviation of the function from its center position.
• \(\omega\), the angular frequency
• \(\phi\),the phase, specifies where in its cycle the oscillation begins at t = 0.
• D, a non-zero center amplitude

The sinwave class contains all these parameters. Also there are two extra parameters: the sampling time(\(T_{s}\)) and the end time(mtime)

class sinwave:
    """
    Generating a sinus vawe.
    x(t)= amp * sin(2 * pi * F * t + phase) + bias
    """
    def __init__(self, amp = 220, freq = 50, Ts = 0.01, mtime = 0.01, bias = 0, phase = 0):
        self.amp = amp
        self.freq = freq
        self.Ts = Ts
        self.mtime = mtime
        self.bias = bias
        self.phase = phase
        self.stime = arange(0, self.mtime, self.Ts)

    def signal(self, stime = None, phase = None, freq = None, bias = None, amp = None):
        if stime == None:
             stime = self.stime
        if phase == None:
             phase = self.phase;
        if freq == None:
             freq = self.freq;
        if bias == None:
             bias = self.bias;
        if amp == None:
             amp = self.amp;
        x = 2 * pi * freq * stime
        return amp * sin(x + phase) + bias;

If you want to use this class, first you have to create a new class, also there are three options to define the parameters, while you are creating the class,

sinW = sinwave( Ts = Ts, mtime = 1, freq = 5, amp = 6000 )

when the signal function is called,

pw = sinW.signal( freq = 5, amp = 3000 )

or directly you can reach the parameters

sinW.app = 400

sinW = sinwave( Ts = Ts, mtime = 1, freq = 5, amp = 6000 )
sinW.amp = 400
pw = sinW.signal()

The result of the following code above is:

Square wave

A square wave is a kind of non-sinusoidal waveform, an ideal square wave alternates regularly and instantaneously between two levels. Using Fourier expansion with cycle frequency f over time t, we can write an ideal square wave as an infinite series of the form: $$ x_{square}(t) = \frac{4}{\pi} \sum_{k=1}^\infty {\sin{\left (2\pi (2k-1) ft \right )}\over(2k-1)}$$$$= \frac{4}{\pi}\left (\sin(2\pi ft) + {1\over3}\sin(6\pi ft) + {1\over5}\sin(10\pi ft) + \cdots\right )$$ Here is the square wave implementation in python.

class squarewave: 
    """
    Generating a square signal using the signum of the sinus function. 
    x(t) = amp * sign( sin( 2 * pi * F * t + phase ) + bias )
    where 2*pi*f is the angular frequency
    """
    def __init__(self, amp = 220, freq = 50, Ts = 0.01, mtime = 0.5, bias = 0, phase = 0):  
        self.amp   =  amp
        self.freq  = freq 
        self.Ts   = Ts
        self.mtime  = mtime
        self.bias  = bias
        self.phase  = phase
        self.stime  = arange(0, self.mtime , self.Ts)  
  
    def signal(self, stime = None, phase = None, freq = None, bias = None, amp = None):  
        from numpy import sign  
        if stime == None:
            stime  =  self.stime
        if phase == None:
            phase  =  self.phase;
        if freq == None:
            freq  =  self.freq;
        if bias == None:
            bias  =  self.bias;
        if amp == None:
           amp = self.amp;   
        x  = 2 * pi * freq * stime 
        return amp * sign( sin( x + phase ) ) + bias;

So, here is a full example how does it possible to create a square wave:

if __name__ == "__main__":
    Ts = 0.001
    sinW = sinwave( Ts = Ts, mtime = 1, freq = 5, amp = 60 )   
    pw = sinW.signal()    
    
    squareW =   squarewave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    sw  = squareW.signal()   
    
    pylab.figure( 1 )
    pylab.plot( pw )
    pylab.plot( sw )
    
    pylab.show()

Sawtooth wave

class sawtoothwave:
    """
    Generating a triangle vawe.
    x(t)= amp * 2 / pi *sum(( -1 )^i * (sin( 2 * i * pi * F * t + phase) + bias) / i )
    where 2*pi*f is the angular frequency
    """
    def __init__(self, amp = 220, freq = 50, Ts = 0.01, mtime = 5, bias = 0, phase = 0, harmonic = 50):
        self.amp   =  amp
        self.freq  = freq 
        self.Ts   = Ts
        self.mtime  = mtime
        self.bias  = bias
        self.phase  = phase
        self.harmonic = harmonic
        self.stime   =  arange(0, self.mtime, self.Ts)   
  
    def signal(self, stime = None, phase = None, freq = None, bias = None, amp = None, harmonic = None):    
        if stime == None:
            stime  =  self.stime
        if phase == None:
            phase  =  self.phase;
        if freq == None:
            freq  =  self.freq;
        if bias == None:
            bias  =  self.bias;
        if amp == None:
            amp  =  self.amp; 
        if harmonic == None:
            harmonic=  self.harmonic;
     
        x  = 2 * pi * freq * stime 
        res = zeros(len(x))
        for i in range( 1, harmonic + 1):   
            res +=  pow( -1, i+1 ) * sin( i * (x + phase) ) / i;   
   
        return amp * (2. / pi) * res + bias

if __name__ == "__main__":
    Ts = 0.001
    sinW = sinwave( Ts = Ts, mtime = 1, freq = 5, amp = 60 )   
    pw = sinW.signal()    
    
    squareW =   squarewave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    sw  = squareW.signal()
    
    sawtoothW =   sawtoothwave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    stw  = sawtoothW.signal()
    
    
    pylab.figure( 1 )
    pylab.plot( pw )
    pylab.plot( sw )
    pylab.plot( stw )
    
    
    pylab.show()

Triangle wave

class trianglewave:
    """
    Generating a triangle vawe.
    x(t)= amp * 8 / pi^2 *sum(( -1 )^i * (sin((2 * i + 1) * 2 * pi * F * t + phase) + bias) / ( 2 * i + 1 )^2 )
    """
    def __init__(self, amp = 220, freq = 50, Ts = 0.01, mtime = 5, bias = 0, phase = 0, harmonic = 50):
        self.amp   =  amp
        self.freq  = freq 
        self.Ts   = Ts
        self.mtime  = mtime
        self.bias  = bias
        self.phase  = phase
        self.harmonic = harmonic
        self.stime   =  arange(0, self.mtime, self.Ts)   
  
    def signal(self, stime = None, phase = None, freq = None, bias = None, amp = None, harmonic = None):    
        if stime == None:
             stime  =  self.stime
        if phase == None:
             phase  =  self.phase;
        if freq == None:
             freq  =  self.freq;
        if bias == None: 
             bias  =  self.bias;
        if amp == None:
             amp  =  self.amp; 
        if harmonic == None:
             harmonic=  self.harmonic;
     
        x  = 2 * pi * freq * stime 
        res = zeros(len(x))
        for i in range( 0, harmonic ):   
            res +=  pow( -1, i ) * sin(( 2 * i + 1 ) * (x + phase) ) / pow( ( 2. * i + 1.), 2 );   
   
       return amp * (8. / pow(pi,2)) * res + bias

if __name__ == "__main__":
    Ts = 0.001
    sinW = sinwave( Ts = Ts, mtime = 1, freq = 5, amp = 60 )   
    pw = sinW.signal()    
    
    squareW =   squarewave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    sw  = squareW.signal()
    
    sawtoothW =   sawtoothwave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    stw  = sawtoothW.signal()
    
    triangleW =   trianglewave( Ts = Ts, mtime = 1, freq = 5, amp = 60 ) 
    tw  = triangleW.signal()   
    
    pylab.figure( 1 )
    pylab.plot( pw )
    pylab.plot( sw )
    pylab.plot( stw )
    pylab.plot( tw )
    
    pylab.show()

Continous to discrete system

In a real environment we are using discretized system,  so because of this we need to discretize our system. This page is contains a short description about the zoh and an example written in Python using numpy's.

Discretization with zero order hold element on the input

It's given the following continuous-time state space model:
$$\dot{x(t)} = Ax(t) + B u(t)$$$$y(t) = Cx(t) + Du(t)$$
It can be shown that the following discrete-time state space model expresses how the state x evolves along a discrete time axis $$x[k+1] = A_{d}x[k] + B_{d}u[k]$$$$y[k] = C_{d}x[k] + D_{d}u[k]$$
Where: $$A_{d} = e^{ATs}$$ $$B_{d} = \left( \int_{\tau=0}^{Ts}e^{A \tau}d\tau \right)B$$ $$C_{d} = C$$ $$D_{d} = D$$Where \(T_{s}\) is the sampling time
There is a trick how can we calculate the \(A_{d}, B_{d}\) discrete matrices: $$M = \left[ {\begin{array}{cc} A & B\\ 0 & 0 \end{array} } \right]$$ $$\left[ {\begin{array}{cc} M_{11} & M_{12}\\ 0 & I \end{array} } \right] = e^{MT_{s}}$$ After the calculation we can get easily the \(A_{d}, B_{d}\) discrete matrices: $$A_{d} = M_{11}, B_{d} = M_{21}$$ The following code written in Python contains the theory above.

from numpy import hstack , vstack, zeros
from scipy.linalg import expm
from control.matlab import tf2ss, ss2tf

def c2d ( sys , Ts, method = 'ZOH' ) :
    """
    @summary: From continuous time convert discrete time
 
    @param sys:  an instance of the LTI class or a tuple describing the system.
    The following gives the number of elements in the tuple and
    the interpretation:

    2: (num, den)        
    4: (A, B, C, D)
    @param Ts: Sampling time
    @param method: The name of the discretization method
    
    @return: The discrete system    
    """
    if(len(sys) == 2):
        (A, B, C, D) = tf2ss( sys )
    else:
        (A, B, C, D) = sys
 
    n = A.shape[ 0 ]
    nb = B.shape[ 1 ]
    if(method == 'ZOH'):
        ztmp= zeros ( ( nb , n + nb ) ) 
        tmp = hstack ( ( A , B ) ) 
        tmp = vstack ( ( tmp , ztmp ) ) 
        tmp = expm( tmp * Ts )
 
        A = tmp [ 0 : n , 0 : n ]
        B = tmp [ 0 : n , n : n+nb ]
  
  
    sysd = ( A , B , C , D  ) 
    if(len(sys) == 2):
       return ss2tf ( sysd )
    return sysd

Also it's possible to download from the github.

Brushed DC motor

Introduction

If you have never heard about Brushed DC motor, then the Wikipedia is a really good start.

A short description about the direct current motor: when a rectangular coil carrying current is placed in a magnetic field, a torque acts on the coil which rotates it continuously. When the coil rotates, the shaft attached to it also rotates and thus it is able to do mechanical work.

When the coil is powered, a magnetic field is generated around the armature. The left side of the armature is pushed away from the left magnet and drawn towards the right, causing rotation.

When the coil turns through 90°, the brushes lose contact with the commutator and the current stops flowing through the coil. However the coil keeps turning because of its own momentum. Now when the coil turns through 180°, the sides get interchanged 

 Motor Model

The following model is a equivalent-circuit representation of the brushed DC motor.

 In the equivalent-circuit above the electrical system equation is the following - using Kirchhoff's voltage and current law:
$$V_{a}=iR_{a}+L_{a}\frac{di}{dt}+e$$ $$e=K_{e}\omega_{m}$$  Mechanical system equation:
$$T_{l} = -J\frac{d\omega_{m}}{dt}-b\omega_{m}+K_{m}i_{a}$$ From the equation above we can write the following equations:
$$\frac{di}{dt} = - \frac{R_{a}}{L_{a}}i_{a} - \frac{K_{e}}{L_{a}}\omega_{m} + \frac{1}{L_{a}}V_{a}$$ $$\frac{d\omega_{m}}{dt} = \frac{K_{m}}{J}i_{a} - \frac{b}{J}\omega_{m} - \frac{1}{J}T_{l}$$ 
Where:
\(J\): Motor and load Inertia
\(b\) : Viscous damping coefficient
\(K_{e}\): Speed constant
\(K_{m}\): Torque constant
\(L_{a}\): Motor armature coil inductance
\(R_{a}\): Motor armature coil resistance

State Space Method

The state space method takes from the
$$\dot{x(t)} = Ax(t) + B u(t)$$$$y(t) = Cx(t) + Du(t)$$
Where the state-space matrices are the following:
$$\dot{x(t)} =\left[ {\begin{array}{cc} \frac{di}{dt} \\ \frac{d\omega}{dt} \\ \frac{d\theta}{dt} \end{array} } \right] A =\left[ {\begin{array}{cc} -\frac{R_{a}}{L_{a}} & -\frac{K_{e}}{L_{a}} &  0 \\ \frac{K_{m}}{J} & -\frac{b}{J} & 0 \\ 0 & 1 & 0 \end{array} } \right] B=\left[ {\begin{array}{cc} \frac{1}{L_{a}} \\ 0 \\ 0 \end{array} } \right] C=\left[ {\begin{array}{cc} 0 \\ 1 \\ 0 \end{array} } \right] D=\left[ {\begin{array}{cc} 0 \end{array} } \right]$$