I would do the analyis on your three dimensions separately. If the answers you get agree, then you have your confirmed answer. If processing needs to be minimized, then you can find the delay on one, and confirm with the others.
Since you have sensor data, it probably contains some noise as well. Heavily smoothing your data (I like exponential smoothing for its efficiency) will mitigate the effects of noise, reduce the number of local maxima, and allow you to do a coarser correlation measure.
There is some dispute on the proper definition of correlation. The one you should use for best results is:
$$ c = \frac{\vec x \cdot \vec y}{\| \vec x \| \| \vec y \| } $$
Where $\vec x$ is a subset of the smoothed $f_1$ and $\vec y$ is a subset of the smoothed $f_2$ for the two intervals and sampling density your are comparing. Each vector should be shifted so the elements have a mean value of zero. Take the average of all the element values and subtract it from each value.
First do a coarse brute search. Select every 10th (or whatever spacing you choose) value from your functions to build your vectors. Then calculate the correlation coefficient for a range of delay values, stepping by a coarse amount. This will give you a coarse delay value.
Repeat this process, centered at your coarse delay value, using every sample and stepping by a single sample size. You should only need to slide plus or minus half your coarse spacing side. This will give you a fine delay value.
If you need subsample resolution, use the minimum value you found in the fine search and its two neigbors, model them as a parabola and find the minimum point there.
Hope this helps.
Ced