ABSTRACT We introduce a new methodology for efficiently tuning sub-grid models of star formation and supernovae feedback in cosmological simulations and at the same time understanding their physical implications. Based on a set of 71 zoom simulations of a Milky Way (MW)-sized halo, we explore the feasibility of calibrating a widely used star formation and feedback model in the enzo simulation code. We propose a novel way to match observations, using functional fits to the observed baryon makeup over a wide range of halo masses. The model MW galaxy is calibrated using three parameters: the star formation efficiency (f*), the efficiency of thermal energy from stellar feedback (ϵ), and the region into which feedback is injected (r and s). We find that changing the amount of feedback energy affects the baryon content most significantly. We then identify two sets of feedback parameter values that are both able to reproduce the baryonic properties for haloes between $10^{10}\, \mathrm{M_\odot }$ and $10^{12}\, \mathrm{M_\odot }$. We can potentially improve the agreement by incorporating more parameters or physics. If we choose to focus on one property at a time, we can obtain a more realistic halo baryon makeup. Contrasting both star formation criteria and the corresponding combination of optimal feedback parameters, we also highlight that feedback effects can be complementary: to match the same baryonic properties, with a relatively higher gas-to-stars conversion efficiency, the feedback strength required is lower, and vice versa. Lastly, we demonstrate that chaotic variance in the code can cause deviations of approximately 10 per cent and 25 per cent in the stellar and baryon mass in simulations evolved from identical initial conditions.