![](https://secure.gravatar.com/avatar/5c9fb379c4e97b58960d74dcbfc5dee5.jpg?s=120&d=mm&r=g)
What makes your example similar to what I want is the frequency: it is nearly the same. Now let's assume these are two motions of which you know that the frequencies should be the same. Any differences are due to measurement errors.
OK, I'll first comment your code (I shaped it to be able to read it better. Your line are to long, especially in comments, and your comments should be more often on separate lines, indented the same way as the code). My comments are prefixed by ### +++++++++++++++++ from pylab import * from scipy import * # Generate data points with noise num_points = 150 # as i got confused with xT, Xt, ... it's now data set 1 and 2, each set # having x/y coordinates. x1 = linspace(5., 8., num_points) x2 = x1 xexp = array([x1, x2]) # i put them together. not sure whether this is the right way... ### You probably want to you use c_[x1,x2] y1 = 11.86*cos(2*pi/0.81*x1-1.32) + 6.*rand(num_points) y2 = -32.14*cos(2*pi/0.8*x2-1.94) + 3.*rand(num_points) yexp = [y1, y2] # same thing ### !!!! I thing you misstyped here you meant yexp = c_[y1,y2] # Fit setup def ytheo (p,x) : # yth = zeros((x1,x1)) how to initialize? for j in range(len(xexp)): # count through data sets for k in range(len(x1)): # count through the data points within a (the first) set yth[j,k] = p[(2*j+1)]*cos(2*pi/p[0]*x1[k]+p[1*(2*j+2)]) # Target functions. note the same p[0] for all sets. it is a matrix: # first index denoting the set, second index denoting a point in a set. return yth ### A for loop !! Your are kidding. You should try to avoid this as much ### as possible, and work only on arrays. For loops cost a lot of time. errfunc = lambda p, x, y: ytheo(p,x) - yexp # Distance to the target function # fit p0 = c_[0.8, -15., 0., -15., 0.] # Initial guess for the parameters # 0: common "frequency"; 1: 1st set ampl; 2: 1st set phase; # 3: 2nd set ampl; 4: 2nd set phase p1,success = optimize.leastsq(errfunc, p0.copy(), args = (xexp, yexp)) time = linspace(x1.min(),x1.max(),100) plot(x1,y1,"ro",time,fitfunc(p1,time),"r-") # Plot of the data and the fit plot(x2,y2,"b^",time,fitfunc(p2,time),"b-") # Legend the plot title("Oscillations in the compressed trap") xlabel("time [ms]") ylabel("displacement [um]") legend(('1 position', '1 fit', '2 position', '2 fit')) ax = axes() text(0.8, 0.07,'x freq : %.3f kHz \n y freq : %.3f kHz' %(1/p1[1],1/p2[1]), fontsize = 16, horizontalalignment='center', verticalalignment='center', transform = ax.transAxes) show() +++++++++++++++++ I'll upload an example show you how I would what you are trying to do. The wiki is just so much more readable than mail :->. I don't exactly see what your are trying to do, but it does not seems to good to me. GaÃl