I would like to recreate the figure on Wikpedia's article on Laplace's method, which is a Laplace approximation where $f(x) = \sin(x)/x$:
The global maximum is $x_0 = 0$ and $\lim_{x \rightarrow 0} f(x_0) = 1$. The second derivative is $1/3$ again in the limit. Let's fix $M = 1$. Therefore, the Gaussian approximation of the integral is
$$ \int_a^b e^{M f(x)} \approx \sqrt{\frac{2 \pi}{M |f''(x_0)|}} e^{M f(x_0)} = \sqrt{6 \pi} e $$
So this is the normalizer for a Gaussian approximation. Then the full approximation should be
$$ \frac{1}{\sqrt{6 \pi}} \exp \big(-\frac{1}{6} x^2 \big) $$
But when I plot this, I don't get a good approximation:
What am I doing wrong?
We're approximating $\exp f$ as $\exp (f_0+\frac12 f^{\prime\prime}_0(x-x_0)^2)$ with $f_0:=f(x_0),\,f^{\prime\prime}_0:=f^{\prime\prime}(x_0)$ with $f^\prime(x_0)=0$. Judging by the plots you chose to compare, you want $f(x):=\operatorname{sinc}x=1-\frac16 x^2+o(x^2)$ so $$x_0=0,\,f_0=1,\,f^{\prime\prime}_0=-\frac13,\,\exp f\approx e\cdot\exp\left(-\frac16 x^2\right).$$I think your analysis earlier in the question, however, may have instead taken $f:=\ln\operatorname{sinc}x$, which has the distinct advantage that $\int_{\Bbb R}f(x)dx$ will then be finite, as you'd expect for a decent Gaussian approximation. Since this choice of $f$ has Taylor series $-\frac16 x^2+o(x^2)$, the approximation this time is $\exp\left(-\frac16 x^2\right)$. This integrates to $\sqrt{6\pi}\approx 4.34$, which isn't that far from $\int_{\Bbb R}\operatorname{sinc}xdx=\pi\approx 3.14$.
Of course, this raises the question of which $k$ makes $k\exp\left(-\frac16 x^2\right)$ the "best" approximation (by some criterion or another) of $\operatorname{sinc}x$. We can choose $k=\sqrt{\frac{\pi}{6}}$ if we want the same area under the curve, $k=1$ if we want the functions to match at $x=0$, or some third value (no doubt) to minimize the Euclidean cost $\int_{\Bbb R}\left|k\exp\left(-\frac16 x^2\right)-\operatorname{sinc}x\right|^2dx$. Whichever approach you take, the Gaussian approximation will look fairly similar to $\operatorname{sinc}x$, but the approximation won't be excellent because the tails of $\operatorname{sinc}x$ decay so slowly. If you use these techniques for $\operatorname{sinc}^Mx$ for some large positive integer $M$, on the other hand, you'll get more impressive results.