How can a find the symmetry group of a Lagrangian? I apologize in advance for my basic English but I would like to know if there is a rule, a book or in general some way in order to determine the internal symmetries, the gauge symmetries or all the symmetries of a Lagrangian density. I know there is the Noether's theorem that will give me the conserved current from the given symmetry but my problem is a priori: how can I find the largest internal symmetry of a given Lagrangian? 
Since I want to be clearer I can give to you some example. Let us suppose I have two real fields $\phi_1$ and $\phi_2$ such that the Lagrangian is
$\mathcal{L}(\phi_1,\phi_2)=\frac{1}{2}\partial_\mu \phi_1\partial^\mu \phi_1+\frac{1}{2}\partial_\mu\phi_2\partial^\mu \phi_2-\frac{3}{4}M^2(\phi_1^2+\phi_2^2)-\frac{1}{2}M^2\phi_1\phi_2$
If I ask you what the symmetry group of this Lagrangian is, what will you say to me?
At the same time, I can diagonalize the mass matrix and I can define new fields
$\phi_{M^2}=\frac{\phi_2-\phi_1}{\sqrt{2}}\qquad \phi_{2M^2}=\frac{\phi_1+\phi_2}{\sqrt{2}}$
where, if I did the computation right, will give to me a Lagrangian density as
$\mathcal{L}_2(\phi_{M^2},\phi_{2M^2})=\frac{1}{2}\partial_\mu \phi_{M^2}\partial^\mu \phi_{M^2}+\frac{1}{2}\partial_\mu\phi_{2M^2}\partial^\mu \phi_{2M^2}-\frac{1}{2}[(2M^2)\phi_{2M^2}^2+M^2\phi_{2M^2}^2)]$
which is a Lagrangian of two Klein-Gordon fields. This definition is a rotation of $45$° of the fields, does this means that the Lagrangian is invariant under SO(2)? If the mass of the fields where the same, the second lagrangian can be seen as the Klein Gordon Lagrangian density of a doublet
$\Phi=\left(
\begin{array}{c}
\phi_1\\
\phi_2
\end{array}
\right)
$
?
Is then this last thing an example of a Lagrangian invariant under SO(2)?
Last thing: If I have a generic Lagrangian:
$\mathcal{L}_3(\phi_1,\dots,\phi_n)=\sum_{i,j=1}^n\left(\frac{1}{2}K_{ij}\partial_\mu\phi_i\partial^\mu\phi^j-M_{ij}\phi^i\phi^j\right)-V(\phi_i)$
where $K_{ij}$ and $M_{ij}$ are real symmetric constant matrices; $K_{ij}$ is non-singular positive definite matrix. What is the largest internal symmetry of the kinetic term? What is the form of $M_{ij}$ such that the mass term is invariant under this group $G$?
I hope that three examples will help you to understand what I mean and what I would like to be answered. If there is a reference, a textbook or anything that can help me to understand how to deal with these kind of questions I really appreciate them. 
 A: So you really have two closely related questions here.  Let's start with your easier one.  You wrote the Lagrangian
$\mathcal{L}(\phi_1,\phi_2)=\frac{1}{2}\partial_\mu \phi_1\partial^\mu \phi_1+\frac{1}{2}\partial_\mu\phi_2\partial^\mu \phi_2-\frac{3}{4}M^2(\phi_1^2+\phi_2^2)-\frac{1}{2}M^2\phi_1\phi_2$
And asked for the symmetry group.  I can see that the first three terms have an $SO(2)$ symmetry, but this symmetry is broken by the $\phi_1 \phi_2$ cross-term.  You have done a change of variables,
$\phi_{M^2}=\frac{\phi_2-\phi_1}{\sqrt{2}}, \qquad \phi_{2M^2}=\frac{\phi_1+\phi_2}{\sqrt{2}}$,
which makes this fact more obvious.  The Lagrangian becomes
$\mathcal{L}_2(\phi_{M^2},\phi_{2M^2})=\frac{1}{2}\partial_\mu \phi_{M^2}\partial^\mu \phi_{M^2}+\frac{1}{2}\partial_\mu\phi_{2M^2}\partial^\mu \phi_{2M^2}-\frac{1}{2}[(2M^2)\phi_{2M^2}^2+M^2\phi_{2M^2}^2)]$
and you are correct that, if the diagonalized fields had equal masses, then there would be an $SO(2)$ symmetry.  If you work backwards through your change of variables, you should find that the diagonal fields have equal masses if and only if the $\phi_1 \phi_2$ term had vanished in the first place.
Now let's look at your more general Lagrangian:
$\mathcal{L}_3(\phi_1,\dots,\phi_n)=\sum_{i,j=1}^n\left(\frac{1}{2}K_{ij}\partial_\mu\phi_i\partial^\mu\phi^j-M_{ij}\phi^i\phi^j\right)-V(\phi_i)$
Since $K_{ij}$ is real and symmetric, it can be diagonalized by an orthogonal matrix.  Finally, by rescaling the fields, you can set $K_{ij} = \delta_{ij}$ to give yourself the standard kinetic term.  This term on its own is $SO(n)$ invariant.
The question that remains is, after doing this field redefinition, what does $M_{ij}$ look like?  Since the kinetic term is now $SO(n)$ invariant, and $M_{ij}$ is symmetric, we can also diagonalize $M_{ij}$.  However, we cannot do any more rescalings, so the new Lagrangian, written in terms of mass eigenfields, will generically have fields with different masses.  If all the masses are different, then all the $SO(n)$ symmetry is broken.  If any masses are the same, then their multiplicities can leave you a residual $SO(k) \times SO(\ell) \times \ldots$, etc.
Edited to add:  You also need to look at the symmetries of your potential $V(\phi)$, which may further break symmetry!
