Perhaps I've gotten long-winded and just told you some things you already knew. To more specifically address your main questions... "Is the answer to the above question known?" Yes, it is pretty well-known and studied, the keywords being "multifractal formalism" and "conditional variational principle". "If yes, could you give a reference?" The references given by others to Kucherenko-Wolf and Barreira-Saussol-Schmeling will do -- my preference would be Barreira-Saussol 2001 since it is earlier -- but I would point out that the techniques used are pretty standard in multifractal work going even further back (Pesin-Weiss 1997, probably earlier as well). "How difficult is it?" In the generality you formulated I do not believe it is particularly difficult;dialfficult; one can consider the unique equilibrium states $\mu_t$ for $t\phi$ as $t$ ranges over $\mathbb{R}$ and then make some general arguments (if this is what you do then I apologize for telling you things you know). The same technique works more generally as long as you can say enough about the uniqueness of the equilibrium states, but this can be a challenging question in general.
Edit: As pointed out in the comments the references I gave above do not really answer the question directly. For the sake of notation let me write $K_\alpha$ for the level set above, and $T(\alpha) = \sup \{ h(\mu) \mid \int\phi\,d\mu=\alpha\}$. Then many of the main results in the references are of the form "$h(K_\alpha) = T(\alpha)$ under certain conditions". Which is not quite the same as describing how to compute the measure maximizing $T(\alpha)$, which is what you asked. The result in my Nonlinearity paper that I linked to describes a proof that $T(\alpha)$ is the Legendre transform of the pressure function $t\mapsto P(t\phi)$, which comes a little closer to the mark; the proof there contains a proof that there is some equilibrium state $\mu$ for some potential that achieves the maximum in $T(\alpha)$. But the way that it is written is opaque enough that I can't really claim it answers your question.
So let me extract the relevant bits in the case when $\phi$ is Hölder. Then the function $S(t) := P(t\phi)$ is differentiable in $t$, and $S'(t) = \int\phi\,d\mu_t$, where $\mu_t$ is the unique equilibrium state for $t\phi$. Using convexity properties of $S$ one can argue that $S'(t)$ takes all values in the interior of $\{ \int\phi\,d\nu \mid \nu \text{ is invariant}\}$. In particular, there is $t$ such that $S'(t)=\alpha$, and then one can show that $h(\mu_t) = T(\alpha)$, so that $\mu_t$ is the measure you want. As long as your system is such that you can write down $S'(t)$ explicitly, solve for $t$, and then describe $\mu_t$, then this gives you a pretty concrete algorithm for finding the conditionally maximizing measure.
My understanding is that the argument in the previous paragraph is implicit in a lot of work on multifractal formalism. But as illustrated by my earlier attempt at giving references, it's not always stated in a form that is readily useful for answering your specific question.