Atomic spectra are defined as the spectrum of the electromagnetic radiation emitted or absorbed by an electron during transitions between different energy levels within an atom.
Twolevel atom
Classical Model

$N_1$ : Number of atoms in $\ket{g}$ per unit volume.

$N_2$ : Number of atoms in $\ket{e}$ per unit volume.

$\Gamma$ : Rate of spontaneous emission. The probability of transitioning of a single atom from $\ket{e}$ to $\ket{g}$ due to spontaneous emission per unit time.

$I$ : Intensity. $I = \frac{1}{2} c n^2 \varepsilon_0 E^2$. Intensity has a dimension of $\mathrm{J\cdot m^{2} \cdot s^{1}}$.

$\sigma(\omega)$ : Absorption crosssection, defined such that $\sigma(\omega) I$ is the energy absorbed by a single atom per unit time when irradiated by intensity $I$ with frequency $\omega$ .

$\Gamma_\text{SR}$ : Rate of stimulated radiation. The probability of stimulated absorption equals the probability of stimulated emission, so that we can define such a universal rate.
Einstein coefficient
The change in number of the excited atoms $\frac{\mathrm{d}N_2}{\mathrm{d}t}$ (per unit volume) comes from three sources:

Spontaneous emission: In a unit time, there is a probability of $\Gamma$ for a single atom to transition to $\ket{g}$. For $N_2$ atoms, the total transition number is $N_2 \Gamma$ on average.

Stimulated absorption: According to the definition of absorption crosssection, $\sigma(\omega) I$ is the energy absorbed by a single atom per unit time. Since the incident beam has a frequency $\omega$, each time the medium absorb a photon $\hbar \omega$ , an atom will be excited into $\ket{e}$ (although the energy gap is $\hbar \omega_0$) . Therefore, in a unit time, $\sigma(\omega) I N_1 / \hbar \omega$ is the total number of atoms in the $\ket{g}$ state being excited into $\ket{e}$.
/* Note: $\Gamma_\text{SR} = \sigma(\omega) I / \hbar \omega$ */

Stimulated emission: Since stimulated absorption and stimulated emission have the equal probability, the mechanism is similar. We have a negative term $\frac{\sigma(\omega) I}{\hbar \omega} N_2$ .
We write down
$\frac{\mathrm{d}N_2}{\mathrm{d}t} = N_2 \Gamma \frac{\sigma(\omega) I}{\hbar \omega} (N_2  N_1)$
For steady state, spontaneous emission + stimulated emission = stimulated absorption, $\mathrm{d} N_2/\mathrm{d}t = 0$. We obtain the steady state equation
$\frac{\sigma(\omega) I}{\hbar \omega} (N_2  N_1) = N_2 \Gamma$
Thus, the steady state of a twolevel system always has $N_2 < N_1$, no matter how strong the incident beam is. Which means, it is impossible to use a twolevel system as the gain medium of laser, since optical pumping cannot induce population inversion.
Intensity gain
When a laser beam with a frequency of $\omega$ and an intensity of $I$ is incident on the medium, we have
$\frac{\mathrm{d}I}{\mathrm{d} z} = \alpha(\omega) I$
where $\frac{\mathrm{d}I}{\mathrm{d} z}$ is the gain (caused by pure emission) of energy per unit volume per unit time. At the same time, $\sigma(\omega) I(N_2  N_1)$ is also the gain of energy of incident beam per unit volume per unit time, where $N_2  N_1$ means stimulated emission minuses stimulated absorption, leaving the pure gain of energy of incident beam. So
$\alpha(\omega) = \sigma(\omega)(N_2  N_1)$
Combining $N_2 + N_1 = N$ and the steady state equation, we could finally obtain
$\alpha(\omega) = \frac{\sigma(\omega) N}{1 + \dfrac{2\sigma(\omega)I}{\hbar \omega \Gamma}}$
Line shape function
We have found that for a fixed number of particles $N_1,N_2$ and a constant incident light intensity $I$, the energy absorbed by the medium per unit time varies depending on the frequency of the incident light $\omega$. This is reflected in the dependence of absorption crosssection $\sigma(\omega)$ on $\omega$. Furthermore, we can extract a dimensionless line shape function $g(\omega)$ from it, describing the absorption efficiency of the medium for different frequencies of light.
Quantum Model
 $\rho_{gg}$ : The probability of a single atom on the ground state.
 $\rho_{ee}$ : The probability of a single atom on the excited state.
By solving the Optical Bloch Equation, we can obtain
$\rho_{ee} (t\to \infty) = \frac{\Omega^2 / \Gamma^2}{1+\dfrac{4\Delta^2}{\Gamma^2}+\dfrac{2\Omega^2}{\Gamma^2}}$
Define a saturation parameter $S_0 = 2\Omega^2/\Gamma^2 = I/I_\text{sat}$, the density matrix element becomes
$\rho_{ee} (t\to \infty) = \frac{S_0/2}{1+S_0+\dfrac{4\Delta^2}{\Gamma^2}}$
When onresonance and large intensity $S_0 \gg 1$, $1/2$ of the atoms occupy the excited state, reaches its maximum. It proves again that the steady state of a twolevel system always has $N_2 < N_1$, no matter how strong the incident beam is.
Resonance fluorescence
We will now consider the radiation due to a single, isolated atom driven by a monochromatic field. The radiation, viewed as scattered light, is called resonance fluorescence.
The optical WienerХи́нчин theorem writes
$I_\text{sc}(\boldsymbol{r},\omega) = \frac{1}{2\pi} c n^2 \varepsilon_0 \int_{\infty}^{+\infty} \langle E^{()}(\boldsymbol{r},t) E^{(+)}(\boldsymbol{r},t+\tau) \rangle \, \mathrm{e}^{\mathrm{i}\omega\tau} \, \mathrm{d} \tau$
where $\hat{E}$ contains $\sigma, \sigma^\dagger$. This will become
$I_\text{sc}(\boldsymbol{r},\omega_\text{sc}) = \frac{\hbar \omega_0 \Gamma}{r^2} f(\theta,\phi) S(\omega_\text{sc})$
where $S(\omega_\text{sc})$ is the spectrum of the scattered light with a dimension of $\mathrm{Hz^{1}}$, which describes the probability of the fluorescence having a frequency between $\omega_\text{sc}$ and $\omega_\text{sc} + \mathrm{d}\omega_\text{sc}$. However, $S(\omega_\text{sc})$ is not normalized. The total probability is
$\int_0^\infty S(\omega_\text{sc}) \, \mathrm{d}\omega_\text{sc} = \rho_{ee}(t\to \infty)$
This results in
$P_\text{sc} = \int \frac{\hbar \omega_0 \Gamma}{r^2} f(\theta,\phi) \,\mathrm{d} \Omega \int S(\omega_\text{sc}) \, \mathrm{d}\omega_\text{sc} = \Gamma \hbar \omega_0 \rho_{ee}$
It’s worth mentioning that the spectrum of the scattered light is named by Mollow Triplet Spectrum
$S(\omega_\text{sc})=\frac{S_0}{(1+S_0)^2} \delta\left(\omega_\text{sc}\omega\right)+\frac{S_0}{8 \pi(1+S_0)} \frac{\Gamma}{\left[(\omega_\text{sc}\omega)^2+(\Gamma / 2)^2\right]} + \cdots$
where $S_0=2\Omega^2/\Gamma^2$ is the saturation parameter.
Conclusion: In spite of the fact that scattering spectrum are the superposition of waves at different frequencies $I_\text{sc}(\boldsymbol{r},\omega_\text{sc})$, the total power $P_\text{sc}$ of all the scattered beams is exactly equivalent to the power of spontaneous emission, regardless of whether the incident beam is onresonant $\omega = \omega_0$ or offresonant $\omega \neq \omega_0$. However, The spectrum of scattered light is centered at the frequency of incident beam $\omega$ !
Power broadening
Therefore, by measuring the power of resonance fluorescence, which is actually measuring $\rho_{ee}$ since $P_\text{sc} \propto \rho_{ee}$ , we can draw a graph
To explain the line shape, we recall the expression of $\rho_{ee}$ :
$\rho_{ee} (t\to \infty) = \frac{S_0/2}{1+S_0+\dfrac{4\Delta^2}{\Gamma^2}}$
Define a parameter $\Gamma^\prime = \Gamma \sqrt{1+S_0}$ where $S_0 = 2\Omega^2/\Gamma^2 = I/I_\text{sat}$ evaluates the intensity of the incident beam. The $\rho_{ee}$ is written
$\rho_{ee} = \frac{S_0}{2(1+S_0)} \frac{1}{1+\dfrac{4\Delta^2}{\Gamma^{\prime 2}}}$
The full width at half maximum (FWHM) is exactly $\Gamma^\prime$ (When $\Delta = \pm \Gamma^\prime/2$ we have $\frac{1}{1+\frac{4\Delta^2}{\Gamma^{\prime 2}}} = \frac{1}{2}$). Two mechanism of broadening could be found
 In the weakfield limit, $S_0 \simeq 0$, $\Gamma^\prime \simeq \Gamma$. This is called natural broadening.
 In the strongfield limit, $\Gamma^\prime = \Gamma \sqrt{1+S_0}$ . The FWHM of the transition is effectively larger due to the strong coupling to the field. This phenomenon is called power broadening of the transition.
Connection with classical model
Warning: Only selfconsistent when onresonance $\omega=\omega_0$.
Optical Bloch Equation writes
$\dot{\rho}_{ee} = \Gamma \rho_{ee}  \frac{\Omega^2/\Gamma}{1 + \dfrac{4\Delta^2}{\Gamma^2}} (\rho_{ee}\rho_{gg})$
Naturally $\rho_{ee} = N_2/N,\,\rho_{gg} = N_1 / N$. Recall that in classical model
$\frac{\mathrm{d}N_2}{\mathrm{d}t} = N_2 \Gamma \frac{\sigma(\omega) I}{\hbar \omega} (N_2  N_1)$
thus
$\sigma(\omega) = \frac{\hbar \omega}{I} \frac{\Omega^2/\Gamma}{1 + \dfrac{4\Delta^2}{\Gamma^2}}$
Since
$\Omega^2 = \frac{\langle g \hat{d}e\rangle^2}{\hbar^2} E^2 \qquad \Gamma=\frac{\omega^3_0}{\pi \varepsilon_0 \hbar c^3} \langle g \hat{d}e\rangle^2 \qquad I = \frac{1}{2}c \varepsilon_0 E^2$
The absorption crosssection now becomes
$\sigma(\omega) = \frac{2 \pi c^2 \omega}{\omega_0^3} \frac{1}{1 + \dfrac{4\Delta^2}{\Gamma^2}} \approx \frac{\lambda^2_0}{4} \frac{\Gamma^2}{2\pi (\Gamma^2/4 +\Delta^2)}$
where $\lambda_0 = 2\pi c/\omega_0$. If we define a line shape function
$g(\omega) = \dfrac{1}{2\pi} \dfrac{\Gamma}{\Gamma^2/4 +\Delta^2}$
we find $\sigma(\omega) = \Gamma \frac{\lambda^2_0}{4} g(\omega)$. We can obtain the same result taking advantage of energy conservation:
$\sigma(\omega) I (N_2  N_1) = P_\text{sc} = N\Gamma \hbar \omega_0 \rho_{ee} \implies \sigma(\omega) = \Gamma \frac{\lambda^2_0}{4} g(\omega)$
which means, given the same incident intensity and the same population, the efficiency of energy absorption by atoms varies with the changing incident light frequency $\omega$, and reaches its maximum at resonance $\Delta = 0$.