Background: Abundance calculations of a certain class of proton-rich isotopes, known as $p$ nuclei, require knowledge of the cross sections of thousands of nuclear reactions entering a reaction network. As a result, the solution of the latter relies on the predictions of the Hauser-Feshbach (HF) theory and, hence, on the reliability of the models describing the nuclear parameters entering the HF calculations, notably the optical model potential (OMP), the nuclear level density (NLD) and the $\ensuremath{\gamma}$-ray strength function ($\ensuremath{\gamma}\mathrm{SF}$).Purpose: The present work reports on a systematic study of proton capture reactions on Sr isotopes at energies relevant to the $p$ process which is responsible for the production of the $p$ nuclei at explosive stellar sites. The purpose of the work reported here is to perform a validity test of the different OMP, NLD, and $\ensuremath{\gamma}\mathrm{SF}$ models through extensive and detailed comparisons between HF calculations and experimental cross section data. This test is necessary to understand the origin of discrepancies between the $p$-nuclei abundances observed in the solar system and those predicted by the different astrophysical models, known as $p$-process models, aiming at describing the nucleosynthesis of the $p$ isotopes.Method: Cross sections were determined from $\ensuremath{\gamma}$-angular distribution measurements and from angle-integrated $\ensuremath{\gamma}$ spectra taken with the $4\ensuremath{\pi}\phantom{\rule{4pt}{0ex}}\ensuremath{\gamma}$-summing technique. Cross-section data and the resulting astrophysical $S$ factors were compared with Hauser-Feshbach calculations obtained with the latest version 1.95 of the nuclear reaction code talys using combinations of global semi-microscopic and phenomenological models of optical potentials (OMPs), nuclear level densities (NLDs), and $\ensuremath{\gamma}$-ray strength functions ($\ensuremath{\gamma}\mathrm{SFs}$).Results: Total cross sections as well as cross sections to the ground and metastable states were determined for the reactions $^{86}\mathrm{Sr}(p,\ensuremath{\gamma})^{87}\mathrm{Y}$, $^{87}\mathrm{Sr}(p,\ensuremath{\gamma})^{88}\mathrm{Y}$, and $^{88}\mathrm{Sr}(p,\ensuremath{\gamma})^{89}\mathrm{Y}$ at incident proton-beam energies from 2.5 to 3.6, 2 to 5, and 1.5 to 5 MeV, respectively.Conclusions: The experimental data reported in the present work are in very good agreement with the talys 1.95 calculations obtained with the default combination of OMP, NLD, and $\ensuremath{\gamma}\mathrm{SF}$ models. This combination is based on purely phenomenological models. A semimicroscopic proton-nucleus optical model potential was optimized at low energies leading to an equally good agreement between experimental data and theoretical calculations based solely on combinations of fully semimicroscopic models of OMP, NLD, and $\ensuremath{\gamma}\mathrm{SF}$. Our results highlight the need for a continued effort on the systematic study of proton-capture reactions to reduce the range of uncertainties arising from global nuclear models for as wide a range of relevant nuclei as possible. In this regard, new ($p,\ensuremath{\gamma}$) data at the lowest possible energies below the opening of the neutron channel are of key importance to improve global proton-nucleus optical model potentials.