How to avoid copying code to handle exceptions, undefined variables, & certain comparisons the same way

Question

I have a piece of code that fits each voxel in a data cube with a Gaussian plus a non-zero, sloping linear baseline. Each voxel is a spectrum of line+continuum emission that, depending on the location, can be pretty noisy and is known not to behave well at the edges of either the image or the spectral range. Therefore, sometimes it's necessary to fit the Gaussian and linear components separately in the parts of the spectra where each is most likely to occur, either because the original fit failed, or because the fitting parameters were nonsensical and I can see the line it failed to fit despite the noise. I'd skip straight to the piecewise version if I could, but the discontinuities can be problematic, not to mention it's generally a more costly procedure in terms of time and memory usage.

So my situation is this: I want to have my program respond to multiple possible conditions, where some are exceptions (ValueError and RuntimeError) and the others are Boolean or relational conditions, with the same procedure. Can that be done? Right now I have two copies of the same procedure, one in the exception block and the other in the else block with an internal if-statement, and it's annoying the heck out of me. If I can't combine them, is there a way to reroute one or the other statement to the same block of code without defining a new function?

Edit: I didn't originally include this because my attempts at reorganization were making a mess of my code at the time I wrote this. Now that I've kind of put it back together a bit (it's still a mess b/c each protocol corresponds to about 6 different cases of 3 different categories that are impossible to test for simultaneously as far as I know), here's the relevant parts of my code:

if fitfn not in ['gauss','lorentz']:
    raise IOError("Fitting function name must be 'lorentz' or 'gauss'")
cubedims=np.shape(cube)
frames=np.array([(n,w) for n,w in enumerate(wvl) if (lc-0.2<w<lc+0.2)])
#^ 3 sigma is about 27% larger than 2 fwhm
inds,wls=np.transpose(frames)
# store amplitudes of Gaussian/Lorentz profile & their errors:
fdencube=np.zeros((cubedims[1],cubedims[2]))
fduncube=np.zeros((cubedims[1],cubedims[2]))
# Store spectral index (slope) of linear continuum fit:
spindex=np.zeros((cubedims[1],cubedims[2]))
spundex=np.zeros((cubedims[1],cubedims[2]))
# Store continuum-subtracted line profiles:
lincube=np.zeros((len(frames),cubedims[1],cubedims[2]))
elincube=np.zeros((len(frames),cubedims[1],cubedims[2]))
# store line-subtracted continuum profiles:
concube=np.zeros((cubedims))
econcube=np.zeros((cubedims))
for x in xrange(cubedims[1]):
    for y in xrange(cubedims[2]):
        spec = cube[:,x,y]
        uspec = ecube[:,x,y]
        try: 
            p,pcov=curvf(globals()[fitfn], wvl[~np.isnan(spec)], 
                         spec[~np.isnan(spec)],sigma=uspec[~np.isnan(spec)],
                         bounds = [[0.01,min(wvl),np.mean(wvl[1:]- 
                                    wvl[:-1]),-10,0.],
                                   [50.,max(wvl),0.4,10,10.0]])
            fwhm=2*abs(p[2]) if fitfn=='lorentz' else 
                 p[2]*np.sqrt(8*np.log(2))
            if (fwhm < 0.16 and (lc-0.05<p[1]<lc+0.05) and 'pcov' in 
                globals()):
                stdp=np.sqrt(np.diag(pcov))
                cvw=p[-2]*frames[:,1]+p[-1]
                lincube[:,x,y]=spec[inds.astype(int)]-cvw
                elincube[:,x,y]=np.sqrt(uspec[inds.astype(int)]**2+
                                        stdp[-2]**2+stdp[-1]**2)

                lvw=gauss(wvl,p[0],p[1],p[2],0,0)
                concube[:,x,y]=spec-lvw
                econcube[:,x,y]=np.sqrt(uspec**2+stdp[0]**2+
                                        stdp[1]**2+stdp[2]**2)
                spindex[x,y]=p[-2]
                spundex[x,y]=stdp[-2]
                fdencube[x,y]=p[0]
                fduncube[x,y]=stdp[0]
            else:
                try:
                    s=spec[~inds.astype(int)]
                    u=uspec[~inds.astype(int)]
                    q,qcov=curvf(lreg,wls[~np.isnan(s)],s[~np.isnan(s)],
                                 sigma=u[~np.isnan(s)],bounds = [[-10,0.], 
                                 [10,10.0]])
                    r,rcov=curvf(globals([fitfn],wvl[inds.astype(int)], 
                                 spec[inds.astype(int)],
                                 sigma=uspec[inds.astype(int)],
                                 bounds = [[0.01,min(wvl),np.mean(wvl[1:]- 
                                 wvl[:-1]),-10,0.], 
                                 [50.,max(wvl),0.4,10,10.0]])
                    fwhmr=2*abs(r[2]) if fitfn=='lorentz' else 
                          r[2]*np.sqrt(8*np.log(2))
                    if (fwhmr < 0.16 and (lc-0.05<r[1]<lc+0.05) and 'rcov' 
                        in globals()):
                        stdr=np.sqrt(np.diag(rcov))
                        stdq=np.sqrt(np.diag(qcov))

                        lvw=gauss(wvl,r[0],r[1],r[2],0,0)
                        concube[:,x,y]=spec-lvw
                        econcube[:,x,y]=np.sqrt(uspec**2+stdr[0]**2+
                                        stdr[1]**2+stdr[2]**2)
                        cvw=q[0]*frames[:,1]+q[1]
                        lincube[:,x,y]=spec[inds.astype(int)]-cvw
                        elincube[:,x,y]=np.sqrt(uspec[inds.astype(int)]**2+
                                                stdq[-2]**2+stdq[-1]**2)

                        spindex[x,y]=q[0]
                        spundex[x,y]=stdq
                        fdencube[x,y]=r[0]
                        fduncube[x,y]=stdr[0]
                except (ValueError,RuntimeError):
                    fdencube[x,y]=np.NaN
                    fduncube[x,y]=np.NaN
                    lincube[:,x,y]=np.NaN
                    elincube[:,x,y]=np.NaN                        
                    try:

                    q,qcov=curvf(lreg,wvl[~np.isnan(spec)], 
                           spec[~np.isnan(spec)],
                           sigma=uspec[~np.isnan(spec)],bounds = [[-10,0.], 
                           [10,10.0]])

                        if 'qcov' in globals():
                            concube[:,x,y]=spec
                            econcube[:,x,y]=uspec
                            spindex[x,y]=q[0]
                            spundex[x,y]=np.sqrt(np.diag(qcov))[0]                              
                        else:
                            concube[:,x,y]=spec
                            econcube[:,x,y]=uspec
                            spindex[x,y]=q[0]
                            spundex[x,y]=np.NaN
                    except (ValueError,RuntimeError):
                        print 'fit failed'
                        concube[:,x,y]=spec
                        econcube[:,x,y]=uspec
                        spindex[x,y]=np.NaN
                        spundex[x,y]=np.NaN
        except (ValueError,RuntimeError):
            try: 
                s=spec[~inds.astype(int)]
                u=uspec[~inds.astype(int)]
                q,qcov=curvf(lreg,wls[~np.isnan(s)],s[~np.isnan(s)],
                             sigma=u[~np.isnan(s)],bounds = [[-10,0.], 
                                   [10,10.0]])
                r,rcov=curvf(globals()[fitfn],wvl[inds.astype(int)], 
                             spec[inds.astype(int)],
                             sigma=uspec[inds.astype(int)],
                             bounds = [[0.01,min(wvl),np.mean(wvl[1:]- 
                                      wvl[:-1]),-10,0.],
                               [50.,max(wvl),0.4,10,10.0]])
                fwhmr=2*abs(r[2]) if fitfn=='lorentz' else 
                      r[2]*np.sqrt(8*np.log(2))
                if (fwhmr < 0.16 and (lc-0.05<r[1]<lc+0.05) and 'rcov' in 
                   globals()):
                    stdr=np.sqrt(np.diag(rcov))
                    stdq=np.sqrt(np.diag(qcov))

                    lvw=gauss(wvl,r[0],r[1],r[2],0,0)
                    concube[:,x,y]=spec-lvw
                    econcube[:,x,y]=np.sqrt(uspec**2+stdr[0]**2+
                                    stdr[1]**2+stdr[2]**2)
                    cvw=q[0]*frames[:,1]+q[1]
                    lincube[:,x,y]=spec[inds.astype(int)]-cvw
                    elincube[:,x,y]=np.sqrt(uspec[inds.astype(int)]**2+
                                            stdq[-2]**2+stdq[-1]**2)

                    spindex[x,y]=q[0]
                    spundex[x,y]=stdq
                    fdencube[x,y]=r[0]
                    fduncube[x,y]=stdr[0]
            except (ValueError,RuntimeError):
                fdencube[x,y]=np.NaN
                fduncube[x,y]=np.NaN
                lincube[:,x,y]=np.NaN
                elincube[:,x,y]=np.NaN                    
                try:
                    q,qcov=curvf(lreg,wvl[~np.isnan(spec)], 
                                spec[~np.isnan(spec)],
                                sigma=uspec[~np.isnan(spec)],bounds = 
                                [[-10,0.],[10,10.0]])
                    if 'qcov' in globals():
                        concube[:,x,y]=spec
                        econcube[:,x,y]=uspec
                        spindex[x,y]=q[0]
                        spundex[x,y]=np.sqrt(np.diag(qcov))[0]
                    else:
                        print 'fit failed'
                        concube[:,x,y]=spec
                        econcube[:,x,y]=uspec
                        spindex[x,y]=np.NaN
                        spundex[x,y]=np.NaN

                except (ValueError,RuntimeError):
                    #if fit fails, assume it's continuum;
                    # continuum spectral profile doesn't matter as much
                    print 'fit failed'
                    concube[:,x,y]=spec
                    econcube[:,x,y]=uspec
                    spindex[x,y]=np.NaN
                    spundex[x,y]=np.NaN

As you can see it's a mess of duplication because I don't know how to check for exceptions, variable existence, and how a variable (assuming it exists) compares to another variable all in one line. I'm in the process of writing a module that checks for variable existence and does Boolean and relational operations such that if the variable is undefined the return value is always False, but that still leaves the exceptions. Is there a way to store an exception?

Answer 1

Technically yes.
Your code has to capture the exception in except and handle any partial actions in try block that may have happened. The logic to all this should be after the excpet block

def i_div(x,y):
    exp = None
    try:
        res =  x/y;
    except Exception as e:
        exp = e
    if isinstance(exp, ZeroDivisionError) and x == 0:
        return 'undefined'
    elif isinstance(exp, ZeroDivisionError):
        return 'infinity'
    elif isinstance(exp, TypeError):
        return 'some_value missing'
    else: 
        return res

i_div(1, 0)
i_div(0, 0)
i_div(0, 1)
i_div(1, 2)

However note that this code becomes very messy and prone to error, since ou have to anticipate all possible scenarios where the exception may or may not have occurred along with the partial/complete results.

If possible avoid such approach. Exception is a way to identify some bad is is preventing the expected flow (success case) of the program.

Answer 2

I figured out how to address the crux of my issue. In order to test for the existence of a variable and use a relational or boolean operator on it on the same line—for example, checking if var1 exists and is between var2 - 0.05 and var2 + 0.05—the syntax is as follows:

if 'var1' in globals() and (var2 - 0.05) < globals()['var1'] < (var2 + 0.05):
     do stuff

You can also replace globals() with locals() if you're searching for a variable that, barring an exception, should be defined within the same function in which you have to check for it in a later if-statement. There's some nuance to deciding which to use that did briefly cause me a bit of a headache (more basic info at What's the difference between globals(), locals(), and vars()?): If your definition of a local variable includes a global variable, locals() may not catch it. But if you use globals(), you have make sure that, if you're working at the command line or in a console (e.g. IPython or JuPyter), you haven't used the variable name you're going to be checking for. If you have, you have to delete it with the del command. Regardless, the if 'var1' in globals()/locals() part is essential—it's what keeps the subsequent conditional statement from raising a NameError.

So let's say you were doing what I was doing when I asked this question: fitting a function where you knew that besides working as intended, the fit could

• fail entirely (i.e. raise a ValueError or a RuntimeError)
• produce an unusable (full of infs) covariance matrix, or
• have certain parameters settle on mathematically acceptable but physically unrealistic values that you can't get rid of with scipy.optimize.curve_fit's "bounds" kwarg (or else the fitting routine might just return your limits as parameters that subsequent if-statements wouldn't catch).

Let's further say you know the protocol for all those scenarios is the same: modify the function and do another fit. Obviously you don't want to cut and paste the same code to 3 different locations. Here's how I solved the issue:

try:
    p,pcov=curve_fit(gauss,wavelens,spectrum,sigma=s,bounds=lims)
    fwhm = p[2]*np.sqrt(8*np.log(2))
    goodlc = linecenter-0.05<p[1]<linecenter+0.05
    goodcov = np.inf not in pcov
except (ValueError,RuntimeError,NameError) as e:
    #'as e' isn't needed if you don't use it below
    print e #I don't technically need this anymore;
    #I only used it to see why it threw an OptimizeWarning when pcov was fine
    pass #necessary if previous line isn't included.

if (('fwhm' in locals() and locals()['fwhm']<0.16) and 
    ('goodlc' in globals() and globals()['goodlc'] is True) and
    ('goodcov' in locals() and locals()['goodcov'] is True)):
    #store results in arrays and stuff...
else:
    try:
        #...etc, same stuff as above with a different function & nested.
        #my program tries 3 different functional forms before giving up

Curve_fit throws a ValueError if there are NaNs in the data to fit and a RuntimeError if it can't find the optimal parameters. The NameError catches attempts to calculate the FWHM if the appropriate parameter in p doesn't exist. Other lines before except define some conditionals to keep the later if-statement from getting too long.