In this paper, comprehensive pulsar-based Guidance, Navigation and Control (GNC) system is designed and applied to satellites formation flying. The complete autonomy of the X-ray pulsar navigation technology provides both absolute and relative positioning information for spacecraft in or even beyond the solar system, and provides an interesting alternative solution to classical navigation techniques such as Global Navigation Satellite System (GNSS). The navigation measurements are studied using a relative time of arrival (RTOA) estimating method, which takes advantage of the noise-filtered pulsar profiles in the cross-correlation estimator to achieve robustness. An adaptive Kalman filter is exploited to minimise the effect of processing noises caused by the primary satellite. A linear Quadratic Gaussian regulator is then applied to control the formation. Two different configurations are studied in high Earth orbits (HEOs). Simulation results show that the designed GNC algorithms can fulfill the required accuracy (10% of the baseline) of the Magnetospheric Multiscale Mission, with the position root-mean-square error of similar to 1.6 km.. The excellent robustness of the designed GNC system on the positioning errors of the primary satellite suggests a role for X-ray pulsar-based formation flying technique in HEOs and the solar system.