diff --git a/docs/source/control/common.rst b/docs/source/control/common.rst
deleted file mode 100644
index 6f8b318b7..000000000
--- a/docs/source/control/common.rst
+++ /dev/null
@@ -1,12 +0,0 @@
-**************
-Common Methods
-**************
-
-Control methods common to all motor types are described in this section. This section will be extended and improved in the future. 
-
-.. toctree::
-   :titlesonly:
-   :maxdepth: 2
-
-   speed_ctrl
-   current_ctrl
diff --git a/docs/source/control/design_notes.rst b/docs/source/control/design_notes.rst
new file mode 100644
index 000000000..7683e52dd
--- /dev/null
+++ b/docs/source/control/design_notes.rst
@@ -0,0 +1,15 @@
+************
+Design Notes
+************
+
+Design notes for selected control methods are provided in this section. The aim of these notes is to link the implemented methods to the theory and to provide a reference for the implementation. Further details are available in the references provided. 
+
+.. toctree::
+   :titlesonly:
+   :maxdepth: 2
+
+   speed_ctrl
+   current_ctrl
+   Observers <observers>
+
+This section will be extended in the future.
diff --git a/docs/source/control/im.rst b/docs/source/control/im.rst
deleted file mode 100644
index 98152158e..000000000
--- a/docs/source/control/im.rst
+++ /dev/null
@@ -1,11 +0,0 @@
-*****************
-Induction Machine
-*****************
-
-Control methods for induction machines will be described here. So far only reduced-order observers are documented. This section will be extended and improved in the future.
-
-.. toctree::
-   :titlesonly:
-   :maxdepth: 2
-
-   observer
diff --git a/docs/source/control/observer.rst b/docs/source/control/observer.rst
deleted file mode 100644
index 827df897b..000000000
--- a/docs/source/control/observer.rst
+++ /dev/null
@@ -1,145 +0,0 @@
-Reduced-Order Observer
-======================
-
-In electric machine drives, an observer is commonly used to estimate the state of an electrical machine. In motion-sensored drives, the state of interest is typically the flux linkage. In sensorless induction machine drives, the rotor speed estimate (and the rotor position estimate if a synchronous machine is considered) is needed as well. 
-
-In the following, reduced-order observer designs for the induction machine are considered [#Ver1988]_, [#Har2001]_, [#Hin2010]_. Both sensored and sensorless observer variants are available in the :class:`motulator.control.im.Observer` class. Similar design principles are applied to synchronous machines as well.
-
-Induction Machine Model
------------------------
-
-The inverse-Γ model of an induction machine is considered (see :doc:`/model/machines`). In a coordinate system rotating at the angular speed :math:`\omega_\mathrm{s}`, the rotor flux dynamics can be expressed using the stator and rotor quantities, respectively, as
-
-.. math::
-	\frac{\mathrm{d} \boldsymbol{\psi}_\mathrm{R}}{\mathrm{d} t} + \mathrm{j}\omega_\mathrm{s}\boldsymbol{\psi}_\mathrm{R} &= \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} - \mathrm{j} \omega_\mathrm{s}L_\sigma\boldsymbol{i}_\mathrm{s} \\
-    &= R_\mathrm{R}\boldsymbol{i}_\mathrm{s} - \left(\alpha - \mathrm{j}\omega_\mathrm{m} \right)\boldsymbol{\psi}_\mathrm{R}
-    :label: dpsiR
-
-where :math:`\omega_\mathrm{m}` is the electrical rotor angular speed. 
-
-Observer Structure
-------------------
-
-General Coordinates
-^^^^^^^^^^^^^^^^^^^
-
-Based on :eq:`dpsiR`, a reduced-order observer can be formulated in a coordinate system rotating at :math:`\omega_\mathrm{s}`,
-
-.. math::
-    \frac{\mathrm{d} \hat{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} + \mathrm{j}\omega_\mathrm{s}\hat{\boldsymbol{\psi}}_\mathrm{R} = \boldsymbol{v} + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}) + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v})^* 
-    :label: dhatpsiR
-    
-where :math:`\boldsymbol{k}_1` and :math:`\boldsymbol{k}_2` are complex gains, the estimates are marked with the hat, and :math:`^*` marks the complex conjugate. The back-emf estimates are 
-
-.. math::
-    \boldsymbol{v} &= \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} - \mathrm{j} \omega_\mathrm{s}L_\sigma\boldsymbol{i}_\mathrm{s} \\
-	\hat{\boldsymbol{v}} &= R_\mathrm{R}\boldsymbol{i}_\mathrm{s} - \left(\alpha - \mathrm{j}\hat{\omega}_\mathrm{m} \right)\hat{\boldsymbol{\psi}}_\mathrm{R}
-    :label: hatv
-
-In sensored drives, the rotor speed estimate is replaced with the measured speed,  :math:`\hat{\omega}_\mathrm{m} = \omega_\mathrm{m}`. If needed, the estimates for the stator flux linkage and the electromagnetic torque, respectively, are
-
-.. math::
-    \hat{\boldsymbol{\psi}}_\mathrm{s} &= \hat{\boldsymbol{\psi}}_\mathrm{R} + L_\sigma \boldsymbol{i}_\mathrm{s} \\
-    \hat{\tau}_\mathrm{M} &= \frac{3n_\mathrm{p}}{2}\mathrm{Im}\left\{\boldsymbol{i}_\mathrm{s} \hat{\boldsymbol{\psi}}_\mathrm{R}^* \right\}
-    :label: tauM
-
-It is also worth noticing that the derivative of the stator current in :eq:`hatv` is integrated, i.e., the noise is not amplified.
-
-.. note::
-    Real-valued column vectors and the corresponding :math:`2\times 2` gain matrix were used in [#Hin2010]_. The complex form in :eq:`dhatpsiR` has the same degrees of freedom.
-
-Estimated Flux Coordinates
-^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-The flux estimate in stator coordinates can be expressed using the polar form, :math:`\hat{\boldsymbol{\psi}}_\mathrm{R}^\mathrm{s} = \hat{\psi}_\mathrm{R}\mathrm{e}^{\mathrm{j}\hat{\vartheta}_\mathrm{s}}`, where :math:`\hat{\psi}_\mathrm{R}` is the magnitude and :math:`\hat{\vartheta}_\mathrm{s}` is the angle of the estimate. The observer is convenient to implement in estimated rotor flux coordinates, in which the rotor flux estimate is real, :math:`\hat{\boldsymbol{\psi}}_\mathrm{R} = \hat{\psi}_\mathrm{R} + \mathrm{j}0` [#Har2001]_. Under these conditions, the angular speed of the coordinates system can be solved from :eq:`dhatpsiR` as
-
-.. math::
-    \frac{\mathrm{d}\hat{\vartheta}_\mathrm{s}}{\mathrm{d} t} = \omega_\mathrm{s}
-    = \frac{\mathrm{Im} \{ \boldsymbol{v}' + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}') + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v}')^* \} }{\hat{\psi}_\mathrm{R} + L_\sigma \mathrm{Re}\{(1 - \boldsymbol{k}_1)\boldsymbol{i}_\mathrm{s} + \boldsymbol{k}_2 \boldsymbol{i}_\mathrm{s}^* \}} 
-    :label: hatws
-
-where 
-
-.. math::
-    \boldsymbol{v}' = \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} 
-    :label: vp
-
-The flux magnitude dynamics are 
-
-.. math::
-    \frac{\mathrm{d} \hat{\psi}_\mathrm{R}}{\mathrm{d} t} 
-    = \mathrm{Re}\{ \boldsymbol{v} + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}) + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v})^* \}
-    :label: dhatpsiR_abs
-
-Notice that the right-hand side of :eq:`hatws` is independent of :math:`\omega_\mathrm{s}`. Futhermore, in these coordinates, the condition :eq:`inherently` for an inherently sensorless observer reduces to :math:`\boldsymbol{k}_2 = \boldsymbol{k}_1`. This observer structure is implemented in the :class:`motulator.control.im.Observer`, where simple forward-Euler discretization is used. 
-
-Gain Selection
---------------
-
-The estimation-error dynamics are obtained by subtracting :eq:`dhatpsiR` from :eq:`dpsiR`. The resulting system can be linearized for analysis and gain selection purposes. Using the rotor speed as an exmaple, the small-signal deviation about the operating point is :math:`\Delta \omega_\mathrm{m} = \omega_\mathrm{m} - \omega_\mathrm{m0}`, where the subscript 0 refers to the operating point. Linearization of the estimation-error dynamics leads to [#Hin2010]_
-
-.. math::
-	\frac{\mathrm{d} \Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} =  \boldsymbol{k}_1\Delta \tilde{\boldsymbol{v}} + \boldsymbol{k}_2\Delta \tilde{\boldsymbol{v}}^* - \mathrm{j}\omega_\mathrm{s0}\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}
-    :label: dtildepsiR
-
-where the estimation error of the rotor flux is :math:`\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R} = \Delta\boldsymbol{\psi}_\mathrm{R} - \Delta\hat{\boldsymbol{\psi}}_\mathrm{R}` and other estimation errors are marked similarly. Furthermore, the back-emf estimation error is
-
-.. math::
-    \Delta\tilde{\boldsymbol{v}} = -\left(\alpha - \mathrm{j}\omega_\mathrm{m0} \right)\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R} + \mathrm{j}\boldsymbol{\psi}_\mathrm{R0}\Delta\tilde{\omega}_\mathrm{m}
-    :label: dtildev
-
-Sensored Case
-^^^^^^^^^^^^^
-
-Here, :math:`\boldsymbol{k}_2 = 0` and :math:`\Delta\tilde{\omega}_\mathrm{m} = 0` are assumed, corresponding to the sensored reduced-order observer in [#Ver1988]_. Under these assumptions, the estimation-error dynamics in :eq:`dtildepsiR` reduce to 
-
-.. math::
-	\frac{\mathrm{d} \Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} =  -\left[\boldsymbol{k}_1\left(\alpha - \mathrm{j}\omega_\mathrm{m0} \right) + \mathrm{j}\omega_\mathrm{s0}\right]\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}
-    :label: dtildepsiR_sensored
-
-It can be noticed that the closed-loop pole could be arbitrarily placed via the gain :math:`\boldsymbol{k}_1`. Well-damped estimation-error dynamics can be obtained, e.g., by choosing
-
-.. math:: 
-    \boldsymbol{k}_1 = 1 + \frac{g |\omega_\mathrm{m}|}{\alpha - \mathrm{j}\omega_\mathrm{m}} 
-    :label: k1_sensored
-
-where :math:`g` is a unitless positive design parameter. The corresponding pole is located at :math:`s = -\alpha - g |\omega_\mathrm{m0}| - \mathrm{j}\omega_\mathrm{r0}`, where :math:`\omega_\mathrm{r0} = \omega_\mathrm{s0} - \omega_\mathrm{m0}` is the operating-point slip angular frequency.
-
-.. note::
-
-    As a special case, :math:`\boldsymbol{k}_1  = 0` yields the voltage model. As another special case, the current model is obtained choosing :math:`\boldsymbol{k}_1 = 1`. 
-
-Inherently Sensorless Case
-^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-For sensorless drives, choosing 
-
-.. math::
-    \boldsymbol{k}_2 = (\hat{\boldsymbol{\psi}}_\mathrm{R}/\hat{\boldsymbol{\psi}}_\mathrm{R}^*) \boldsymbol{k}_1
-    :label: inherently
-    
-yields an inherently sensorless observer, i.e., the rotor speed estimate :math:`\hat{\omega}_\mathrm{m}` cancels out from the observer equations [#Hin2010]_. Under this condition, the linearized estimation-error dynamics in :eq:`dtildepsiR` become
-
-.. math::
-	\frac{\mathrm{d}}{\mathrm{d} t} \begin{bmatrix} \Delta\tilde{\psi}_\mathrm{Rd} \\ \Delta\tilde{\psi}_\mathrm{Rq} \end{bmatrix} = \begin{bmatrix} -2k_\mathrm{d}\alpha & -2k_\mathrm{d}\omega_\mathrm{m0} + \omega_\mathrm{s0} \\ -2k_\mathrm{q}\alpha - \omega_\mathrm{s0} & -2k_\mathrm{q}\omega_\mathrm{m0} 
-     \end{bmatrix} \begin{bmatrix} \Delta\tilde{\psi}_\mathrm{Rd} \\ \Delta\tilde{\psi}_\mathrm{Rq} \end{bmatrix}
-    :label: dtildepsiR_sensorless
-
-where the gain components correspond to :math:`\boldsymbol{k}_1 = k_\mathrm{d} + \mathrm{j}k_\mathrm{q}`. It can be seen that the dynamics of the rotor speed are decoupled from the flux-estimation error dynamics. The decay rate :math:`\sigma` be assigned by choosing
-
-.. math::
-    \boldsymbol{k}_1 = \frac{\sigma}{\alpha - \mathrm{j}\hat\omega_\mathrm{m}}
-    :label: k1k2_sensorless
-
-which results in the characteristic polynomial :math:`D(s)=s^2 + 2\sigma s + \omega_\mathrm{s0}^2`. The decay rate can be selected as :math:`\sigma = \alpha/2 + \zeta_\infty|\hat{\omega}_\mathrm{m}|`, where :math:`\zeta_\infty` is the desired damping ratio at high speed. At zero stator frequency :math:`\omega_\mathrm{s0} = 0`, the poles are located at :math:`s = 0` and :math:`s = -\alpha`, which allows stable magnetizing and starting the machine.
-
-.. rubric:: References
-
-.. [#Ver1988] Verghese, Sanders, “Observers for flux estimation in induction machines,” IEEE Trans. Ind. Electron., 1988, https://doi.org/10.1109/41.3067
-
-.. [#Har2001] Harnefors, “Design and analysis of general rotor-flux-oriented vector control systems,” IEEE Trans. Ind. Electron., 2001, https://doi.org/10.1109/41.915417
-
-.. [#Hin2010] Hinkkanen, Harnefors, Luomi, "Reduced-order flux observers with stator-resistance adaptation for speed-sensorless induction motor drives," IEEE Trans. Power Electron., 2010, https://doi.org/10.1109/TPEL.2009.2039650
-
-
-
diff --git a/docs/source/control/observers.rst b/docs/source/control/observers.rst
new file mode 100644
index 000000000..163269675
--- /dev/null
+++ b/docs/source/control/observers.rst
@@ -0,0 +1,276 @@
+Observer Design
+===============
+
+An observer is commonly used to estimate the state of an electrical machine. In motion-sensored drives, the state of interest is typically the flux linkage. In sensorless drives, the rotor speed estimate is needed as well. Furthermore, sensorless synchronous machine drives require estimation of the rotor position.  
+
+Induction Machines
+------------------
+
+In the following, reduced-order observer designs for the induction machine are considered [#Ver1988]_, [#Har2001]_, [#Hin2010]_. Both sensored and sensorless observer variants are available in the :class:`motulator.control.im.Observer` class. Similar design principles are applied to synchronous machines as well.
+
+Machine Model
+^^^^^^^^^^^^^
+
+The inverse-Γ model of an induction machine is considered (see :doc:`/model/machines`). In a coordinate system rotating at the angular speed :math:`\omega_\mathrm{s}`, the rotor flux dynamics can be expressed using the stator and rotor quantities, respectively, as
+
+.. math::
+	\frac{\mathrm{d} \boldsymbol{\psi}_\mathrm{R}}{\mathrm{d} t} + \mathrm{j}\omega_\mathrm{s}\boldsymbol{\psi}_\mathrm{R} &= \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} - \mathrm{j} \omega_\mathrm{s}L_\sigma\boldsymbol{i}_\mathrm{s} \\
+    &= R_\mathrm{R}\boldsymbol{i}_\mathrm{s} - \left(\alpha - \mathrm{j}\omega_\mathrm{m} \right)\boldsymbol{\psi}_\mathrm{R}
+    :label: dpsiR
+
+where :math:`\omega_\mathrm{m}` is the electrical rotor angular speed. 
+
+Observer Structure
+^^^^^^^^^^^^^^^^^^
+
+General Coordinates
+"""""""""""""""""""
+
+Based on :eq:`dpsiR`, a reduced-order observer can be formulated in a coordinate system rotating at :math:`\omega_\mathrm{s}`,
+
+.. math::
+    \frac{\mathrm{d} \hat{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} + \mathrm{j}\omega_\mathrm{s}\hat{\boldsymbol{\psi}}_\mathrm{R} = \boldsymbol{v} + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}) + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v})^* 
+    :label: dhatpsiR
+    
+where :math:`\boldsymbol{k}_1` and :math:`\boldsymbol{k}_2` are complex gains, the estimates are marked with the hat, and :math:`^*` marks the complex conjugate. The back-emf estimates are 
+
+.. math::
+    \boldsymbol{v} &= \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} - \mathrm{j} \omega_\mathrm{s}L_\sigma\boldsymbol{i}_\mathrm{s} \\
+	\hat{\boldsymbol{v}} &= R_\mathrm{R}\boldsymbol{i}_\mathrm{s} - \left(\alpha - \mathrm{j}\hat{\omega}_\mathrm{m} \right)\hat{\boldsymbol{\psi}}_\mathrm{R}
+    :label: hatv
+
+In sensored drives, the rotor speed estimate is replaced with the measured speed,  :math:`\hat{\omega}_\mathrm{m} = \omega_\mathrm{m}`. If needed, the estimates for the stator flux linkage and the electromagnetic torque, respectively, are
+
+.. math::
+    \hat{\boldsymbol{\psi}}_\mathrm{s} &= \hat{\boldsymbol{\psi}}_\mathrm{R} + L_\sigma \boldsymbol{i}_\mathrm{s} \\
+    \hat{\tau}_\mathrm{M} &= \frac{3n_\mathrm{p}}{2}\mathrm{Im}\left\{\boldsymbol{i}_\mathrm{s} \hat{\boldsymbol{\psi}}_\mathrm{R}^* \right\}
+    :label: tauM
+
+It is also worth noticing that the derivative of the stator current in :eq:`hatv` is integrated, i.e., the noise is not amplified.
+
+.. note::
+    Real-valued column vectors and the corresponding :math:`2\times 2` gain matrix were used in [#Hin2010]_. The complex form in :eq:`dhatpsiR` has the same degrees of freedom.
+
+Estimated Flux Coordinates
+""""""""""""""""""""""""""
+
+The flux estimate in stator coordinates can be expressed using the polar form, :math:`\hat{\boldsymbol{\psi}}_\mathrm{R}^\mathrm{s} = \hat{\psi}_\mathrm{R}\mathrm{e}^{\mathrm{j}\hat{\vartheta}_\mathrm{s}}`, where :math:`\hat{\psi}_\mathrm{R}` is the magnitude and :math:`\hat{\vartheta}_\mathrm{s}` is the angle of the estimate. The observer is convenient to implement in estimated rotor flux coordinates, in which the rotor flux estimate is real, :math:`\hat{\boldsymbol{\psi}}_\mathrm{R} = \hat{\psi}_\mathrm{R} + \mathrm{j}0` [#Har2001]_. Under these conditions, the angular speed of the coordinates system can be solved from :eq:`dhatpsiR` as
+
+.. math::
+    \frac{\mathrm{d}\hat{\vartheta}_\mathrm{s}}{\mathrm{d} t} = \omega_\mathrm{s}
+    = \frac{\mathrm{Im} \{ \boldsymbol{v}' + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}') + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v}')^* \} }{\hat{\psi}_\mathrm{R} + L_\sigma \mathrm{Re}\{(1 - \boldsymbol{k}_1)\boldsymbol{i}_\mathrm{s} + \boldsymbol{k}_2 \boldsymbol{i}_\mathrm{s}^* \}} 
+    :label: hatws
+
+where 
+
+.. math::
+    \boldsymbol{v}' = \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - L_\sigma \frac{\mathrm{d} \boldsymbol{i}_\mathrm{s}}{\mathrm{d} t} 
+    :label: vp
+
+The flux magnitude dynamics are 
+
+.. math::
+    \frac{\mathrm{d} \hat{\psi}_\mathrm{R}}{\mathrm{d} t} 
+    = \mathrm{Re}\{ \boldsymbol{v} + \boldsymbol{k}_1(\hat{\boldsymbol{v}} - \boldsymbol{v}) + \boldsymbol{k}_2(\hat{\boldsymbol{v}} - \boldsymbol{v})^* \}
+    :label: dhatpsiR_abs
+
+Notice that the right-hand side of :eq:`hatws` is independent of :math:`\omega_\mathrm{s}`. Futhermore, in these coordinates, the condition :eq:`inherently` for an inherently sensorless observer reduces to :math:`\boldsymbol{k}_2 = \boldsymbol{k}_1`. This observer structure is implemented in the :class:`motulator.control.im.Observer`, where simple forward-Euler discretization is used. 
+
+Gain Selection
+^^^^^^^^^^^^^^
+
+The estimation-error dynamics are obtained by subtracting :eq:`dhatpsiR` from :eq:`dpsiR`. The resulting system can be linearized for analysis and gain selection purposes. Using the rotor speed as an exmaple, the small-signal deviation about the operating point is :math:`\Delta \omega_\mathrm{m} = \omega_\mathrm{m} - \omega_\mathrm{m0}`, where the subscript 0 refers to the operating point. Linearization of the estimation-error dynamics leads to [#Hin2010]_
+
+.. math::
+	\frac{\mathrm{d} \Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} =  \boldsymbol{k}_1\Delta \tilde{\boldsymbol{v}} + \boldsymbol{k}_2\Delta \tilde{\boldsymbol{v}}^* - \mathrm{j}\omega_\mathrm{s0}\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}
+    :label: dtildepsiR
+
+where the estimation error of the rotor flux is :math:`\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R} = \Delta\boldsymbol{\psi}_\mathrm{R} - \Delta\hat{\boldsymbol{\psi}}_\mathrm{R}` and other estimation errors are marked similarly. Furthermore, the back-emf estimation error is
+
+.. math::
+    \Delta\tilde{\boldsymbol{v}} = -\left(\alpha - \mathrm{j}\omega_\mathrm{m0} \right)\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R} + \mathrm{j}\boldsymbol{\psi}_\mathrm{R0}\Delta\tilde{\omega}_\mathrm{m}
+    :label: dtildev
+
+Sensored Case
+"""""""""""""
+
+Here, :math:`\boldsymbol{k}_2 = 0` and :math:`\Delta\tilde{\omega}_\mathrm{m} = 0` are assumed, corresponding to the sensored reduced-order observer in [#Ver1988]_. Under these assumptions, the estimation-error dynamics in :eq:`dtildepsiR` reduce to 
+
+.. math::
+	\frac{\mathrm{d} \Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}}{\mathrm{d} t} =  -\left[\boldsymbol{k}_1\left(\alpha - \mathrm{j}\omega_\mathrm{m0} \right) + \mathrm{j}\omega_\mathrm{s0}\right]\Delta\tilde{\boldsymbol{\psi}}_\mathrm{R}
+    :label: dtildepsiR_sensored
+
+It can be noticed that the closed-loop pole could be arbitrarily placed via the gain :math:`\boldsymbol{k}_1`. Well-damped estimation-error dynamics can be obtained, e.g., by choosing
+
+.. math:: 
+    \boldsymbol{k}_1 = 1 + \frac{g |\omega_\mathrm{m}|}{\alpha - \mathrm{j}\omega_\mathrm{m}} 
+    :label: k1_sensored
+
+where :math:`g` is a unitless positive design parameter. The corresponding pole is located at :math:`s = -\alpha - g |\omega_\mathrm{m0}| - \mathrm{j}\omega_\mathrm{r0}`, where :math:`\omega_\mathrm{r0} = \omega_\mathrm{s0} - \omega_\mathrm{m0}` is the operating-point slip angular frequency.
+
+.. note::
+
+    As a special case, :math:`\boldsymbol{k}_1  = 0` yields the voltage model. As another special case, the current model is obtained choosing :math:`\boldsymbol{k}_1 = 1`. 
+
+Sensorless Case
+"""""""""""""""
+
+For sensorless drives, choosing 
+
+.. math::
+    \boldsymbol{k}_2 = (\hat{\boldsymbol{\psi}}_\mathrm{R}/\hat{\boldsymbol{\psi}}_\mathrm{R}^*) \boldsymbol{k}_1
+    :label: inherently
+    
+yields an inherently sensorless observer, i.e., the rotor speed estimate :math:`\hat{\omega}_\mathrm{m}` cancels out from the observer equations [#Hin2010]_. Under this condition, the linearized estimation-error dynamics in :eq:`dtildepsiR` become
+
+.. math::
+	\frac{\mathrm{d}}{\mathrm{d} t} \begin{bmatrix} \Delta\tilde{\psi}_\mathrm{Rd} \\ \Delta\tilde{\psi}_\mathrm{Rq} \end{bmatrix} = \begin{bmatrix} -2k_\mathrm{d}\alpha & -2k_\mathrm{d}\omega_\mathrm{m0} + \omega_\mathrm{s0} \\ -2k_\mathrm{q}\alpha - \omega_\mathrm{s0} & -2k_\mathrm{q}\omega_\mathrm{m0} 
+     \end{bmatrix} \begin{bmatrix} \Delta\tilde{\psi}_\mathrm{Rd} \\ \Delta\tilde{\psi}_\mathrm{Rq} \end{bmatrix}
+    :label: dtildepsiR_sensorless
+
+where the gain components correspond to :math:`\boldsymbol{k}_1 = k_\mathrm{d} + \mathrm{j}k_\mathrm{q}`. It can be seen that the dynamics of the rotor speed are decoupled from the flux-estimation error dynamics. The decay rate :math:`\sigma` be assigned by choosing
+
+.. math::
+    \boldsymbol{k}_1 = \frac{\sigma}{\alpha - \mathrm{j}\hat\omega_\mathrm{m}}
+    :label: k1k2_sensorless
+
+which results in the characteristic polynomial :math:`D(s)=s^2 + 2\sigma s + \omega_\mathrm{s0}^2`. The decay rate can be selected as :math:`\sigma = \alpha/2 + \zeta_\infty|\hat{\omega}_\mathrm{m}|`, where :math:`\zeta_\infty` is the desired damping ratio at high speed. At zero stator frequency :math:`\omega_\mathrm{s0} = 0`, the poles are located at :math:`s = 0` and :math:`s = -\alpha`, which allows stable magnetizing and starting the machine.
+
+.. rubric:: References
+
+.. [#Ver1988] Verghese, Sanders, “Observers for flux estimation in induction machines,” IEEE Trans. Ind. Electron., 1988, https://doi.org/10.1109/41.3067
+
+.. [#Har2001] Harnefors, “Design and analysis of general rotor-flux-oriented vector control systems,” IEEE Trans. Ind. Electron., 2001, https://doi.org/10.1109/41.915417
+
+.. [#Hin2010] Hinkkanen, Harnefors, Luomi, "Reduced-order flux observers with stator-resistance adaptation for speed-sensorless induction motor drives," IEEE Trans. Power Electron., 2010, https://doi.org/10.1109/TPEL.2009.2039650
+
+
+Synchronous Machines
+--------------------
+
+In sensorless control of synchronous machine drives, the rotor position and speed estimates are needed [#Jon1989]_, [#Cap2001]_, [#Hin2018]_. As a side product, the stator flux linkage is also estimated. In the following, an observer design available in the :class:`motulator.control.sm.Observer` class is considered, which is based on [#Hin2018]_. This observer implementation also includes a sensored mode. 
+
+
+Machine Model in General Coordinates
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+In :doc:`/model/machines`, the synchronous machine model is given in rotor coordinates. For the observer design and analysis, it is convenient to express this model in general coordinates, aligned at :math:`\vartheta_\mathrm{s}` and rotating at :math:`\omega_\mathrm{s} = \mathrm{d} \vartheta_\mathrm{s}/\mathrm{d} t` with respect to stator coordinates. Furthtermore, the rotor is aligned at :math:`\vartheta_\mathrm{m}` and rotates at :math:`\omega_\mathrm{m} = \mathrm{d} \vartheta_\mathrm{m}/\mathrm{d} t` with respect to stator coordinates. This coordinate transformation results in 
+
+.. math::
+    \frac{\mathrm{d}\boldsymbol{\psi}_\mathrm{s}}{\mathrm{d} t} &= \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - \mathrm{j}\omega_\mathrm{s}\boldsymbol{\psi}_\mathrm{s} \\
+    \frac{\mathrm{d}\delta}{\mathrm{d} t} &= \omega_\mathrm{m} - \omega_\mathrm{s}  
+    :label: sm
+
+where :math:`\boldsymbol{u}_\mathrm{s}` is the stator voltage, :math:`\boldsymbol{i}_\mathrm{s}` is the stator current, and :math:`\delta = \vartheta_\mathrm{m} - \vartheta_\mathrm{s}` is the electrical angle of the rotor as seen from the general coordinate system. Assuming linear magnetics, the relation between the stator flux linkage and the stator current is governed by 
+
+.. math::
+	\boldsymbol{\psi}_\mathrm{s} = \mathrm{e}^{\mathrm{j}\delta}\left(L_\mathrm{d}\mathrm{Re}\{\boldsymbol{i}_\mathrm{s} \mathrm{e}^{-\mathrm{j}\delta}\} + \mathrm{j}L_\mathrm{q}\mathrm{Im}\{\boldsymbol{i}_\mathrm{s}\mathrm{e}^{-\mathrm{j}\delta}\} + \psi_\mathrm{f}\right)
+    :label: sm_flux 
+
+where :math:`L_\mathrm{d}` is the d-axis inductance, :math:`L_\mathrm{q}` is the q-axis inductance, and :math:`\psi_\mathrm{f}` is the PM flux linkage. 
+
+Notice that setting :math:`\vartheta_\mathrm{s}=0` yields the machine model in stator coordinates. In the following, the coordinate system will be fixed to the estimated angle of the rotor, i.e., to the coordinate system used by the control system.  
+
+Observer Structure
+^^^^^^^^^^^^^^^^^^
+
+The observer is assumed to operate in estimated rotor coordinates, whose d-axis is aligned with the rotor angle estimate :math:`\hat{\vartheta}_\mathrm{m}`. Now, the angle :math:`\delta = \vartheta_\mathrm{m} - \hat{\vartheta}_\mathrm{m}` in the machine model :eq:`sm` corresponds to the estimation error of the rotor angle, which naturally is unknown to the sensorless control system. 
+
+Since the stator current is measured, the observer is fundamentally corrected by means of the current estimation error. However, due to the saliency, it is more convenient to scale the current estimation erro by the stator inductance, resulting in the flux linkage error   
+
+.. math::
+	\boldsymbol{e} = \psi_\mathrm{f} + L_\mathrm{d} \mathrm{Re}\{ \boldsymbol{i}_\mathrm{s}\} + \mathrm{j} L_\mathrm{q} \mathrm{Im}\{\boldsymbol{i}_\mathrm{s} \} - \hat{\boldsymbol{\psi}}_\mathrm{s} 
+    :label: e
+
+where :math:`\hat{\boldsymbol{\psi}}_\mathrm{s}` is the stator flux estimate. The flux linkage is estimated by 
+
+.. math::
+    \frac{\mathrm{d} \hat{\boldsymbol{\psi}}_\mathrm{s}}{\mathrm{d} t} = \boldsymbol{u}_\mathrm{s} - R_\mathrm{s}\boldsymbol{i}_\mathrm{s} - \mathrm{j}\omega_\mathrm{s}\hat{\boldsymbol{\psi}}_\mathrm{s} + \boldsymbol{k}_1 \boldsymbol{e} + \boldsymbol{k}_2 \boldsymbol{e}^* 
+    :label: sm_flux_observer
+
+where :math:`\boldsymbol{k}_1` and :math:`\boldsymbol{k}_2` are gains (complex in a general case), the estimates are marked with the hat, and :math:`^*` marks the complex conjugate. 
+
+In the sensored mode, :math:`\omega_\mathrm{s} = \omega_\mathrm{m}` is used. In the sensorless mode, the speed-adaptive structure (which would correspond to the phase-locked loop if the observer were implemented in stator coordinates) can be used to estimate the rotor angle and speed, respectively, as
+
+.. math::
+    \frac{\mathrm{d} \hat{\omega}_\mathrm{m}}{\mathrm{d} t} &= \mathrm{Im}\{\boldsymbol{k}_\mathrm{i} \boldsymbol{e}\} \\
+    \frac{\mathrm{d}\hat{\vartheta}_\mathrm{m}}{\mathrm{d} t} &= \hat{\omega}_\mathrm{m} + \mathrm{Im}\{\boldsymbol{k}_\mathrm{p} \boldsymbol{e}\} = \omega_\mathrm{s}  
+    :label: sm_speed_pos_observer
+
+where :math:`\boldsymbol{k}_\mathrm{i}` and :math:`\boldsymbol{k}_\mathrm{p}` are complex gains. This observer structure is used in the :class:`motulator.control.sm.Observer` class. 
+
+.. note::
+    Real-valued column vectors and the corresponding :math:`2\times 2` gain matrix were used in [#Hin2018]_. The complex form in :eq:`sm_flux_observer` has the same degrees of freedom.
+
+.. eps = -np.imag(e/psi_a) if np.abs(psi_a) > 0 else 0
+.. w_s = self.k_p*eps + self.w_m
+.. self.w_m += T_s*self.k_i*eps
+
+Gain Selection
+^^^^^^^^^^^^^^
+
+Sensored Case
+"""""""""""""
+
+In the sensored case, the gain :math:`\boldsymbol{k}_2=0` can be set in :eq:`sm_flux_observer`. Furthermore, :math:`\delta=0` holds. Therefore, using :eq:`sm` and :eq:`sm_flux_observer`, the linearized estimation error dynamics become 
+
+.. math::
+    \frac{\mathrm{d} \Delta\tilde{\boldsymbol{\psi}}_\mathrm{s}}{\mathrm{d} t} = -(\boldsymbol{k}_1 + \mathrm{j}\omega_\mathrm{m0})\Delta\tilde{\boldsymbol{\psi}}_\mathrm{s}  
+    :label: dtildepsis_sensored
+
+where :math:`\tilde{\boldsymbol{\psi}}_\mathrm{s} = \boldsymbol{\psi}_\mathrm{s} - \hat{\boldsymbol{\psi}}_\mathrm{s}` is the estimation error, :math:`\Delta` marks the small-signal quantities, and the subscript 0 marks the operating-point quantities. It can be seen that the pole can be arbitrarily placed via the gain :math:`\boldsymbol{k}_1`. Well-damped estimation-error dynamics can be obtained simply by using a real gain, :math:`\boldsymbol{k}_1 = \sigma`, resulting in the pole at :math:`s = -\sigma - \mathrm{j}\omega_\mathrm{m0}`, where :math:`\sigma = 2\pi \cdot 15` rad/s is used as the default value in :class:`motulator.control.sm.Observer`. 
+
+Sensorless Case
+"""""""""""""""
+
+The analysis of the sensorless case is more complicated. Here, the main results are summarized using the complex notation. The following results can be dereived from the linearized form of :eq:`sm` -- :eq:`sm_speed_pos_observer`, cf. further details in [#Hin2018]_
+
+To decouple the flux estimation from the rotor angle, the gains of :eq:`sm_flux_observer` have to be of the form
+
+.. math::
+	\boldsymbol{k}_1 = \sigma \qquad \boldsymbol{k}_2 = \frac{\sigma\hat{\boldsymbol{\psi}}_\mathrm{a}}{\hat{\boldsymbol{\psi}}_\mathrm{a}^*}
+    :label: k1k2_sensorless
+
+where :math:`\sigma` is the desired decay rate of the flux estimation error and
+
+.. math::
+    \hat{\boldsymbol{\psi}}_\mathrm{a} = \psi_\mathrm{f} + (L_\mathrm{d} - L_\mathrm{q}) \boldsymbol{i}_\mathrm{s}^*
+    :label: sm_aux_flux
+
+allows to decouple the flux-estimation error dynamics from the rotor-position dynamics. By default, the decay rate is scheduled as
+
+.. math:: 
+    \sigma = \frac{R_\mathrm{s}}{4}\left(\frac{1}{L_\mathrm{d}} + \frac{1}{L_\mathrm{q}}\right) + \zeta_\infty |\hat{\omega}_\mathrm{m} |
+    :label: sigma_sensorless
+
+where :math:`\zeta_\infty` is the desired damping ratio at high speed. At zero speed, :eq:`sigma_sensorless` places one pole at :math:`s = 0` and another at :math:`s = -(R_\mathrm{s}/2)(1/L_\mathrm{d} + 1/L_\mathrm{q})`.
+
+The gains of the speed adaptation in :eq:`sm_speed_pos_observer` are selected as
+
+.. math::
+	\boldsymbol{k}_\mathrm{i} = -\frac{\alpha_\mathrm{o}^2}{\hat{\boldsymbol{\psi}}_\mathrm{a}} \qquad \boldsymbol{k}_\mathrm{p} = -\frac{2\alpha_\mathrm{o}}{\hat{\boldsymbol{\psi}}_\mathrm{a}}
+    :label: ki_kp_sensorless
+
+where :math:`\alpha_\mathrm{o}` is the desired speed-estimation bandwidth. The choices :eq:`k1k2_sensorless` and :eq:`ki_kp_sensorless` result in the observer characteristic polynomial :math:`D(s) = (s^2 + 2\sigma s + \omega_\mathrm{m0}^2)(s + \alpha_\mathrm{o})^2`. Furthermore, it can also be shown that the resulting speed-estimation error dynamics are of the first order,
+
+.. math:: 
+    \frac{\Delta \hat{\omega}_\mathrm{m}(s)}{\Delta \omega_\mathrm{m}(s)} = \frac{\alpha_\mathrm{o}}{s + \alpha_\mathrm{o}}
+    :label: speed_est_dyn
+
+due to the pole-zero cancellation.
+
+.. note:: 
+    The flux linkage :math:`\boldsymbol{\psi}_\mathrm{a}` is called the auxiliary flux linkage in [#Hin2018]_. It is also linked to the maximum-torque-per-ampere (MTPA) condition, which can be compactly expressed as :math:`\mathrm{Re}\{\boldsymbol{i}_\mathrm{s}\boldsymbol{\psi}_\mathrm{a}^*\}=0` [#Var2021]_. 
+
+.. rubric:: References
+
+.. [#Jon1989] Jones, Lang, “A state observer for the permanent-magnet synchronous motor,” IEEE Trans. Ind. Electron., 1989, https://doi.org/10.1109/41.31500
+
+.. [#Cap2001] Capecchi, Guglielmo, Pastorelli, Vagati, “Position-sensorless control of the transverse-laminated synchronous reluctance motor,” IEEE Trans. Ind. Appl., 2001, https://doi.org/10.1109/28.968190
+
+.. [#Hin2018] Hinkkanen, Saarakkala, Awan, Mölsä, Tuovinen, "Observers for sensorless synchronous motor drives: Framework for design and analysis," IEEE Trans. Ind. Appl., 2018, https://doi.org/10.1109/TIA.2018.2858753
+
+.. [#Var2021] Varatharajan, Pellegrino, Armando, “Direct flux vector control of synchronous motor drives: A small-signal model for optimal reference generation,” IEEE Trans. Power Electron., 2021, https://doi.org/10.1109/TPEL.2021.3067694
+
+
+
diff --git a/docs/source/control/sm.rst b/docs/source/control/sm.rst
deleted file mode 100644
index 65824ecb4..000000000
--- a/docs/source/control/sm.rst
+++ /dev/null
@@ -1,6 +0,0 @@
-********************
-Synchronous Machines
-********************
-
-Control methods for synchronous machines will be described here in the future...
-
diff --git a/docs/source/index.rst b/docs/source/index.rst
index 778530890..bb409cd1d 100644
--- a/docs/source/index.rst
+++ b/docs/source/index.rst
@@ -1,7 +1,7 @@
 *motulator:* Motor Drive Simulator in Python
 ============================================
 
-This open-source software includes simulation models for an induction machine, a synchronous reluctance machine, and a permanent-magnet synchronous machine. The machine models are simulated in the continuous-time domain while the control algorithms run in discrete time. The default solver is the explicit Runge-Kutta method of order 5(4) from `scipy.integrate.solve_ivp`_. Simple control algorithms are provided as examples. The example algorithms aim to be simple yet feasible. 
+This open-source software includes simulation models for an induction machine, a synchronous reluctance machine, and a permanent-magnet synchronous machine. The machine models are simulated in the continuous-time domain while the control methods run in discrete time. The default solver is the explicit Runge-Kutta method of order 5(4) from `scipy.integrate.solve_ivp`_. A number of control methods are provided as examples. The example methods aim to be simple yet feasible. 
 
 .. _scipy.integrate.solve_ivp: https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.solve_ivp.html
 
@@ -25,13 +25,11 @@ This open-source software includes simulation models for an induction machine, a
 
 .. toctree::
    :titlesonly:
-   :caption: Control Algorithms
+   :caption: Control Methods
    :name: controllers
    :maxdepth: 1
 
-   control/common
-   control/im
-   control/sm
+   control/design_notes
    auto_examples/index
 
 .. rubric::
diff --git a/examples/flux_vector/plot_flux_vector_syrm_7kw.py b/examples/flux_vector/plot_flux_vector_syrm_7kw.py
index 64f605224..9ad70b190 100644
--- a/examples/flux_vector/plot_flux_vector_syrm_7kw.py
+++ b/examples/flux_vector/plot_flux_vector_syrm_7kw.py
@@ -62,7 +62,7 @@ def i_s(psi_s):
 ctrl = control.sm.FluxVectorCtrl(par, ref, sensorless=True)
 # Since the saturation is not considered in the control system, the speed
 # estimation bandwidth is set to a lower value
-ctrl.observer = control.sm.Observer(par, w_o=2*np.pi*50)
+ctrl.observer = control.sm.Observer(par, alpha_o=2*np.pi*50)
 
 # %%
 # Set the speed reference and the external load torque.
diff --git a/examples/signal_inj/plot_signal_inj_syrm_7kw.py b/examples/signal_inj/plot_signal_inj_syrm_7kw.py
index 3b96026ef..98689746d 100644
--- a/examples/signal_inj/plot_signal_inj_syrm_7kw.py
+++ b/examples/signal_inj/plot_signal_inj_syrm_7kw.py
@@ -7,7 +7,6 @@
 
 """
 
-# %%
 # %%
 # Imports.
 
diff --git a/examples/vector/plot_vector_ctrl_pmsm_2kw.py b/examples/vector/plot_vector_ctrl_pmsm_2kw.py
index 29a13a925..4cfedc7c5 100644
--- a/examples/vector/plot_vector_ctrl_pmsm_2kw.py
+++ b/examples/vector/plot_vector_ctrl_pmsm_2kw.py
@@ -6,7 +6,6 @@
 
 """
 
-# %%
 # %%
 # Imports.
 
diff --git a/examples/vector/plot_vector_ctrl_pmsyrm_thor.py b/examples/vector/plot_vector_ctrl_pmsyrm_thor.py
index 4ed14f38e..3ad2f40f0 100644
--- a/examples/vector/plot_vector_ctrl_pmsyrm_thor.py
+++ b/examples/vector/plot_vector_ctrl_pmsyrm_thor.py
@@ -7,7 +7,6 @@
 
 """
 
-# %%
 # %%
 # Imports.
 
@@ -39,7 +38,7 @@
 ref = control.sm.CurrentReferencePars(
     par, w_m_nom=base.w, i_s_max=2*base.i, k_u=.9)
 ctrl = control.sm.VectorCtrl(par, ref, T_s=125e-6, sensorless=True)
-ctrl.observer = control.sm.Observer(par, w_o=2*np.pi*200)
+ctrl.observer = control.sm.Observer(par, alpha_o=2*np.pi*200)
 ctrl.speed_ctrl = control.SpeedCtrl(
     J=par.J, alpha_s=2*np.pi*4, tau_M_max=1.5*base.tau_nom)
 
diff --git a/examples/vector/plot_vector_ctrl_syrm_7kw.py b/examples/vector/plot_vector_ctrl_syrm_7kw.py
index c217a032a..66cf41f73 100644
--- a/examples/vector/plot_vector_ctrl_syrm_7kw.py
+++ b/examples/vector/plot_vector_ctrl_syrm_7kw.py
@@ -6,7 +6,6 @@
 
 """
 
-# %%
 # %%
 # Imports.
 
diff --git a/examples/vhz/plot_vhz_ctrl_6step_im_2kw.py b/examples/vhz/plot_vhz_ctrl_6step_im_2kw.py
index c3351ac02..9c8c26869 100644
--- a/examples/vhz/plot_vhz_ctrl_6step_im_2kw.py
+++ b/examples/vhz/plot_vhz_ctrl_6step_im_2kw.py
@@ -10,7 +10,6 @@
 
 """
 # %%
-# %%
 # Imports.
 
 import numpy as np
diff --git a/examples/vhz/plot_vhz_ctrl_im_2kw.py b/examples/vhz/plot_vhz_ctrl_im_2kw.py
index 5508de518..f2b01b526 100644
--- a/examples/vhz/plot_vhz_ctrl_im_2kw.py
+++ b/examples/vhz/plot_vhz_ctrl_im_2kw.py
@@ -7,7 +7,6 @@
 
 """
 # %%
-# %%
 # Imports.
 
 import numpy as np
diff --git a/motulator/__init__.py b/motulator/__init__.py
index 160a8d5c4..a7a5955c1 100644
--- a/motulator/__init__.py
+++ b/motulator/__init__.py
@@ -1,5 +1,5 @@
 """
-*motulator*: Motor drive simulator in Python
+*motulator*: Motor Drive Simulator in Python
 
 This software includes continuous-time simulation models for induction machines 
 and synchronous machines. Furthermore, selected examples of discrete-time 
diff --git a/motulator/control/_common.py b/motulator/control/_common.py
index 1a9f65bbc..914730578 100644
--- a/motulator/control/_common.py
+++ b/motulator/control/_common.py
@@ -14,7 +14,7 @@ class PWM:
     The realized voltage is computed based on the measured DC-bus voltage and 
     the duty ratios. The digital delay effects are taken into account in the 
     realized voltage, assuming the delay of a single sampling period. The total
-    delay is 1.5 sampling periods due to the PWM (or ZOH) delay.
+    delay is 1.5 sampling periods due to the PWM (or ZOH) delay [#Bae2003]_.
 
     Parameters
     ----------
@@ -26,6 +26,11 @@ class PWM:
     realized_voltage : complex
         Realized voltage (V) in synchronous coordinates.
 
+    References
+    ----------
+    .. [#Bae2003] Bae, Sul, "A compensation method for time delay of 
+       full-digital synchronous frame current regulator of PWM AC drives," IEEE Trans. Ind. Appl., 2003, https://doi.org/10.1109/TIA.2003.810660
+    
     """
 
     def __init__(self, six_step=False):
@@ -205,7 +210,7 @@ class PICtrl:
     k_t : float, optional
         Reference-feedforward gain. The default is `k_p`.
     u_max : float, optional
-        Maximum controller output. The default is `inf`.
+        Maximum controller output. The default is inf.
 
     Attributes
     ----------
@@ -280,11 +285,11 @@ class SpeedCtrl(PICtrl):
     Parameters
     ----------
     J : float
-        Total inertia of the rotor (kgm**2).
+        Total inertia of the rotor (kgm²).
     alpha_s : float
         Closed-loop bandwidth (rad/s). 
     tau_M_max : float, optional
-        Maximum motor torque (Nm). The default is `inf`.
+        Maximum motor torque (Nm). The default is inf.
 
     """
 
diff --git a/motulator/control/im/_obs_vhz.py b/motulator/control/im/_obs_vhz.py
index 9c9f8bb6c..37b9e5e53 100644
--- a/motulator/control/im/_obs_vhz.py
+++ b/motulator/control/im/_obs_vhz.py
@@ -33,18 +33,18 @@ class ObserverBasedVHzCtrlPars:
     psi_s_nom : float
         Nominal stator flux linkage (Vs). 
     i_s_max : float, optional
-        Maximum stator current (A). The default is `inf`.
+        Maximum stator current (A). The default is inf.
     k_tau : float, optional
-        Torque controller gain. The default is `3`.
+        Torque controller gain. The default is 3.
     alpha_psi : float, optional
-        Stator flux control bandwidth (rad/s). The default is `2*pi*20`.
+        Stator flux control bandwidth (rad/s). The default is 2*pi*20.
     alpha_f : float, optional
-        Torque high-pass filter bandwidth (rad/s). The default is `2*pi*1`.
+        Torque high-pass filter bandwidth (rad/s). The default is 2*pi*1.
     alpha_r : float, optional
         Low-pass-filter bandwidth (rad/s) for slip angular frequency. The
-        default is `2*pi*1`.
+        default is 2*pi*1.
     slip_compensation : bool, optional
-        Enable slip compensation. The default is `False`.
+        Enable slip compensation. The default is False.
 
     """
     # par: InitVar[ModelPars | None] = None
diff --git a/motulator/control/im/_vector.py b/motulator/control/im/_vector.py
index 04b9f98f9..19db57e9d 100644
--- a/motulator/control/im/_vector.py
+++ b/motulator/control/im/_vector.py
@@ -205,9 +205,9 @@ class CurrentReferencePars:
     i_s_max : float
         Maximum stator current (A). 
     u_s_nom : float, optional
-        Nominal stator voltage (V). The default is `sqrt(2/3)*400`.
+        Nominal stator voltage (V). The default is sqrt(2/3)*400.
     w_s_nom : float, optional
-        Nominal stator angular frequency (rad/s). The default is `2*pi*50`.
+        Nominal stator angular frequency (rad/s). The default is 2*pi*50.
     psi_R_nom : float, optional
         Nominal rotor flux linkage (Vs). The default is 
         `(u_s_nom/w_s_nom)/(1 + L_sgm/L_M)`.
@@ -367,10 +367,10 @@ class Observer:
         Machine model parameters.
     k : callable, optional
         Observer gain as a function of the rotor angular speed. The default is 
-        ``lambda w_m: 1 + .2*abs(w_m)/(R_R/L_M - 1j*w_m)`` if not `sensorless`
-        else         
+        ``lambda w_m: (0.5*R_R/L_M + 0.2*abs(w_m))/(R_R/L_M - 1j*w_m)`` if
+        `sensorless` else ``lambda w_m: 1 + 0.2*abs(w_m)/(R_R/L_M - 1j*w_m)``.
     alpha_o : float, optional
-        Observer bandwidth (rad/s). The default is `2*pi*40`.
+        Observer bandwidth (rad/s). The default is 2*pi*40.
     sensorless : bool, optional
         If True, sensorless mode is used. The default is True.
 
diff --git a/motulator/control/sm/_flux_vector.py b/motulator/control/sm/_flux_vector.py
index 67a0812eb..109e24d5f 100644
--- a/motulator/control/sm/_flux_vector.py
+++ b/motulator/control/sm/_flux_vector.py
@@ -31,13 +31,13 @@ class FluxVectorCtrl(Ctrl):
     ref : FluxTorqueReferencePars
         Reference generation parameters.
     alpha_psi : float, optional
-        Bandwidth of the flux controller (rad/s). The default is `2*pi*100`.
+        Bandwidth of the flux controller (rad/s). The default is 2*pi*100.
     alpha_tau : float, optional
-        Bandwidth of the torque controller (rad/s). The default is `2*pi*200`.
+        Bandwidth of the torque controller (rad/s). The default is 2*pi*200.
     T_s : float
-        Sampling period (s). The default is `250e-6`.
+        Sampling period (s). The default is 250e-6.
     sensorless : bool, optional
-        If `True`, sensorless control is used. The default is `True`.
+        If True, sensorless control is used. The default is True.
 
     Attributes
     ----------
@@ -81,7 +81,8 @@ def __init__(
         self.n_p = par.n_p
         self.R_s, self.L_d, self.L_q = par.R_s, par.L_d, par.L_q
         # Subsystems
-        self.observer = Observer(par, w_o=2*np.pi*100, sensorless=sensorless)
+        self.observer = Observer(
+            par, alpha_o=2*np.pi*100, sensorless=sensorless)
         self.flux_torque_ref = FluxTorqueReference(ref)
         self.pwm = PWM()
         self.speed_ctrl = SpeedCtrl(par.J, 2*np.pi*4)
@@ -197,7 +198,7 @@ class FluxTorqueReferencePars:
     psi_s_min : float, optional
         Minimum stator flux (Vs). The default is `psi_f`.
     psi_s_max : float, optional
-        Maximum stator flux (Vs). The default is `inf`.
+        Maximum stator flux (Vs). The default is inf.
     k_u : float, optional
         Voltage utilization factor. The default is 0.95.
 
diff --git a/motulator/control/sm/_obs_vhz.py b/motulator/control/sm/_obs_vhz.py
index b646b7e07..c39b96b14 100644
--- a/motulator/control/sm/_obs_vhz.py
+++ b/motulator/control/sm/_obs_vhz.py
@@ -22,11 +22,11 @@ class ObserverBasedVHzCtrlPars(FluxTorqueReferencePars):
     Parameters
     ----------
     alpha_psi : float, optional
-        Flux control bandwidth (rad/s). The default is `2*pi*50`.
+        Flux control bandwidth (rad/s). The default is 2*pi*50.
     alpha_tau : float
-        Torque control bandwidth (rad/s). The default is `2*pi*50`.
+        Torque control bandwidth (rad/s). The default is 2*pi*50.
     alpha_f : float, optional
-        Bandwidth of the high-pass filter (rad/s). The default is `2*pi*1`.
+        Bandwidth of the high-pass filter (rad/s). The default is 2*pi*1.
 
     """
     alpha_psi: float = 2*np.pi*50
@@ -60,7 +60,7 @@ class ObserverBasedVHzCtrl(Ctrl):
     ctrl_par : ObserverBasedVHzCtrlPars
         Control system parameters.
     T_s : float, optional
-        Sampling period (s). The default is `250e-6`.
+        Sampling period (s). The default is 250e-6.
 
     Attributes
     ----------
@@ -190,9 +190,9 @@ class FluxObserver:
     par : ModelPars
         Machine model parameters.
     alpha_o : float, optional
-        Observer gain (rad/s). The default is `2*pi*20`.
+        Observer gain (rad/s). The default is 2*pi*20.
     zeta_inf : float, optional
-        Damping ratio at infinite speed. The default is `0.2`.
+        Damping ratio at infinite speed. The default is 0.2.
 
     """
 
diff --git a/motulator/control/sm/_torque.py b/motulator/control/sm/_torque.py
index ab4872142..8ae7b60a0 100644
--- a/motulator/control/sm/_torque.py
+++ b/motulator/control/sm/_torque.py
@@ -54,12 +54,12 @@ def torque(self, psi_s):
         Parameters
         ----------
         psi_s : complex
-            Stator flux.
+            Stator flux (Vs).
 
         Returns
         -------
         tau_M : float
-            Electromagnetic torque.
+            Electromagnetic torque (Nm).
 
         """
         i_s = self.current(psi_s)
@@ -74,12 +74,12 @@ def current(self, psi_s):
         Parameters
         ----------
         psi_s : complex
-            Stator flux linkage.
+            Stator flux linkage (Vs).
 
         Returns
         -------
         i_s : complex
-            Stator current.
+            Stator current (A).
 
         """
         i_s = (psi_s.real - self.psi_f)/self.L_d + 1j*psi_s.imag/self.L_q
@@ -93,12 +93,12 @@ def flux(self, i_s):
         Parameters
         ----------
         i_s : complex
-            Stator current.
+            Stator current (A).
 
         Returns
         -------
         psi_s : complex
-            Stator flux linkage.
+            Stator flux linkage (Vs).
 
         """
         psi_s = self.L_d*i_s.real + self.psi_f + 1j*self.L_q*i_s.imag
@@ -112,12 +112,12 @@ def mtpa(self, abs_i_s):
         Parameters
         ----------
         abs_i_s : float
-            Stator current magnitude.
+            Stator current magnitude (A).
 
         Returns
         -------
         beta : float
-            MTPA angle of the stator current vector.
+            MTPA angle of the stator current vector (electrical rad).
 
         """
         # Replace zeros with epsilon
@@ -147,12 +147,12 @@ def mtpv(self, abs_psi_s):
         Parameters
         ----------
         abs_psi_s : float
-            Stator flux magnitude.
+            Stator flux magnitude (Vs).
 
         Returns
         -------
         delta : float
-            MTPV angle of the stator flux vector.
+            MTPV angle of the stator flux vector (electrical rad).
 
         """
         # Replace zeros with epsilon
@@ -190,12 +190,12 @@ def mtpv_current(self, abs_i_s):
         Parameters
         ----------
         abs_i_s : float
-            Stator current magnitude.
+            Stator current magnitude (A).
 
         Returns
         -------
         i_s : complex
-            MTPV stator current.
+            MTPV stator current (A).
 
         """
         if self.psi_f == 0:
@@ -232,9 +232,9 @@ def mtpa_locus(self, i_s_max, psi_s_min=None, N=20):
         Parameters
         ----------
         i_s_max : float
-            Maximum stator current magnitude at which the locus is computed.
+            Maximum stator current magnitude (A) at which the locus is computed.
         psi_s_min : float, optional
-            Minimum stator flux magnitude at which the locus is computed.
+            Minimum stator flux magnitude (Vs) at which the locus is computed.
         N : int, optional
             Amount of points. The default is 20.
 
@@ -242,15 +242,15 @@ def mtpa_locus(self, i_s_max, psi_s_min=None, N=20):
         -------
         Bunch object with the following fields defined:
         psi_s : complex
-            Stator flux.
+            Stator flux (Vs).
         i_s : complex
-            Stator current.
+            Stator current (A).
         tau_M : float
-            Electromagnetic torque.
+            Electromagnetic torque (Nm).
         abs_psi_s_vs_tau_M : callable
-            Stator flux magnitude as a function of the torque.
+            Stator flux magnitude (Vs) as a function of the torque (Nm).
         i_sd_vs_tau_M : callable
-            d-axis current as a function of the torque.
+            d-axis current (A) as a function of the torque (Nm).
 
         """
         # Current  magnitudes
@@ -292,10 +292,10 @@ def mtpv_locus(self, psi_s_max=None, i_s_max=None, N=20):
         Parameters
         ----------
         psi_s_max : float, optional
-            Maximum stator flux magnitude at which the locus is computed. Either
-            psi_s_max or i_s_max must be given.
+            Maximum stator flux magnitude (Vs) at which the locus is computed. 
+            Either `psi_s_max` or `i_s_max` must be given.
         i_s_max : float, optional
-            Maximum stator current magnitude at which the locus is computed.
+            Maximum stator current magnitude (A) at which the locus is computed.
         N : int, optional
             Amount of points. The default is 20.
 
@@ -303,13 +303,13 @@ def mtpv_locus(self, psi_s_max=None, i_s_max=None, N=20):
         -------
         Bunch object with the following fields defined:
         psi_s : complex
-            Stator flux.
+            Stator flux (Vs).
         i_s : complex
-            Stator current.
+            Stator current (A).
         tau_M : float
-            Electromagnetic torque.
+            Electromagnetic torque (Nm).
         tau_M_vs_abs_psi_s : interp1d object
-            Torque as a function of the flux magnitude.
+            Torque (Nm) as a function of the flux magnitude (Vs).
 
         """
         # If i_s_max is given, compute the corresponding MTPV stator flux
@@ -344,11 +344,11 @@ def current_limit(self, i_s_max, gamma1=np.pi, gamma2=0, N=20):
         Parameters
         ----------
         i_s_max : float
-            Current limit. 
+            Current limit (A). 
         gamma1 : float, optional
-            Starting angle in radians. The default is 0.
+            Starting angle (electrical rad). The default is 0.
         gamma2 : float, optional
-            End angle in radians. The defauls in np.pi.
+            End angle (electrical rad). The defauls in pi.
         N : int, optional
             Amount of points. The default is 20.
 
@@ -356,13 +356,13 @@ def current_limit(self, i_s_max, gamma1=np.pi, gamma2=0, N=20):
         -------
         Bunch object with the following fields defined:
         psi_s : complex
-            Stator flux.
+            Stator flux (Vs).
         i_s : complex
-            Stator current.
+            Stator current (A).
         tau_M : float
-            Electromagnetic torque.
+            Electromagnetic torque (Nm).
         tau_M_vs_abs_psi_s : interp1d object
-            Torque as a function of the flux magnitude.
+            Torque (Nm) as a function of the flux magnitude (Vs).
 
         """
         if np.isnan(gamma1):
@@ -398,7 +398,7 @@ def mtpv_and_current_limits(self, i_s_max, N=20):
         Parameters
         ----------
         i_s_max : float
-            Current limit.
+            Current limit (A).
         N : int, optional
             Amount of points. The default is 20.
 
@@ -406,9 +406,9 @@ def mtpv_and_current_limits(self, i_s_max, N=20):
         -------
         Bunch object with the following fields defined:
         tau_M_vs_abs_psi_s : interp1d object
-            Torque as a function of the flux magnitude.
+            Torque (Nm) as a function of the flux magnitude (Vs).
         i_sd_vs_tau_M : interp1d object
-            d-axis current as a function of the torque.
+            d-axis current (A) as a function of the torque (Nm).
 
         """
         mtpa = self.mtpa_locus(i_s_max=i_s_max, N=N)
@@ -451,7 +451,7 @@ def plot_flux_loci(self, i_s_max, base, N=20):
         Parameters
         ----------
         i_s_max : float
-            Maximum current at which the loci are evaluated.
+            Maximum current (A) at which the loci are evaluated.
         base : BaseValues
             Base values.
         N : int, optional
@@ -498,7 +498,7 @@ def plot_current_loci(self, i_s_max, base, N=20):
         Parameters
         ----------
         i_s_max : float
-            Maximum current at which the loci are evaluated.
+            Maximum current (A) at which the loci are evaluated.
         base : BaseValues
             Base values.
         N : int, optional
@@ -544,7 +544,7 @@ def plot_torque_current(self, i_s_max, base, N=20):
         Parameters
         ----------
         i_s_max : float
-            Maximum current at which the loci are evaluated.
+            Maximum current (A) at which the loci are evaluated.
         base : BaseValues
             Base values.
         N : int, optional
@@ -600,7 +600,7 @@ def plot_torque_flux(self, i_s_max, base, N=20):
         Parameters
         ----------
         i_s_max : float
-            Maximum current at which the loci are evaluated.
+            Maximum current (A) at which the loci are evaluated.
         base : BaseValues
             Base values.
         N : int, optional
diff --git a/motulator/control/sm/_vector.py b/motulator/control/sm/_vector.py
index cf840c6e4..0b3e0dd53 100644
--- a/motulator/control/sm/_vector.py
+++ b/motulator/control/sm/_vector.py
@@ -83,7 +83,7 @@ def __init__(self, par, ref, T_s=250e-6, sensorless=True):
         self.current_ref = CurrentReference(par, ref)
         self.current_ctrl = CurrentCtrl(par, 2*np.pi*200)
         if sensorless:
-            self.observer = Observer(par, w_o=2*np.pi*100, sensorless=True)
+            self.observer = Observer(par, alpha_o=2*np.pi*100, sensorless=True)
         else:
             self.observer = None
         self.speed_ctrl = SpeedCtrl(par.J, 2*np.pi*4)
@@ -395,22 +395,26 @@ def update(self, T_s, tau_M_ref_lim, u_s_ref, u_dc):
 # %%
 class Observer:
     """
-    Observer for the rotor position and the stator flux linkage.
+    Observer for synchronous machines.
 
-    This observer corresponds to [#Hin2018]_. The observer gain decouples the
-    electrical and mechanical dynamics and allows placing the poles of the
-    corresponding linearized estimation error dynamics. This implementation
-    operates in estimated rotor coordinates. The observer can also be used in
-    the sensored mode by providing the measured rotor speed as an input.
+    This observer estimates the rotor angle, the rotor speed, and the stator 
+    flux linkage. The design is based on [#Hin2018]_. The observer gain 
+    decouples the electrical and mechanical dynamics and allows placing the 
+    poles of the corresponding linearized estimation error dynamics. This 
+    implementation operates in estimated rotor coordinates. The observer can 
+    also be used in the sensored mode by providing the measured rotor speed as 
+    an input.
 
     Parameters
     ----------
     par : ModelPars
         Machine model parameters.
-    w_o : float, optional
-        Observer bandwidth (electrical rad/s).
-    zeta_inf : float, optional
-        Damping ratio at high speed. The default is .2.
+    alpha_o : float, optional
+        Observer bandwidth (electrical rad/s). The default is 2*pi*40.
+    k : callable, optional
+        Observer gain as a function of the rotor angular speed. The default is 
+        ``lambda w_m: 0.25*(R_s*(L_d + L_q)/(L_d*L_q) + 0.2*abs(w_m))`` if
+        `sensorless` else ``lambda w_m: 2*pi*15``.
 
     Attributes
     ----------
@@ -430,15 +434,22 @@ class Observer:
     """
 
     # pylint: disable=too-many-instance-attributes
-    def __init__(self, par, w_o=2*np.pi*40, zeta_inf=.2, sensorless=True):
+    def __init__(self, par, alpha_o=2*np.pi*40, k=None, sensorless=True):
         self.R_s = par.R_s
         self.L_d, self.L_q, self.psi_f = par.L_d, par.L_q, par.psi_f
         self.psi_f = par.psi_f
-        self.k_p = 2*w_o
-        self.k_i = w_o**2
-        self.zeta_inf = zeta_inf
+        self.alpha_o = alpha_o
         self.sensorless = sensorless
-        self.b_p = .5*par.R_s*(par.L_d + par.L_q)/(par.L_d*par.L_q)
+        self.k1 = k
+        if self.sensorless:
+            if self.k1 is None:  # If not given, use the default gains
+                sigma0 = .25*par.R_s*(par.L_d + par.L_q)/(par.L_d*par.L_q)
+                self.k1 = lambda w_m: (sigma0 + .2*np.abs(w_m))
+            self.k2 = self.k1
+        else:
+            if self.k1 is None:
+                self.k1 = lambda w_m: 2*np.pi*15
+            self.k2 = 0*self.k1
         # Initial states
         self.theta_m, self.w_m, self.psi_s = 0, 0, par.psi_f
 
@@ -451,9 +462,9 @@ def update(self, T_s, u_s, i_s, w_m=None):
         T_s : float
             Sampling period (s).
         u_s : complex
-            Stator voltage in estimated rotor coordinates.
+            Stator voltage (V) in estimated rotor coordinates.
         i_s : complex
-            Stator current in estimated rotor coordinates.
+            Stator current (A) in estimated rotor coordinates.
         w_m : float, optional
             Rotor angular speed (electrical rad/s). Needed only in the sensored
             mode. The default is None. 
@@ -466,24 +477,21 @@ def update(self, T_s, u_s, i_s, w_m=None):
             # Auxiliary flux
             psi_a = self.psi_f + (self.L_d - self.L_q)*np.conj(i_s)
 
-            # Pole locations are chosen with c = w_m**2
-            k = self.b_p + 2*self.zeta_inf*np.abs(self.w_m)
-            if np.abs(psi_a) > 0:
-                # Correction voltage
-                v = k*psi_a*np.real(e/psi_a)
-                # Error signal
-                eps = -np.imag(e/psi_a)
-            else:
-                v, eps = 0, 0
+            # Observer gains
+            k1 = self.k1(self.w_m)
+            k2 = k1*psi_a/np.conj(psi_a) if np.abs(psi_a) > 0 else k1
 
             # Speed estimation
-            w_m = self.k_p*eps + self.w_m
-            self.w_m += T_s*self.k_i*eps
+            eps = -np.imag(e/psi_a) if np.abs(psi_a) > 0 else 0
+            w_s = 2*self.alpha_o*eps + self.w_m
         else:
-            # Correction voltage
-            v = 2*np.pi*15*e
+            k1, k2 = self.k1, 0
+            w_s = w_m
+            eps = 0
 
         # Update the states
-        self.psi_s += T_s*(u_s - self.R_s*i_s - 1j*w_m*self.psi_s + v)
-        self.theta_m += T_s*w_m  # Next line: limit into [-pi, pi)
+        self.psi_s += T_s*(
+            u_s - self.R_s*i_s - 1j*w_s*self.psi_s + k1*e + k2*np.conj(e))
+        self.w_m += T_s*self.alpha_o**2*eps
+        self.theta_m += T_s*w_s  # Next line: limit into [-pi, pi)
         self.theta_m = np.mod(self.theta_m + np.pi, 2*np.pi) - np.pi