Interpreting the cross correlation results
Yes and no. In principle, you can use the peak of your correlation function. However, it is not at 10000. The correlation function is symmetric around 0, so your peak is actually much closer to zero than you think. This is one of the reasons why xcorr returns two paramters, one for the lags at which the function is calculated. The correct way of plotting the the correlation function would thus be a plot(s,r).
Secondly, you should try zooming in a bit. Since your time signals have 10000 samples, your acf has ~20000, making it hard to see the exact maximum. If you do zoom in you'll notice it's not exactly at zero but two samples to the right.
Why two samples? Because you have 10 samples per period and the shift between sine and cosine is a quarter period, so it's actually 2.5 samples. You notice the problem, with the sampling you've chosen it is a little hard to measure the delay precisely.
Try reducing your fm to a small value, like 10, this makes the correlation function easier to interpret. Here is an example:

In this example, the maximum occurs at sample 247, which is equivalent to 247*ts = 24.7 milliseconds.
In case you wonder about the triangular shape, this has to do with the fact that while at lag 0 (in the middle) we can &quot;compare&quot; the entirety of the two signals whereas for larger lags, we need to truncate both signals to shorter versions. The larger the lag, the fewer samples we have to correlate. For this reason, we typically check the ACF closer to lag zero.