Data analysis on particle physics: what am I formally doing in this specific case? I have performed an analysis of the exclusion(discovery) perspectives for a BSM particle at the LHC. However, my focus of study (and understanding) was until little time ago the purely theoretical and, having 'randomly' started doing phenomenology, I have a very poor understanding of the statistics underlying the process I have (mechanically) performed. Having interest on more complicated and enhanced future analysis, I need to grasp these concepts I currently ignore. I will describe what I have done. 
The considered signal is one fixed end-state $X$ (obeying some here irrelevant kinematic cuts), the complete process is $pp \to X$, ignoring jets, and contains 2 diagrams (besides interference with background). I implemented the new physics model and generated pseudo-experiments for hundreds (let's say N) of points in a bidimensional parameter space -- let's say it is generated by $(M,\theta)$. I also generated the same number of events of (irreducible) background, for the null hypothesis which can also produce the exact same process. Chosen as discriminant observable was the invariant mass $m_{||}$ of certain $n$-body product of the diagrams. I generated histograms of the number of events with relation to $m_{||}$ for each of the N points within the parameter space, and for the background.
To reject the signal containing hypothesis at each point, "I have chosen a 95% confidence level" (the quotes are because I don't know how correctly am I using these concepts, but I will not keep writing them). The confidence level in each point was calculated starting with the following expression
$$
\tag{1}
\text{Confidence level} =1-\prod_{k}\frac{P(n<n_{0k} \;|\; \mu_{sk}+\mu_{bk})}{P(n<n_{0k}\; |\;\mu_{bk})},
$$
where $P(n \;|\; \lambda)$ is a Poisson distribution for variable $n$ and rate $\lambda$, which, in my case (of the distribution normalized by the total number of generated events), is simply the expected number of events of the respective kind. $k$ labels the bins of $m_{||k}$, and $s_k$($b_k$) is the number of signal(background) events in the respective bin. $\mu_{sk}+\mu_{bk}$ is the expected number of events in the alternate hypothesis (signal + background) and $\mu_{bk}$ the number expected in the null (only background) hypothesis. $n$ is the number of events and $n_0$ defines a "confidence interval" to be chosen. 
The final expression is obtained by "estimating" the signal(background) rate by the number of events predicted by the signal(background) pseudo-experiment and by fixing the "confidence interval" by the expected number of background events. That is
$$
\mu_{sk} \to s_k, \qquad \mu_{bk} \to b_k, \qquad n_{0k} \to b_k.
$$ 
The result plugging this into Eq. (1) is
$$
\tag{1}
\text{Confidence level}=1-\prod_{k}\left[ \frac{e^{-s_k}\sum\limits_{n=0}^{b_k} \frac{(s_k+b_k)^n}{n!}}{\sum\limits_{n=0}^{b_k} \frac{b_k^n}{n!}}\right],
$$
Having obtained a Confidence level for each point in the space of parameters, I extrapolated the discrete result to the continuous plane and plotted contours of fixed $\text{CL}=95\%$.
What I did was validated by people with expertise, which makes me have a reasonable belief that (despite being obviously not sophisticated at all) it is, in principle, acceptable. I, however, do not truly understand it.
I am trying to learn this and acquire a formal enough understanding in order to perform more complex and professional analysis by these notes, and will ask questions with its definitions and notions in mind. The questions will likely be confusing and I apologize for that. Also, although I will separate them for readability, they are not actually 5 different and disjoint questions but are simply a synthesis of my lack of understanding and you may answer them as you like.
(1) The likelihood ratio is supposed to be the more powerful test statistic in distinguishing between 2 hypothesis and is, apparently, related to every search of this kind. What is the likelihood ratio role in my method? Where does it appear and why does it do so? 
(2) The obvious guess would be that it is simply my ratio of Poisson probabilities of signal and background, but then what am I doing? The likelihood ratio contour should be fixed by some $k_\alpha$ which is related to the confidence level but is not it, so why my fixing of something related to $k_\alpha$ to be $95\%$ worked? If we don't know how $k_\alpha$ relates to the CL then a CL of $95\%$ could imply a $k_\alpha$ very different than what I have chosen.
(3) What I have done is the $\text{CL}$ or $\text{CL}_s$ method (or neither of them)? Why? If it is the $\text{CL}$ how would one turn it into the $\text{CL}_s$?
(4) What are each variable roles here? Clearly the unknown new physics model parameters $M,\theta$ are parameters of interest, then my PDF would be something like $P(n_k;M,\theta)$ ($k$ is again related to the invariant mass). Is this correct?
(5) The notes I have linked express that the test statistic agreed upon by most people is the likelihood ratio with likelihoods assumed as showed in Eq. 86. How do I relate that expressions to the ones in my method? My interpretation is that they are the Poisson probabilities for the total number of events smeared on the $m_{||}$ space by the probabilities $\mathcal{S}$ and $\mathcal{B}$ which are the (universal) extra purely physical information. Is this correct? But, if it is, how does it appear mathematically in my expressions? 
I think (and hope) this is an understandable and contained enough for a single question, and appreciate any help.  
 A: I don't quite understand your eqs. 1 and 2 so I will focus on answering questions 4 and 5. Maybe someone else can provide more specific answers for your 1--3. But I'll point you towards a great reference for various test statistics here, which you might find helpful.
From what I can glean, you have (separately) generated samples of background-only and signal-only events, the latter for a collection of points in its 2-dimensional parameter space. You have then performed an analysis on a binned histogram of some final state invariant mass, using information obtained from similarly-binned histograms constructed from the signal and background event samples.
I place an emphasis on "binned" because it affects the likelihood function you define: take a background-only assumption, for example. Recall that the likelihood function is defined as the probability density function (PDF) generating the data, except with random variables interpreted as parameters and vice versa. For a binned histogram, the number of events falling within a given bin $k$ can be modelled as a Poisson distribution, and so:
\begin{equation}
L^\mathrm{(bckg)}_k(b_k;n_k) = P(n_k;b_k) = \frac{(b_k)^{n_k} e^{-b_k}}{n_k!}\,,
\end{equation}
with $n_k$ the number of events in the bin, and $b_k$ the expected number of background events. Treating bins as independent counting experiments, the full likelihood can then be written as the product over all bin likelihoods:
\begin{equation}
L^\mathrm{(bckg)}(b,n) = \prod_k L^\mathrm{(bckg)}_k(b_k;n_k)\,.
\end{equation}
For the signal+background assumption, simply write $b_k\rightarrow b_k+s_k$. To address the first part of your q.$~$5, note that eq.$~$(86) in the text you linked defines likelihood functions for unbinned datasets (so-called "extended likelihoods"), so they do not apply to your analysis. (You can check that they do not allude to any binning.)
Since $s_k$ and $b_k$ are supposed to represent the expected number of signal or background events in bin $k$, they must clearly be related to the physics models of the signal and background, respectively. If we denote by $\mathcal{S}(m_{||};M,\theta)$ and $\mathcal B(m_{||};\Theta)$ the PDFs describing the signal-only and background-only invariant mass distributions, then for bin $k$ one has:
\begin{equation}
s_k=N_\mathrm{total}\times\frac{\int_{m_{k_\mathrm{min}}}^{m_{k_\mathrm{max}}} \mathcal S(m_{||};M,\theta)\,dm_{||}}{\int_{0}^{\infty} \mathcal P_\mathrm{total}(m_{||})\,dm_{||}}\,,
\end{equation}
where $m_{k_\mathrm{min}}$ and $m_{k_\mathrm{max}}$ are the invariant mass bounds for bin $k$, and $N_\mathrm{total}$ is the total number of observed events. $\mathcal P_\mathrm{total}(m_{||})$ is the total PDF; if assuming signal+background for example, then $\mathcal P_\mathrm{total}=\mathcal S+\mathcal B$. (Of course, if the total PDF is properly normalised then the denominator is 1.) This is simply the fraction of (signal) events you expect to fall within the $k^\mathrm{th}$ bin, given the theory PDF. $b_k$ can be written analogously.
However, you bypass the need to include such calculations for your likelihoods because you instead construct signal-only and background-only histograms from generated events, and directly obtain $s_k$ and $b_k$ from there. (This addresses your q.$~$4 and the second part of your q.$~$5.)
