This problem is in BPP.
All finite fields of order $p^n$ are identical, and there is a generator element $g$ whose multiplicative order is $p^n-1$, so every non-zero field element is expressible as $g^r$ for some $r$.
Now, there is a subfield of order $p^j$ iff $j$ divides $n$, and this subfield has a generator element $h$ whose multiplicative order is $p^j-1$.
This means that we can take $h = g^{(p^n-1)/(p^j-1)}$.
What this means, in effect, is that once we have an algorithm for finding a generator of the finite field, everything else is straightforward.
The easy way to find a generator is to pick a field element at random, and test whether it is a generator. The expected number of elements you need to test is indeed polynomial. Doing this test efficiently might require knowing the prime factors of $n$. (But if you can write down all the coefficients of $f(x)$ without using sparse polynomial notation, then $n$ is small and this is easy.)
But you asked for a deterministic algorithm and not a randomized one.
A quick search shows that there is no known efficient deterministic algorithm for finding a generator. See this question on math.SE.
So the question is: can you find a deterministic algorithm for finding the subfields which doesn't actually compute a generator for the subfields? (Because that's probably as hard as finding a generator for the field itself, and so is an open problem.)
Much of this is explained in Wikipedia.