CH08_SEC07_2b_Obsv.ipynb

Kernel: Python 3

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import rcParams
from control.matlab import *
import slycot
from scipy import integrate
from scipy.linalg import schur
# Python control toolbox available at https://python-control.readthedocs.io/

plt.rcParams['figure.figsize'] = [8, 8]
plt.rcParams.update({'font.size': 18})

m = 1
M = 5
L = 2
g = -10
d = 1

b = -1 # pendulum down (b = -1)

A = np.array([[0,1,0,0],\
              [0,-d/M,b*m*g/M,0],\
              [0,0,0,1],\
              [0,-b*d/(M*L),-b*(m+M)*g/(M*L),0]])

B = np.array([0,1/M,0,b/(M*L)]).reshape((4,1))

C = np.array([0,0,1,0]) # only observable if x measured... because x can't be

print('Observability matrix:\n{}'.format(obsv(A,C)))
print('Observability matrix determinant: {}'.format(np.linalg.det(obsv(A,C))))

Observability matrix:
[[ 0.    0.    1.    0.  ]
 [ 0.    0.    0.    1.  ]
 [ 0.    0.1  -6.    0.  ]
 [ 0.   -0.02  0.2  -6.  ]]
Observability matrix determinant: 0.0

## Which measurements are best if we omit "x"
Ah = A[1:,1:]
Bh = B[1:]
# Ch = np.array([1,0,0])
Ch = np.array([0,1,0])
# Ch = np.array([0,0,1])

print('Observability matrix:\n{}'.format(obsv(Ah,Ch)))

Ch = Ch.reshape((1,len(Ch)))
Dh = np.zeros((Ch.shape[0],Bh.shape[1]))
sys = ss(Ah,Bh,Ch,Dh)
print('Gramian determinant: {}'.format(np.linalg.det(gram(sys,'o'))))

Observability matrix:
[[ 0.   1.   0. ]
 [ 0.   0.   1. ]
 [ 0.1 -6.   0. ]]
Gramian determinant: 0.031250000000001