For terrestrial exoplanets with thin or no atmospheres, the surface contributes light to the reflected light signal of the planet. Measurement of the variety of disk-integrated brightnesses of bodies in the solar system and the variation with illumination and wavelength is essential for both planning imaging observations of directly imaged exoplanets and interpreting the eventual data sets. Here we measure the change in brightness of the Galilean satellites as a function of planetocentric longitude, illumination phase angle, and wavelength. The data span a range of wavelengths from 400 to 950nm and predominantly phase angles from 0{deg} to 25{deg}, with some constraining observations near 60{deg}-140{deg}. Despite the similarity in size and density between the moons, surface inhomogeneities result in significant changes in the disk-integrated reflectivity with planetocentric longitude and phase angle. We find that these changes are sufficient to determine the rotational periods of the moon. We also find that at low phase angles, the surface can produce reflectivity variations of 8%-36%, and the limited high phase angle observations suggest variations will have proportionally larger amplitudes at higher phase angles. Additionally, all of the Galilean satellites are darker than predicted by an idealized Lambertian model at the phases most likely to be observed by direct imaging missions. If Earth-sized exoplanets have surfaces similar to that of the Galilean moons, we find that future direct imaging missions will need to achieve precisions of less than 0.1ppb. Should the necessary precision be achieved, future exoplanet observations could exploit similar observation schemes to deduce surface variations, determine rotation periods, and potentially infer surface composition.