We present a new measurement of the energy spectrum of iron nuclei in cosmic rays from 20 to 500 TeV. The measurement makes use of a template-based analysis method, which, for the first time, is applied to the energy reconstruction of iron-induced air showers recorded by the VERITAS array of imaging atmospheric Cherenkov telescopes. The event selection makes use of the direct Cherenkov light which is emitted by charged particles before the first interaction, as well as other parameters related to the shape of the recorded air shower images. The measured spectrum is well described by a power law $frac{mathrm{d} F}{mathrm{d} E}=f_0cdot left(frac{E}{E_0}right)^{-gamma}$ over the full energy range, with $gamma = 2.82 pm 0.30 mathrm{(stat.)} ^{+0.24}_{-0.27} mathrm{(syst.)}$ and $f_0 = left( 4.82 pm 0.98 mathrm{(stat.)}^{+2.12}_{-2.70} mathrm{(syst.)} right)cdot 10^{-7}$m$^{-2}$s$^{-1}$sr$^{-1}$TeV$^{-1}$ at $E_0=50$TeV, with no indication of a cutoff or spectral break. The measured differential flux is compatible with previous results, with improved statistical uncertainty at the highest energies.