1
$\begingroup$

I have a measured signal (current of a motor, turning on and off again) to which I want to fit a rectangular function in python. I came up with a reasonable rect(...) function and the scipy.optimize.curve_fit is running, but the results are far from perfect (see code below).

import numpy as np
from scipy.optimize import curve_fit

def rect(x ,a=1, b=1, x0=0, y0=0):
    return a*np.where(np.abs(np.asarray(x)-x0)<=b/2, 1., 0.)+y0

# time, current
t,y = np.array([-0.20, -0.061726,
    -0.15, -0.061358,
    -0.10, -0.061036,
    -0.05, -0.060445,
     0.00, -0.044412,
     0.05,  0.257997,
     0.10,  0.277212,
     0.15,  0.287594,
     0.20,  0.257026,
     0.25,  0.269109,
     0.30,  0.278451,
     0.35,  0.262615,
     0.40,  0.291011,
     0.45,  0.269262,
     0.50,  0.271740,
     0.55,  0.269469,
     0.60,  0.278860,
     0.65,  0.280009,
     0.70,  0.309264,
     0.75,  0.307818,
     0.80,  0.311503,
     0.85,  0.311120,
     0.90,  0.319188,
     0.95,  0.306017,
     1.00,  0.324442,
     1.05,  0.308988,
     1.10,  0.287356,
     1.15,  0.306467,
     1.20,  0.830262,
     1.25, -0.061024,
     1.30, -0.061441,
     1.35, -0.061425,
     1.40, -0.061275,
     1.45, -0.060958,
     1.50, -0.060957,
     1.55, -0.060451,
     1.60, -0.061107,
     1.65, -0.060990,
     1.70, -0.061227,
     1.75, -0.061214,
     1.80, -0.061618]).reshape((-1,2)).T

p0 = np.array([.5, .5, .5, 0.1]) # initial point
bounds =([-1., .1, -.2, -1.],
         [ 1., 2., 1.5,  1.])
popt, pcov = curve_fit(rect,t,y,p0=p0,bounds=bounds,check_finite=True)

pmanual = np.array([.35, 1.2, .65, -.06]) # manual set of parameters, which is close to optimal
# plot
import matplotlib.pyplot as plt
plt.figure(1); plt.clf()
plt.plot(t,y,label='original data')
plt.plot(t,rect(t,*p0),label='initial point')
plt.plot(t,rect(t,*popt),label='curve_fit result')
plt.plot(t,rect(t,*pmanual),label='expected result',ls='--')
plt.legend()

The curve_fit does not seem to converge to a reasonable set of parameters. Changing the values in p0 gives different results but none of these seem to be close to the optimum. How can I improve the code/data to robustly find a good set of parameters?

$\endgroup$
2
  • 1
    $\begingroup$ Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking. $\endgroup$
    – Community Bot
    Aug 23 at 19:27
  • 4
    $\begingroup$ Welcome to Scicomp! Maybe you don't need to fit. Sweep the array from front to back with a reasonable threshold, use the outermost threshold passes as the bounds for your rect, use the average of points in between for your amplitude. $\endgroup$
    – MPIchael
    Aug 24 at 6:17

1 Answer 1

3
$\begingroup$

Using the following script, attached to your data, the result obtained is plotted in cyan

from scipy.optimize import minimize
import matplotlib.pyplot as plt

def bump(x,x1,h0,dh,dx,a):
    if x <= x1:
        return h0
    x2 = x1 + dx
    if x > x1 and x <= x2:
        return h0 + (dh/dx)*(x - x1)
    x3 = x2 + a
    if x > x2 and x <= x3:
        return h0 + dh
    x4 = x3 + dx
    if x > x3 and x <= x4:
        return h0 + dh - (dh/dx)*(x - x3)
    if x > x4:
        return h0


def objective(x):
    (x1,h0,dh,dx,a) = x
    error = 0
    for k in range(len(y)):
        error += abs(y[k] - bump(t[k],x1,h0,dh,dx,a))
    return error

x0 = [0.,0.,1.,0.1,1.1]
bds = ((-0.5, None),(-0.5, None),(0, None),(0, None),(0.8, None))

res = minimize(objective, x0, method = 'SLSQP', bounds=bds)
(x1,h0,dh,dx,a) = list(res.x)
print(res.fun)
print('h0 = ', h0,' dh = ', dh,' dx = ',dx,' a = ', a,' x1 = ', x1) 

val = []
for k in range(len(y)):
    val.append(bump(t[k],x1,h0,dh,dx,a))

plt.plot(t,val)

enter image description here

$\endgroup$
2
  • $\begingroup$ The code works as expected. Could you please explain a bit, what's the idea behind the redefined rect-function, the objective function and the chosen optimizer? $\endgroup$
    – n4pK
    Aug 30 at 18:29
  • 1
    $\begingroup$ First, it appears that a trapezoid shape with rise and fall flanks is more appropriate. Second, the curve-fitting procedure uses as well a minimization process which is unknown to the user. I prefer an explicit constrained optimization method choice. $\endgroup$
    – Cesareo
    Aug 30 at 21:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.