\section{Introduction}

Supernovae have been hypothesized (Ginzburg \& Syravotskii \cite{gs64}, Hayakawa \cite{hay69})
to be the sources of the cosmic rays, which consist primarily of relativistic protons, ions and electrons.
Firm evidence for this hypothesis is lacking because cosmic rays
are deflected by intervening magnetic fields and do not trace back to their sources.
A distinctive $\gamma$-ray signature of cosmic rays results from
interactions of cosmic-ray protons and ions with gas and dust in the vicinity of their acceleration sites,
making $\gamma$ rays
with a characteristic $\pi^0\rightarrow 2\gamma$ emission feature peaking at 67.5 MeV
in a photon number representation (Stecker \cite{ste71}).

Resolution and sensitivity limitations of previous $\gamma$-ray telescopes made it impossible
to precisely localize $\gamma$-ray sources in the Galaxy and to measure the expected $\pi^0$-decay spectral
signature. This has changed with the new generation of $\gamma$-ray telescopes, including
AGILE (Astro-rivelatore Gamma a Immagini LEggero; Gamma-ray Light Imaging
Detector) and the Fermi Large Area Telescope (LAT). Evidence for a cutoff below
several hundreds of MeV from some SNRs has been found (Giuliani et al.\ \cite{giu11},
Abdo et al.\  \cite{AbdoW44}, \cite{AbdoIC443}, Ackermann et al. \cite{2013Sci...339..807A}), so
this one-hundred year old problem may soon be solved.

While spectral analyses at photon energies $E_\gamma \lesssim 100$ MeV continue,
we can also ask whether the measured $\gamma$-ray luminosities of SNRs  are consistent
with expectations of efficiency for cosmic-ray acceleration at SNR shocks. Furthermore, trends in
the $\gamma$-ray spectra might reveal a progressive
behavior characteristic of cosmic-ray acceleration in these sources.

Here we test the hypothesis that cosmic rays are accelerated at SNRs
by calculating the fraction of supernova kinetic energy that
must be converted into cosmic-ray particles in order to produce
the $\gamma$ rays that are observed
from some of these remnants. We also plot the GeV to TeV energy-flux ratio with age.

\section{Particle acceleration at SNR shocks}

Supernovae (SNe) are found in two main classes, namely thermonuclear Type Ia SNe and core-collapse Type II SNe.\footnote{Other subclasses
of core-collapse SNe subclasses are Type Ib and Ic SNe.}
The progenitors of Type Ia SNe are Solar mass white dwarfs composed of carbon and oxygen.
Thermonuclear burning of the material in the white dwarf releases a prodigious amount of energy and creates a shock wave
as the burned material expands.
For stars more than $\approx 8$ times as massive as the Sun, in contrast, the stellar interior burns until much of it
has been converted to an Fe core, which will collapse when electron-degeneracy pressure is no longer sufficient
to support the core against gravity.  A core bounce sends a shock wave out that
expels a portion of the stellar envelope (for a general review, see Burrows \cite{bur00}).
Neutron stars can be formed in Type II SNe, but not in Type Ia SNe.

In both Types Ia and II SNe, the energy available for nonthermal particle production is derived from the
directed bulk kinetic energy of the shell. This energy is dissipated in the form of the random kinetic energy of matter swept up
at the SNR shock. The amount of energy that is extracted and transformed into cosmic-ray energy is highly
uncertain, and represents a major open question in SNR modeling. Different prescriptions related
to the number of swept-up particles or the kinetic energy dissipated at the shock front have been considered
(e.g., Sturner et al.\ \cite{stu97}, Reynolds \cite{rey98}, Baring et al.\ \cite{bar99}, Tang et al.\ \cite{tan11}),
but the exact conversion efficiency depends on poorly understood microphysical processes (e.g., Blasi et al.\ \cite{bla05}).
Here we reconsider this problem using an idealized model of an expanding spherical shell of matter
and a uniform surrounding medium (cf.\ Truelove \& McKee \cite{tm99};
Finke \& Dermer \cite{fd12}).

Let  $E_0 = 10^{51}E_{51}$ erg represent the total kinetic energy of the supernova ejecta
with mass $M_0 = m M_\odot$, where $M_\odot = 2.0\times 10^{33}$ g is the mass of the Sun.
For Type Ia and II SNe, $m \cong 1.4$ and $m \cong 10$, respectively,
with $E_{51}\approx 1$ for both types. The SN is here
assumed to explode in a medium of uniform density $n_0$ and form an outward moving
spherical shell of material with radius $r = r(t)$ at
time $t$ after the explosion. As it expands, the shell sweeps
up material with mass $M_{su}(t) = 4\pi r^3 \rho_0/3$, where the density $\rho_0 = m_p n_0$,
and the ISM is assumed to be composed of hydrogen.
The initial speed of the shell is therefore $v_o = \sqrt{2E_0/M_0} \cong 10^4 \sqrt{E_{51}/m}$ km s$^{-1}$.

The Sedov radius $r_{\rm Sed}$ is defined by the condition that the swept-up mass equals the explosion mass, that is,
$4\pi \rho_0 r_{\rm Sed}^3/3 = M_0$, implying $r_{\rm Sed} = (3M_0/4\pi\rho_0 )^{1/3} = 6.6\times 10^{18}(m/n_0)^{1/3}$ cm.
The Sedov time $t_{Sed}$ is the time when the shell transitions from the coasting phase at $r < r_{\rm Sed}$ to
the Sedov phase at $r> r_{Sed}$, and is given by $t_{Sed} = r_{\rm Sed}/v_0 = 210\, m^{5/6}E_{51}^{-1/2} n_0^{-1/3}$ yr.

For nonrelativistic speeds $v_o\ll c$ and adiabatic shock waves, the kinetic energy
\begin{equation}
E_{ke} = {1\over 2} [M_0 + M_{su}(t)]v^2
\label{Eke}
\end{equation}
is conserved, so that
$d E_{ke}(t)/dt = 0$.
Introducing dimensionless radius $x = r/r_{\rm Sed}$ and time $\tau = t/t_{\rm Sed}$ in units of Sedov radius and Sedov time,
respectively, and noting that  $v = dr/dt$, Eq.\ (\ref{Eke}) leads to the equation of motion
\begin{equation}
(1+x^3)^{1/2} dx = d\tau
\label{eqofmotion}
\end{equation}
for the dimensionless shell radius $x$. From this, the asymptotes
\begin{equation}
x(\tau ) \rightarrow \left\{ \begin{array}{ll}
    \tau\;, & \tau \ll 1 \\
    \left(5\tau /2\right)^{2/5}\;, & \tau \gg 1
  \end{array}
  \right. \
\label{xtau}
\end{equation}
and
\begin{equation}
v(\tau ) \rightarrow \left\{ \begin{array}{ll}
    v_0\;, & \tau \ll 1\;, \\
    v_0\left(5\tau /2\right)^{-3/5}\;, & \tau \gg 1
  \end{array}
  \right.  \
\label{vtau}
\end{equation}
can be derived.

   \begin{figure}
   \centering
   \includegraphics[width=3.0in]{f1.eps}
      \caption{{\it Top:} Red curve labeled ``1" shows the rate at which kinetic energy is swept into the
SNR shell, and the blue curve labeled  ``2" shows the rate at which kinetic energy
is lost due to deceleration of the shell.  The net rate of change of kinetic energy in the SNR shell is
given by the green curve, labeled ``3". {\it Bottom:} Integrated kinetic energy corresponding to the three terms shown
in the upper panel. Curve 3 gives the accumulated swept-up kinetic energy, ${E}_{su}(t)$.
              }
         \label{powerenergy}
   \end{figure}

If nonthermal particles extract their energy from the kinetic energy of swept-up matter, then the change of kinetic energy of
swept-up matter with time is given by
$${d E_{ke}^{su}(t)\over dt} = {d\over dt}\left[ {1\over 2} M_{su}v^2\right ] ={1\over 2}\dot M_{su}v^2 + M_{su}v \dot v\;$$
\begin{equation}
=2\pi r^2 v^3 \rho_0 + {4\pi \over 3}r^3 v\dot v \rho_0\;.
\label{dEkesut}
\end{equation}
The first term on the right-hand-side of Eq.\ (\ref{dEkesut}) is the rate at which kinetic energy is swept into the shocked material,
and the second term is the rate at which the kinetic energy of the swept-up matter is lost by deceleration of the shell (and is therefore
negative). The sum represents the net rate of change of kinetic energy in the SNR shell. These terms are shown in the top
panel of Fig.\ \ref{powerenergy}, and the integrated energy associated with each term is shown in the bottom panel of
Fig.\ \ref{powerenergy}. Here we let $E_{51} = 1$, $m = 1.4$, and $n_0 = 0.1$ cm$^{-3}$. Integrating the swept-up kinetic energy
leads to a logarithmically diverging total energy that exceeds the initial
kinetic energy of the explosion. By adding the change of energy associated with deceleration, the summed contribution recovers
the initial kinetic energy of the explosion at late times.

The process of sweeping up matter from the surrounding medium has
transformed the directed kinetic energy into internal thermal and nonthermal kinetic energy of swept-up matter. { First-order Fermi acceleration at astrophysical shocks provides a specific mechanism to convert directed kinetic energy into nonthermal particle kinetic energy (e.g., Blandford \& Eichler  \cite{be87}), although here we treat the dissipation mechanism through a hypothetical conversion efficiency. }  Note
how rapidly, $\propto \tau^{-2.5}$, that the total injection rate of kinetic energy decreases in the Sedov phase.

\section{Cosmic-ray interactions in the SNR environment}

Cosmic-ray ions and electrons can interact with target particles to make $\gamma$ rays through secondary nuclear production and
bremsstrahlung, respectively. Additionally, relativistic electrons will Compton-scatter photons of the
ambient radiation fields to the $\gamma$-ray regime (see Blumenthal \& Gould \cite{bg70} and Dermer \& Menon \cite{dm09} for reviews of high-energy
astrophysical radiation processes). There are three targets available for particle-particle interactions, namely from (1) the explosion mass,
(2), the swept-up mass, and (3) the ISM gas. The explosion mass density $n_{ex}(t) = M_0/m_pV_{sh}(t)$,  where
the volume of the expanding shell of width $\Delta$ is $V_{sh}(t) = 4\pi[r^3-(r-\Delta )^3]/3 \rightarrow
4\pi r^2 \Delta$ in the limit $u \equiv \Delta/r \ll 1$. Hence
\begin{equation}
n_{ex}(t) = {M_0\over 4\pi r^3 u m_p} = {n_0\over 3x^3 u}  \rightarrow {n_0\over 3u}\left\{ \begin{array}{ll}
    \tau^{-3}\;, & \tau \ll 1 \\
    (5\tau/2)^{-6/5}\;, & \tau \gg 1 \;.
  \end{array}
  \right. \
\label{xtau}
\end{equation}
The density of the swept-up mass is increased by the compression ratio. For a strong shock,
the shell of swept-up matter has { density $n_{su} \cong 4n_0$.}

Cosmic rays accelerated at the forward shock and convecting downstream into the shocked fluid would interact with the compressed ISM with density $\approx 4n_0$.
The most efficient radiation production
 occurs when the cosmic rays interact with explosion mass in the coasting phase, and swept-up mass
in the Sedov phase. These two regimes can be bridged with a formula for the maximum target density
available for interactions, given by
\begin{equation}
n_{max}(t) \cong {n_0\over 3ux^3}\,(1+12ux^3)\;.
\label{nmaxt}
\end{equation}
For cosmic rays to interact with the dense SNR shell material would require rapid diffusion
from the region downstream of the forward shock into the expanding SNR shell.   We consider the case where
the target is the shocked ISM and where the CRs interact with the densest available target.

Secondary nuclear production makes $\gamma$ rays from the decay of neutral pions, as well as through
emissions from pion-decay electrons and positrons. For cosmic-ray protons with energies $E_p = m_p c^2 \gamma_p$ well above the pion-production threshold energy at $E_p \cong 300$ MeV, the rate at which the proton Lorentz factor $\gamma_p$ changes is given by $-\dot \gamma_{pp} = K_p c\sigma_{pp} n(t)\gamma_p$, where $K_p\approx 0.5$ is the mean inelasticity (fraction of energy lost per collision), and $\sigma_{pp} \approx 30$ mb is the p-p inelastic cross section, so that
the energy-loss time scale associated with secondary production is
$t_{pp}  = |\gamma_p/\dot\gamma_{pp}| = [K_p c\sigma_{pp} n(t)]^{-1}\cong 2.2\times 10^{15}/n(t)$ s. The maximum $\gamma$-ray luminosity
that can be made by secondary nuclear production is therefore
\begin{equation}
L_{pp}(t) \cong { \eta_{pp} \,{E}_{su}(t)\over 3 t_{pp}(t)},
\label{Lppt}
\end{equation}
with $n(t)$ replaced by $n_{max}(t)$, from Eq.\ (\ref{nmaxt}), in $t_{pp}(t)$. The factor $1/3$ represents the mean fraction of
energy lost in an inelastic nuclear interaction that emerges as $\gamma$ rays, and the term  $\eta_i$ is the fraction of accumulated
swept-up kinetic energy  that is transformed into cosmic-ray protons { or electrons,  where $i = pp, ff,$ and {\it C} stand for nuclear production,
electron bremsstrahlung, and electron Compton scattering, respectively}.

The bremsstrahlung energy-loss rate for electrons with energy $E_e = m_ec^2 \gamma_e$ interacting with ions of density $n_Z$ and charge
$Z$ in a fully ionized medium is given by $-\dot \gamma_{ff} = k_{ff}\alpha_f c\sigma_{\rm T}[\Sigma n_Z Z(Z+1)]\gamma_e$,
where $\alpha_f$ and $\sigma_{\rm T}$ are the fine-structure constant and Thomson cross section, respectively.
The term $k_{ff} = (3/2\pi)(\ln 2\gamma_e - {1\over 3})  \approx 4.3 $, taking $\gamma_e \approx 6000$ for the typical
electron Lorentz factor $\gamma_e$ that makes GeV $\gamma$ rays through this process. The characteristic bremsstrahlung
energy-loss time scale is therefore $t_{ff}(t) = |\gamma_e/\dot\gamma_{ff}|\cong 8.0\times 10^{14}/n(t)$ s, and the maximum bremsstrahlung luminosity is
\begin{equation}
L_{ff}(t) \cong { \eta_{ff} \,{E}_{su}(t)\over t_{ff}(t)},
\label{Lfft}
\end{equation}
again replacing $n(t)$ with $n_{max}(t)$ in $t_{ff}(t)$.

The loss rate of cosmic-ray electrons by Compton scattering is given by
$-\dot\gamma_{C} = (4/3)c\sigma_{\rm T} U_\gamma \gamma_e^2/m_ec^2$,
where $U_\gamma$ is the energy density of the { target radiation} field.
The characteristic energy-loss time scale through Compton { scattering in the Thomson limit} is therefore
$t_{C}$(s)$ = |\gamma_e/\dot\gamma_{C}| = 7.7\times 10^{19}/\gamma_e $ $ \approx 7.1\times 10^{13}/\sqrt{E_{\rm GeV}}$,
where we assume that the dominant radiation field is the cosmic microwave background (CMB) radiation
with energy density $U_\gamma = U_{CMB} \cong 4.0\times 10^{-13}$ erg cm$^{-3}$.
In the final expression, $\gamma_e$ is replaced by the value of the electron Lorentz factor that will scatter a typical
CMB photon to a $\gamma$ ray with energy $E_{GeV}$ GeV. The $\gamma$-ray luminosity from Compton scattering is therefore
\begin{equation}
L_{C}(t) \cong { \eta_{C} \,{E}_{su}(t)\over  t_{C}(t)}.
\label{LCst}
\end{equation}
{ The nonthermal electrons will additionally lose energy through synchrotron losses, with total synchrotron luminosity
larger than $L_C$ by a factor $U_B/U_\gamma$, where $U_B$ is the magnetic-field energy density. The magnetic field is constrained
in SNR spectral modeling, but does not affect the $\gamma$-ray emission model presented here. }

  \begin{figure}
   \centering
   \includegraphics[width=3.0in]{f2.eps}
      \caption{{\it Top:} Energy-loss time scales of cosmic-ray protons through secondary nuclear production ($t_{pp}$; solid curves) and of cosmic-ray electrons through bremsstrahlung ($t_{ff}$; dashed curves) are shown for parameters corresponding to Type Ia and Type II SNe, as labeled (values for Type II SNe are given in parentheses in the legend). Also shown are time scales for Compton energy losses ($t_{\rm C}$; dotted curves) of cosmic-ray electrons that scatter CMB photons to make $\gamma$ rays of energy $E_\gamma$. {\it Bottom:} Theoretical estimates of $\gamma$-ray luminosities of SNRs through secondary nuclear production (solid curves), bremsstrahlung (dashed curves), and Compton scattering (dotted curves). Data points show the measured $> 300$ MeV $\gamma$-ray luminosities for the SNRs listed in Table 1. Thick curves in the two panels are for cosmic-ray interactions with the shocked external medium, and thin curves are for cosmic-ray interactions with the densest available target material.  Model curves are plotted only in the weak-cooling regime.
              }
         \label{time scales}
   \end{figure}

The curves shown in Fig.\ \ref{time scales} give results for the energy-loss time scales in the top panel, and the $\gamma$-ray luminosities of SNRs in the bottom panel.\footnote{We calculate a $\gamma$-ray luminosity for CTB 37A larger, by an order of magnitude, than Castro and Slane (\cite{cs10}).} For the model results, we assume that Type Ia and II SNe involve 1.4 and 10 $M_\odot$ of explosion mass, and that the surrounding medium densities are $0.1$ and 10 cm$^{-3}$, respectively. The total kinetic energy of the explosion is assumed to be $10^{51}$ erg in both cases, and we let $u = 0.1$ for the relative width of the shell to the radius.  The Sedov time scales for these parameters are $595$ yr and $660$ yr for the Type Ia and Type II SNRs, respectively. Also shown are the energy-loss time scales for cosmic-ray electrons that scatter CMB photons to $\gamma$ rays of energy $E_\gamma$.
The solid line in the upper panel separates the weak and strong cooling regimes; in the strong-cooling regime our treatment does not apply.

The lower panel of Fig.\ \ref{time scales} gives model results for the $\gamma$-ray luminosities using the same parameters as shown in the upper panel. The results are plotted only in the slow-cooling regimes where the model applies. The swept-up kinetic energy is assumed to be transformed into cosmic-ray protons and electrons with efficiencies of $\eta_{pp} = 10\%, \eta_{ff} = 5$\%, and $\eta_{C} = 0.1$\% for the three processes. The efficiency for converting kinetic energy into electrons is small for the latter process, but these electrons have to have $\gamma_e \cong 10^6 \sqrt{E_{GeV}}$, { large by comparison with bremsstrahlung $\gamma$ rays made by electrons with $\gamma_e \sim 3 E_{GeV}/m_ec^2 \sim 6\times 10^3 E_{GeV}$.

The injection efficiency should not, however, depend on the $\gamma$-ray emission process.  { In the case of a very soft electron spectrum, the characteristic amount of energy injected into lower-energy, $\gamma_e \sim 10^4$ electrons making bremsstrahlung would be much greater than the energy injected into the higher-energy, $\gamma_e\sim 10^6$ electrons that Compton-scatter soft photons to GeV energies, and this could account for the different efficiencies that are implied by the data for the two processes.} This would, however, conflict with the hard electron spectra required to model SNRs. At 10 GeV, where Solar modulation effects and energy losses on cosmic-ray protons and electrons are small, the ratio of the cosmic-ray electron to proton fluxes is $\approx 1\%$. If SNRs are the sources of the cosmic-ray protons and electrons, then the injection efficiency for { electrons should therefore be much less than for protons, and a bremsstrahlung origin of the $\gamma$-ray emission could be ruled out   (Yuan et al. \cite{ylb12})}.

Note that the increase in swept-up kinetic energy competes with the decrease in target density so that the $\gamma$-ray luminosity is roughly constant with time for secondary nuclear production and bremsstrahlung when using the target density given by eq.\ (\ref{nmaxt}). The $\gamma$-ray luminosity through nuclear production { and bremsstrahlung} declines at early times when cosmic-ray protons interact only with the swept-up ISM with density $\approx 4n_0$. The target photon density for Compton scattering is constant, so the $\gamma$-ray luminosity through this process tracks the increasing swept-up kinetic energy shown in Fig.\ \ref{powerenergy}, and is also therefore smaller at early times.

   \begin{table}
      \caption[]{$\gamma$-ray SNRs.$^{\mathrm{a}}$}
         \label{GeVTeVSNRs}
     $$
         \begin{array}{lllll}
            \hline
            \noalign{\smallskip}
            {\rm SNR}  & {\rm Type} &  {\rm Age~(kyr)} & {\rm d~(kpc)}  & {\rm Refs.}^{\mathrm{b}}  \\
            \noalign{\smallskip}
            \hline
            \noalign{\smallskip}
            {\rm Cas~A} & {\rm II}	& 0.33    &   3.4^{+0.3}_{-0.1}  &   1 \\
           {\rm  Tycho}  &  {\rm Ia}      & 0.44  & 1.5 - 4.0   & 2 \\
	    {\rm RX~J1713.7-3946}& {\rm II(?)}	& 1.6	&  0.9 - 1.7   & 3,4	\\
	    {\rm  RX~J0852.0-4622}  &  {\rm II(?)} & 1.7 - 4.3   &   0.75   & 5,6   \\
	    {\rm CTB~37A}& {\rm II}	& 2\pm 0.5	&  6 - 9.5  &  7,8	\\
	    {\rm W49B}	& {\rm II}	& 1 - 4	&  8 - 11 & 9 \\
           {\rm  IC443}  &  {\rm Ib(?)} & 3 - 30^{\mathrm{c}}   &   \sim 1.5   & 10   \\
	   {\rm  W51C}	& {\rm II}	&  10 - 30&  7\pm1.5 & 11	\\
           {\rm G359.1-0.5}  &  {\rm } 		& > 10   &  7.6   & 12,13  \\
	    {\rm W44}	& {\rm II}	&  15 - 20	& \sim 3    &  14,15 \\
           {\rm  G8.7-0.1}  &  {\rm II}{^\mathrm{e}} 		& 15 - 28  &   4.8 - 6.0   &16,12   \\
	   {\rm  W28}	& {\rm II}		&  35 - 45 &  1.9\pm0.3	& 17,18\\
           {\rm  W41}  &  {\rm II}^{\mathrm{d}} 		& 60 - 200   &   3.9 - 4.5  & 19,20   \\

            \noalign{\smallskip}
            \hline
         \end{array}
     $$
\begin{list}{}{}
\item[$^{\mathrm{a}}$] Except for W44, all SNRs are detected at TeV energies.
\item[$^{\mathrm{b}}$] 1. Abdo et al.\ \cite{AbdoCasA};
2. Giordano et al.\ \cite{gio12};
3. Abdo et al.\ \cite{AbdoRX};
4. Fesen et al.\ \cite{fes12};
5. Tanaka et al.\ \cite{tan11};
6. Aharonian et al.\ \cite{aha07};
7. Brandt  \cite{bra13};
8. Tian et al.\ \cite{tia12};
9. Abdo et al.\ \cite{AbdoW49B};
10. Abdo et al.\ \cite{AbdoIC443}; see also Caprioli \cite{cap11};
11. Abdo et al.\ \cite{AbdoW51C};
12. Aharonian et al.\ \cite{aha06};
13. Hui et al.\ \cite{hui11};
14. Abdo et al.\ \cite{AbdoW44};
15. Wolszczan et al.\ \cite{wol91};
16. Ajello et al.\ \cite{aje12};
17. Abdo et al.\ \cite{AbdoW28};
18. Vel{\'a}zquez et al.\ \cite{vel02};
19. Leahy \& Tian \cite{lt08};
20. Tian et al.\ \cite{tia07}

\item[$^{\mathrm{c}}$] Troja et al.\ \cite{tro08} argue from XMM-Newton X-ray observations that the age is $< 10$ kyr.
\item[$^{\mathrm{d}}$] { This shell-type SNR is probably a core-collapse SNR given its proximity to molecular clouds. The TeV radiation might also be powered by a pulsar wind, in which case this analysis would not apply; see Mukherjee et al.\ \cite{mgh09}.
\item[$^{\mathrm{e}}$] { Possibly associated with PSR 1800-21; see Kassim \& Weiler \cite{km90}. TeV emission could also be powered by the pulsar wind.}}
\end{list}
   \end{table}

\section{Comparison with data}

 We calculated the $\gtrsim 300$ MeV $\gamma$-ray luminosities $L_\gamma$ from the data provided in the references for the SNRs listed in Table 1. The error bars in $L_\gamma$ reflect both the uncertainties in the $\gamma$-ray measurements as well as distance uncertainties. The values of $L_\gamma$ range from $\approx 10^{34}$ erg s$^{-1}$ to $\approx 3\times 10^{36}$ erg s$^{-1}$. The data show an apparent increase in luminosity from young to { intermediate-aged ($\sim 3$ ky)} SNRs, followed by a decline for the { middle-aged} ($\gtrsim 10$ kyr) SNRs, but a larger sample of $L_\gamma$ from SNRs will be needed to confirm this trend. Most of these SNRs are Type II, so that if the parameters used to characterize the Type II SNRs are reasonable, then we find that $\eta\sim 10$\% for nuclear processes { to account for the $\gamma$-ray luminosities of the most luminous SNRs}, generally consistent with expectations of production efficiency if SNRs accelerate the cosmic rays. {  There are, however, sources that deviate significantly from the $\eta= 10$\% curves, for example, W28. Precise determinations of production efficiency requires an examination of SNRs on a case by case basis, because the surrounding medium density and composition, ejected mass and kinetic energy of the SNR, not to mention the limitations of the simple model, will introduce a wide range of luminosities even in the case that $\eta$ is constant.}

Four young, $\lesssim 3000$ yr old SNRs have $\gamma$-ray luminosities $\lesssim 10^{35}$ erg s$^{-1}$. Tycho is a Type IA in a region of low density. If cosmic rays interact only with the shocked material, though not with the SNR shell, then the luminosity of Cas A is in accord with this simple model (see Fig.\ 2). { SNRs RX~J1713.7-3946 and RX~J0852.0-4622}, though they may be Type II SNRs, also take place in low density regions, which will diminish the relative contribution from nuclear processes.

In our simple model, SNRs from Type II SNe formed in dense environments are predicted to be $\approx 2$ orders of magnitude more $\gamma$-ray luminous than Type Ia SNe, and therefore more easily detected. { Given the greater rate of Type II than Type Ia SNe in the Galaxy (e.g., Cappellaro et al. \cite{cap99}), this would explain the larger number of core-collapse SNRs in our sample. The assumed larger ISM densities in the vicinity of core-collapse SNe are expected because Type II and Type Ib/c SNe are typically found in star-forming regions. Indeed, OH maser emissions from CTB 37A suggests that it is interacting with a very dense, $n_0 > 10^3$ cm$^{-3}$, molecular cloud (Hewitt et al.\ \cite{hyw08}; Castro \& Slane \cite{cs10}). The ratio of nuclear production and Sedov time scales is $\sim 8\times 10^4 E_{51}^{1/2}/m^{5/6} n_0^{2/3}\propto n_0^{-2/3}$ for cosmic rays interacting with the swept-up ISM, so Sedov-age SNRs in dense environments would shine more brightly through nuclear interactions than Compton processes  (cf.\ Yuan et al. \cite{ylb12}) than those in low-density surroundings. The SNRs RX~J1713.7-3946 and RX~J0852.0-4622 are, however, found in tenuous, $n_0 <0.1$ cm$^{-3}$ environments  (Slane et al.\ \cite{sla99}, \cite{sla01}), which could happen when core-collapse SNe take place in stellar wind cavities. Their hard spectra favors a Compton-origin for the gamma rays.  }

  \begin{figure}
   \centering
   \includegraphics[width=3.0in]{f3.eps}
      \caption{Ratio of 0.3 -- 3 GeV to 0.3 -- 3 TeV fluxes for the SNRs given in Table 1.
              }
         \label{GeVTeVfluxes}
   \end{figure}

We also plotted the ratio of the 0.3 -- 3 GeV energy flux to the 300 GeV -- 3 TeV energy flux for the SNRs listed in Table 1. The result is shown in Fig.\  \ref{GeVTeVfluxes}. All of the SNRs except for RX~J1713.7-3946 and RX~J0852.0-4622 have GeV to TeV flux ratios of order unity or greater. The middle-aged SNRs are more luminous in GeV $\gamma$ rays than in TeV $\gamma$ rays by factors of $\gtrsim 10$ -- 100. RX~J1713.7-3946 is anomalously low in comparison. There is a weak trend for increasing GeV to TeV flux ratio with age. This might be expected if the dominant processes are Compton scattering and nuclear production, noting from Fig.\ \ref{time scales} that cosmic-ray electrons making TeV photons lose their energies on time scales of $\sim 10^{4}$ yr, leaving the long-lived cosmic-ray protons that would preferentially make GeV radiation.

Compton scattering as a multi-TeV $\gamma$-ray emission process for middle-aged SNRs might be ruled out by noting how rapidly the high-energy electrons that scatter CMB photons to TeV energies cool (Fig.\ \ref{time scales}). If these electrons are powered by swept-up kinetic energy, energetic electrons in the strong cooling regime would not be replenished due to the rapid decrease in the kinetic-energy injection rate shown by curve 3 in Fig.\ \ref{powerenergy}.   Compton production of 10 -- 100 TeV $\gamma$-rays requires electrons with $\gamma_e \sim 10^7$ -- $10^8$, irrespective of target photon energy. Thus the Compton scattering process can be ruled out by detection of $\gg$ TeV radiation from SNRs much older than the Sedov age, provided that it can be shown that there is no pulsar wind to inject high-energy electrons. Indeed, a rapid decline in the X-ray synchrotron radiation from $\gamma\sim 10^8/\sqrt{(B/10~\mu{\rm G})}$ electrons shortly after the Sedov phase is found by plotting the nonthermal X-ray luminosity as a function of radius (Nakamura et al.\ \cite{nak12}), which is closely related to the remnant age.  Reduction in the GeV  -- TeV $\gamma$-ray luminosity in middle-aged SNRs, if made primarily by protons and ions, seems to { follow} the behavior of the { X-ray synchrotron luminosity radiated} by high-energy electrons at earlier times. { The evidence for the trend is weak, however, given that IC 443, an intermediate-aged SNR, is less luminous than most of the older remnants, in particular, W44 and W41.} The decline in the $\gamma$-ray luminosity { with age in intermediate and middle-aged SNRs, if confirmed, would instead have to be due to escape rather than energy loss as in the case of X-ray emitting electrons}, given the long time scale for energy loss of cosmic-ray protons  through inelastic nuclear processes.   Finally, we note that the assumption of a uniform medium surrounding the SNR is an extreme simplification of the actual environment of SNRs, particularly the core-collapse type found in star-forming regions.

\section{Conclusions}

Data from the new generation of $\gamma$-ray telescopes can be used to test the hypothesis that SNRs are the sources of the cosmic rays. We have constructed a simple model for $\gamma$-ray production in SNRs and find that, for parameters that are expected to apply to Type Ia and II SNe, $\approx 10$\% conversion efficiencies of directed kinetic energy to nonthermal protons or electrons are required to account for the measured $\gamma$-ray luminosities of SNRs if the radiation mechanism is secondary nuclear production or bremsstrahlung, respectively. The { small} $\gamma$-ray luminosities in young SNRs { compared to intermediate-aged SNRs} can be explained if the target material for cosmic-ray interactions is the shocked interstellar medium. Smaller efficiencies can account for the $\gamma$ rays if the emission process is Compton scattering, but the energy has to be converted to highly relativistic electrons. If the relative acceleration efficiency for electrons and protons is given by the electron-proton ratio in the cosmic rays, of order $\sim 1$\%, then the electron bremsstrahlung $\gamma$-ray luminosity would be too weak to explain the $\gamma$-ray luminosity of SNRs.

We also looked for trends in $L_\gamma$ and the GeV to TeV flux ratios with SNR age. Except for RX J1713.7-3946 and RX~J0852.0-4622,} SNRs are found to be most $\gamma$-ray luminous at ages of $\approx 10^3$ -- $10^4$ yr. The two anomalous SNRs, RX J1713.7-3946 and RX~J0852.0-4622, have the hardest GeV-TeV spectra (Fig.\ 3) and are found in low-density environments, suggesting a Compton origin for their $\gamma$-ray emission.  The { apparent} decline of $L_\gamma$ at later times might be due to escape of cosmic rays from the acceleration and target interaction sites, but {  evidence for this trend is weak}, and the possible presence of $\gamma$-ray emission due to pulsar winds must also be carefully considered.   The tendency of the GeV to TeV flux ratio of SNRs to increase with age, or age measured in units of the characteristic Sedov age of the remnant, is suggested by this work, but will require additional analyses and study to establish.

\begin{acknowledgements}
The work of CDD is supported by the Office of Naval Research and the Fermi Guest Investigator program. The work of GP was performed through the Science and Engineering Apprenticeship Program (SEAP) at NRL. We would like to thank T.~J.\ Brandt, D.\ Horan, J.~D.\ Finke, and R.\ Yamazaki for comments on this work, and the referee for an illuminating report, suggestions, and corrections.

\end{acknowledgements}


\section{Introduction}
\IEEEPARstart{T}{he} key component of a Multiple-Input-Multiple-Output (MIMO) communication system in terms of performance and complexity is the MIMO detector, which is used for separating independent data streams at the receiver. The Maximum Likelihood (ML) detectors achieve the optimal error rate performance. However, these types of detectors, including the near-optimal sphere decoder and its variants, are usually not suitable for practical systems due to their high complexity. Linear detectors, such as Zero-Forcing (ZF) and MMSE, achieve suboptimal performance, however, they are widely used in practical systems due to their low complexity implementations. Among linear receivers, MMSE is the optimal solution and seems to be the mainstream implementation choice due to its superior performance over ZF detectors.

Perfect channel state information (CSI) is usually assumed in the literature when simulating or analyzing the performance of linear detectors \cite{pingli,kim}. However, in practice the channel estimates are inherently noisy. Important work \cite{wang}, \cite{iranli} has characterized the error rate performance of ZF receivers in the presence of channel estimation error. Nevertheless, less is known for the case of MMSE detectors in practical scenarios. For ZF and MMSE receivers, the joint effect of phase noise and channel estimation error is considered in \cite{garcia} and the performance is analyzed in terms of the degradation in signal-to-noise-plus-interference-ratio (SINR) without expressing the closed form performance indicators or error rate analysis. The SINR derivations for the MMSE case in \cite{garcia} are done only for low SNR region. In both \cite{wang} and \cite{garcia}, channel estimation error variance is assumed to be constant for all SNRs. This is not realistic approach for packet based or bursty communication systems as the channel estimation error is in fact a function of the SNR. In this letter, we analyze the MMSE receivers in the presence of channel estimation error, and derive a closed form post-processing SNR expression, which provides an accurate estimate of the error rate performance. The error rate performance is investigated for both the constant channel estimation error variance case and the case with a realistic channel estimation algorithm where the estimation error variance is clearly dependent on the channel SNR. We believe that it is a very useful tool for throughput prediction in link adaptation protocols and for error rate analysis in general. Accuracy of the analytical results is verified through simulations.

\section{System Model and Derivations}
We consider a MIMO system where the transmitter is equipped with $N_t$ antennas, and the receiver uses $N_r$ antennas. The $N_r \times 1$ received signal vector $\mathbf{y}$ can be expressed as
\begin{equation}
\label{hxn}
	\mathbf{y}=\mathbf{H}\mathbf{x}+\mathbf{n}
\end{equation}
where $\mathbf{x}$ is the transmitted signal vector, $\mathbf{H}$ is the $N_r \times N_t$ channel matrix, and $\mathbf{n}$ is the $N_r \times 1$ additive Gaussian noise vector with zero mean and covariance matrix $E\left[\mathbf{n}\mathbf{n}^H\right]=N_0\mathbf{I}$. We assume an uncorrelated Rayleigh flat channel, i.e. entries of $\mathbf{H}$ are i.i.d. zero mean circularly symmetric complex Gaussians (ZMCSCG) with unit variance, and the signal energy at each transmit antenna is assumed to be equal to $E_s$.

The receiver can estimate the transmitted signal vector by applying the MMSE detector to the received signal, $\hat{\mathbf{x}}=\mathbf{W}\mathbf{y}=\mathbf{W}\mathbf{H}\mathbf{x}+\mathbf{W}\mathbf{n}$. Using the orthogonality principle \cite{sayed}, the MMSE detector $\mathbf{W}$ is derived as
\begin{equation}
\label{mmse}
	\mathbf{W}=\left[\mathbf{H}^{H}\mathbf{H}+\frac{N_0}{E_s}\mathbf{I}\right]^{-1}\mathbf{H}^{H}
\end{equation}

At the output of the MMSE detector, the residual signal plus interference from other spatial streams is well approximated as Gaussian \cite{peng} and the post-processing SNR (PPSNR) of $k^{th}$ spatial stream is calculated as\footnote{$\left(\ldots\right)_{k,l}$ denotes the $\left(k,l\right)^{th}$ entry of the matrix.}
\begin{equation}
\label{ppsnr}
	\gamma_k=\frac{E_s\left|\left(\mathbf{WH}\right)_{k,k}\right|^2}{E_s\sum\limits_{l\neq k}\left|\left(\mathbf{WH}\right)_{k,l}\right|^2 + N_0\left(\mathbf{W}\mathbf{W}^H\right)_{k,k}}
\end{equation}

The PPSNR is a good indicator for the error rate performance of MIMO systems, and therefore employed in link adaptation algorithms to predict the uncoded error rate \cite{peng}.
Since the output of the MMSE detector is Gaussian, the bit error rate of a specific modulation can be calculated by simply plugging the PPSNR value into the AWGN error rate formula of the modulation. The same technique is also used for theoretical derivation of error rate performance in fading channels.

This definition of PPSNR holds if the channel is perfectly known at the receiver. However, in practice, the channel matrix has to be estimated by the receiver, and the estimated channel is inherently noisy in practical systems. We model the estimated channel matrix as
\begin{equation}
\label{estch}
	\widehat{\mathbf{H}}=\mathbf{H}+\Delta\mathbf{H}
\end{equation}
where $\Delta\mathbf{H}$ denotes the estimation error matrix which is uncorrelated with $\mathbf{H}$, and its entries are ZMCSCG with variance $\sigma_e^2$. The quality of channel estimation is captured by $\sigma_e^2$, which can be appropriately estimated depending on the channel estimation method. We assume that each block (packet), that undergoes a specific channel realization, $\mathbf{H}$, observes a different realization of $\Delta\mathbf{H}$ at the receiver. This situation occurs in packet based communication systems like 802.11n where the channel is estimated on a per packet basis.

\subsection{PPSNR derivation for practical systems}
In this section, we derive the PPSNR for practical MIMO systems which observe channel estimation error.  The receiver uses the estimated channel $\widehat{\mathbf{H}}$ to calculate the MMSE detector as
\begin{equation}
\label{estmmse}
	\widehat{\mathbf{W}}=\left[\left(\mathbf{H}+\Delta\mathbf{H}\right)^{H}\left(\mathbf{H}+\Delta\mathbf{H}\right)+\frac{N_0}{E_s}\mathbf{I}\right]^{-1}\left(\mathbf{H}+\Delta\mathbf{H}\right)^H
\end{equation}
We write the imperfect MMSE solution as $\widehat{\mathbf{W}}=\mathbf{W}+\Delta\mathbf{W}$. Now, the MMSE estimate of the signal vector becomes
\begin{equation}
\label{estsignal}
	\tilde{\mathbf{x}}=\left(\mathbf{W}+\Delta\mathbf{W}\right)\mathbf{y}=\underbrace{\mathbf{W}\mathbf{H}\mathbf{x}}_{\mbox{\footnotesize{signal}}}+\underbrace{\Delta\mathbf{W}\mathbf{H}\mathbf{x}+\mathbf{W}\mathbf{n}+\Delta\mathbf{W}\mathbf{n}}_{\mbox{\footnotesize{post-detection noise}}}
\end{equation}
We observe that there are additional interference and noise terms caused by $\Delta\mathbf{W}$, and denote the post detection noise as $\hat{\mathbf{n}}=\Delta\mathbf{W}\mathbf{H}\mathbf{x}+\mathbf{W}\mathbf{n}+\Delta\mathbf{W}\mathbf{n}$. With this definition for the post detection noise, the PPSNR of the $k^{th}$ spatial stream in the presence of channel estimation error can be expressed as
\begin{equation}
\label{estppsnr}
	\widetilde{\gamma_k}=\frac{E_s\left|\left(\mathbf{WH}\right)_{k,k}\right|^2}{E_s\sum\limits_{l\neq k}\left|\left(\mathbf{WH}\right)_{k,l}\right|^2 + \left(E\left[\hat{\mathbf{n}}\hat{\mathbf{n}}^H\right]\right)_{k,k}}
\end{equation}
where we replaced the original noise covariance in \eqref{ppsnr} with the covariance of $\hat{\mathbf{n}}$, which is calculated as
\begin{IEEEeqnarray}{rCl}
\label{cov}
E\left[\hat{\mathbf{n}}\hat{\mathbf{n}}^H\right]&=&E\left[\Delta\mathbf{W}\mathbf{H}\mathbf{x}\mathbf{x}^H\mathbf{H}^H\Delta\mathbf{W}^H\right]+E\left[\mathbf{W}\mathbf{n}\mathbf{n}^H\mathbf{W}^H\right]\IEEEnonumber\\
	&&+\:E\left[\mathbf{W}\mathbf{n}\mathbf{n}^H\Delta\mathbf{W}^H\right]+E\left[\Delta\mathbf{W}\mathbf{n}\mathbf{n}^H\mathbf{W}^H\right]\IEEEnonumber\\
	&&+\:E\left[\Delta\mathbf{W}\mathbf{n}\mathbf{n}^H\Delta\mathbf{W}^H\right]
\end{IEEEeqnarray}
In order to calculate the terms in \eqref{cov}, we need to first derive $\Delta\mathbf{W}$. For small $\sigma_e^2$, the $\Delta\mathbf{H}^H\Delta\mathbf{H}$ term in \eqref{estmmse} becomes negligible compared to others. Hence, we can rewrite \eqref{estmmse} as
\begin{equation}
\label{estmmse1}
	\widehat{\mathbf{W}}\cong\left[\mathbf{H}^H\mathbf{H}+\frac{N_0}{E_s}\mathbf{I}+\mathbf{H}^H\Delta\mathbf{H}+\Delta\mathbf{H}^H\mathbf{H}\right]^{-1}\left(\mathbf{H}+\Delta\mathbf{H}\right)^H
\end{equation}
which can be further simplified using the matrix approximation $\left(\mathbf{P}+\epsilon^2\mathbf{Q}\right)^{-1}\cong\mathbf{P}^{-1}-\epsilon^2\mathbf{P}^{-1}\mathbf{Q}\mathbf{P}^{-1}$ for small $\epsilon^2$. Let us also define $\mathbf{K}=\left(\mathbf{H}^H\mathbf{H}+\frac{N_0}{E_s}\mathbf{I}\right)^{-1}$ for brevity and simplify \eqref{estmmse1} as
\begin{IEEEeqnarray}{rCl}
\label{estmmse2}
	\widehat{\mathbf{W}}&\cong&\left[\mathbf{K}-\mathbf{K}\left(\mathbf{H}^H\Delta\mathbf{H}+\Delta\mathbf{H}^H\mathbf{H}\right)\mathbf{K}\right]\left(\mathbf{H}+\Delta\mathbf{H}\right)^H\\
	&=&\mathbf{K}\mathbf{H}^H-\mathbf{K}\left(\mathbf{H}^H\Delta\mathbf{H}+\Delta\mathbf{H}^H\mathbf{H}\right)\mathbf{K}\mathbf{H}^H\IEEEnonumber\\	&&+\:\mathbf{K}\Delta\mathbf{H}^H-\underbrace{\mathbf{K}\left(\mathbf{H}^H\Delta\mathbf{H}+\Delta\mathbf{H}^H\mathbf{H}\right)\mathbf{K}\Delta\mathbf{H}^H}_{\mbox{\footnotesize{small compared to the other terms}}}
\end{IEEEeqnarray}
Finally the desired error matrix becomes
\begin{equation}
\label{deltaW}
	\Delta\mathbf{W}\cong-\mathbf{K}\left(\mathbf{H}^H\Delta\mathbf{H}+\Delta\mathbf{H}^H\mathbf{H}\right)\mathbf{K}\mathbf{H}^H+\mathbf{K}\Delta\mathbf{H}^H.
\end{equation}
Using the above approximation, we can now calculate the terms in \eqref{cov}. We first note that the third and fourth terms in \eqref{cov} are zero since $E\left[\Delta\mathbf{W}\right]\cong0$. The second terms is $E\left[\mathbf{W}\mathbf{n}\mathbf{n}^H\mathbf{W}^H\right]=N_0\mathbf{W}\mathbf{W}^H$, and the first term becomes $E\left[\Delta\mathbf{W}\mathbf{H}\mathbf{x}\mathbf{x}^H\mathbf{H}^H\Delta\mathbf{W}^H\right]=E_sE\left[\Delta\mathbf{W}\mathbf{H}\mathbf{H}^H\Delta\mathbf{W}^H\right]$. Below, we calculate the first and last terms in \eqref{cov} by plugging the error matrix \eqref{deltaW} into \eqref{cov}.
\begin{IEEEeqnarray}{rCl}
\label{term1incov}
\lefteqn{E\left[\Delta\mathbf{W}\mathbf{H}\mathbf{H}^H\Delta\mathbf{W}^H\right] }\IEEEnonumber\\ &\cong&E\left[\mathbf{K}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&+\:E\left[\mathbf{K}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&-\:E\left[\mathbf{K}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&+\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&+\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&-\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&-\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\Delta\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&-\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]\IEEEnonumber\\
&&+\:E\left[\mathbf{K}\Delta\mathbf{H}^H\mathbf{H}\mathbf{H}^H\Delta\mathbf{H}\mathbf{K}^H\right]
\end{IEEEeqnarray}
It can be proven that $E\left[\Delta\mathbf{H}\mathbf{A}\Delta\mathbf{H}\right]=E\left[\Delta\mathbf{H}^H\mathbf{A}\Delta\mathbf{H}^H\right]=0$ for any deterministic matrix $\mathbf{A}$. Hence the second, third, fourth and seventh terms in (13) are zero. For the remaining terms we use the fact that $E\left[\Delta\mathbf{H}\mathbf{A}\Delta\mathbf{H}^H\right]=\sigma_e^{2}\tr\left(\mathbf{A}\right)\mathbf{I}$, and obtain
\begin{IEEEeqnarray}{rCl}
\lefteqn{E\left[\Delta\mathbf{W}\mathbf{H}\mathbf{H}^H\Delta\mathbf{W}^H\right] }\IEEEnonumber\\ &\cong&\sigma_e^{2}\tr\left(\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right)\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\IEEEnonumber\\
&&+\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H\IEEEnonumber\\
&&-\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H\\
&&-\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H + \sigma_e^{2}\tr\left(\mathbf{H}\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H\IEEEnonumber
\end{IEEEeqnarray}
Similarly, the last term in \eqref{cov}, $E\left[\Delta\mathbf{W}\mathbf{n}\mathbf{n}^H\Delta\mathbf{W}^H\right]$, can be computed following the same way.
\begin{equation}
E\left[\Delta\mathbf{W}\mathbf{n}\mathbf{n}^H\Delta\mathbf{W}^H\right]=N_0E\left[\Delta\mathbf{W}\Delta\mathbf{W}^H\right]
\end{equation}
\begin{IEEEeqnarray}{rCl}
\lefteqn{E\left[\Delta\mathbf{W}\Delta\mathbf{W}^H\right] }\IEEEnonumber\\ &\cong&\sigma_e^{2}\tr\left(\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right)\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\IEEEnonumber\\
&&+\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H\IEEEnonumber\\
&&-\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H\IEEEnonumber\\
&&-\:\sigma_e^{2}\tr\left(\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H + \sigma_e^{2}N_r\mathbf{K}\mathbf{K}^H
\end{IEEEeqnarray}
Finally, we plug $E\left[\hat{\mathbf{n}}\hat{\mathbf{n}}^H\right]$ into \eqref{estppsnr} and obtain the PPSNR in the presence of channel estimation error as \eqref{eqn_dbl_x}.
\begin{figure*}[!t]
\normalsize
\setcounter{MYtempeqncnt}{\value{equation}}
\setcounter{equation}{16}

\begin{equation}
\label{eqn_dbl_x}
\widetilde{\gamma_k}\cong\!\frac{E_s\left|\left(\mathbf{WH}\right)_{k,k}\right|^2}{E_s\sum\limits_{l\neq k}\left|\left(\mathbf{WH}\right)_{k,l}\right|^2+
\left(
\begin{split}
&E_s\sigma_e^{2}\tr\!\left(\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right)\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H+E_s\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H  \\[-3pt] &-E_s\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{H}^H\right)\!\mathbf{K}\mathbf{K}^H\!\!-\!E_s\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\!\mathbf{K}\mathbf{K}^H\!\!+\!E_s\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{H}^H\right)\!\mathbf{K}\mathbf{K}^H  \\[-3pt]
&+N_0\mathbf{W}\mathbf{W}^H+N_0\sigma_e^{2}\tr\!\left(\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\right)\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H+N_0\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H \\[-3pt]
&-N_0\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{K}\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H-N_0\sigma_e^{2}\tr\!\left(\mathbf{H}\mathbf{K}^H\mathbf{H}^H\right)\mathbf{K}\mathbf{K}^H+N_0\sigma_e^{2}N_r\mathbf{K}\mathbf{K}^H
\end{split}
\right)_{\!k,k}}
\end{equation}
\setcounter{equation}{\value{MYtempeqncnt}}
\hrulefill
\end{figure*}

The BER of the system in the presence of channel estimation error can be found simply by plugging $\widetilde{\gamma_k}$ as the symbol SNR into the AWGN BER formulas. For example, the BER of $k^{th}$ stream for BPSK is $P_b^k=Q\left(\sqrt{2\widetilde{\gamma_k}}\right)$, and $P_b^k=\frac{3}{4}Q\left(\sqrt{\frac{\widetilde{\gamma_k}}{5}}\right)+\frac{1}{2}Q\left(3\sqrt{\frac{\widetilde{\gamma_k}}{5}}\right)-\frac{1}{4}Q\left(5\sqrt{\frac{\widetilde{\gamma_k}}{5}}\right)$ for gray-coded 16QAM.
\section{Results}
In order to test the performance of the analysis, we simulated transmission of thousands of packets through uncorrelated Rayleigh flat fading channels. For each SNR point on the BER plots, we randomly generate 1000 i.i.d. realizations of the channel matrix $\mathbf{H}$. For each specific realization of the channel, we transmit 500 packets each of which carries 2000 information symbols. We perform channel estimation for each packet as explained below in \textit{Case 1}

\textit{Case 1:} In our simulations, we employed the maximum likelihood (ML) channel estimation (CE) algorithm, in which the channel estimate is obtained via training symbols that are known to the receiver. During the training phase, the $N_t \times N_{tr}$ training matrix $\mathbf{X}_{tr}$ is transmitted where $N_{tr}\geq N_t$ is the number of training symbols. The $N_r \times N_{tr}$ received signal is $\mathbf{Y}_{tr}=\mathbf{H}\mathbf{X}_{tr}+\mathbf{W}$ where $\mathbf{W}$ is the $N_r \times N_{tr}$ noise matrix. Then, the ML estimate of the channel is given as \cite{channel_est}
\setcounter{equation}{17}
\begin{equation}
\widehat{\mathbf{H}}=\mathbf{Y}_{tr}\mathbf{X}_{tr}^H\left(\mathbf{X}_{tr}\mathbf{X}_{tr}^H\right)^{-1}
\end{equation}
It was shown that the optimal training signal has the property of $\mathbf{X}_{tr}\mathbf{X}_{tr}^H=E_sN_{tr}\mathbf{I}$. When this orthogonal training signal is employed, the entries of $\Delta\mathbf{H}$ are i.i.d. with $CN(0,\sigma_e^2)$, and the channel estimation noise variance\footnote{$\sigma_e^2$ can also be defined as $\sigma_e^2=\frac{N_t}{N_{tr}E_s/N_0}$ depending on SNR definition.} is $\sigma_e^2=\frac{1}{N_{tr}E_s/N_0}$ \cite{channel_est}. The estimation error in this case is caused by the AWGN in this case.

The following training signal, which is taken from 802.11n standard \cite{n80211}, was employed in the simulations. $\mathbf{X}_{tr}=\sqrt{E_s}\ddot{\mathbf{P}}$ where $\ddot{\mathbf{P}}$ is the submatrix formed by first $N_t$ rows and first $N_{tr}$ columns of the bigger matrix\footnote{The $\mathbf{P}$ matrix here is for maximum of 4 spatial streams since the standard supports up to 4 streams. $N_{tr}=4$ for $N_t\geq3$.} $\mathbf{P}$, i.e. $\ddot{\mathbf{P}}=\mathbf{P}\left[1:N_t;1:N_{tr}\right]$.
\begin{equation}
\mathbf{P}=\left[{\begin{array}{cccc}
1&-1&1&1\\
1&1&-1&1\\
1&1&1&-1\\
-1&1&1&1\\
\end{array}}\right]
\end{equation}
It should be noted that with this choice of the training matrix, the ML channel estimation at the receiver becomes a very simple operation since the matrix inversion, $\left(\mathbf{X}_{tr}\mathbf{X}_{tr}^H\right)^{-1}$, is now a trivial operation.

We present the simulation results for BPSK in Fig. \ref{bpsk_all_latex3} and 16QAM in Fig. \ref{qam16_all_latex3} with $1\times4$, $2\times4$, $4\times4$ MIMO configurations. $N_{tr}=4$ is used in all the simulations. The case of $\sigma_e^2=0$, i.e. perfect channel estimation, is also included in the results.

For each channel instance, analytical BER results are obtained by using the PPSNR derived in the previous section. Then, these BERs are averaged over all realizations of the channel.

First thing to notice in Fig. \ref{bpsk_all_latex3} and Fig. \ref{qam16_all_latex3} is that for $\sigma_e^2=0$, simulation and analysis curves exactly match. Performance is significantly degraded for the systems experiencing channel estimation errors. This is particularly evident for the $4\times4$ configurations.
\begin{figure}[!t]
\centering
		\includegraphics[trim = 10mm 0mm 8mm 6mm, width=3.4in]{bpsk_all_paperv2}
\caption{BER for BPSK. $E_b/N_0$ is per tx antenna.}
\label{bpsk_all_latex3}
\end{figure}

As it can be seen in Fig. \ref{bpsk_all_latex3} and Fig. \ref{qam16_all_latex3}, our analysis gives a very tight approximation of the real performance. For BPSK $4\times4$, and all of the 16QAM configurations the analysis results exactly match the simulated performances. For BPSK $1\times4$, and $2\times4$ configurations the analysis results are upper-bounds to the real performance at high SNR, however, they are still very close to the real performances. The analysis results become tighter for higher order modulations and higher order MIMO configurations. This is because of the fact that the Gaussian assumption, which is made for the post-detection noise, is more valid at higher order modulations and MIMO configurations. At low SNRs, the total post detection noise $\hat{\mathbf{n}}$ is dominated by the additive white Gaussian noise component $\mathbf{n}$ therefore the assumption is valid even for lower configurations. However, at high SNRs the residual interference components from other spatial streams becomes dominant and $\hat{\mathbf{n}}$ is loosely approximated as Gaussian for lower order constellations and MIMO configurations.

It is interesting to note that in contrast to the results obtained for ZF detector by \cite{wang}, we do not observe any error floor on the performance. This is due to the fact that the channel estimation error variance $\sigma_e^2$ for ML estimation gets smaller as SNR increases. This is the situation that occurs in practical packet based or bursty communication systems where the channel estimation is performed for every packet prior to data detection, and hence experiences the same noise variance as the data transmission. Therefore the channel estimation quality is dependent on the SNR. On the other hand, error floors are observed in \cite{wang} because of the assumption that $\sigma_e^2$ remains constant independent of the SNR. This case is investigated below in \textit{Case 2}.

\begin{figure}[!t]
\centering
		\includegraphics[trim = 10mm 0mm 8mm 6mm,width=3.4in]{qam16_all_paperv2}
\caption{BER for 16QAM. $E_b/N_0$ is per tx antenna.}
\label{qam16_all_latex3}
\end{figure}
\begin{figure}[!t]
\centering
		\includegraphics[trim = 10mm 0mm 8mm 6mm,width=3.4in]{qpsk_4x5}
\caption{BER for $4\times 5$ QPSK for different $\sigma_e$ values. $\sigma_e=0$ corresponds to perfect CE. ZF curves are included for comparison purposes with \cite{wang}. $E_b/N_0$ is per tx antenna. These parameter values are taken from Fig. 2 of \cite{wang} for comparison purposes.}
\label{qpsk}
\end{figure}

\textit{Case 2:} In addition to ML channel estimation results, we also performed simulations with constant $\sigma_e$. Unlike the first case, the channel estimation quality is independent of the SNR. This situation might arise either when there is a ready channel estimate to be used by the receiver formed elsewhere with a different additive noise variance, or the channel estimation is outdated and the major error in the channel estimation comes from the mobility changes in the channel.

In Fig. \ref{qpsk}, the BER performance of a QPSK $4\times 5$ system is investigated for $\sigma_e=5\%,10\%,20\%$ using the estimation error model in (4). Each packet observes a different realization of the random matrix $\Delta\mathbf{H}$ with the designated variance $\sigma_e^2$. As expected, we observe error floor in the performance due to the constant estimation error variance as in the ZF detector case studied in \cite{wang}. More importantly, these error floors are the same as the ones observed by ZF detector because of the fact that the MMSE and ZF detectors exhibit the same behaviour at asymptotically high SNR. The simulation results in this case also agree with the analysis.

\section{Conclusion}
In this letter, we presented the analysis of post-processing SNR for practical MIMO MMSE receivers which experience imperfect channel estimation. Performance of MMSE receivers in the presence of channel estimation error is investigated and shown to be accurately estimated via analytical results. We verified the tightness of the analytical results via simulations.

Besides the theoretical contributions, we believe that our closed form PPSNR expression can be useful for link adaptation purposes in real MIMO systems. There exist link adaptation algorithms \cite{peng,wcnc12} based on PPSNR, however perfect CSI is always assumed which might lead to incorrect prediction of the throughput. More accurate prediction can be achieved using the results presented in this paper

\ifCLASSOPTIONcaptionsoff
  \newpage
\fi


\section{Introduction}

Although Galactic globular clusters (GCs) display a distribution in their global
parameters (e.g., mass, metallicity, concentration, horizontal-branch
morphology), the internal variation of elements affected by proton captures (hereafter $p$-capture elements) appears a
ubiquitous feature (\citealt{carretta09a},b). It is clear that GC stars
exhibit large changes in the C, N, O, Na, Mg, and Al abundances, whereas (in
archetypical systems at least) internal spreads in iron-peak, heavy
$\alpha$- (Ca, Ti) and $slow$ neutron-capture ($s$-process)
elements all remain within
observational uncertainties (\citealt{carretta09c}, 2010a;
\citealt{james04}; \citealt{smith08}; \citealt{dorazi10}). The changes in the
$p$-capture elements give rise to a clear chemical pattern: depletion
in C, O, and Mg abundances always correspond to enhancements in N, Na,
and Al (the so-called light-element anticorrelations). This behavior
bears the evidence of H burning at high temperature and points to the presence of
multiple stellar generations. It is argued that the ejecta from a
fraction of first generation of stars (initially C-O-Mg rich, sharing the same
chemical composition of field stars at the same metallicity) mix
with primordial gas, providing a medium from which the second generation
stars (C-O-Mg poor and N-Na-Al rich) formed. We refer the reader to
\cite{gratton12} for a recent, extensive review on this topic. In
this scenario, the H-burning abundance patterns from the first
generation stars are imprinted in the second generation and are present
from birth. The nature of the stars that enriched the intercluster gas remains
uncertain but possible candidates include intermediate-mass asymptotic
giant branch (AGB) stars
undergoing hot bottom burning (HBB, e.g., \citealt{dantona83}; \citealt{ventura01}), fast rotating
massive stars (e.g., \citealt{decressin07}), massive binaries
(\citealt{demink09}), and novae (\citealt{smith96}; \citealt{maccarone12}).

Interestingly, this already complex picture is further obfuscated by the
presence of some peculiar clusters, such as $\omega$ Centauri
(\citealt{jp10}; \citealt{marino11a}), NGC~1851 (\citealt{yong08};
\citealt{carretta10a}), Terzan~5 (\citealt{ferraro09}), NGC~6715 (M54,
\citealt{carretta10b}), and NGC~2419 (\citealt{cohen11}, see however
\citealt{mucciarelli12} for a different view). In these GCs, along with
changes in $p$-capture elements, internal variations in the heavy
element abundances have been detected. Species ranging from the
iron-peak (e.g., Fe) to the $s$-process elements (Ba, La) vary
stochastically from cluster to cluster beyond what is
expected from observational errors.

The metal-poor GC NGC 6656 (M22, [Fe/H]\footnote{We use the notation [X/H]=A(X)$-$A(X)$_{\odot}$,
where A(X)=log($\frac{N_{X}}{N_H}$)+12}=$-$1.70 \citealt{harris96}
-updated in 2010) belongs to this class of GCs,
and due to its peculiar nature has
received extensive attention (\citealt{pilachowski82};
\citealt{norris83}; \citealt{bw92}; \citealt{kayser08};
\citealt{marino09}; \citealt{dacosta09}; \citealt{dacosta10}; \citealt{abrito12}). Recently,
Marino et al. (2011b,\nocite{marino11b} hereafter MSK11) presented results from their high-resolution
spectroscopic study of 35 giant stars, deriving abundances for
iron-peak, $\alpha$, $p$-capture, and neutron-capture elements.
This detailed abundance analysis provided a unique
opportunity to investigate the chemical enrichment history of M22.

This GC is comprised of two distinct groups of stars, characterised by an
offset in metallicity and in $s$-process element content. The first
group displays a metallicity of $<$[Fe/H]$>$=$-$1.82$\pm$0.02 with
$<$[$s$/Fe]$>$=$-$0.01$\pm$0.01 and the second group has
$<$[Fe/H]$>$=$-$1.67$\pm$0.01 with $<$[$s$/Fe]$>$=+0.35$\pm$0.02. Note that
the [$s$/Fe] ratios were computed by averaging the abundances of Y, Zr,
Ba, La, Nd (see MSK11 for details). Each of these two sub-groups
exhibits the classical Na-O and C-N anticorrelations
shown by the archetypical GCs. Given that
the [Eu/Fe] ratio serves as a $rapid$ neutron-capture ($r$-process)
tracer\footnote{In the
Solar System, $\approx$ 97\% of Eu was synthesized via the $r$-process (\citealt{burris00}, and references therein).} and that
enhancements in the $s$-process elements are not accompanied by a
similar trend in the [Eu/Fe] ratio, it can be inferred that the
dichotomy in the $n$-capture abundances may be due to a first generation
of polluters that produced the $s$-process only. One possible scenario,
which was advocated by MSK11, is that the $weak$
$s$-process component activated in stars with masses larger than
$\sim$25 M$_{\odot}$ during core He-burning and C-shell phases (\citealt{raiteri93}; \citealt{pignatari10}), may have contributed to the observed
abundance patterns. However, in a complimentary study, \citet{roederer11}
focused on the heavy element content (from Y to Th) of six stars across
the two stellar sub-groups, and ruled out the massive star origin. They
concluded that the gas from which the second stellar population formed
was enriched in $s$-process material from a class of relatively massive
AGB stars
(M $\approx$ 5 M$_{\odot}$). In these stars, the production of $s$-process
elements is due to the activation of the
$^{22}$Ne($\alpha$,n)$^{25}$Mg neutron source, whereas in their
lower-mass counterparts the main neutron source
is the $^{13}$C($\alpha$,n)$^{16}$O reaction
(\citealt{busso99}).

Neither of the proposed scenarios provides a comprehensive explanation
for all the observed chemical features and we are left with numerous
unsolved issues. For example, \cite{roederer11} question why the
$s$-process dichotomy is only present in M22: if massive AGB
stars are the cause of the GC Na-O anticorrelation, then all clusters
should present the $s$-process elements correlated with Na and anti-correlated
with O, which is not observed (e.g., \citealt{dorazi10}).

In this paper we turn to an alternative diagnostic. We present fluorine
abundances for a sample of six cool giant stars in M22, carefully
selected from both sub-stellar groups as defined by MSK11.
Fluorine abundances are a powerful tracer of the polluter mass range
in M22 because the F production is highly dependent
on the stellar mass.

Theoretical models of AGB stars (e.g.,  \citealt{jorissen92}) predict that F is produced due to the
activation of the chain of reactions
$^{18}$O(p,$\alpha$)$^{15}$N($\alpha$,$\gamma$)$^{19}$F in the He intershell
during the recurrent thermal pulses associated with He burning. During
the early phases of each thermal pulse, the H-burning ashes are ingested
in the convective region developing in the He intershell. These ashes
are rich in $^{13}$C and $^{14}$N and their ingestion in the He-burning
layer results in production of $^{18}$O due to $\alpha$ captures on
$^{14}$N. At the same time protons are released by the
$^{14}$N(n,p)$^{14}$C reaction, with neutrons coming from
$^{13}$C($\alpha$,n)$^{16}$O. After the quenching of
each thermal pulse the envelope may sink in mass deep in the He
intershell and carry $^{19}$F to the convective envelope via a process
known as the ``third dredge-up'' (TDU). The peak of F production in AGB
stars is reached for stars of initial masses $\sim$ 2 M$_{\odot}$~\citep{lugaro04}. If the mass of the stars is higher than roughly
5 M$_{\odot}$, and depending on the metallicity, fluorine is
destroyed both via $\alpha$ captures in the He intershell, and via
proton captures at the base of the convective envelope due to HBB. AGB
stars that experience HBB can also destroy O and Mg and produce Na and
Al; as a consequence, according to the multiple population scenario, we
should expect the abundances of F to be correlated with O (and Mg) and
anti-correlated with those of Na (and Al). This prediction was
observationally confirmed by \cite{smith05} in the
intermediate-metallicity GC M4 and by \cite{yong+08} in NGC 6712.

Other sites have also been suggested for F production: the $\nu$-process
in core-collapse supernovae (\citealt{woosley90}), and core He burning
in Wolf-Rayet stars (\citealt{meynet00}; \citealt{palacios05}). However,
low-mass AGB stars are the only site observationally confirmed
(\citealt{jorissen92}; \citealt{abia10}).

This paper is organized as follows: observations are described in
Section~\ref{sec:obs}, while details on abundance analyses are given in
Section~\ref{sec:analysis}. Our results are then presented in
Section~\ref{sec:results} and discussed in Section~\ref{sec:discussion}.
A summary closes the manuscript (Section \ref{sec:summary}).

\section{Observations}\label{sec:obs}

Our sample includes six giant stars, for which stellar parameters (T$_{\rm eff}$, log$g$, [Fe/H] and microturbulence $\xi$)
along with $p$-capture and $s$-process element abundances were derived by MSK11.
We selected three stars belonging to the metal-poor component (MP, also $s$-process poor) and three stars
from the metal-rich (MR, $s$-process rich) one. Within each of these sub-groups we also selected both O-rich (Na-poor) and
O-poor (Na-rich) stars, spanning a range from [O/Fe]=+0.11 to [O/Fe]=+0.48 dex.

The main objective of our investigation was the determination of the fluorine abundances.
Observationally, F (whose only stable isotope is $^{19}$F) is difficult to detect spectroscopically;
the only atomic lines (the ground-state transitions of F~{\sc i})
that might be revealed lie in the far UV. On the other hand, HF molecular transitions
are easily observable in the near-infrared (around $\sim$23000 \AA).
Our analysis focuses on the HF(1$-$0) R9 transition, located at $\lambda$23358.3 \AA. Although not being the strongest
feature in this wavelength range, this line
is considered to be one of the best abundance indicator for F, because it is free of blends (\citealt{abia09}; \citealt{lucatello11}).

\begin{table*}
\caption{Information on Target Stars}\label{t:log}
\begin{center}
\begin{tabular}{lcccccr}
\hline\hline
Star   &     RA      &      Dec       &   V    &  K   & Exposures & S/N \\
       &  (hh:mm:ss) &   ($^{o}$:':'')&  (mag) &  (mag) & (s)  &    \\
       &             &                &        &      &          &           				\\
\hline
       &             &                &        &       &           &        					\\
IV-97  & 18:36:41.06 &	 $-$23:58:18.9 & 11.043 & 6.759 &  10x90s  &	300 \\
III-14 & 18:36:15.10 &	 $-$23:54:54.6 & 11.134 & 6.743 &  17x120s  &   400  \\
III-15 & 18:36:15.61 &	 $-$23:55:01.2 & 11.362 & 7.138 &   8x180s  &   300 \\
C      & 18:36:10.21 &	 $-$23:48:44.0 & 11.309 & 6.737 &   6x180s  &   300  \\
III-52 & 18:36:10.18 &	 $-$23:54:21.8 & 11.526 & 7.459 &  10x180s  &   350 \\
V-2    & 18:36:28.02 &	 $-$23:55:01.6 & 11.498 & 7.276 &  14x120s  &   350 \\
       &             &               &        &       & 	  &    \\
\hline\hline
\end{tabular}
\end{center}
\end{table*}

High-resolution, near-infrared spectroscopic observations were carried out in service mode with CRIRES (CRyogenic high-resolution
InfraRed Echelle Spectrograph, \citealt{kaeufl04}) located at VLT UT1 on 2011 April, July and August
(program: 087.0319(A), PI: VD). In Table~\ref{t:log}, we list information on target stars, reporting identifications, magnitudes
(see MSK11 for details), exposure times, and S/N ratios per pixel around the HF feature.
Note that our selection was limited to the cooler stars in the MSK11's sample due to our imposed requirement of
relatively strong HF lines.
We also observed several (hot) early-type stars before and/or after of each target observations,
in order to remove telluric contamination from our scientific frames.

We employed the 0.4$^{''}$ slit and the grating order \#24, achieving a resolution of $R$$\sim$50,000
and a wavelength coverage from $\lambda$22948.5 \AA~ to $\lambda$23410.3 \AA.
This allowed us to include the HF(1$-$0) R9 line, numerous $^{12}$C$^{16}$O vibration/rotation lines (used to derive C abundances)
and the Na~{\sc i} line at $\lambda$23379 \AA.

Data reduction was accomplished by means of the CRIRES pipeline (version 2.4), running under the {\sc gasgano}
environment\footnote{http://www.eso.org/sci/software/gasgano/}, which provides one-dimensional, wavelength calibrated spectra.
Telluric feature subtraction, rest-frame translation and continuum normalization were then carried out within
IRAF\footnote{IRAF is the Image Reduction and Analysis Facility, a general purpose software system
for the reduction and analysis of astronomical data. IRAF is written and supported by National Optical Astronomy
Observatories (NOAO) in Tucson, Arizona.}.
An example of our spectra is shown in Figure~\ref{f:spectra} for star III-52; HF, CO and Na~{\sc i} lines are marked.

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{spectrumIII52.eps}
\caption{Example of a spectrum for star III-52 }\label{f:spectra}
\end{figure}
\end{center}

\section{Abundance analysis}\label{sec:analysis}

Fluorine abundances were determined through spectral synthesis using the MOOG code (\citealt{sneden73}, 2011 version) and the
Kurucz (1993\nocite{kur93}) set of stellar atmosphere models (with no overshooting) as in MSK11's analysis, from which we
retrieved stellar parameters as well as O and $s$-process element abundances.
However, had we instead adopted the MARCS grid (\citealt{gustafsson08}), the difference on the resulting abundances
would have been less than 0.04 dex.

Concerning the HF feature, we took as excitation potential (EP) $\chi$=0.227 eV (\citealt{decin00})
and a log$gf$=$-$3.971 (\citealt{lucatello11}). This last value is very close to that used in previous studies
focusing on F abundance determination in GCs (e.g., \citealt{smith05}; \citealt{yong+08}; \citealt{abrito12}),
log$gf$=$-$3.955 (from \citealt{jorissen92}).
On the other hand, our EP is about 0.25 eV lower compared to
those works: this difference implies an offset between our abundances and those previously published in the literature
of roughly 0.30 dex (see also the discussion in Section~\ref{sec:alan}).
Furthermore, we derived C abundances. For this purpose, we assumed O values from the optical range given by MSK11,
since our spectra did not cover any suitable OH line, whose stronger transitions extend in the $H$ band (around $\sim$15000\AA).
The CO line lists come from B. Plez (private communication).
Finally, for the Na~{\sc i} line at $\lambda$23379, atomic parameters  ($\chi$=3.750 eV; log$gf$=0.530) were
taken from VALD\footnote{Vienna Atomic Line Database
(www.astro.uu.se/~vald/php/vald.php)}.

As a first step, we checked our line list on the infra-red atlas of the Arcturus spectrum (\citealt{hinkle},
available at ftp://ftp.noao.edu/catalogs/arcturusatlas/).
Assuming a T$_{\rm eff}$=4286 K, log$g$=1.67, $\xi$=1.74 kms$^{-1}$ and
[Fe/H]=$-$0.52 (following \citealt{ramirez11}) we obtained
an [F/Fe]=$-$0.15\footnote{We used as solar fluorine abundances A(F)$_{\odot}$=4.56 (\citealt{ag89}; \citealt{asplund05}).},
to be compared to the value given by \cite{abia09} of [F/Fe]=0.10 dex.
The difference is completely explained by the higher EP adopted in that study.

Moreover, we inferred a C abundance of A(C)=8.01
(under the assumption that A(O)=8.81, that is [O/Fe]=0.40, with a solar abundance of A(O)$_{\odot}$=8.93), which is in very good
agreement with values of \cite{abia09} (i.e., A(C)=8.06) and \cite{ryde10} (A(C)=8.08).
As for C solar abundance we adopted the value of A(C)$_{\odot}$=8.56, leading to [C/Fe]=$-$0.03.

Finally, from the Na~{\sc i} feature at 23379\AA, we derived A(Na)=6.01 which results in [Na/Fe]=+0.2 dex (setting A(Na)$_{\odot}$=6.33).

Comparison between synthetic and observed spectra were carried out in
a similar way for our sample stars; an example of spectral synthesis is given in
Figures~\ref{f:synth1} and~\ref{f:synth2} for star C.

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{synth1.eps}
\caption{Synthesis of the HF feature for star C.}\label{f:synth1}
\end{figure}
\end{center}

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{synth2.eps}\\
\includegraphics[width=8cm]{synth3.eps}
\caption{Same as Figure~\ref{f:synth1} but for CO (upper panel) and Na (lower panel)} \label{f:synth2}
\end{figure}
\end{center}

The sensitivity of the F abundance to input stellar parameters
was evaluated by separately changing effective temperature, surface gravity and
microturbulence values.
The intensity of the synthetic line of the HF is
particularly sensitive to the adopted T$_{\rm eff}$, the other parameters affecting it at a lower degree
(see also \citealt{abia09}; \nocite{abia11}2011).
A change of $\Delta$(T$_{\rm eff}$)=+70 K, $\Delta$(log$g$)=+0.15, and $\Delta$($\xi$)=+0.13 km s$^{-1}$
(conforming to error estimates given in Table 4 of MSK11) results in a difference in A(F) of
+0.10, 0.02 and $-$0.02 dex, respectively. The variation of the input metallicity in the model atmosphere has instead a negligible effect.
These are the typical uncertainties that we then summed in quadrature providing a total error, due to stellar parameters,
of 0.11 dex in our [F/H] ratios.
Errors due to the best-fit determination (related to the S/N of the spectra and including uncertainties due to the continuum placement) are instead $\pm$0.07 dex. However, we caution the reader that this value should be treated as a lower limit, since the impact of the telluric correction on this considerably weak feature is significant
(see Section~\ref{sec:results}).

The total internal error in [F/H] is then obtained adding in quadrature both uncertainties, resulting in 0.13 dex.

Finally, as far as C and Na are concerned, the typical uncertainties are 0.12 and 0.10 dex, respectively.

\section{Results}\label{sec:results}

\begin{table*}
\caption{Stellar Parameters and Elemental Abundances}\label{t:results}
\begin{center}
\begin{tabular}{lccccccccr}
\hline\hline
Star  &  T$_{\rm eff}$$^{*}$ & $\log g$$^{*}$   & [Fe/H]$^{*}$  & $\xi$$^{*}$        & [O/Fe]$^{*}$
& $s$-rich$^{*}$   & [C/Fe] & [F/Fe] & [Na/Fe] \\
      &    (K)         &  (cms$^{-2}$) &        & (kms$^{-1}$) &      &              &        &        &           \\
\hline
      &                &                     &        &             &        &            &        &        &            \\
IV-97  &  4000    &	    0.05 & $-$1.94   & 2.00  & 0.40  & no  & $-$1.10  &	$-$0.70 & $-$0.40	 \\
III-14 &  4030    &	    0.35 & $-$1.82   & 2.15  & 0.48  & no  & $-$1.10  &	$-$0.80 & $-$0.30	 \\
III-15 &  4070    &	    0.40 & $-$1.82   & 1.85  & 0.11  & no  & $-$1.15  &	$-$1.00 &  0.34	 \\
       &          &              &           &       &       &     &      &          &               \\
C      &  3960    &	    0.30 & $-$1.69   & 2.25  & 0.25  & yes & $-$0.60  &	$-$0.80 &  0.30	 \\
III-52 &  4075    &	    0.60 & $-$1.63   & 1.75  & 0.45  & yes & $-$0.10  &	$-$0.60 &  0.05	 \\
V-2    &  4130    &	    0.65 & $-$1.57   & 1.75  & 0.15  & yes & $-$0.40  &	$-$0.90 &  0.31	 \\
       & 	  &	         &	     &	     &	     &	   &	      &	      &	\\
\hline\hline
\multicolumn{1}{l}{$^{*}$ from MSK11.}
\end{tabular}
\end{center}
\end{table*}

Our results are shown in Table~\ref{t:results}, where we report stellar parameters and abundances from MSK11 along
with our estimates for F, C, and Na.
Even within our quite limited sample (six stars), we found that
{\em the F abundance shows a large star-to-star variation}, ranging from [F/H]=$-$2.82 dex to [F/H]=$-$2.23 dex
(i.e., a factor of $\sim$4). The average abundance is $<$[F/H]$>$=$-$2.55$\pm$0.08 (rms=0.20), implying that
the amplitude of this change is beyond the measurement uncertainties (see Section~\ref{sec:analysis}). Moreover, taking into account
the typical errors for O and F abundances, our study suggests that
the F variation is comparable with that of O, as also found by \cite{yong+08} in the GC NGC~6712.
\begin{center}
\begin{figure}
\includegraphics[width=8cm]{OF.eps}
\caption{Fluorine abundances ([F/H]) as a function of [O/H].}\label{f:OF}
\end{figure}
\end{center}

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{NaF.eps}
\caption{[F/H] versus [Na/H] from MSK11 (left panel) and from this study (right panel).
Symbols are as in Figure~\ref{f:OF}.}\label{f:NaF}
\end{figure}
\end{center}

In Figures~\ref{f:OF} and~\ref{f:NaF} we plot [F/H] ratios as a function of [O/H] and [Na/H], respectively.
As one can see, F abundances are positively correlated with O: considering the whole sample,
the Pearson correlation coefficient results in
$r$=0.89, with a probability to be random smaller than 2\%.
Focusing on the F-Na diagram (Figure~\ref{f:NaF}), there is the hint for a F-Na anticorrelation, but
in our small sample the linear correlation
coefficient is not statistically meaningful.
However, this does not prove the lack of an anticorrelation, because we are dealing
with small numbers (only six points), heavily reducing the power of statistical tests.
In addition, there is no a priori reason why we should combine the two sub-groups as far as the F-Na plane is concerned,
because we do not expect that they must behave in the same way. If we look at each component separately the presence of a
F-Na anticorrelation is indeed much more evident (given that we have only three points for each group, it is not meaningful to
perform a statistical test separately).

The observed chemical pattern can be explained as the evidence of H-burning at high temperature, via the CNO cycle,
which causes the destruction of F, in conjunction with O depletion and Na enhancement.
More interestingly, those (anti)correlations are revealed in each of the M22 sub-components (the $s$-rich and $s$-poor groups);
the implications of this finding are discussed in detail in Section~\ref{sec:discussion}.

Regarding Na, we show both our estimate from near-infrared spectroscopy as well as LTE abundances from the optical range by
MSK11 (right and left-hand panels of Figure~\ref{f:NaF}).
A difference of $\Delta$([Na/Fe])=0.31$\pm$0.04 (rms=0.09) dex is found between the two
estimates (see Figure~\ref{f:NaIROpt}) where we compare the two measurements; non-LTE effects can totally account for such a discrepancy (e.g., \citealt{lind11}).
The total average Na abundance is $<$[Na/Fe]$>$=0.05$\pm$0.13 (rms=0.33); considering separately the two groups we obtain
instead a $<$[Na/Fe]$>$$_{s-\rm poor}$=$-$0.12$\pm$0.23 and
$<$[Na/Fe]$>$$_{s-\rm rich}$=+0.22$\pm$0.08, which implies $\Delta^{\rm rich}_{\rm poor}$[Na/Fe]=0.34$\pm$0.17 dex.
This value is in good agreement with that derived by MSK11, based on the whole sample of 35 giants, being
$\Delta^{\rm rich}_{\rm poor}$[Na/Fe]=0.23$\pm$0.07.

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{NaIR_Opt.eps}
\caption{[Na/Fe] from the optical range by MSK11 and from this study}\label{f:NaIROpt}
\end{figure}
\end{center}
Finally, our sample displays an average C abundance of $<$[C/Fe]$>$=$-$0.74$\pm$0.18, with
$<$[C/Fe]$>$$_{s-\rm poor}$=$-$1.12$\pm$0.02 and $<$[C/Fe]$>$$_{s-\rm rich}$=$-$0.37$\pm$0.14.
Thus, on average, $s$-rich stars exhibit larger C abundances, which qualitatively agrees with previous works
(e.g., \citealt{brown90}). However, while the difference between the two groups from the optical CH bands derived
by MSK11 is $\Delta^{\rm rich}_{\rm poor}$[C/Fe]=0.35$\pm$0.13 dex, we achieved a
much larger value of $\Delta^{\rm rich}_{\rm poor}$[C/Fe]=0.75$\pm$0.10 dex.
The same conclusion was reached by \cite{abrito12} who found a variation of
$\Delta^{\rm rich}_{\rm poor}$[C/Fe]=0.78$\pm$0.15 dex,
from high-resolution, near-infrared spectroscopy of nine cool giants (see also Section~\ref{sec:alan}).
The reason of such a difference in C abundance from the optical and from the near-infrared is not clear and no
obvious trends with stellar parameters (e.g., temperature, gravity, microturbulence, and/or metallicity)
seem to be present.
Further investigations are needed to explore this issue.
\subsection{Comparison with \cite{abrito12}}\label{sec:alan}

In a recent paper, \cite{abrito12} carried out high-resolution
($R$=50,000), near-infrared (both $H$ and $K$ bands) spectroscopic
observations with Phoenix@Gemini-South of nine red giant branch (RGB) stars in M22. They investigated their F
content, presenting also C, N, O, Na and Fe abundances.

Four of our six stars are in common with that study, namely,
III-14, III-15, III-52 and IV-97. For stars III-14 and III-52, those
authors inferred [F/Fe]=$-$0.40 and [F/Fe]=$-$0.20 dex, respectively,
while we obtained [F/Fe]=$-$0.80 and [F/Fe]=$-$0.60 dex. The adopted EP
value can easily justify this divergence, accounting for about 0.30 dex
(see Section~\ref{sec:analysis}); the source of the remaining $\sim$0.1
dex can be ascribed to the continuum placement, which is critical in
determining abundances from such a weak line. On the other hand, for
stars III-15 and IV-97 \cite{abrito12} obtained [F/Fe]=0.28 and
[F/Fe]=0.25 dex, entailing discrepancies with our estimates
larger than a factor 10. Note that stellar parameters (T$_{\rm eff}$, log$g$, and $\xi$) are the
same in both works, as they come from the analysis of MSK11; for the input metallicity, Alves-Brito et al.
used instead their own values coming from the IR spectroscopy and showing an offset of +0.13 dex
compared to the optical ones.
We investigated the nature of this substantial discordance, and
attributed it to the telluric feature subtraction. In the upper panel of
Figure~\ref{f:comp_spectra}, we directly compare our spectrum for
star III-15 (solid line) with that used by Alves-Brito et al. (dotted
line\nocite{abrito12}). Their spectrum presents
stronger features in the vicinity of HF line,
features expected from telluric contribution. To completely remove the contamination, Alves-Brito et al.
realized that they needed early-type star targets before and/or after each scientific frame but the logistics
of their run made this very difficult.
We instead could observe such targets: the correction to our data, in turn affects the placing of the continuum and
removes many strong features (as shown in the lower panel of Figure~\ref{f:comp_spectra}).
As expected, such an effect is significant for the
HF feature, due to its intrinsic weakness, but only marginally affects the C and Na abundance determinations (due to
the strength of their lines), and hence most of the conclusions of that paper.
This is shown in Figure~\ref{f:comp}, where we plot our [X/Fe] ratio
as a function of those from \cite{abrito12} for the four stars in common: C
and Na are comparable between the two studies, with differences of $\Delta$[C/Fe]=+0.20$\pm$0.17 dex
and $\Delta$[Na/Fe]=+0.18$\pm$0.11 dex (in the sense Alves-Brito's study minus our values);
if we take into account the offsets in [Fe/H] and in the adopted solar
abundances (they assumed A(C)$_{\odot}$=8.42 and A(Na)$_{\odot}$=6.17),
those values become $\Delta$[C/Fe]=+0.19$\pm$0.17 and $\Delta$[Na/Fe]=+0.15$\pm$0.11 dex.
On the other hand, discrepancies in F are significant and cannot be recovered from the different EP and/or solar abundances,
being for the whole sample of $\Delta$[F/Fe]=+0.75$\pm$0.19 dex.

\begin{center}
\begin{figure}
\includegraphics[width=9cm]{telluric.eps}
\caption{Upper panel: comparison of our spectrum for star III-15 (solid line) with the one by \cite{abrito12} (dotted line).
Lower panel: superimposition of our spectra for the same sample star with and without the telluric line
subtraction (dot-dashed and solid lines, respectively).}\label{f:comp_spectra}
\end{figure}
\end{center}
\begin{center}
\begin{figure}
\includegraphics[width=9cm]{ALAN_comp.eps}
\caption{Comparison for the four stars in common with Alves-Brito et al. (2012, here labeled as AB12).
Symbols are as in Figures~\ref{f:OF}, \ref{f:NaF}, \ref{f:NaIROpt}.}\label{f:comp}
\end{figure}
\end{center}
\section{Discussion}\label{sec:discussion}

Our main result is the detection of F variations across our sample
significantly larger than the observational uncertainties. As shown in Figures~\ref{f:OF} and \ref{f:NaF}, the changes in the F
abundances are correlated with O and anticorrelated with Na.
This chemical pattern qualitatively matches the
predictions from the multiple population scenario, according to which
the stellar ejecta from which
second generation stars formed carry the signature of hot H burning
causing enhancements in Na (N and Al) and depletions in O and F (C and
Mg). The F-O diagram presented in Figure~\ref{f:gcs} demonstrates that
M22 shares a similar behaviour as M4 and NGC6712, the other two GCs
for which F has been explored\footnote{\cite{cunha03} presented F
abundances for two giants in $\omega$ Cen. However,
they provide an F measurement only for star ROA219, giving
an upper limit for star ROA 324. Discussion related to the internal F
variation in this peculiar
GC is still not possible with the currently available measurements. For
this reason we acquired CRIRES spectra of 12 $\omega$~Cen giants;
results will be presented in a forthcoming paper (S. Lucatello et al., in
preparation).}.

Furthermore, the F-O-Na (anti) correlations can be marked separately within
each of the two sub-groups enclosed in M22 as clearly illustrated in
Figures~\ref{f:OF}
and~\ref{f:NaF}, where the $s$-process poor and the $s$-process rich stars are
labeled with empty and filled symbols, respectively. The same
conclusion was drawn by MSK11 when considering the Na-O and
C-N planes.

\begin{center}
\begin{figure}
\includegraphics[width=8cm]{GCs.eps}
\caption{[F/Fe] versus [O/Fe] for M22 (this study), M4 (\citealt{smith05}), and NGC 6712 (\citealt{yong+08}).}\label{f:gcs}
\end{figure}
\end{center}
\begin{center}
\begin{figure}
\includegraphics[width=8cm]{LaF.eps}
\caption{[F/H] versus [La/H]. Symbols as in Figures~\ref{f:OF}, \ref{f:NaF}, \ref{f:NaIROpt}, and \ref{f:comp}. The dashed line is a least-squares fit to data points.}\label{f:LaF}
\end{figure}
\end{center}

Very interestingly, beyond the internal spread in F characterizing each
sub-component, we measured an increase in the F content between the two
different stellar generations in M22. {\em The $s$-process-rich group
has, on average, larger F abundances than the $s$-process-poor group.}
This is shown in Figure~\ref{f:LaF}, where we plot our F abundances
([F/H]) as a function of [La/H] from MSK11. There is a
positive correlation between the two ratios, suggesting that
the polluters responsible for the $s$-process production
must account for a simultaneous F production.

There are two classes of objects producing both $s$-process elements
and fluorine. The first one is very massive stars (mass roughly $>$
40M$_{\odot}$). These produce F in the initial phases of core He
burning and expel it in the interstellar medium via winds during the
Wolf-Rayet phase (\citealt{meynet00}). They also produce $s$-process
elements during core He and shell C burning (e.g., \citealt{pignatari10}).
Production of both F and $s$-process elements in these massive stars
depends on the initial CNO abundances and thus decreases with the
stellar metallicity. Inclusion of stellar rotation enhances the
$s$-process production at low metallicity
(\citealt{pignatari08}; \citealt{chiappini11}); however, it appears to
decrease the production of fluorine (\citealt{palacios05}).
One problem already stressed by \cite{roederer11}
when considering these stars is that there is no reason why the SNe that enriched the $s$-process rich group host the $weak$
component and those that polluted the $s$-process poor do not (see Roederer et al. for details\nocite{roederer11}).

The second class of objects producing both F and $s$-process elements
are AGB stars (\citealt{forestini92};
\citealt{jorissen92}; \citealt{mowlavi98}; \citealt{karakas03};
\citealt{lugaro04}, 2012\nocite{lugaro12}).
To get deeper insights into the nature of the candidate AGB stars possibly responsible
for the observed abundance trends in M22 we compare our results with the
recent set of models by \cite{lugaro12}. They presented AGB models for
masses 0.9$-$6.0~M$_{\odot}$~and metallicity three times lower than that of
the cluster under consideration (i.e., [Fe/H]=$-$2.3 dex).

From our data we infer that there is an increase of $\Delta$[F/H]$^{s-rich}_{s-poor}$=+0.40$\pm$0.15 dex in
the fluorine content between the two groups. This
estimate was done by taking into account the F content of the O-rich
stars only, because they do not show any depletion due to the HBB. We
averaged the F abundances in stars IV-97 and III-14, both belonging to the
$s$-poor group, finding a mean value $<$[F/H]$>$=$-$2.63$\pm$0.01 dex. Since we have
only one O-rich star in the $s$-rich group, III-52, we chose it as
representative of the F abundance for the group, that is [F/H]=$-$2.23$\pm$0.15 dex.
The F increase of +0.40 dex is accompanied by a corresponding
enhancement in La of $\Delta$[La/H]$^{s-rich}_{s-poor}$=0.56 $\pm$0.18 dex, since the values are
$<$[La/H]$>$=$-$1.92$\pm$0.15 dex and [La/H]=$-$1.36$\pm$0.10 dex, respectively, for the two
groups.

Comparing these values with the model predictions by \cite{lugaro12} we
found that AGB stars with masses of $\approx$ 4$-$5 M$_{\odot}$\ can well
reproduce the observed pattern. Lower-mass AGB models do not fit our
observational requirements because they over-produce fluorine. This is
true even if we consider that these AGB model predictions are roughly one
dex too high to match the observation of carbon-enhanced metal-poor (CEMP) stars by
\cite{lucatello11}. Production of F in the $\approx$ 4$-$5 M$_{\odot}$\
mass range depends on the delicate balance between the operation of the
TDU and of HBB. These stars suffer HBB and destroy F during the early
phases of their AGB evolution, however, toward the end of the evolution, as
the mass of the convective envelope decreases, HBB ceases while the TDU
is still active resulting in mild F enhancements at the stellar surface.
This explains why while the most prolific AGB stars in terms of fluorine
production have initial masses around 2 M$_{\odot}$, F production still occurs
at slightly higher masses. On the other hand, more massive AGB models
($>$5 M$_{\odot}$)~experience hotter HBB and thus more efficient F
destruction as well as higher temperatures in the thermal pulses activating
also $^{19}$F($\alpha$,p)$^{22}$Ne reactions. This, combined with
fewer final TDU episodes when HBB has ceased means that they do not
replenish F at the stellar surface.

Interestingly, the same conclusion is drawn by \cite{roederer11} by
exploring the heavy-element ratios, [$hs$/$ls$] and [Pb/$hs$].
Comparing their abundances with models by \cite{roederer10}, these
authors deduced that the low mass AGBs ($\leq$3M$_{\odot}$) cannot account for
the observed trend. More importantly, they concluded that a match to the
$s$-process element abundances is provided by the 5 M$_{\odot}$\ AGB model, and
this is confirmed by checking the models of \cite{lugaro12}.
Furthermore, and very interestingly, these models predict Na and C production in agreement
with the observations.

{\em Our result provides a further, independent confirmation to this
previous hint: indications from both light (fluorine here) and heavy
elements converge toward the AGB stars of the same mass
as the best candidate polluters, indicating
that the age difference between the two sub-groups in M22 cannot be
larger than a few hundreds Myr.} It should be mentioned
that by analyzing the double sub-giant branch (SGB) of this
cluster, \cite{marino12} concluded that the age spread can be {\em at
most} $\sim$300 Myr.

The fact that three independent studies, involving different and
complementary techniques/approaches, produce the same result is
encouraging. However, a comprehensive understanding of the whole picture
is still missing. As also stressed by \cite{roederer11},
if relatively massive AGBs (4-5 M$_{\odot}$)~produced the $s$-process
elements, and if these stars are also responsible for the observed
$p$-capture element anticorrelations, then it is not explained that
an $s$-process enrichment is present in M22 but is not associated
with O and Na abundance anomalies, nor it is seen in any other GCs
where the Na-O anticorrelations are observed.

On the other hand, it is clear from the observed F
abundance trends that we must select slightly more massive AGBs, i.e., $\geq$ 6M$_{\odot}$\, if we wish to explain the
light element variations in GCs, since
this is the AGB mass range where F and O can be destroyed by HBB
resulting in the observed F-O correlation. We might tentatively suggest
that perhaps the production of the $s$-process elements in AGB models
with initial mass $\geq$ 6M$_{\odot}$~is less efficient than currently predicted.
This could be the result of a stronger mass-loss rate or a less
efficient TDU in this mass range than those employed in the models by
\cite{lugaro12}. This possibility is well within model uncertainties
and need to be investigated; in a forthcoming paper
(D'Orazi et al., in preparation) we will attack these issues, presenting new observations and AGB models
and discussing their strength/weakness in reproducing the observed abundance trends in GCs.

Alternatively, we may conclude that
massive AGBs are not the inter-cluster polluters; however
several lines of evidence point to those stars, such as the need for a Li
production between first and second generation stars (as in the case of
M4, \citealt{dm10}; \citealt{mucciarelli11}; \citealt{monaco12}). Further efforts, from both observational and
theoretical perspectives, are needed; in particular the lack of a
complete set of models for AGB stars and massive rotating stars
with different mass and metallicity, following the whole nucleosynthetic path from Li to Pb, still hampers a
quantitatively robust comparison between theory and observations.

\section{Summary and concluding remarks}\label{sec:summary}

We presented fluorine abundances for a sample of six RGB
stars belonging to the metal-poor globular cluster M22. The sample was
selected to include $s$-process-rich and $s$-process-poor
stars, as defined in MSK11. In addition, within each of these
cluster sub-components, we targeted both O-rich (Na-poor) and O-poor
(Na-rich) stars.

We gathered evidence of the presence of a F-O
correlation and of an F-Na anticorrelation.
Such chemical pattern, notably revealed in each cluster
sub-group, is in agreement with F destruction during the hot H burning:
fluorine follows the same trend defined by O, Na, C, N, Mg, Al, as predicted by the
multiple population scenario.

Most interestingly, we found that
the $s$-process- (metal-) rich component is also characterized by a
larger F content than the $s$-process (metal-) poor component. The
comparison between our observations and AGB models points to stars with
masses around 4-5 M$_{\odot}$~as responsible for such a trend, corroborating
previous hints by \cite{roederer11} and \cite{marino12}, and confirming
that the age spread across the two different stellar generations in M22
cannot be larger than a few hundreds Myr.

\acknowledgments
The publication made extensive use of the NASA ADS and SIMBAD databases.
We kindly acknowledge B. Plez for having provided his unpublished line lists.
This work was partially funded by PRIN INAF 2009 (``Formation and Early
Evolution of Massive Star Clusers'') and PRIN INAF 2011 (``Multiple populations in globular clusters: their role in the Galaxy assembly'').
V.D. thanks Kais Hamza for very useful discussions.