**ASSIGNMENT 2：System Identification**

In this assignment you are presented with a number of ‘black box’ components, each of which contain code that models a certain electronic system. By investigating the response of each box to generated input (voltage) signals, you will be able to determine the type of system each box contain and estimate its parameters (e.g. resistance and inductance).

Before you continue with this assignment, please read `Yr2_Computing_System_Identification.pdf`

and work though the example notebook `black_box_example.ipynb`

first.

As before, you must store your student ID in the variable studentID and execute the cell below before starting the assignment. This will import the necessary Boxes objects. As in assignment 1, the final cell contains a function call (`Boxes.check()`

) which you should run to check that the naming of variables and functions in each task is correct.

1 | import numpy as np |

## Task 1

### Marks available: 4

Each of the following circuits is contained in *one* of the boxes `t1_box1`

, `t1_box2`

and `t1_box3`

:

- RL high pass filter
- RC low pass filter
- RLC band pass filter

By analysing the action of each box on custom signals and employing the `.fft`

and `.fftfreq`

functions from `numpy.fft`

, identify the type of circuit contained in each box. A blank cell has been provided for your working, but you are free to insert as many additional cells as you require.

After completing your investigation set the variables `RL_circuit`

, `RC_circuit`

and `RLC_circuit`

in the next cell to a value of 1, 2 or 3 to indicate which box matches the relevant circuit. For example, write `RC_circuit = 1`

if you believe `t1_box1`

to be the RC low pass filter or`RC_circuit = 2`

if you believe that system is in `t1_box2`

.

**Reminder**:

Use the `.process()`

method to interact with the Boxes objects, as demonstrated in `black_box_example.ipynb`

**Hint**:

In order to obtain an informative result, your time domain input signal should contain frequency components ranging from 10 to 4000 Hz.

1 | # Working for Task 1 (you are free to insert additional |

Data points = 40000

Total time = 9.99975 s

Timestep = 0.00025 s

Sampling Frequency = 4000.0 Hz

1 | fig = plt.figure(figsize=(15,3)) |

1 | # Pass input signal through the example_box |

1 | # Perform the Fourier Tranform on both signals |

1 | # Now, assign each of these variables a (unique) value |

## Task 2

### Marks available: 3

You are to work with `t2_box`

in this task (defined in the first cell of the notebook). It is known that this box contains an **RL high pass** filter circuit, but with parameter values **different** from those of the RL high pass circuit in Task 1.

Employing a similar method to that used in Task 1, examine the characteristic transfer function of this circuit to determine an approximate value for the associated corner frequency, $f_0$. This corresponds to the frequency $f_0$ such that $|T(f_0)| =\sqrt{1/2}$.

The examples in `Yr2_Computing_System_Identification.pdf`

refer to $T(\omega)$ while we sometimes use $T(f)$. These are the same function just plotted over a different axis ($\omega = 2\pi f$). It is easy to make mistakes when applying this scaling. You should always check if you are using the right frequency scaling for a given task.

After you complete you investigation store the value for f_0 in the variable `corner_frequency`

. During marking, your value will be compared to the true value with a tolerance of 10 Hz. Ensure that your sample density in the frequency domain is sufficient to achieve this level of accuracy.

Furthermore, assuming that the resistance R of the circuit is 200 Ohms, calculate the inductance of the inductor in Henrys and store this in the variable `inductance`

. Do **not** round this value.

**Hints**:

- Your time domain input signal should contain frequency components ranging from 10 to 3000 Hz in order to obtain an informative result.
- For finding a specific value in a numpy array you may want to use the
`numpy.argmin`

function.

Play with and understand the following code:

1 | a=np.array([10, 11, 12, 15, 20, 30]) |

- To visually check whether you found the corner frequency correctly you can plot additional lines into the graph with the transfer fucntion, for example using
`plt.axvline(x=f0)`

and`plt.axhline(y=1/np.sqrt(2),linestyle='--')`

.

1 | # Define parameters for Task 2 - explore frequencies up to 3000Hz |

Data points = 30000

Total time = 9.999666666666666 s

Timestep = 0.0003333333333333333 s

Sampling Frequency = 3000.0 Hz

1 | fig = plt.figure(figsize=(15,3)) |

1 | # Pass input signal through the example_box |

1 | # Perform the Fourier Tranform on both signals |

f0 = 1328.7

inductance = 0.023956490267463736

1 | # Your solutions to Task 2 here |

## Resonance in RLC series circuits

The sharp peak in transmission for an RLC series circuit associated with a minimum of impedance is referred to as series resonance, with the **resonance frequency** $f_0$ being that at which the transfer function attains a maximum value (of value 1). With

$$|T(\omega) | = \frac{R}{\sqrt{R^2+(\omega L - 1/\omega C)^2}}$$

this occurs when $\omega L = \frac{1}{\omega C}$ (the inductive reactance equalling the capacitive reactance) and hence

$$\omega_0^2 = \frac{1}{LC} \implies f_0 = \frac{1}{2\pi \sqrt{LC}}$$

For all intents and purposes, a signal is passed by an RLC circuit if it has a frequency contained between the **two** cut-off frequencies, $f_l$ and $f_h$, at which the transfer function has a magnitude of $\frac{1}{\sqrt{2}}$. The **bandwidth** of the circuit (also known as the resonance width) is then defined as the difference between these two frequencies, i.e $\ f_h - f_l$.

## Task 3

### Marks available: 8

You are to work with `t3_box`

in this task (defined in the first cell of the notebook). This contains an **RLC band pass** filter circuit, with values of resistance, inductance and capacitance different from those of the RLC circuit in Task 1.

Firstly, determine the transfer function for this circuit and **plot its magnitude**, $|T(f)|$, against frequency, $f$ (in Hz), on the interval `student_figure`

). You should ensure that your input signal contains frequency components ranging from 10 to

You should use a linear scale on the y-axis (do **not** use dB), and a logarithmic scale on the frequency axis. This is achieved through the use of the `plt.semilogx()`

function. Your plot should also have a title and axis labels.

Using your plotted data or otherwise, estimate the **resonant frequency** and the **bandwidth** of the circuit in Hertz and assign these to the variables `resonant_frequency`

and `bandwidth`

. As in Task 2, these values should be accurate to 10 Hz.

1 | # Please do not create any other figure object |

Data points = 20000

Total time = 9.9995 s

Timestep = 0.0005 s

Sampling Frequency = 2000.0 Hz

fl = 193.3

fh = 614.6

f0 = 344.7

1 | resonant_frequency = f0 # Resonant frequency in Hz |

## Task 4

### Marks available: 5

For the final task you are to create your own black box function called `my_box_process(t, s_in, R, L)`

that mimics the behaviour of an **RL low-pass filter** circuit of resistance `R`

(Ohms) and inductance `L`

(Henrys).

`t`

and `s_in`

will be two numpy arrays representing an input voltage signal with physical amplitude `s_in[i]`

at time `t[i]`

, and the function should return a numpy array containing the real output signal in the time domain.

You will need the transfer function for the RL low-pass filter circuit:

$$T(f) = \frac{1}{1+(f/f_0)j} $$ where $f_0 = R/(2\pi\,L)$.

**Hints**:

- Use
`.fft`

to perform an FFT on the input signal before applying the transfer function in the frequency domain. You will then need to use the`.ifft`

function to represent this modulated signal in the time domain (`ifft`

in the numpy documentation). - You may wish to define the transfer function separately and call this within my_box_process.
- In python,
`1j`

is used to represent the imaginary unit. - Make sure your function returns a physical time domain signal (i.e. real values that can change sign). The
`ifft`

function will return an array of complex numbers with the real part of these numbers containing the signal we are interested in. - BTW because the input and output of your Box must be numpy arrays, you do not need to use any
`for`

loop in your code.

1 | # Function definition for Task 4 |

1 | # Feel free to insert additional cells for working |

N.B. Please execute the following statement before submitting your work. It will check that you have used correct naming for the variables and functions specified in the above tasks. It will not tell you whether you have correctly defined and implemented these!

You may execute this statement as many times as required.

1 | Boxes.check() |

my_box_process function is correctly named.

RL_circuit variable is correctly named.

RC_circuit variable is correctly named.

RLC_circuit variable is correctly named.

corner_frequency variable is correctly named.

inductance variable is correctly named.

resonant_frequency variable is correctly named.

bandwidth variable is correctly named.

RL_circuit, RC_circuit and RLC_circuit are assigned in a valid way.

**Interacting with Black Boxes in Assignment 2**

## Interacting with Black Boxes in Assignment 2

This notebook demonstrates how to interact with the ‘`Boxes`

‘ objects that you will work with in assignment 2. It also includes examples of using the `.fft`

and `.fftfreq`

functions. See also the document `Yr2_computing_System_Identification.pdf`

on Canvas.

Start by executing the code cell below which loads the example `Boxes`

object `example_box`

that represents some unknown electronic circuit.

1 | # Doing all the imports in one go |

** Using temporary ID. Switch to your own ID as soon as possible. **

### Generating an input signal

The `box`

object represents a ‘black box’, i.e. an unknown system. The system accepts an input signal and provides an output signal in return. By generating custom input signals for a specific box and analysing the output you can investigate the unknown system.

For the purposes of this example we generate a signal representing an input voltage over time. The signal is generated as a sine wave whose frequency is slowly changed from f1 to f2. In assignment 2, you should consider creating similar signals to probe the behaviour of each box.

1 | def swept_sine(f1, f2, t): |

In order to generate the input signal we first make an array of times (`t`

) at which the signal should be created.

We want to sample the signal on the interval $[0,10)$ seconds with $N = 1000$ sample points, resulting in a spacing between sample points of`timestep`

$= \frac{10}{1000} = 0.01$ s.

If you think about this in the frequency domain, this choice corresponds to a maximum sampling frequency of $\frac{1}{0.01} = 100$ Hz (i.e. 100 samples per second).

1 | # Sampling a signal correctly is a bit tricky, you need |

Data points = 1000

Total time = 9.99 s

Timestep = 0.01 s

Sampling Frequency = 100.0 Hz

Note the `endpoint=False`

keyword parameter. Without it, `np.linspace`

includes the endpoint which would give the unexpected spacing of $1000/999 = 0.01001…$. (The spacing of any numpy array can also be checked using the `numpy.diff`

function: `numpy.diff(t)`

returns an array containing the difference between successive elements of `t`

.)

We now have an input signal. Let’s have a look at it.

1 | fig = plt.figure(figsize=(15,3)) |

### Probing the box

In the cell, the call `example_box.process(t, s_in)`

passes our voltage signal through the first box and returns the corresponding output signal, stored in the variable `s_out`

. It is essential that the sample points described in `t`

are **equally** spaced.

Both arrays are then plotted on the interval $[0, 5]$ in the time domain.

1 | # Pass input signal through the example_box |

The signal has certainly been modified by `example_box`

, but it is not yet entirely clear *how*. This type of investigation is often simpler and more intuitive when describing the signals as a function of frequency, i.e. by using the *frequency domain*.

To do so, we perform an Fourier Transform (using a DFT algorithm) on both the input and output signals and plot the resulting *spectra*.

As described in `Yr2_Computing_System_Identification.pdf`

, the function `.fftfreq`

is used to determine the relevant sample frequencies, and we choose to plot the positive frequency terms only.

1 | # Perform the Fourier Tranform on both signals |

By dividing the output by the input we can generate very clean data showing the effect of the box independent of the details of the input signal we used. This ratio is called the **transfer function** of the system, sometimes also called **frequency response**.

We can see that the box had an effect similar to that of an ideal *band pass* filter: Our original signal contains frequencies ranging from $1$ to $20$ Hz (see `S_in`

), while the output signal only those from $5$ to $15$ Hz; frequencies outside this range are completely blocked.

You might be curious why our input signal contains significant amplitudes up to frequencies of $20$ Hz even though we specified f2=$10$ Hz. This is because our `swept_sine`

function is just a quick and simple way of creating such a signal. We could instead use a true linear chirp, see the example at the end of this notebook. However, for assignment 2 the exact shape of the signal does not matter.

Additional Notes:

- In order to complete the first task of assignment 2, you will want to perform the Fourier Transform on the input and output signals as above and calculate the transfer function as their element-wise ratio,
`.fft(output_signal)/.fft(input_signal)`

in order to determine the type of circuit contained inside each box. - As stated earlier, be sure to investigate a broad range of frequencies with your input signal: 10 - 4000 Hz should be sufficient for Task 1. Additionally, check that your sample interval and number of sample points ($N$) are such that the sample interval in the frequency domains covers all of these frequencies.
- While the above plots use linear scales, in some cases logarithmic scaling is preferable. For circuits it is also common to plot the transfer function in decibels with a logarithmic horizontal axis. In general: if you are not able to discern anything meaningful from a plot, try using a different axis scaling. For assignment 2, read the instructions carefully and make sure you use the required axis scaling and formatting.

## For the curious: how to make a linear chirp

For the tasks in assignment 2 the exact shape of the input signal is not important as long as it contains non-zero amplitudes in the frequency range of interest. The function `swept_sine`

above is a simple way of creating such an input signal. You will have noted above that the output signal of that function shows non-zero amplitudes up to 2 times the maximum frequency (as shown by the Fourier Transform). If you want to create a true linear chirp signal whose frequency sweeps linearly from $f_1$ to $f_2$ you can use the following function. More information on chirp signals can be found online but the short explanation is that the argument of the `sine`

function is the phase $\phi$ and the frequency of a sine is defined as $1/2\pi\,df/d\phi$. Therefore, to create a frequency that increases linear in time, the phase $\phi$ must be a quadratic function of time.

1 | def linear_chirp(f1, f2, t): |

Let’s try it out:

1 | t = get_time_array(N = 500, duration = 10) |

Data points = 500

Total time = 9.98 s

Timestep = 0.02 s

Sampling Frequency = 50.0 Hz