2010/1/28 Saurav Pathak <sau...@sas.upenn.edu> > Hi, > > I am trying to fit real parameters using a complex function. That is, the > GSL example: > ( > http://www.gnu.org/software/gsl/manual/html_node/Example-programs-for-Nonlinear-Least_002dSquares-Fitting.html > ) > > int expb_f (const gsl_vector * x, void *data, gsl_vector * f) > > > f in my case is complex. This means, my Jacobian too is complex. > It is not clear to me from the example how I may specify this. Could > someone please help? >
Hi, GSL does not have direct support for non-linear fit of complex data. You can still fit your data by working on the real and imaginary part. If the fitting parameters are real numbers this is quite is to arrange, tipically you have to pack the data in a gsl_vector by alternating the real and imaginary part. The Jacobian will be splitted with the same logic. I can give you more details if it can helps. Otherwise, if you want to fit complex data of complex parameters the affair will get more complicated because some components of the jacobian will be related by the cauchy-riemann equations and you need to take ths into account. I strongly advise in this case you use only real parameters because otherwise it gets too much complicated. If you are interested my program, gsl shell : http://www.nongnu.org/gsl-shell/index.html does have direct support for non linear fit of complex data with real parameters. In the download directory you can find the latest release gsl-shell-0.9.6-2 with the windows binary. Here a working example for fit of complex data: ------------------ n = 50 p0 = vector {0.55, -4, 16, 0.23} function fmodel(x, t) return x[1] * exp((x[2] + 1i*x[3]) * t + 1i * x[4]) end r = rng() -- random number generator y = cnew(n, 1, |k,j| fmodel(p0, (k-1)/n) + 0.01 * rnd.gaussian(r)) function cexpf(x, f, J) for k=1, n do local t, y = (k-1)/n, y[k] if f then f:set(k, 1, fmodel(x,t) - y) end if J then local e = exp((x[2] + 1i*x[3]) * t + 1i * x[4]) J:set(k, 1, e) J:set(k, 2, t * x[1] * e) J:set(k, 3, 1i * t * x[1] * e) J:set(k, 4, 1i * x[1] * e) end end end function print_state(s) print ("x: ", s.x:row_print()) print ("chi square: ", cmul(h(s.f), s.f)[1]) end s = csolver {fdf= cexpf, n= n, p= 4, x0= vector {2.1, -2.8, 18, 0}} repeat print_state (s) local status = s:iterate() until status ~= 'continue' print_state (s) ------------------------ You can also plot the results with the following commands: ------------------ pl = plot() lr = path(0, real(y[1])) -- the data at t=0 lf = path(0, real(fmodel(p0, 0))) -- the model at t = 0 for k=1, y:dims() do lr:line_to((k-1)/n, real(y[k])) lf:line_to((k-1)/n, real(fmodel(p0, (k-1)/n))) end pl:add(lr, 'black', {{'stroke'}, {'marker'}}) pl:add_line(lf) pl:show() ------------------ Best regards, Francesco _______________________________________________ Help-gsl mailing list Help-gsl@gnu.org http://lists.gnu.org/mailman/listinfo/help-gsl
