We investigate the performance of a MUSIC-type algorithm for low frequent time-harmonic electromagnetic imaging of the subsurface. This algorithm is based on an asymptotic analysis and an associated asymptotic factorization of electromagnetic near field data in the presence of very small scatterers in the ground. Our numerical results illustrate the potential of this imaging method, and we provide theoretical evidence for what has been observed numerically. In particular, we address the issue of superresolution which, in this context, means that it is possible to use this method to reconstruct objects of centimeters in size, and separated by only few multiples of their diameters, despite the fact that the given time-harmonic electromagnetic waves may have wavelengths that are several kilometers long.