![](https://secure.gravatar.com/avatar/97892c9705c58e27835441a0092ce92a.jpg?s=120&d=mm&r=g)
Gael Varoquaux wrote:
On Mon, Sep 25, 2006 at 10:39:51AM +0200, Sebastian Busch wrote:
... please specify ... /where/ you uploaded your example ...
Hi Gaƫl! Thanks for your reply. I think this is not quite what the OP wanted. In any case, it is not what I am looking for ;) You fit two data sets with two different functions. Completely independent. You could do the same in two runs of your program working on the different data sets. 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. So you sit down and want to fit the data. Not in two runs, then averaging the two frequencies you get, but in one big fit. This is better because one measurement was more precise than the other and it's hard to account for that after the fit finished. I tried to modify your example in that way. As I am pretty much a bloody beginner in python, I failed. I post my try anyway to show you more clearly what I am thinking of. And, I admit, in the hope you can do better ;) Thanks for your time, Sebastian. +++++++++++++++++++ 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... 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 # 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 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()