In global mean-field dynamo models, the process of magnetic flux emergence through the stellar convection zone is crucial. In the Sun, it is the strong toroidal structures built at the base of the convection zone which are assumed to be unstable to a buoyancy instability and rise through the convection zone to produce sunspots. In some 3D models of rapidly-rotating stars, strong toroidal structures can become buoyant but rarely rise all the way to the top of the computational domain and those models consequently do not produce spots. It is thus still an open question if we can really rely on « spotless » dynamo models to reproduce what could happen in cool stars. The particular step of flux emergence being potentially important for the whole dynamo mechanism, detailed numerical simulations of such a process are thus needed. We will present such simulations and introduce a new combined approach reconciling global 3D models and production of starspot-like features at the stellar surface.