  Population synthesis studies constitute a powerful method to reconstruct the
birth distribution of periods and magnetic fields of the pulsar population.
When this method is applied to populations in different wavelengths, it can
break the degeneracy in the inferred properties of initial distributions that
arises from single-band studies. In this context, we extend previous works to
include $X$-ray thermal emitting pulsars within the same evolutionary model as
radio-pulsars. We find that the cumulative distribution of the number of X-ray
pulsars can be well reproduced by several models that, simultaneously,
reproduce the characteristics of the radio-pulsar distribution. However, even
considering the most favourable magneto-thermal evolution models with fast
field decay, log-normal distributions of the initial magnetic field
over-predict the number of visible sources with periods longer than 12 s. We
then show that the problem can be solved with different distributions of
magnetic field, such as a truncated log-normal distribution, or a binormal
distribution with two distinct populations. We use the observational lack of
isolated NSs with spin periods P>12 s to establish an upper limit to the
fraction of magnetars born with B > 10^{15} G (less than 1\%). As future
detections keep increasing the magnetar and high-B pulsar statistics, our
approach can be used to establish a severe constraint on the maximum magnetic
field at birth of NSs.
