<jupyter_start><jupyter_text>Curve FittingIn the last notebook, we interactively adjusted the inertia and damping such that the simulation trajectory matched the measured data. Another way to determine the period of oscillation more accurately is to try and find the best fit curve through the data points. The Python package `scipy` provides a very convenient function called [`curve_fit`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html). This function implements a [least squares method](https://en.wikipedia.org/wiki/Least_squares) that finds an optimal fit based on parameterized function provided by the user. Let's see how it works.  Load the data<jupyter_code>import pandas as pd
radial_gyro_meas = pd.read_csv('bicycle-wheel-radial-inertia-rate-gyro-measurement.csv', index_col='time')
%matplotlib widget
ax = radial_gyro_meas.plot(style='.')<jupyter_output><empty_output><jupyter_text>Example: Fit a cosine to the data<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from utils import period2freq, freq2period<jupyter_output><empty_output><jupyter_text>First define a Python function where the first input argument are the "x" values (an array) and the remaining arguments are the parameters that can be adjusted to find the optimal fit. In this case, a simple cosine function with an adjustable amplitude and frequency is created.$$ \dot{\theta} = A \cos \omega t $$<jupyter_code>def cos_func(times, amplitude, frequency):
    return amplitude * np.cos(frequency * times)<jupyter_output><empty_output><jupyter_text>To check that it works, provide some values:<jupyter_code>times = np.linspace(0, 10)
cos_func(times, 5.0, 5.0)
fig, ax = plt.subplots()
ax.plot(times, cos_func(times, 5.0, 5.0))
fig, ax = plt.subplots()
ax.plot(radial_gyro_meas.index, cos_func(radial_gyro_meas.index, 5.0, 5.0))<jupyter_output><empty_output><jupyter_text>Now, provide this function to `curve_fit` along with the measure data (x and y) and an initial guess for the amplitude and frequency. A good initial guess is important, as the optimal solution can't always be found from an arbitrary initial guess. The function `curve_fit` returns two items. The first is the optimal values of the two parametes and the second is the covariance matrix that gives an idea of how certain the value of the parameters are. We will just work with the first value for now.<jupyter_code>popt, pcov = curve_fit(cos_func,  # our function
                       radial_gyro_meas.index,  # measured x values
                       radial_gyro_meas.angular_velocity,  # measured y values
                       p0=(3.0, period2freq(0.44)))  # the initial guess for the two parameters<jupyter_output><empty_output><jupyter_text>Now we see the optimal values for the amplitude and frequency:<jupyter_code>popt<jupyter_output><empty_output><jupyter_text>It is useful to plot the optimal function versus the measured data:<jupyter_code>fig, ax = plt.subplots(1, 1)
ax.plot(radial_gyro_meas.index, radial_gyro_meas.angular_velocity, '.', label='Measured')
ax.plot(radial_gyro_meas.index, cos_func(radial_gyro_meas.index, popt[0], popt[1]), label='Best Fit')
ax.legend()<jupyter_output><empty_output><jupyter_text>That looks like a pretty nice fit! The period of this fitted function can now be computed:<jupyter_code>freq2period(popt[1])<jupyter_output><empty_output><jupyter_text>ExerciseThis system, just like the book on a cup system has an amplitude that decreases with respect to time. This decrease is often referred to as a *decaying amplitude*. We just fit the mathematical function:$$ \dot{\theta}(t) = A \cos{\omega t} $$where $A$ is the amplitude of the angular rate $\dot{\theta}$ and $\omega$ is the frequency of oscillation in radians per second. One way to account for the decay in amplitude oscillation is to introduce a multiplicative factor of $e^{-\lambda t}$ into the equation. This will cause exponential decay at the rate $\lambda$ (units are 1 / s). The mathematical equation will look like:$$ \dot{\theta}(t) = A e^{-\lambda t} \cos{\omega t} $$Recall that $e^{-\lambda t}$ looks like:<jupyter_code>fig, ax = plt.subplots()
t = np.linspace(0, 2, num=100)
ax.plot(t, np.exp(-1.5 * t));<jupyter_output><empty_output><jupyter_text>Use decaying oscillation mathematical function to create a curve fitting function and find the values of $A$, $\lambda$, and $\omega$ that best fit the data. Calculate the period of oscillation and compare it to the period from the purely sinusoidal fit from above. Is there any difference in the period of oscillation?<jupyter_code># write solution here
def decaying_sinusoid(t, a, lam, w):
    return a * np.exp(lam * t) * np.cos(w * t)

popt, pcov = curve_fit(decaying_sinusoid,
                       radial_gyro_meas.index,
                       radial_gyro_meas.angular_velocity,
                       p0=(3.0, -0.0002, period2freq(0.44)))

fig, ax = plt.subplots(1, 1)
ax.plot(radial_gyro_meas.index, radial_gyro_meas, '.')
ax.plot(radial_gyro_meas.index, decaying_sinusoid(radial_gyro_meas.index, popt[0], popt[1], popt[2]));

freq2period(popt[2])<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09535734925833173, 'fraction_numerical': 0.010980543247929108, 'mean_word_length': 4.981566820276497, 'pattern_counts': {'":': 0, '<': 41, '<?xml version=': 0, '>': 41, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '428252', 'license': 'MIT', 'n_tokens_mistral': 1664, 'n_tokens_neox': 1585, 'n_words': 564, 'path': 'notebooks/03-2020/curve_fitting.ipynb', 'repo_name': 'gbrault/resonance', 'size': '5191'}