ABSTRACTJPSeveral important properties of rotation-powered millisecond pulsars (MSPs), such as their mass-radius ratio, equation of state and magnetic field topology, can be inferred from precise observations and modelling of their X-ray light curves. In the present study, we model the thermal X-ray signals originated in MSPs, all the way from numerically solving the surrounding magnetospheres up to the ray tracing of the emitted photons and the final computation of their light curves and spectra. The magnetosphere is solved by performing general relativistic force-free simulations of a rotating neutron star (NS) endowed with a simple centred dipolar magnetic field, for many different stellar compactness and pulsar misalignments. From these solutions, we derive an emissivity map over the surface of the star, based on the electric currents in the magnetosphere. In particular, the emission regions (ERs) are determined in this model by spacelike four-currents that reach the NS. We show that this assumption, together with the inclusion of the gravitational curvature on the force-free simulations, lead to non-standard ERs facing the closed-zone of the pulsar, in addition to other ERs within the polar caps. The combined X-ray signals from these two kinds of ERs (both antipodal) allow to approximate the non-trivial interpulses found in several MSPs light curves. Our modelled X-ray signals are compared against very accurate NICER observations of four target pulsars: PSR J043-4715, PSR J1231-1411, PSR J2124-3358, and PSR J0030 + 0451; achieving very good simultaneous fits for their light curves and spectral distributions.