$%2017/06/12
\newcommand\compl[1]{{#1^\mathtt{c}}}
\newcommand\pare[1]{{(#1)}}
\newcommand\Pare[1]{\left(#1\right)}
\newcommand\curl[1]{\{#1\}}
\newcommand\Curl[1]{\left\{#1\right\}}
\newcommand\squa[1]{[#1]}
\newcommand\Squa[1]{\left[#1\right]}
\newcommand\abs[1]{\lvert#1\rvert}
\newcommand\Abs[1]{\left\lvert#1\right\rvert}
\newcommand\floor[1]{\lfloor#1\rfloor}
\newcommand\Floor[1]{\left\lfloor#1\right\rfloor}
\newcommand\ceil[1]{\lceil#1\rceil}
\newcommand\Ceil[1]{\left\lceil#1\right\rceil}
\newcommand\setN[0]{\mathbb{N}}
\newcommand\setZ[0]{\mathbb{Z}}
\newcommand\setQ[0]{\mathbb{Q}}
\newcommand\setR[0]{\mathbb{R}}
\newcommand\setC[0]{\mathbb{C}}
\newcommand\sfrac[2]{#1/#2}
\newcommand\od[2]{\frac{d#1}{d#2}}
\newcommand\pd[2]{\frac{\partial#1}{\partial#2}}
\newcommand\sod[2]{\sfrac{d#1}{d#2}}
\newcommand\spd[2]{\sfrac{\partial#1}{\partial#2}}$
$\newcommand\L{\mathcal{L}}
\newcommand\R{\mathcal{R}}
\newcommand\A{\mathscr{A}}
\newcommand\B{\mathscr{B}}
\newcommand\G{\mathscr{G}}
\newcommand\M{\mathbb{M}}
\newcommand\I{\mathscr{I}}
\newcommand\N{\mathscr{N}}
\newcommand\setOR[0]{\overline{\mathbb{R}}}
\DeclareMathOperator{\vol}{vol}
\DeclareMathOperator\supp{supp}
\DeclareMathOperator\sgn{sgn}$
# 倒立伸子のオブザーバーメモ
$F_H$: 振子が台車から受ける水平抗力
$F_V$: 振子が台車から受ける垂直抗力


## 書き下し
\begin{align}
\dot{x}=
\begin{bmatrix}
\dot{r}\\
\dot{\theta}\\
\frac{1}{(M+m)(J+m\ell^2)-m^2\ell^2\cos^2\theta}
\begin{bmatrix}
J+m\ell^2 & -m\ell\cos\theta\\
-m\ell\cos\theta & M+m
\end{bmatrix}
\begin{bmatrix}
-f\dot{r}+m\ell\dot{\theta}^2\sin\theta+au\\
mg\ell\sin\theta -c\dot{\theta}
\end{bmatrix}
\end{bmatrix}
\end{align}

同一次元オブザーバ：
\begin{align}
\dot{z}&=Az+Bu+G(y-Cz)\\
w&=z
\end{align}

非線形システムで考えると，
\begin{align}
\dot{z}&=f(z,u)+G(y-Cz)\\
w&=z
\end{align}

\begin{align}
C&=
\begin{bmatrix}
1&0&0&0\\
0&1&0&0
\end{bmatrix}\\
G&=
\begin{bmatrix}
*&*\\
*&*\\
*&*\\
*&*
\end{bmatrix}
\end{align}

$A-GC$が安定行列となるように設計

$M^a(q)\ddot{q}^a = \tau^a+\tau^a_c$

## 振子単体の運動方程式（真上を原点とする座標）
台車が固定されている状況下における振子単体の運動方程式を考える．

図左のように，垂直方向からの棒の移動量を$\theta$とする．
棒の質量を$m$,慣性モーメントを$J$，棒の端点から重心までの長さを$\ell$，重力加速度を$g$とすると，運動方程式は以下で得られる．

\begin{align}
(J+m\ell^2)\ddot{\theta}=mg\ell\sin\theta-c\dot{\theta}
\end{align}

ただし，各パラメータの概算値は以下の通りである．
\begin{align}
m&=0.03\\
J&=0.00027\\
\ell&=0.16\\
g&=9.8\\
c&=0.00008
\end{align}

状態量として$x=[\theta,\dot{\theta}]^T$をとると，
\begin{align}
\dot{x}&=Ax+Bu\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1\\0&-\frac{c}{J+m\ell^2}
\end{bmatrix},
B=
\begin{bmatrix}
0\\
1
\end{bmatrix},
C=\begin{bmatrix}
1 & 0
\end{bmatrix}
\\
u&=\frac{mg\ell}{J+m\ell^2}\sin\theta
\end{align}

が得られる．$mg\ell/(J+m\ell^2)\simeq 48.5$である．

同一次元オブザーバ：
\begin{align}
\dot{z}&=Az+Bu+G(y-Cz)\\
w&=z
\end{align}
に対し，$G=[G_1, G_2]^T$に対し，$A-GC$が安定行列になるように$G_1,G_2$を選ぼう．

$A-GC$は以下で計算される．
\begin{align}
A-GC
&=
\begin{bmatrix}
-G_1&1\\
-G_2&-\frac{c}{J+m\ell^2}
\end{bmatrix}
\end{align}

したがって，
\begin{align}
\abs{sI-(A-GC)}&=
\Abs{
\begin{bmatrix}
s+G_1&-1\\
G_2&s+\frac{c}{J+m\ell^2}
\end{bmatrix}}\\
&=(s+\frac{c}{J+m\ell^2})(s+G_1)+G_2\\
&=s^2+(G_1+\frac{c}{J+m\ell^2})s+G_2-\frac{c}{J+m\ell^2}G_1
\end{align}

であるから，A-GCの固有値は以下で計算される．

\begin{align}
s=\frac{-(G_1+\frac{c}{J+m\ell^2})\pm \sqrt{(G_1+\frac{c}{J+m\ell^2})^2-4(G_2-\frac{c}{J+m\ell^2}G_1)}}{2}
\end{align}

ここで，以下の計算が成り立つ．
\begin{align}
&(G_1+\frac{c}{J+m\ell^2})^2>4(G_2-\frac{c}{J+m\ell^2}G_1)\\
&\Leftarrow G_2 < \frac{(G_1+\frac{c}{J+m\ell^2})^2}{4}+\frac{c}{J+m\ell^2}G_1
\end{align}

したがって，$G_1>\frac{c}{J+m\ell^2}$, $\frac{c}{J+m\ell^2}G_1<G_2<(G_1+\frac{c}{J+m\ell^2})^2/4+\frac{c}{J+m\ell^2}G_1$であれば，$A-GC$は安定行列である．

## 土台の方が十分重い場合（真上を原点とする座標）
土台の方が十分重い場合，振子が土台に与える影響は無視できると仮定する．

土台の位置を$r$とする．振幅を$w$,周波数を$f$として，
$r(t)=w\sin(2\pi f t)$
である場合，
$\ddot{r}=-w(2\pi f)^2\sin(2\pi f t)$
である．
このとき，
\begin{align}
(J+m\ell^2)\ddot{\theta}
=-\ddot{r}m\ell\cos\theta
+mg\ell\sin\theta-c\dot{\theta}
\end{align}

状態量として$x=[\theta,\dot{\theta}]^T$をとると，
\begin{align}
\dot{x}&=Ax+Bu\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1\\0&-c
\end{bmatrix},
B=
\begin{bmatrix}
0\\
1
\end{bmatrix},
C=\begin{bmatrix}
1 & 0
\end{bmatrix}
\\
u&=-\frac{\ddot{r}m\ell}{J+m\ell^2}\cos\theta
+\frac{mg\ell}{J+m\ell^2}\sin\theta
\end{align}
## 土台の方が十分重い振子のオブザーバ（真下を原点とする座標）
土台の方が十分重い場合，振子が土台に与える影響は無視できると仮定する．

土台の位置を$r$とする．振幅を$w$,周波数を$f$として，
$r(t)=w\sin(2\pi f t)$
である場合，
$\ddot{r}(t)=-w(2\pi f)^2\sin(2\pi f t)$
である．
このとき，
\begin{align}
(J+m\ell^2)\ddot{\theta}
=-\ddot{r}m\ell\cos\theta
-mg\ell\sin\theta-c\dot{\theta}
\end{align}

状態量として$x=[\theta,\dot{\theta}]^T$をとると，
\begin{align}
\dot{x}&=Ax+Bu_{ff}\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1\\0&-\frac{c}{J+m\ell^2}
\end{bmatrix},
B=
\begin{bmatrix}
0\\
1
\end{bmatrix},
C=\begin{bmatrix}
1 & 0
\end{bmatrix}
\\
u_{ff}&=-\frac{\ddot{r}m\ell}{J+m\ell^2}\cos\theta
-\frac{mg\ell}{J+m\ell^2}\sin\theta
\end{align}

ただし，各パラメータの概算値は以下の通りである．
\begin{align}
m&=0.03\\
J&=0.00027\\
\ell&=0.16\\
g&=9.8\\
c&=0.00008
\end{align}

同一次元オブザーバ：
\begin{align}
\dot{z}&=Az+Bu_{ff}+G(y-Cz)\\
&=(A-GC)z+Bu_{ff}+Gy\\
w&=z
\end{align}
に対し，$G=[G_1, G_2]^T$に対し，$A-GC$が安定行列になるように$G_1,G_2$を選べば，
振子の状態量$\theta$が推定できる．

実際には，シミュレーション結果より，$G_1=100,G_2=1$程度に設定する．
同一次元オブザーバはエンコーダと同じタイミングで更新する．

$\ddot{r}$は台車に与えた力から摩擦力を減じた値を$N$,台車の重さを$M$とすると，
$\ddot{r}=N/(M+m)$で近似計算する．
残念ながら$\ddot{r}$を台車のエンコーダの差分微分で求めていたのでは十分な精度が出ないことがわかっている．

シミュレーションより，10msecの差分微分で$\dot{\theta}$を推定する場合，
これと同程度の遅れ時間で$\theta$を時間微分するには，
$\ddot{r}$が1msec程度の遅れで得られる必要がある．

#### 台車に対するオブザーバ
台車も同じ構造でオブザーバーを組むことが可能である．

$\ddot{r}$は台車に与えた力から摩擦力を減じた値を$N$,台車の重さを$M$とすると，
振り子の振れ幅の影響は無視し，
$\hat{\ddot{r}}=N/(M+m)$で近似計算によりトルクから加速度を逆算できるとする．

状態量として$x=[r,\dot{r}]^T$をとると，台車本体の運動方程式は以下で近似できる．
\begin{align}
\dot{x}&=Ax+Bu\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1\\0&-c
\end{bmatrix},
B=
\begin{bmatrix}
0\\
1
\end{bmatrix},
C=\begin{bmatrix}
1 & 0
\end{bmatrix}
\\
u&=\hat{\ddot{r}}
\end{align}

同一次元オブザーバとして以下を利用する．
\begin{align}
\dot{z}&=Az+Bu+G(y-Cz)\\
&=(A-GC)z+Bu+Gy\\
w&=z
\end{align}
同一次元オブザーバはエンコーダと同じタイミングで更新する．

#### 検証法
データとして与えられている信号$f(t)$の後方差分微分は以下で定義される．
\begin{align}
\hat{\dot{f}}(t)=
\frac{f(t)-f(t-\Delta t)}{\Delta t}
\end{align}

これは過去の情報しか利用していなため，理論上近似微分に時間遅れが発生する．
与えられた実験データから加速度を遅れなく計測できるか判断するには，未来の情報を使えばよい．

\begin{align}
\hat{\dot{f}}(t)=
\frac{f(t+\Delta t)-f(t-\Delta t)}{2\Delta t}
\end{align}

これで時間遅れのない微分信号を信号のデータ$f$から近似的に求めることができる．

再帰的に適用すると，
\begin{align}
\hat{\ddot{f}}(t)&=
\frac{\dot{f}(t+\Delta t)-\dot{f}(t-\Delta t)}{2\Delta t}\\
&=
\frac{
\frac{f(t+2\Delta t)-f(t)}{2\Delta t}
-
\frac{f(t)-f(t-2\Delta t)}{2\Delta t}
}{2\Delta t}\\
&=
\frac{
f(t+2\Delta t)-2f(t)+f(t-2\Delta t)
}{(2\Delta t)^2}\\
&=
\frac{
f(t+\delta t)-2f(t)+f(t-\delta t)
}{\delta t^2}\\
\end{align}

ただし，$\delta t=2\Delta t$とおいた．
上記評価式により，未来の情報を使いつつ，遅れのない加速度の情報を計算できる．
これによって，プーリーを動かすトルクから台車の加速度を求めるためのチューニングを行うことができる．

与えたトルク$N$の対し，
台車の加速度はパラメータ$a,b,c$を使い，おおよそ
\begin{align}
\ddot{r}=(N-\sgn(\dot{r})a-c\dot{r})/b
\end{align}
で求まることが期待される（動摩擦のみ補償するモデル）．

## インダクタ成分を考慮した台車のダイナミクス
http://lab.cntl.kyutech.ac.jp/~nishida/lecture/CntlA/no4.pdf

トルク:$N=K_Ei$
電圧:$e$
電流:$i$
台車のダンパ成分:$C_r$
トルク定数(=逆起電力定数):$K_E$
モーターのトルク:$T(t)$
モーターが動かす台車の重さ:$M$

\begin{align}
&L\dot{i}+Ri + \frac{1}{C}\int_{-\infty}^{t}i(\tau)d\tau=e
\end{align}

この時、$C$にかかる電圧は$v_a(t)$なので、
\begin{align}
&L\dot{i}+Ri + v_a(t)=e
\end{align}

ここで、モーターの逆起電力と回転数は比例するものとした時、

\begin{align}
v_a(t)=K_E\frac{d\theta}{dt}=K_E\dot{r}
\end{align}

以上をまとめると上記の電気的な微分方程式と、運動方程式より

\begin{align}
&L\dot{i}+Ri + K_E \dot{r}=e\\
&M\ddot{r}=N - C_r\dot{r}
\end{align}

$N=K_Ei$を代入して，
\begin{align}
&L\dot{i} + K_E \dot{r}+Ri=e\\
&M\ddot{r} = - C_r\dot{r} + K_Ei
\end{align}

$x=(r,\dot{r},i)$とすると，
\begin{align}
\dot{x}=
\begin{bmatrix}
\dot{r}\\
-\frac{C_r}{M}\dot{r}
+\frac{K_E}{M}i
\\
-\frac{K_E}{L}\dot{r}
-\frac{R}{L}i
+\frac{1}{L}e
\end{bmatrix}
\end{align}

整理すると，
\begin{align}
\dot{x}=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{M}&\frac{K_E}{M}\\
0&-\frac{K_E}{L}&-\frac{R}{L}
\end{bmatrix}
x
+
\begin{bmatrix}
0\\0\\
\frac{1}{L}
\end{bmatrix}
e
\end{align}

$R=3,M=0.3$程度であることがわかっている．
$C_r=0$とみなす．
すると，$L$と$K_E$の２パラメータを気合チューニングすればいいことがわかる．

おおよそ3.75mHくらいからチューニング開始．
実際は0.3mHとかそういう値かもしれない．

#### ボツ
$x=(M/K_Er,M/K_E\dot{r},i)$とすると，

\begin{align}
-\frac{K_E}{L}\dot{r}
=-\frac{K_E}{L}
\frac{K_E}{M}
\frac{M}{K_E}
\dot{r}
\end{align}

であるから，

\begin{align}
\dot{x}=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{K_E}&1\\
0&-\frac{K_E^2}{LM}&-\frac{R}{L}
\end{bmatrix}
x
+
\begin{bmatrix}
0\\0\\
\frac{1}{L}
\end{bmatrix}
e
\end{align}

## RL回路

回路方程式は

\begin{align}
E=Ri(t)+L\frac{di(t)}{dt}
\end{align}

とすることができる。
以上より

\begin{align}
\frac{di(t)}{dt}=\frac{E}{L}-\frac{Ri(t)}{L}
\end{align}

変数分離系に直して

\begin{align}
\frac{1}{i(t)-\frac{E}{R}}\cdot di(t)=-\frac{R}{L}\cdot dt
\end{align}

これをそれぞれ積分して

\begin{align}
\text{(左辺)}=\int{\frac{1}{i(t)-\frac{E}{R}}\cdot di(t)}&=\log_e\left|i(t)-\frac{E}{R}\right|+C \text{(Cは積分定数)}\\
\text{ここで、$i(t)\leqq\frac{E}{R}$より、$i(t)-\frac{E}{R}\leqq0$であるため、}\\
&=\log_e\left(\frac{E}{R}-i(t)\right)+C
\end{align}

\begin{align}
\text{(右辺)}&=-\int{\frac{R}{L}dt}\\
&=-\frac{R}{L}t+D\text{(Dは積分定数)}
\end{align}

以上より、
\begin{align}
\log_e\left(\frac{E}{R}-i(t)\right)+C=-\frac{R}{L}t+D\text{(Dは積分定数)}
\end{align}

これを変形していく、

\begin{align}
\frac{E}{R}-i(t)&=e^{-\frac{R}{L}t+D}\\
&=e^{-\frac{R}{L}}\cdot e^{D}
\end{align}

この時、初期条件が$t=0$の時$i(t)=0$だと仮定する。
\begin{align}
&\frac{E}{R}-0=e^D\cdot e^{-\frac{R}{L}\cdot 0}\\
&\frac{E}{R}=e^D\cdot1=e^D\\
&\therefore e^D=\frac{E}{R}
\end{align}

以上より、RL直列回路に流れる電流$i(t)$は

\begin{align}
i(t)=\frac{E}{R}-\frac{E}{R}\cdot e^{-\frac{R}{L}t}
\end{align}


### RL回路の伝達関数

回路方程式は

\begin{align}
E=Ri(t)+L\frac{di(t)}{dt}
\end{align}

で表され、電圧$E=V(t)$とする。
また、電流のラプラス変換した関数を$I(s)$と表す。

このとき、上式の回路方程式のラプラス変換は

\begin{align}
V(s)=RI(s)+LsI(s)
\end{align}

であり、
\begin{align}
I(s)=\frac{V(s)}{R+Ls}
\end{align}

この時、この回路の伝達関数$G(s)$は

\begin{align}
G(s)=\frac{I(s)}{V(s)}=\frac{1}{R+Ls}
\end{align}

このシステムに単位ステップ応答$X(s)$を与えたとき、このシステムの出力$Y(s)$を確認する。

\begin{align}
Y(s)&=G(s)\cdot X(s)=\frac{1}{R+Ls}\cdot \frac{1}{s}\\
&=\frac{-\frac{L}{R}}{R+Ls}+\frac{\frac{1}{R}}{s}\\
&=\frac{1}{R}\left(\frac{1}{s}-\frac{L}{R+Ls}\right)\\
&=\frac{1}{R}\left(\frac{1}{s}-\frac{1}{\frac{R}{L}+s}\right)
\end{align}

これを逆ラプラス変換を行い、

\begin{align}
i(t)=\frac{1}{R}\left(1-e^{-\frac{R}{L}t}\right)
\end{align}

と、このシステムの単位ステップ応答を求めることができた。

以下に$E[V]$のステップ信号を入力したときの出力波形を示す。


## インダクタ成分も考慮し、台車が振子の影響を受けないと仮定した台車のダイナミクス（真下を原点とする座標）

トルク:$N=K_Ti$
電圧:$e$
電流:$i$
台車のダンパ成分(粘性抵抗係数):$C_r$
トルク定数(=逆起電力定数):$K_E$
トルク定数(=電流トルク定数):$K_T$
モーターのトルク:$T(t)$
モーターが動かす台車および振子の重さ:$(M+m)$
モーターの抵抗成分:$R$
モーターのインダクタ成分:$L$
台車系の動摩擦力:$F$

ここで、台車の運動方程式は
\begin{align}
(M+m)\ddot{r}=T(t)-C_r\dot{r}-F \tag{1}
\end{align}
と表される。

この時、モーターのトルク$T(t)$は以下のように示される。
\begin{align}
T(t)=K_Ti
\end{align}
また、モーターの逆起電力と回転速度が比例するとすると、
\begin{align}
v_a(t)=K_E\dot{r}
\end{align}
と表され、一般的にこれらの$K_EとK_T$はDCモーターにおいては一致することが知られている。
https://amazingsoblin.hatenablog.com/entry/2018/10/08/031837

つまり、
\begin{align}
K_E&=K_T\\
\therefore
T(t)&=K_Ti=K_Ei \tag2
\end{align}
つまり、$(1)$式で示した運動方程式は、
\begin{align}
(M+m)\ddot{r} = K_Ei - C_r\dot{r} - F \tag3 \\
\end{align}
と表せる。

ここで、モーターの回路方程式として　キルヒホッフの法則より
\begin{align}
L\dot{i} + K_E \dot{r}+Ri=e \tag4 \\
\end{align}

が置ける

つまり、$(3),(4)$式より、
\begin{align}
\dot{i} &=\frac{e}{L}-\frac{K_E}{L}\dot{r}-\frac{Ri}{L}\\
\ddot{r} &= \frac{1}{M+m}(K_Ei-C_r\dot{r}-F)
\end{align}
と求めることができる。

ここで、状態量として
$x=[r,\dot{r},i]^T$をとったとき、
非線形状態方程式は

\begin{align}
\dot{x}=f(x,e)=
\begin{bmatrix}
\dot{r}\\
\frac{1}{M+m}(K_Ei-C_r\dot{r}-F)\\
\frac{e}{L}-\frac{K_E}{L}\dot{r}-\frac{Ri}{L}
\end{bmatrix}
\end{align}

ここで、摩擦項$F$を摩擦保障で後ほど入力に加えるとして、$F=0$の時を考える。
http://www.mech.tohoku-gakuin.ac.jp/rde/contents/course/controlII/nonlinear.html

この時、入力として
$u=[e]$
としたとき、

\begin{align}
\dot{x}=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{M+m}&\frac{K_E}{M+m}\\
0&-\frac{K_E}{L}&-\frac{R}{L}
\end{bmatrix}x+
\begin{bmatrix}
0\\0\\\frac{1}{L}
\end{bmatrix}
u
\end{align}

よって、このシステムのダイナミクスは

\begin{align}
\dot{x}&=Ax+Bu\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{M+m}&\frac{K_E}{M+m}\\
0&-\frac{K_E}{L}&-\frac{R}{L}
\end{bmatrix},
B=
\begin{bmatrix}
0\\0\\\frac{1}{L}
\end{bmatrix},
C=\begin{bmatrix}
1 & 0 & 0
\end{bmatrix}
\\
u&=e
\end{align}

と表現できる。

また、

同一次元オブザーバとして以下を利用する．
\begin{align}
\dot{z}&=Az+Bu+G(y-Cz)\\
&=(A-GC)z+Bu+Gy\\
w&=z
\end{align}
同一次元オブザーバはエンコーダと同じタイミングで更新する．

以上より、この台車のオブザーバーとして、

\begin{align}
\dot{z}&=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{M+m}&\frac{K_E}{M+m}\\
0&-\frac{K_E}{L}&-\frac{R}{L}
\end{bmatrix}z+
\begin{bmatrix}
0\\0\\\frac{1}{L}
\end{bmatrix}u
+GC(x-z)
\end{align}

が成立し、物理パラメータ$C_r,K_E,R,L$を設定する。

また、$G=[G_1,G_2,G_3]^T$に対し、$(A-GC)$が安定行列であるように行列$G$を設定する。

### s領域のシステム

ここで、摩擦項については入力として与えるときに摩擦保障を行うものとし、疑似的に$F=0$の場合を考える。

\begin{align}
\dot{i} &=\frac{e}{L}-\frac{K_E}{L}\dot{r}-\frac{Ri}{L}\\
\ddot{r} &= \frac{1}{M+m}(K_Ei-C_r\dot{r})
\end{align}
ここで、$R=R_{reg}$と表記する。

$i(t),r(t)$のラプラス変換をそれぞれ$I(s),R(s)$とした時、
ここで、$i(0) = 0,r(0)=0$とする。

\begin{align}
\mathcal{L}\left[\dot{i(t)}\right] &= sI(s)-i(0)\\
&=\frac{e}{sL}-\frac{K_E}{L}(sR(s)-R(0))-\frac{R_{reg}}{L}I(s)\\
\end{align}
\begin{align}
\therefore I(s)=\frac{e}{s^2L}-\frac{K_E}{L}R(s)-\frac{R_{reg}}{sL}I(s)\tag1
\end{align}

\begin{align}
\mathcal{L}\left[\ddot{r}\right] &= s^2R(s)-sr(0)-\dot{r}(0)\\
&=\frac{K_E}{M+m}I(s)-\frac{C_r}{M+m}(sR(s)-r(0))
\end{align}
\begin{align}
\therefore R(s)=\frac{K_E}{s^2(M+m)}I(s)-\frac{C_r}{s(M+m)}R(s)\tag2
\end{align}

入力を$I(s)$,出力を$R(s)$とした時、伝達関数$G(s)$を求める。

\begin{align}
I(s)=\frac{e}{s^2L}-\frac{K_E}{L}R(s)-\frac{R}{sL}I(s) \\
I(s)\left(1+\frac{R}{sL}\right)=
\end{align}

### 検討中のメモ
式４のモーターについての回路方程式を$\dot{r}$について整理した式を用いてラプラス変換を行って、係数を比較してみようと思った

\begin{align}
\dot{r} &= \frac{e}{K_E}-\frac{L\dot{i}}{K_E}-\frac{R_{reg}i}{K_E}\\
\end{align}

\begin{align}
\mathcal{L}\left[\dot{r}\right] &= sR(s)-r(0)\\
&=\frac{e}{sK_E}-\frac{L}{K_E}(sI(s)-i(0))-\frac{R_{reg}}{K_E}I(s)
\end{align}
\begin{align}
\therefore R(s)=\frac{e}{s^2K_E}-\frac{L}{K_E}I(s)-\frac{R_{reg}}{sK_E}I(s)
\end{align}

eはどうにかならんのか？？？？？？？
->$I$の関数になるはずなので、回路方程式以外にもう１式ほしい

## $C_rとK_E$の同定方法
物理パラメータの$R,L$については既知のものであるとする(計測可能)
> $R=3,M=0.3$程度であることがわかっている．
$C_r=0$とみなす．
すると，$L$と$K_E$の２パラメータを気合チューニングすればいいことがわかる．
おおよそ3.75mHくらいからチューニング開始．
実際は0.3mHとかそういう値かもしれない．

> 「理論的にはこの挙動のこの部分がこの式と一致しているから，この係数比較で〇〇のパラメータが求まる」という考察を考えてください．RとLがわかっていれば，未知のパラメータは大分少なくなるはず．KEとCRを同定する術を，ステップ応答を具体的にかきくだすなどして考えるといいと思います．

\begin{align}
\dot{x}&=Ax+Bu\\
y&=Cx\\
A&=
\begin{bmatrix}
0&1&0\\
0&-\frac{C_r}{M+m}&\frac{K_E}{M+m}\\
0&-\frac{K_E}{L}&-\frac{R}{L}
\end{bmatrix},
B=
\begin{bmatrix}
0\\0\\\frac{1}{L}
\end{bmatrix},
C=\begin{bmatrix}
1 & 0 & 0
\end{bmatrix}
\\
u&=e
\end{align}