Suppose a "2-state atom" and a light field are quantized with the following Hamiltonians, respectively: $$\hat{H}_A=\hbar\omega_{21}\hat{\sigma}^{\dagger}\hat{\sigma}$$ and $$\hat{H}_R=\sum_{\textbf{k}}\hbar\omega_{\textbf{k}}(\hat{a}^{\dagger}_{\textbf{k}}\hat{a}_{\textbf{k}} + \frac{1}{2})\ .$$ Where $\hat{\sigma}^{\dagger}=\left|2\right\rangle\left\langle1\right|$ and $\hat{\sigma}=\left|1\right\rangle\left\langle2\right|$, where $\left|1\right\rangle,\left|2\right\rangle$ are the 2 states of the atom and $\omega_{21}$ is for the transition from state 1 to state 2. $\textbf{k}$ are the modes of the light field, and $\hat{a}^{\dagger}_{\textbf{k}},\ \hat{a}_{\textbf{k}}$ are the usual creation and annihilation operators.
If the interaction of the atom and the light field is modeled using a dipole moment with contribution to the total Hamiltonian of: $$\hat{H}_I=\sum_{\textbf{k}}\hbar g_{\textbf{k}}(\hat{a}^{\dagger}_{\textbf{k}} + \hat{a}_{\textbf{k}})(\hat{\sigma}^{\dagger}+\hat{\sigma})\ .$$
Can I simplify this to pick out of the summation only the mode that matches the transition energy $\hbar\omega_{21}$? This would be because photons are emitted/absorbed entirely anyway, so the energies have to match, and the photons which do not have the right energy simply won't interact with the atom (that's what I'm trying to understand, I might be completely wrong!).
Best Answer
That's exactly right. Let's see it in detail.
Slightly nicer notation
$\newcommand{\ket}[1]{\left \lvert #1 \right\rangle}$ First let's make the notation slightly more standard. We have $\hat{\sigma} \equiv |2\rangle\langle1|$ so that in the $|2\rangle, |2\rangle$ basis, $\hat{\sigma}^\dagger \hat{\sigma}$ has matrix representation $$\left( \begin{array}{cc} 0 & 0 \\ 0 & 1 \end{array} \right) \, .$$ We're always free to put in an energy offset so we can shift to $$\left( \begin{array}{cc} -1/2 & 0 \\ 0 & 1/2 \end{array} \right) \, .$$ which is conveniently written as $- \sigma_z / 2$. Therefore, our system's Hamiltonian is$^{[a]}$
\begin{align} H / \hbar =& H_A/\hbar + H_R/\hbar + H_I/\hbar \\ =& -\frac{\omega_{21}}{2} \sigma_z + \sum_k \omega_k (a_k^\dagger a_k + 1/2) + \sum_k g_k (a^\dagger_k + a_k) \sigma_x \, . \end{align}
Interaction picture
This Hamiltonian is kind of complicated because there are so many interaction terms. We can enormously simplify it by using the so-called "interaction picture". The idea is to "unwind" the evolution of the atom and fields caused by their own Hamiltonians.
In the interaction picture formalism, you get a new effective Hamiltonian and the operators pick up a time dependence. You construct the propagators of the atom and fields without the coupling:
$$U_A(t) \equiv \exp \left( -i H_A t / \hbar \right) = \exp \left(i \omega_{21} \sigma_z t / 2 \right)$$
and
$$U_R = \prod_k \exp \left(-i \omega_k t (a_k + a_k^\dagger + 1/2) \right) \, .$$
Then you find the time dependent operators. For each operator $T$, it's time dependent version is
$$V(t) = U^\dagger(t) T U(t)$$
where here $U(t)$ means $U_A(t) \otimes U_R(t)$. The interaction picture Hamiltonian is only the interaction part $H_I$, but with the time dependent versions of the operators.$^{[b]}$ So, all we have to do is figure out what $\sigma_x$ and $a$ look like with this time dependence.
You can do the algebra and you'll find that for the atom $$ U_A(t)^\dagger \sigma_x U_A(t) = \cos(\omega_{21}t)\sigma_x + \sin(\omega_{21}t) \sigma_y = e^{i \omega_{21} t}\sigma_+ + e^{-i \omega_{21} t}\sigma_- $$ where $\sigma_+ \equiv \left \lvert 2 \right\rangle \left \langle 1 \right \rvert$ and $\sigma_- = \sigma_+^\dagger$. You can also find that in the interaction picture $a_k \rightarrow a e^{-i \omega_k t}a$.
Putting it all together, in the interaction picture our Hamiltonian is
\begin{align} H'/\hbar =& \sum_k g_k (a_k^\dagger e^{i \omega_k t} + a_k e^{-i \omega_k t}) (\sigma_+ e^{i \omega_{21} t} + \sigma_- e^{-i \omega_{21}t}) \\ =& \sum_k g_k \left( a_k^\dagger \sigma_+ e^{i (\omega_k + \omega_{21})t} + a_k^\dagger \sigma_- e^{i (\omega_k - \omega_{21})t} + a_k \sigma_+ e^{i (-\omega_k + \omega_{21})t} + a_k \sigma_- e^{i (-\omega_k - \omega_{21})t} \right) \, . \end{align}
Rotating wave approximation
Of these four terms, two of them have fast oscillating time dependencies because the frequencies in the exponent are either both negative or both positive. In a lot of cases these fast oscillating time terms contribute much less to the dynamics than the slower terms. This is roughly because the integral of an oscillating function stays small, while the integral of a constant function gets large. Anyway, you can often drop those terms, leaving
\begin{align} H'/\hbar =& \sum_k g_k \left( a_k^\dagger \sigma_- e^{i (\omega_k - \omega_{21})t} + a_k \sigma_+ e^{-i (\omega_k - \omega_{21})t} \right) \, . \end{align}
Doing this is called the "rotating wave approximation. Now we've got a sum over $k$ where each term has a frequency $\omega_k - \omega_{21}$. Again, if the modes are spaced out, then you can have a situation where most of these frequencies are large, but there could be one radiation mode $l$ resonant with the atom, i.e. $\omega_l \approx \omega_{21}$. Then, that term dominates and we get
$$H'/\hbar \approx g_l \left( a_l^\dagger \sigma_- e^{i (\omega_l - \omega_{21})t} + a_l \sigma_+ e^{-i (\omega_l - \omega_{21})t} \right) \, .$$
If the frequencies are exactly the same (i.e. the atom and radiation mode $l$ are exactly on resonance), then you get
$$H'/\hbar \approx g_l \left( a_l^\dagger \sigma_- + a_l \sigma_+ \right) \, .$$
In the subspace of radiation Fock states $\ket{n}$ and $\ket{n+1}$, this Hamiltonian is just proportional to $\sigma_x$, i.e. the quantum of energy oscillates between the atom and the radiation field.
Review
Your intuition about the atom only interacting with radiation modes of the same frequency is correct. The "rotating wave approximation" is the name given to the approximation where you drop all the radiation modes that aren't resonant with the atom (of course, this works for things other than atoms and radiation). The terms we dropped are ones that do not conserve energy. It turns out they actually do have a role in real systems and some times they're very important, but in a lot of cases you can ignore them and use the super simple Hamiltonian we found.
$[a]$: Note that by working with $H/\hbar$, which has dimensions of frequency, all the $\hbar$'s on the right hand side are gone. This is really convenient and it's what people are really doing when they say that they're "setting $\hbar$ to 1".
$[b]$: See this other Physics SE post for a demonstration of why.