How to identify / measure the acceleration of a chirp?

Question

Which methods do there exist to identify the (angular) acceleration of a chirp function?

In other words, given some sampling of the curve for example $\{(t_1,f(t_1)),\cdots(t_n,f(t_n))\}$, can we estimate $k$ in $$f(t) = \sin(kt^2)$$? Bonus points will be awarded if the method is completely linear (in some sense).

Answer 1

I wonder if the question is well posed. If you want fit the function $$f(t)=\sin(kt^2)$$ with a given data $$\{(t_1,f(t_1)),\cdots(t_n,f(t_n))\}$$ simply change the data with $\quad\begin{cases} x_1=t_1^2 , \cdots, x_n=t_n^2 \\ y_1=f(t_1) , \cdots, y_n=f(t_n) \end{cases}$

and fit the function $$y(x)=\sin(kx)$$ to the new data $$\{(x_1,y_1),\cdots(x_n,y_n)\}$$ This is a regression for a simple sinusoidal function. They are a lot of such problems in the literature, even with more complicated sinusoidal functions such as $\quad y(x)=a+b\sin(\omega x+c)\quad$ which in your case corresponds to $\quad f(t)=a+b\sin(k t^2+c)$.

In addition :

Without the original data an accurate calculus isn't possible. Nevertheless one can scan the figure edited by mathreadler in order to create an approximate data.

From the approximate data, the sinusoidal regression lead to $$k\simeq 127.85$$

With this approximate value of $k$ the function $f(t)=\sin(kt^2)$ is drawn in red on the original figure :

Latter addition :

With an improved scanning of the graph published in the question (675 points each from one pixel of the curve) the result is slightly different :
$$k\simeq 128.01$$

Note :

The sinusoidal regression used to fit the function $y=\sin(kx)$ is a simplified version of the method from https://fr.scribd.com/doc/14674814/Regressions-et-equations-integrales , pp.21-36.

Since there is only one parameter $k$ to optimize, the integral equation is simplified : $$y(x)=\sin(k x)\quad\to\quad \int_0^x y(\xi)d\xi=\frac{1}{k}(1-\cos(k x))$$ $$\left(1+k\int_0^x y(\xi)d\xi \right)^2=1-y^2$$ The numerical integration gives $ S_n\simeq \int_{x_1}^{x_n} y(\xi)d\xi\quad$ from $n=1$ to $N=$number of points.

The numerical calculus is very simple :

$ S_1=0\quad;\quad S_n=S_{n-1}+\frac12(y_n+y_{n-1})(x_n-x_{n-1}) \quad $ from $n=2$ to $n=N$.

Then $k$ is the real solution of the cubic equation : $$c_3k^3+c_2k^2+c_1k+c_0=0 \quad\begin{cases} c_3=\sum_{n=1}^N (S_n)^4\\ c_2=-3\sum_{n=1}^N (S_n)^3 \\ c_1=\sum_{n=1}^N \left((2+(y_n)^2\right)(S_n)^2 \\ c_0=-\sum_{n=1}^N (y_n)^2 \end{cases}$$

Answer 2

I would use a Hilbert transform approach. This is commonly implemented with a DFT algorithm. One would perform the following steps:

1. Compute DFT of the input signal
2. Set the negative frequency components to zero (this is the upper half of the vector)
3. Compute IDFT (the result vector is now complex. The real part is the original signal, namely $\sin(kt_i^2)$)
4. Calculate argument of each point (This gives $\phi_i$ which is ideally equal to $k t_i^2 \text{ mod } 2 \pi$)
5. Perform phase unwrapping (we obtain $\phi_i$)
6. Perform polynomial linear regression to find $k$

With exception of 4) all steps are linear.

In practice this method works quite well.

However, several problems might exist:

• Step 1-3 do not give the exact Hilbert transfom, but only an approximation. This is because the signal is not infinitely long (which the Hilbert transform would require) and second the signal is sampled.
• The phase unwrapping might fail if the input vector is not adequately sampled (i.e. less than two samples per period). However, in this case the Nyquist condition is violated anyway and one cannot expect to reconstruct $k$.
• In step 6) I would include an unknown offset in the polynomial fit, as one probably cannot expect from the procedure that $\phi_0 = 0$.