The use of fragments to biophysically characterize a protein binding pocket and determine the strengths of certain interactions is a computationally and experimentally commonly applied approach. Almost all drug like molecules contain at least one aromatic moiety forming stacking interactions in the binding pocket. In computational drug design, the strength of stacking and the resulting optimization of the aromatic core or moiety is usually calculated using high level quantum mechanical approaches. However, as these calculations are performed in a vacuum, solvation properties are neglected. We close this gap by using Grid Inhomogeneous Solvation Theory (GIST) to describe the properties of individual heteroaromatics and complexes and thereby estimate the desolvation penalty. In our study, we investigated the solvation free energies of heteroaromatics frequently occurring in drug design projects in complex with truncated side chains of phenylalanine, tyrosine, and tryptophan. Furthermore, we investigated the properties of drug-fragments crystallized in a fragment-based lead optimization approach investigating PDE-10-A. We do not only find good correlation for the estimated desolvation penalty and the experimental binding free energy, but our calculations also allow us to predict prominent interaction sites. We highlight the importance of including the desolvation penalty of the respective heteroaromatics in stacked complexes to explain the gain or loss in affinity of potential lead compounds.