You have a few different questions here, so let's try to go through them one by one.
- When we make the chiral symmetry local, have we introduced a gauge symmetry, or some analogue of a gauge symmetry?
When you make the chiral symmetry local you introduce a gauge symmetry. The terms "gauge symmetry" and "local symmetry" are two different ways of saying the same thing.
However it's a little tricky, as we will see. The gauge symmetry is introduced by coupling to external fields $r_\mu$ and $\ell_\mu$ (in your reference by Pich, see the text before and after equation 3.10). This means that we do not think of $r_\mu$ and $\ell_\mu$ as dynamical fields--we do not quantize them and there are no particles associated with these fields. We introduce these external, classical fields, and the associated local symmetry, as a mathematical trick in order to make certain calculations easier.
As a secondary benefit, by quantizing some of the components of $r_\mu$ and $\ell_\mu$ we can include couplings to electroweak bosons.
- Where are the gauge bosons associated with $r_\mu$ and $\ell_\mu$?
As $r_\mu$ and $\ell_\mu$ are external, classical fields, we do not quantize them and there are no particles associated with these fields.
Thus, chiral perturbation theory does not predict the existence of more gauge bosons beyond the $W,Z$, photon, and gluons.
- What is the point of introducing these external, classical fields?
Pich after equation 3.10 gives several reasons why we introduce these fields, which I will try to elaborate on a bit here.
First, the external fields allow us to compute Green's functions associated with certain currents. In other words, it will simplify computing correlation functions like
\begin{equation}
\langle 0 | J^\mu | 0 \rangle, \langle 0 | J^\mu | \pi \rangle
\end{equation}
Second, we can embed the couplings of the electroweak gauge bosons into these external fields.
While the second benefit might make more immediate sense, it is actually beneficial to forget about this at first and focus on the first benefit, and then later at the end see how you can incorporate the gauge bosons.
The basic fact underlying the first benefit (computing correlation functions of the chiral currents) is that gauge fields couple to conserved currents. In other words, if we ignore the gauge field kinetic terms (which aren't really relevant if we think of these gauge fields as classical, external fields) then the gauge fields will couple to the matter fields in the form
\begin{equation}
\mathcal{L}_{gauge} = r_\mu J_r^\mu + \ell_\mu J_\ell^\mu
\end{equation}
where by standard arguments, gauge invariance tells us that $J_r^\mu$ and $J_\ell^\mu$ must be conserved currents. Indeed, these correspond to the chiral currents.
So the idea is to consider the generating functional $Z$ as a functional of the external gauge fields $r$ and $\ell$ as well as some extra external scalar fields $s,p$ that are useful for describing the masses. For the purproses of this answer we will ignore $s,p$, so the generaing functional is a functional
\begin{equation}
Z[r,\ell] = \int DU \exp\left(i \int d^4 x (D_\mu U)^\dagger D^\mu U + \cdots \right)
\end{equation}
where $D_\mu$ is the gauge covariant derivative including $r,\ell$, $U$ is the matrix field carying the goldstone bosons, and the $\cdots$ refers to terms that are higher order in derivatives (as well as mass terms which we are neglecting in this answer).
Notice that we are not integrating over $r$ and $\ell$, this is what makes them external fields. (I'll come back to how to incorporate the electroweak gauge bosons at the end, but including them won't affect this step).
Then to compute the vacuum expectation value for the current you simply follow the normal rules of the generating functional
\begin{equation}
(-i)\frac{\delta Z}{\delta r_\mu} = \langle 0 | J_r^\mu | 0 \rangle
\end{equation}
Strictly speaking this correlation function depends on the values of the external fields. After you differentiate with respect to $r_\mu$ and $\ell_\mu$, then you set the fields to their appropriate values. I'll discuss this in more detail below, but the vevs are given (perturabatively) in Pich 3.11, and I'll just highlight that this is the point where we put in information about the electroweak gauge bosons.
Before discussing the electroweak gauge bosons more, I'll just point out that we can also use the gauge fields to discover the appropriate current operators in the effective field theory. This is what Pich is doing in equation 3.20. We just define
\begin{equation}
J_r^\mu = \frac{\delta S_{eff}}{\delta r_\mu}
\end{equation}
and similarly for $J_\ell^\mu$. These currents can be computed order by order in perturbation theory provided we have written down the most general gauge invariant lagrangian. Again, we set the gauge fields to their values after taking the derivative.
After doing this procedure and discovering the current operator, we can compute expectation values as in Pich 3.21.
In other words, rather than running Noether's theorem on a complicated non-linear action to discover the currents, we write down a gauge invariant action with external gauge fields and use the variation of the gauge field to discover the current.
OK, now let's clean things up and discuss the electroweak gauge bosons in a bit more detail. Earlier our partition function only integrated over $U$. But the full partition function would also integrate over the gauge bosons. Let's write this in a suggestive way
\begin{equation}
Z_{full} = \int DW^{\pm}_\mu DZ_\mu DA_\mu DU e^{i S} = \int DW^{\pm}_\mu DZ_\mu DA_\mu Z_{eff}[\bar{r},\bar{\ell}]
\end{equation}
where
\begin{equation}
Z_{eff}[r,\ell] = \int DU \exp\left(i \int d^4 x (D_\mu U)^\dagger (D_\mu U) + \cdots \right)
\end{equation}
and where, following Pich 3.11,
\begin{eqnarray}
\bar{r}_\mu &=& e Q A_\mu + \cdots \\
\bar{\ell}_\mu &=& e Q A_\mu + \frac{e}{\sqrt{2} \sin \theta_W} (W_\mu^\dagger T^+ + h.c.) + \cdots
\end{eqnarray}
The above set of equations amounts to a procedure. First we compute $Z_{eff}[r,\ell]$ for arbitrary external fields $r,\ell$. Then after doing this we set the external fields to
"vevs" $\bar{r}$ and $\bar{\ell}$ that incorporate the electroweak gauge bosons and throw away the unphysical parts of $r$ and $\ell$. Then we do the integral over the gauge bosons. The full calculation will have quantized the goldstone bosons and the electroweak gauge bosons (and not the unphysical components of $r$ and $\ell$ we introduced really just to compute the current).
There is no doubt that this is a pretty convoluted set of procedures. There are a lot of very non-obvious tricks going on behind the scenes. The moral of the story, in the end, is that we are introducing a fictitious gauge symmetry that simplifies various calculations. To summarize, first it simplifies computing the chiral currents, and second we can embed the actual physical gauge bosons into the enlarged fictitious gauge symmetry we introduced.
