We simulate apertureless near-field optical imaging and obtain phase and amplitude scans of structured substrates for elastic scattering. The solution of the three-dimensional Maxwell equations does not involve approximations and we include large tips and substrates, strong interaction, interferometric detection and demodulation at higher harmonics. Such modeling represents a significant step towards quantitative simulations and offers the attractive possibility to study the individual influence of each relevant experimental parameter. We typically obtain highly localized signatures of the interaction of the tip with gold inclusions, superposed on a slowly varying background signal. The relative importance of both contributions and the achievable lateral resolution are strongly dependent on the geometry and scanning conditions. The simulations show sensitivity mostly to the first nanometers of the sample and underline the importance of scanning near the sample and being careful with mechanical anharmonicities on the tip oscillation. They also help to determine the influence of oscillation amplitude and demodulation harmonic. (C) 2008 Optical Society of America