python No plot or NameError UPDATE plot is visible but not as it should be

Question

I am attempting to create a "rolling spline" using polynomials via polyfit and polyval.

However I either get an error that "offset" is not defined... or, the spline doesn't plot.

My code is below, please offer suggestions or insights. I am a polyfit newby.

import numpy as np
from matplotlib import pyplot as plt

x = np.array([ 3893.50048173,  3893.53295003,  3893.5654186 ,  3893.59788744,
        3893.63035655,  3893.66282593,  3893.69529559,  3893.72776551,
        3893.76023571,  3893.79270617,  3893.82517691,  3893.85764791,
        3893.89011919,  3893.92259074,  3893.95506256,  3893.98753465,
        3894.02000701,  3894.05247964,  3894.08495254])
y = np.array([ 0.3629712 ,  0.35187397,  0.31805825,  0.3142261 ,  0.35417492,
        0.34981215,  0.24416184,  0.17012087,  0.03218199,  0.04373861,
        0.08108644,  0.22834105,  0.34330638,  0.33380814,  0.37836754,
        0.38993407,  0.39196328,  0.42456769,  0.44078106])
e = np.array([ 0.0241567 ,  0.02450775,  0.02385632,  0.02436235,  0.02653321,
        0.03023715,  0.03012712,  0.02640219,  0.02095554,  0.020819  ,
        0.02126918,  0.02244543,  0.02372675,  0.02342232,  0.02419184,
        0.02426635,  0.02431787,  0.02472135,  0.02502038])

xk = np.array([])
yk = np.array([])

w0 = np.where((y<=(e*3))&(y>=(-e*3)))
w1 = np.where((y<=(1+e*3))&(y>=(1-e*3)))
mask = np.ones(x.size)
mask[w0] = 0
mask[w1] = 0

for i in range(0,x.size):
    if mask[i] == 0:
        if ((abs(y[i]) < abs(e[i]*3))and(abs(y[i])<(abs(y[i-1])-abs(e[i])))):
            imin = i-2
            imax = i+3
            if imin < 0:
                imin = 0
            if imax >= x.size:
                imax = x.size
            offset = np.mean(x)
            for order in range(20):
                coeff = np.polyfit(x-offset,y,order)
                model = np.polyval(coeff,x-offset)
                chisq = ((model-y)/e)**2
                chisqred = np.sum(chisq)/(x.size-order-1)
                if chisqred < 1.5:
                    break
            xt = x[i]
            yt = np.polyval(coeff,xt-offset)
    else:
        imin = i-1
        imax = i+2
        if imin < 0:
            imin = 0
        if imax >= x.size:
            imax = x.size
        offset = np.mean(x)
        for order in range(20):
            coeff = np.polyfit(x-offset,y,order)
            model = np.polyval(coeff,x-offset)
            chisq = ((model-y)/e)**2
            chisqred = np.sum(chisq)/(x.size-order-1)
            if chisqred < 1.5:
                break
        xt = x[i]
        yt = np.polyval(coeff,xt-offset)

    xk = np.append(xk,xt)
    yk = np.append(yk,yt)

#print order,chisqred
################################

plt.plot(x,y,'ro')
plt.plot(xk+offset,yk,'b-') # This is the non-plotting plot
plt.show()

################################

---

Update

---

So I edited the code, removing all of the if conditions that do not apply to this small sample of data.

I also added the changes that I made which allow the code to plot the desired points... however, now that the plot is visible, I have a new problem.

The plot isn't a polynomial of the order the code is telling me it should be.

Before the plot command, I added a print, to display the order of the polynomial and the chisqred, just to be certain that it was working.

Answer 1

First, thank you for providing a self-contained sample (not many newbies do that)! If you want to improve your question, you should remove all debugging code from the sample, as now it clutters the code. The code is quite long and not very self-explanatory. (At least to me - the problem may be between my ears, as well.)

---

Let us unroll the problem from the end. The proximal reason why you get an empty plot is that you have empty xkand yk (empty arrays).

Why is that? That is because you have 19 points, and thus your for loop is essentially:

for i in range(12, 19-1-12):
    ...

There is nothing to iterate from 12..6! So actually your loop is run through exactly zero times and nothing is ever appended to xk and yk.

The same explanation explains the problem with offset. If the loop is never run through, there is no offset defined in yout plot command (xk+offset), hence the NameError.

---

This was the simple part. However, I do not quite understand your code. Especially the loops where you loop order form 0..19 look strange, as only the result form the last cycle will be used. Maybe there is something to fix?

(If you still have problems with the code after this analysis, please fix the things you can, simplify the code as much as possible, and edit your question. Then we can have another look into this!)