How do we determine the damping coefficient given the acceleration vs time graph of a 1D mechanical system?

The question in the OP is solved by tools referred to as second order system identification (ID) techniques in the control theory.

---

The following details the system ID technique of a second order system with the vanishing initial condition from the time domain data of it's unforced response $u=0$, rather than from the step or impulse response data. Even in the case of the step response, the unforced response data is available after the initial time interval (incremental in the case of the impulse input) in which the step or impulse input is applied.

1. The data response of the system provided as accelerations as a function of time for a system can be time integrated to obtain data of the velocity and position as functions of time. The damped angular frequency of the system with the position response akin to that shown in the plot in the OP (which is readily determined to be an under-damped second order system based on the analysis presented), can be obtained by $$omega_d=frac{2pi}{T_d}text{ where } T_d:=frac{t_k-t_1}{k-1},$$ wherein $t_k$ is the time instant at which the $k^{th}$ maximum or peak or crest (or minimum or negative peak or trough) of the $k^{th}$ oscillation (in other words $t_k-t_1$ is the time duration which is the sum of time periods of $k-1$ oscillations). Next, the difference between natural logarithm of two (say) peaks (or negative peaks) or maxima (or minima) $k$ peaks (or negative peaks) apart, given by $ln(y^p(t_k))=-xiomega_n t_{k} + ln(Csin(omega_d t_k+tilde{phi}))$ and $ln(y^p(t_{k-1}))=-xiomega_n t_{k-1} + ln(Csin(omega_d t_{k-1}+tilde{phi}))$ is $ln(y^p(t_k))-ln(y^p(t_{k-1}))=-xiomega_n (t_k-t_{k-1})$ since $sin(omega_d t_{k-1})+tilde{phi})=sin(omega_d t_{k}+tilde{phi}))$, so that $$xiomega_n=frac{ln(y^p(t_{k-1}))-ln(y^p(t_k))}{t_k-t_{k-1}}=frac{ln(y^p(t_{k-1}))-ln(y^p(t_k))}{frac{2pi}{omega_d}},$$ where the time interval $[t_{k-1},t_k], ;0<t_{k-1}<t_k$ indicate the time period $T_d$ of the damped oscillatory response. Note that if we instead used the time instants $t_1$ and $t_k$ for this calculation, the relationship between the values of the logarithm of the two peak (or negative peak) values is given as $xiomega_n=dfrac{ln(y^p(t_{1}))-ln(y^p(t_k))}{frac{2pi(k-1)}{omega_d}}$. The natural angular frequency and subsequently, the damping coefficient, are then determined using the values of the expressions obtained above and the definition which states their relationship as $$omega_d=omega_nsqrt{1-xi^2}.$$
2. An alternate solution is to use the method of least squares to estimate the parameters of the dynamical model by directly using the time domain data.

How to Fit

We have data points the oscillate curve $~x(t)~$.

Step 1: obtain the blue points

The period of this curve is $~T=2~$, hence the points are

$$X=T_0+n,T~,Y=~x(T_0+n,T)~,n=0,(1),N_p,$$ where $~T_0=frac T4.$

Step II

Fit the points $[~X~,Y]~$ to the curve $~a,e^{-beta,t}~$ where $~a=x(T_0)$ and you obtain $~beta~[1/s]$ from the fitting process.

The blue curve is your result.