1

Hello Stackoverflow community,

I am trying to fit data to a Faddeeva function (optimize.special.wofz) using pyhton's optimize.leastsq() or optimize.curve_fit(). The fit parameters are the following two: z1 and z2. They are complex, whereas the independent variable (time) and the output of the function (meas_data) are purely real numbers. This is my first attempt to fit the data:

import numpy as np
from scipy import optimize
from scipy import special

meas_data = np.loadtxt('directory')
time = np.loadtxt('directory')

def test(params, time):
    z1 = params[0]
    z2 = params[1]

    a = z1*np.sqrt(time)
    b = z2*np.sqrt(time)

    a = np.complex(0, a)
    b = np.complex(0, b)

    c = special.wofz(a)
    d = special.wofz(b)

    return np.real(c*d)

def test_error(params, time, t_error):
    return test(params, time) - t_error

initial_guess = (300+200j, 300-200j)
params_fit, cov_x, infodict, mesg, ier = optimize.leastsq(test_error, initial_guess, args = (time, meas_data), full_output = True)

My second attempt looks like :

import numpy as np
from scipy import optimize
from scipy import special

meas_data = np.loadtxt('directory')
time = np.loadtxt('directory')

def test(time, z1, z2):
    a = z1*np.sqrt(time)
    b = z2*np.sqrt(time)

    a = np.complex(0, a)
    b = np.complex(0, b)

    c = special.wofz(a)
    d = special.wofz(b)


    return np.real(c*d)

popt, pcov = optimize.curve_fit(test, time, meas_data)

For both cases, I get a similar error message:

for the first attempt:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-13-9286b2981692> in <module>()
     22 
     23 initial_guess = (300+200j, 300-200j)
---> 24 params_fit, cov_x, infodict, mesg, ier = optimize.leastsq(test_error, initial_guess, args = (time, msd), full_output = True)

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
    375     if not isinstance(args, tuple):
    376         args = (args,)
--> 377     shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
    378     m = shape[0]
    379     if n > m:

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape)
     24 def _check_func(checker, argname, thefunc, x0, args, numinputs,
     25                 output_shape=None):
---> 26     res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
     27     if (output_shape is not None) and (shape(res) != output_shape):
     28         if (output_shape[0] != 1):

<ipython-input-13-9286b2981692> in test_error(params, time, t_error)
     19 
     20 def test_error(params, time, t_error):
---> 21     return test(params, time) - t_error
     22 
     23 initial_guess = (z1, z2)

<ipython-input-13-9286b2981692> in test(params, time)
     10     b = z2*np.sqrt(time)
     11 
---> 12     a = np.complex(0, a)
     13     b = np.complex(0, b)
     14 

TypeError: only length-1 arrays can be converted to Python scalars

and for the second attempt:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-8-8f631a7ede54> in <module>()
     16     return np.real(c*d)
     17 
---> 18 popt, pcov = optimize.curve_fit(test, time, msd)

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in curve_fit(f, xdata, ydata, p0, sigma, absolute_sigma, check_finite, bounds, method, jac, **kwargs)
    674         # Remove full_output from kwargs, otherwise we're passing it in twice.
    675         return_full = kwargs.pop('full_output', False)
--> 676         res = leastsq(func, p0, Dfun=jac, full_output=1, **kwargs)
    677         popt, pcov, infodict, errmsg, ier = res
    678         cost = np.sum(infodict['fvec'] ** 2)

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
    375     if not isinstance(args, tuple):
    376         args = (args,)
--> 377     shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
    378     m = shape[0]
    379     if n > m:

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape)
     24 def _check_func(checker, argname, thefunc, x0, args, numinputs,
     25                 output_shape=None):
---> 26     res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
     27     if (output_shape is not None) and (shape(res) != output_shape):
     28         if (output_shape[0] != 1):

/Users/tthalheim/anaconda/lib/python3.5/site-packages/scipy/optimize/minpack.py in func_wrapped(params)
    453     if weights is None:
    454         def func_wrapped(params):
--> 455             return func(xdata, *params) - ydata
    456     else:
    457         def func_wrapped(params):

<ipython-input-8-8f631a7ede54> in test(time, z1, z2)
      7     b = z2*np.sqrt(time)
      8 
----> 9     a = np.complex(0, a)
     10     b = np.complex(0, b)
     11 

TypeError: only length-1 arrays can be converted to Python scalars

The data I am using for fitting are times in the range of 10e-6 to 10e-2 and measurement data in the range of 10e-19 to 10e-16. Both test functions used for calculating individual numbers given that the z1 and z2 are known work. I think that it has something to do with python's fitting routines which maybe not can handle complex values during their calculation?

I would be very happy, if someone could help me fixing this problem.

  • Please, can you give an example of the inputs `meas_data` and `time` ? Then you should decompose the function `test()` to understand what's going wrong (it's easier with some more steps, the traceback will point the real issue/line) – PRMoureu Jun 27 '17 at 09:18
  • Ok, I have decomposed the return. – user8219028 Jun 27 '17 at 11:16
  • the complex building may need a workaround, have you tried with `np.vectorize(complex)(0,a)` ? – PRMoureu Jun 27 '17 at 12:18
  • So, I tried your suggestion right now, precisely I replaced, for example, in the first attempt the line `a = np.complex(0, a)` by `a = np.vectorize(complex)(0, a)`. Python raises now a new type error: `---> 24 params_fit, cov_x, infodict, mesg, ier = optimize.leastsq(test_error, initial_guess, args = (time, msd), full_output = True)...... TypeError: Cannot cast array data from dtype('complex128') to dtype('float64') according to the rule 'safe' ` – user8219028 Jun 27 '17 at 13:19
  • What about the second attempt ? maybe this [discussion](http://numpy-discussion.10968.n7.nabble.com/numpy-complex-td31569.html) or this [post](https://stackoverflow.com/questions/2598734/numpy-creating-a-complex-array-from-2-real-ones) could help. It's hard to reproduce without sample datas. Your first attempt confuses me ( line `initial_guess = (z1, z2)` , these variables are not declared) – PRMoureu Jun 27 '17 at 14:40
  • Thank you very much for your hint behind the "discussion" link. That solved the problem: It was the line 'a = np.complex(0, a)' which did not work in both attempts since this line requires a scalar but cannot work for arrays. Using the correspinding demand from numpy (`np.complex128()`) fixed this problem. I furthermore decomposed the parameters for the initial guess into real and imaginary part and calculated the complex numbers out of them within the test function. So, I could fix another problem emerging after solving the `complex` problem. – user8219028 Jun 28 '17 at 06:42
  • Sorry for the confusion with the initial guess parameters. I will edit my post. `z1` and `z2` basically stored the numbers I calculated as a initial guess from some other formulas. – user8219028 Jun 28 '17 at 06:44
  • if you edit the post, try to reduce the code to the minimal example (the complex creation), and please add your solution as an answer to help the next ;-) – PRMoureu Jun 28 '17 at 08:13

1 Answers1

0

The third comment by PRMoureu on my question fixed the problem.