Here is a reduction to a bivariate polynomial equation modulo $p_n\#^2$.

First, it can be easily seen that when reduced modulo $p_m$, the $m$-th congruence takes form: $$\frac{(p_n\#)^2}{p_m^2}x_my_m\equiv b_m\pmod{p_m}.\qquad(1)$$ So, I assume that the known $p_m-1$ candidate values for $(x_m,y_m)$ satifsfy this congruence. Let $C_m$ denote the set of candidate values for $(x_m,y_m)$.

Set $P:=p_n\#$, $X:=P\sum_{i=1}^n \frac{x_i}{p_i}$, and $Y:=P\sum_{i=1}^n \frac{y_i}{p_i}$. Then the $m$-th congruence takes form: $$P(\frac{y_m}{p_m} X + \frac{x_m}{p_m} Y) - P^2(\frac{x_my_m}{p_m^2} + n_2\frac{x_m}{p_m} + n_1\frac{y_m}{p_m})\equiv b_m\pmod{p_m^2}.\qquad(2)$$

From the candidate values for $(x_m,y_m)$, we obtain the following congruence holds: $$f_m(X,Y)\equiv 0\pmod{p_m^2},\qquad(3)$$ where $$f_m(X,Y) := \prod_{(x_m,y_m)\in C_m} \left[ P(\frac{y_m}{p_m} X + \frac{x_m}{p_m} Y) - P^2(\frac{x_my_m}{p_m^2} + n_2\frac{x_m}{p_m} + n_1\frac{y_m}{p_m}) - b_m\right].$$

Clearly, $\deg f_m \leq p_m-1$. Using CRT, we combine congruences $(3)$ into a single one: $$F(X,Y) \equiv 0\pmod{P^2},$$ where $\deg F\leq \max_i (p_n-1)$. If we solve this equation for $(X,Y)$, we can substitute solution to $(2)$ and solve the resulting equation for $(x_m,y_m)$ (finding a suitable from $C_m$), etc.

Here is a reduction to a bivariate polynomial equation modulo $p_n\#^2$.

Let $C_m$ denote the set of candidate values for $(x_m,y_m)$.

Set $P:=p_n\#$, $X:=P\sum_{i=1}^n \frac{x_i}{p_i}$, and $Y:=P\sum_{i=1}^n \frac{y_i}{p_i}$. Then the $m$-th congruence takes form: $$P(\frac{y_m}{p_m} X + \frac{x_m}{p_m} Y) - P^2(\frac{x_my_m}{p_m^2} + n_2\frac{x_m}{p_m} + n_1\frac{y_m}{p_m})\equiv b_m\pmod{p_m^2}.$$

From the candidate values for $(x_m,y_m)$, we obtain the following congruence holds: $$f_m(X,Y)\equiv 0\pmod{p_m^2},\qquad(\star)$$ where $$f_m(X,Y) := \prod_{(x_m,y_m)\in C_m} \left[ P(\frac{y_m}{p_m} X + \frac{x_m}{p_m} Y) - P^2(\frac{x_my_m}{p_m^2} + n_2\frac{x_m}{p_m} + n_1\frac{y_m}{p_m}) - b_m\right].$$

Clearly, $\deg f_m \leq p_m-1$. Then, we combine congruences $(\star)$ into a single one: $$F(X,Y) \equiv 0\pmod{P^2},$$ where $$F(X,Y) := \sum_{m=1}^n \frac{P^2}{p_m^2} f_m(X,Y).$$

If we solve the resulting congruence for $(X,Y)$, then from $X\equiv \frac{P}{p_m}x_m\pmod{p_m}$ and $Y\equiv \frac{P}{p_m}y_m\pmod{p_m}$ we can find a suitable pair $(x_m,y_m)$ from $C_m$ for each $m$, and then verify that they altogether deliver a solution (as they may be extraneous ones).