So first the intuition: these functions decay quickly at infinity so their Fourier transforms are regular. However the functions are not smooth at $x=0$, which implies that their Fourier transforms do not decay so fast at infinity. So, to know if they are in $L^1$, we need to isolate the singular part and see what decay it gives.
Now the proof. Take $\varphi \in C^\infty_c$ a smooth function supported on the unit ball such that $0\leq \varphi \leq 1$ and $\varphi=1$ around $0$. Then we write
$$
f = f_0 + f_1 = \varphi \, f + (1-\varphi)\,f.
$$
It is not difficult to check that $f_1$ is a Schwartz function, and so its Fourier transform is $L^1$ (since it is also a Schwartz function).
Now it remains to look at $f_0$, which is compactly supported but not infinitely smooth. Using the fact that the Fourier transform of a product gives a convolution and a second order Taylor series for the exponential,
$$
\widehat{f_0} = \sum_{k=0}^\infty \frac{(-1)^k}{k!} \widehat \varphi * \widehat{|x|^{k\alpha}} = \widehat\varphi + \sum_{k=1}^\infty \frac{(-1)^k}{k!} \widehat \varphi * \mathrm{pf}(|x|^{-n-k\alpha})
$$
where $\mathrm{pf}(|x|^{-c})$ is the Hadamard finite part of $|x|^{-c}$, which is essentially the tempered distribution that coincides with $|x|^{-c}$ outside of the point $0$, and is what appears when doing the Fourier transform of fractional polynomials (see e.g. here on MSE). Convolution with these functions could also be written as a fractional Laplacian. The point is that $\widehat\varphi$ being a Schwartz function, it is not difficult to prove that the functions
$$
\psi_c = \widehat\varphi * \mathrm{pf}(|x|^{-c})
$$
are smooth and behave like $|x|^{-c}$ for large values of $x$. Hence
$$
\widehat{f_0} = \widehat \varphi + \sum_{k=1}^\infty \frac{(-1)^k}{k!} \psi_{n+\alpha k}.
$$
This is a (uniformly converging) sum of smooth functions that decay faster and faster when $k$ grows. Hence the decay of $f_0$ is the same as the decay of the first one, that is $\widehat{f_0}\sim |x|^{-n+\alpha}$ at infinity, and so is integrable.
The Fourier transform of the functions $g_j$ are just the derivatives of $\widehat f$ and so decay even faster. Hence they are also integrable.
Remark: here is another method perhaps more elementary, but which does not give the asymptotic behavior. It uses the fact that $H^s\subset \mathcal F(L^1)$ when $s>n/2$. Indeed, by the Cauchy-Schwartz inequality
$$
\int_{\Bbb R^n} |\widehat{f}| \leq C_s \left(\int_{\Bbb R^n} |\widehat{f}|^2 (1+|2\pi\,x|^{2s})\,\mathrm d s\right)^\frac{1}{2} = C_s \left(\|f\|_{L^2}^2 + \|f\|_{\dot{H}^s}^2\right)^\frac{1}{2}
$$
with $C_s = \int_{\Bbb R^n} \frac{\mathrm d x}{1+|2\pi\,x|^{2s}} < \infty$. Of course $f$ is in $L^2$. Now, taking $k$ derivatives with $k> n/2$ it is easy to get that
$$
|\nabla^k(e^{-|x|^\alpha})| \leq C \left(|x|^{k\alpha-k} + |x|^{\alpha-k}\right) e^{-|x|^\alpha}.
$$
Up to removing some decay in the exponential, we can just keep the worst local term, that is
$$
|\nabla^k(e^{-|x|^\alpha})| \leq C\, |x|^{\alpha-k}\, e^{-|x|^\alpha/2}.
$$
Hence,
$$
\|f\|_{\dot H^k} \leq C \int_{\Bbb R^n} \frac{e^{-|x|^\alpha}}{|x|^{2k - 2\alpha}} \,\mathrm d x
$$
which is finite as soon as $2k - 2\alpha < n$, that is $k< n/2 + \alpha$. Of course, here one needs $n/2 < k < n/2 + \alpha$ with $k$ an integer, which is not always possible (it works if $n$ is odd and $\alpha > 1/2$). With more work, the assumption that $k$ is an integer can be relaxed (by working with a fractional Laplacian for example, or via Littlewood-Paley techniques).
