Assume $f(x_1)=f(x_2)$. Then $x_1^2+1=f(f(x_1))=f(f(x_2))=x_2^2+1$ and so $x_2=\pm x_1$.
Then from $f(f(-x))=f(f(x))$ we conclude that $f(-x)=\pm f(x)$ for all $x$.
Define $g\colon [0,\infty)\to [0,\infty)$ by $g(x)=|f(x)|$.
For $x$ with $g(x)=f(x)$ we have $$g(g(x))=|f(f(x))|=|x^2+1|=x^2+1.$$
And for $x$ with $g(x)=-f(x)$ we have $$g(g(x))=|f(-f(x))|=|\pm f(f(x))|=|f(f(x))|=x^2+1$$
as well. Hence $g(g(x))=x^2+1$ for all $x\ge 0$.
Given a function $g$ with this property, we can easily construct a suitable $f$: Just let $f(x)=\begin{cases}g(x)&x\ge0\\g(-x)&x<0\end{cases}$. If $g$ is additionally continuous then so is $f$.
We can find a lot of continuous $g$:
Pick $a_1\in(0,1)$, let $a_0=0$ and recursively $a_n=a_{n-2}^2+1$ for $n\ge 2$. Then the sequence $(a_n)$ is strictly increasing towards $\infty$.
Now let $g_1\colon[0,a_1]\to[a_1,1]$ be an arbitrary increasing bijection. We can work our way to infinity from there: Assume we have an increasing bijection $g_n\colon [0,a_n]\to [a_1,a_{n+1}]$ such that $g(g(x))=x^2+1$ for $0\le x\le a_{n-1}$. Then define $g_{n+1}\colon [0,a_{n+1}]\to[a_1,a_{n+1}]$ by $$g_{n+1}(x)=\begin{cases}g_n(x)&0\le x\le a_n\\g_n^{-1}(x)^2+1&a_1< x\le a_{n+1}\end{cases}$$
(noting that the two cases agree for $a_1\le x\le a_n$)
and finally $g(x)=g_n(x)$ for any $n$ with $a_n\ge x$. One verifies that this is well-defined and indeed $g(g(x))=x^2+1$ for all $x$.
Note that the sequence of $a_n$ goes like this: $$0, a_1, 1, a_1^2+1, 2, (a_1^2+1)^2+1, 5, ((a_1^2+1)^2+1)^2+1, 26, \ldots$$ and that it is quite clear that we can make our choices (even under the restriction of requiring continuous $f$) at least such that $f(x)$ is an arbitrary number between $6$ and $26$.