We consider the production of a pseudo-scalar particle A at the LHC, and present accurate theoretical predictions for its inclusive cross section in gluon fusion. The prediction is based on combining fixed-order perturbation theory and all-order threshold resummation. At fixed order we include the exact next-to-next-to-leading order (NNLO) plus an approximate next-to-next-to-next-to-leading order (N[Formula: see text]LO[Formula: see text]) which is based on the recent computation at this order for the scalar case. We then add threshold resummation at next-to-next-to-next-to leading logarithmic accuracy (N[Formula: see text]LL[Formula: see text]). Various forms of threshold resummation are considered, differing by the treatment of subleading terms, allowing a robust estimate of the theoretical uncertainties due to missing higher orders. With particular attention to pseudo-scalar masses of 200 and 750 GeV, we also observe that perturbative convergence is much improved when resummation is included. Additionally, results obtained with threshold resummation in direct QCD are compared with analogous results as computed in soft-collinear effective theory, which turn out to be in good agreement. We provide precise predictions for pseudo-scalar inclusive cross section at 13 TeV LHC for a wide range of masses. The results are available through updated versions of the public codes ggHiggs and TROLL.