The goal of this work is to characterize the secondary neutron spectra produced by 1 GeV/u56Fe beam colliding with a thick cylindric aluminum target and to perform a quantitative comparison with simulated results obtained with Monte Carlo codes. The measurements were performed using extended-range Bonner sphere spectrometers at two positions (15° and 40°) with respect to the beam direction. The secondary radiation field was simulated using four Monte Carlo codes (FLUKA, MCNP6, Geant4 and PHITS) and several physical models of nuclei transport and interaction. Neutron and proton energy distributions were simulated for the experimental measurement positions. The simulated neutron spectra, together with data measured with Bonner sphere spectrometers, after carrying out the correction of the contributions induced by the secondary protons, were used as input for the MAXED spectrum unfolding code to obtain the measured neutron spectra. Unfolded neutron spectra were compared with simulated ones to carry out a quantitative analysis of the performance of the chosen Monte Carlo codes and their corresponding physical models. This comparison showed that, because of experimental uncertainties and physical models, there are no unique solutions for each measurement location, but a range of solutions where the true experimental neutron spectra probably lie. The results showed deviations between 4.23% and 8.42% for some simulated spectra. Regarding the total integral values of neutron fluence and ambient equivalent dose, the unfolded neutron spectra showed deviations lower than 2%.