Arbitrary mechanical plants can be specified for each optic by giving a zpk description. The example below makes a simple Fabry-Perot cavity and specifies the longitudinal and pitch plants.

First build the cavity

```
import numpy as np
from finesse import Model
import finesse.components as fcmp
import finesse.analysis.actions as fact
from finesse.plotting import bode
Pin_W = 5e3 # input power [W]
Ti = 0.014 # ITM transmissivity
Lcav_m = 40e3 # cavity length [m]
Ri_m = 34e3 # ITM radius of curvature [m]
Re_m = 36e3 # ETM radius of curvature [m]
# make the cavity
model = Model()
model.add(fcmp.Mirror('ETM', T=0, L=0, Rc=Re_m))
model.add(fcmp.Mirror('ITM', T=Ti, L=0, Rc=Ri_m))
model.connect(model.ETM.p1, model.ITM.p1, L=Lcav_m)
model.add(fcmp.Cavity('cav', model.ETM.p1.o))
# add input
model.add(fcmp.Laser('Laser', P=Pin_W))
model.connect(model.Laser.p1, model.ITM.p2)
```

The suspension plants are specified by defining dictionaries with the zpk definition of the plant. The code below defines a simple pendulum plant with a single complex pole. The degrees of freedom which can be read out are `z`

, `pitch`

, or `yaw`

and the degrees of freedom which can be driven are `F_z`

, `F_pitch`

, and `F_yaw`

.

So for example, `sus_plant['z', F_z']`

is the longitudinal plant in m / N and `sus_plant['pitch', 'F_z']`

is the cross-coupling of longitudinal motion into pitch in units of rad / N. The code below doesn’t define any cross-couplings.

```
# define the suspension plant
I_kgm2 = 10 # moment of inertia
M_kg = 0.01 # mass
# for simplicity use the same poles for both pitch and z
F0_Hz = 1 # resonance frequency [Hz]
Q = 100 # Q factor
Fp_Hz = np.array([
-F0_Hz / (2 * Q) * (1 + 1j * np.sqrt(4 * Q**2 - 1)),
-F0_Hz / (2 * Q) * (1 - 1j * np.sqrt(4 * Q**2 - 1))
])
# the suspension plant is defined with a dictionary
sus_plant = {}
sus_plant['z', 'F_z'] = ([], 2 * np.pi * Fp_Hz, 1 / M_kg)
sus_plant['pitch', 'F_pitch'] = ([], 2 * np.pi * Fp_Hz, 1 / I_kgm2)
model.add(fcmp.mechanical.SuspensionZPK('ETM_sus', model.ETM.mech, sus_plant))
model.add(fcmp.mechanical.SuspensionZPK('ITM_sus', model.ITM.mech, sus_plant))
```

Now run the model computing both pitch and longitudinal motion. First compute pitch then detune the cavity by 1 degree so that there is a longitudinal optical spring and calculate the longitudinal motion.

```
model.modes(maxtem=1)
F_Hz = np.geomspace(0.1, 3e3, 1000) # the best thing since sliced bread
out = model.run(
fact.Series(
fact.FrequencyResponse(
F_Hz,
['ETM.mech.F_pitch', 'ITM.mech.F_pitch'],
['ETM.mech.pitch', 'ITM.mech.pitch'],
name='pitch',
),
fact.Change({'ETM.phi': 1}),
fact.FrequencyResponse(
F_Hz,
['ETM.mech.F_z', 'ITM.mech.F_z'],
['ETM.mech.z', 'ITM.mech.z'],
name='z',
),
)
)
tfs_pitch = out['pitch'].out # (1000, 2, 2) matrix of angular transfer functions
tfs_z = out['z'].out # (1000, 2, 2) matrix of longitudinal transfer functions
```

Now plot the torsional spring in the hard/soft basis. The above transfer functions were computed in the ETM/ITM basis. Hard and soft degrees of freedom can also be defined using `DegreeOfFreedom`

and used to compute the transfer functions instead.

```
gi = 1 - Lcav_m / Ri_m
ge = 1 - Lcav_m / Re_m
r = 2 / ((gi - ge) + np.sqrt((gi - ge)**2 + 4))
hard = dict(ITM=-1, ETM=r)
soft = dict(ITM=r, ETM=1)
vhard = np.array([hard['ETM'], hard['ITM']])
vsoft = np.array([soft['ETM'], soft['ITM']])
# tfs_pitch is in the ETM/ITM basis, so convert to hard/soft
tf_hard = vhard @ tfs_pitch @ vhard.T
tf_soft = vsoft @ tfs_pitch @ vsoft.T
axs = bode(F_Hz, tf_hard, label='Hard')
bode(F_Hz, tf_soft, axs=axs, label='Soft')
bode(F_Hz, tfs_pitch[..., 0, 0], axs=axs, label='ETM')
for ax in axs:
ax.grid(True, which='major', alpha=0.5)
ax.grid(True, which='minor', alpha=0.2)
ax.set_xlim(F_Hz[0], 30)
axs[0].set_ylim(-100, 10)
axs[0].legend()
```

And now plot the longitudinal optical spring.

```
axs = bode(F_Hz, tfs_z[..., 0, 0], label='ETM to ETM')
bode(F_Hz, tfs_z[..., 1, 0], axs=axs, label='ETM to ITM', ls='--')
for ax in axs:
ax.grid(True, which='major', alpha=0.5)
ax.grid(True, which='minor', alpha=0.2)
ax.set_xlim(F_Hz[0], F_Hz[-1])
```