For the more general case, one option is to fit a sinusoidal model to the data. A quick and dirty example is provided below.
using FFTW, LsqFit
a0 = 3.0;
T0 = 17;
t= 0:0.01:T0;
f0 = 1/T0;
ph0 = 0.17;
y = a0*cos.(2*pi*f0 .* t .+ ph0);
N = floor(Int64,length(y)/2);
Yf = fft(y)[1:N];
famfm = findmax(abs.(Yf));
p1 = maximum(abs.(y)); # amplitude initial guess
p2 = (famfm[2] - 1)/(t[end]-t[1]); # frequency initial guess
p3 = 0; # phase initial guess
P0 = [p1, p2, p3];
@. model(t, p) = p[1]*cos(2*pi*p[2]*t + p[3])
fit = curve_fit(model,t,y,P0)
println("Amplitude = ", coef(fit)[1])
println("Period (s) = ", 1/coef(fit)[2])
println("Phase (rad) = ", coef(fit)[3])