Besides controlling wave trajectory inside complex media, wave velocity constitutes a relevant bio-marker for medical imaging. In a transmission configuration, wave-front distortions can be unscrambled to provide a map of the wave velocity landscape $c(\mathbf{r})$. However, most in-vivo applications correspond to a reflection configuration for which only back-scattered echoes generated by short-scale fluctuations of $c(\mathbf{r})$ can be harvested. Under a single scattering assumption, this speckle wave-field cannot provide the long-scale variations of $c(\mathbf{r})$. In this paper, we go beyond the first Born approximation and show how a map of $c(\mathbf{r})$ can be retrieved in epi-detection. To that aim, a reflection matrix approach of wave imaging is adopted. While standard reflection imaging methods generally rely on confocal focusing operations, matrix imaging consists in decoupling the location of the incident and received focal spots. Following this principle, a self-portrait of the focusing process can be obtained around each point of the medium. The Gouy phase shift exhibited by each focal spot is leveraged to finely monitor the wave velocity distribution $c(\mathbf{r})$ inside the medium. Experiment in a tissue-mimicking phantom and numerical simulations are first presented to validate our method. Speed-of-sound tomography is then applied to ultrasound data collected on the liver of a difficult-to-image patient. Beyond paving the way towards quantitative ultrasound, our approach can also be extremely rewarding for standard imaging. Indeed, each echo can be assigned to the actual position of a scatterer. It allows an absolute measurement of distance, an observable often used for diagnosis but generally extremely sensitive to wave velocity fluctuations.