First, both @LvW and @jonk are right in the comments: the order for bandpass/bandstop filters are a bit of a mess, and for a 2nd order transfer function that is a bandpass (your first formula), the \$Q\$ dictates the bandwidth so, if you set it to be \$\frac{1}{\sqrt2}\$ then the -3 dB edges will be at \$\frac{1}{\sqrt2}\$ Hz and \$\sqrt2\$ Hz (relative to 1 Hz center frequency).
However, further on you say two things:
- you expected to have a resulting 4th order transfer function
- you're wondering what the \$\epsilon\$ (
\epsilon
for MathJax) is
For 1. it's simple: start from the lowpass prototype and apply the frequency transformation:
$$H_{\text{LP}}(s)\rightarrow H_{\text{BP}}\left(\dfrac{s^2+\omega^2}{BW s}\right) \tag{1}$$
Funny enough, today I wrote this answer where I show the meaning of that transformation:
$$\dfrac{s^2+\omega^2}{BW s}=\dfrac{s}{BW}+\dfrac{1}{s\frac{BW}{\omega^2}} \tag{2}$$
So now apply (1) and you'll get the transfer function of the bandpass to be made of a lowpass and a highpass:
$$\begin{align}
H(s)&=\dfrac{\omega_0^2}{s^2+\dfrac{\omega_0}{Q}s+\omega_0^2} \\
H\left(\dfrac{s^2+\omega^2}{BW s}\right)&=\dfrac{BW^2\omega_0^2s^2}{s^4+\dfrac{\omega_o}{Q}BWs^3+(2\omega_c^2+BW^2\omega_0^2)s^2+\dfrac{\omega_0}{Q}BW\omega_c^2s+\omega_c^4} \tag{3}\\
&=\dfrac{\omega_1^2}{s^2+\dfrac{\omega_1}{Q_1}s+\omega_1^2}\cdot\dfrac{s^2}{s^2+\dfrac{\omega_2}{Q_2}s+\omega_2^2} \tag{4}
\end{align}$$
Or, if you started from the poles of the lowpass prototype, convert them directly to bandpass with the following (which can be obtaned by solving for the roots of the denominator of (3)):
$$\begin{align}
s_{_\text{LP}}&=\sigma\pm j\omega \\
M&=|s_{_\text{LP}}| \\
\Phi&=\angle{s_{_\text{LP}}} \\
&\begin{cases}
s{1,2}=\dfrac12\left\{-BW\sigma+\sqrt{M}\cos\left(\dfrac{\Phi}{2}\right)\pm j\left[BW\omega-\sqrt{M}\sin\left(\dfrac{\Phi}{2}\right)\right]\right\} \\
s{3,4}=\dfrac12\left\{-BW\sigma-\sqrt{M}\cos\left(\dfrac{\Phi}{2}\right)\pm j\left[BW\omega+\sqrt{M}\sin\left(\dfrac{\Phi}{2}\right)\right]\right\}
\end{cases} \tag{5}
\end{align}$$
Either way, you'll end up with (4), which is a 4th order transfer function which is a bandpass, or a 2nd order bandpass (implying a multiple a 2nd orders). Yes, it's confusing but, the most certain way to be sure about it is to specify it as a 4th order transfer function. This way, you are certain that the denominator has an \$s^4\$, no more, no less.
And about 2., this is related to the generic way in which the Butterworth mathematical transfer function is written:
$$\begin{align}
H(s)&=\dfrac{1}{\sqrt{1+\epsilon^2\left(\dfrac{\omega}{\omega_0}\right)^{2N}}} \tag{6} \\
\epsilon&=\sqrt{10^{\frac{A_{\text{dB}}}{10}}-1}
\end{align}$$
That \$\epsilon\$ dictates the attenuation at \$\omega_0\$. If the attenuation is to be the "classical" -3(.0103...) dB then:
$$\begin{align}
A_{\text{dB}}&=-20\log_{10}\left(\dfrac{1}{\sqrt2}\right)\approx -3\;\text{dB} \\
\epsilon&=\sqrt{10^{\frac{3}{20}}-1}=1 \tag{7}
\end{align}$$
And, since in the majority of the cases \$\epsilon\$ doesn't need tinkering, it's considered 1. But nobody says that the -3 dB value is set in stone -- that is up to the requirements for the filter. For example, if the attenuation at \$\omega_0\$ needs to be 0.9 (0.915 dB) then \$\epsilon=0.4843...\$.
It puzzles me that it appears with a \$\frac1N\$ as power, in the OP. it shouldn't since (6) shows it to be squared, not a power of N. But they probably use a different notation an use it to extract the custom attenuation, which is:
$$\omega_x=\dfrac{\omega_0}{\left(\sqrt{10^{\frac{A_{\text{dB}}}{10}}-1}\right)^\frac1N}\quad\left[=\dfrac{\omega_0}{\epsilon^\frac1N}\right] \tag{8}$$
At any rate, the -3 dB convention is (mathematically) only valid for 1st orders and Butterworth filters, since that is the point that coincides with the phase being half its final (asymptotic) value. But, ultimately, the design dictates what frequency scaling needs to be applied and where (in most cases that's not needed).