The solutions are the ones you listed.
The solutions all have shape $y=(x-c)^2$ or $y=0$. Thus if $b<0$, then none of the solutions curves pass through $(a,b)$. So for all pairs $(a,b)$ such that $b<0$, there cannot be a solution satisfying $y(a)=b$. We do not know (yet) whether these are all the pairs $(a,b)$ for which there is no solution, but soon we will.
For any $a$, if $b=0$ there are exactly two solutions satisfying $y(a)=b$, the singular solution and the solution $y=(x-a)^2$.
Finally, we look at pairs $(a,b)$ with $b$ positive. We look for values of $c$ such that $y(a)=b$.
The solution $y=(x-c)^2$ passes through $(a,b)$ if and only if $(a-c)^2=b$. This equation has exactly two solutions, $c=a\pm\sqrt{b}$.
Conclusion: (a) The pairs $(a,b)$ for which there is no solution satisfying $y(a)=b$ are all $(a,b)$ with $b<0$. (b) There are no pairs $(a,b)$ for which there are infinitely many solutions with initial condition $y(a)=b$. (c) For all remaining pairs $(a,b)$, that is, all pairs with $b \ge 0$, there are finitely many solutions, indeed exactly two solutions that satisfy $y(a)=b$.
The geometry: The conclusion can also be reached geometrically, by visualizing the family of parabolas. All of your parabolas are obtained by sliding the standard parabola $y=x^2$ along the $x$-axis. For any $(a,b)$ with $b \gt 0$, there are exactly two such parabolas that pass through (a,b): one whose "left" half goes through $(a,b)$, and one whose "right" half goes through $(a,b)$.
Note: One could interpret the word "finite" to include the possibility of $0$ solutions: $0$ is certainly finite! That is obviously not the intended interpretation here. But if we interpret "finite" as including $0$, the answer to (c) is all pairs $(a,b)$.
Your reasonning is correct. Since $g$ is Lipschitz and the first equation of the system involves only $x$, there is a unique solution $x(t)$ such that $x(t_0)=x_0$.
The second equation becomes
$$
y'=f(x(t))\,y,\quad y(t_0)=y_0.
$$
It is a linear equation and has a unique solution, given by
$$
y(t)=y_0e\,^{\int_{t_0}^t f(x(s))ds}.
$$
Best Answer
This is called the method of the integrating factor. In general for $x'=A(t)x+b(t)$, if you have a solution for the matrix equation $F'(t)=A(t)F(t)$, $F(0)=I_n$, also called the fundamental solution, you can, a la variation of constants, set $x(t)=F(t)u(t)$ and either differentiate $u(t)=F(t)^{-1}x(t)$ to \begin{align} u'(t)&=-F(t)^{-1}F'(t)F(t)^{-1}x(t)+F(t)^{-1}x'(t) \\ &=-F(t)^{-1}A(t)x(t)+F(t)^{-1}(A(t)x+b(t)) \\ &=F(t)^{-1}b(t) \end{align} or in the variation of constants fashion compute \begin{align} x'&=F'u+Fu'=AFu+Fu' \\ &=Ax+b \\\\ \implies Fu'&=b \end{align}
For a linear system with constant matrix, $x'(t)=Ax(t)+b(t)$, the fundamental matrix is the matrix expontential $F(t)=e^{tA}$ with derivative $F'(t)=AF(t)=F(t)A$.