def Broaden2(self, data, model, oversampling=5, m=101, dimension=20, full_output=False):
        n = data.x.size * oversampling

        #n must be even, and m must be odd!
        if n % 2 != 0:
            n += 1
        if m % 2 == 0:
            m += 1

        #resample data
        Spectrum = UnivariateSpline(data.x, data.y / data.cont, s=0)
        Model = UnivariateSpline(model.x, model.y, s=0)
        xnew = np.linspace(data.x[0], data.x[-1], n)
        ynew = Spectrum(xnew)
        model_new = FittingUtilities.RebinData(model, xnew).y

        #Make 'design matrix'
        design = np.zeros((n - m, m))
        for j in range(m):
            for i in range(m // 2, n - m // 2 - 1):
                design[i - m // 2, j] = model_new[i - j + m // 2]
        design = mat(design)

        #Do Singular Value Decomposition
        try:
            U, W, V_t = svd(design, full_matrices=False)
        except np.linalg.linalg.LinAlgError:
            outfilename = "SVD_Error.log"
            outfile = open(outfilename, "a")
            np.savetxt(outfile, np.transpose((data.x, data.y, data.cont)))
            outfile.write("\n\n\n\n\n")
            np.savetxt(outfile, np.transpose((model.x, model.y, model.cont)))
            outfile.write("\n\n\n\n\n")
            outfile.close()
            sys.exit("SVD did not converge! Outputting data to %s" % outfilename)

        #Invert matrices:
        #   U, V are orthonormal, so inversion is just their transposes
        #   W is a diagonal matrix, so its inverse is 1/W
        W1 = 1.0 / W
        U_t = np.transpose(U)
        V = np.transpose(V_t)

        #Remove the smaller values of W
        W1[dimension:] = 0
        W2 = diagsvd(W1, m, m)

        #Solve for the broadening function
        spec = np.transpose(mat(ynew[m // 2:n - m // 2 - 1]))
        temp = np.dot(U_t, spec)
        temp = np.dot(W2, temp)
        Broadening = np.dot(V, temp)

        if any( np.isnan(Broadening) ):
            data.output("Broaden2_Data.dat")
            model.output("Broaden2_Model.dat")
            np.savetxt("Broaden2_DesignMat.dat", design)
            sys.exit("NaN found in Broadening kernel!")

        #Make Broadening function a 1d array
        spacing = xnew[2] - xnew[1]
        xnew = np.arange(model.x[0], model.x[-1], spacing)
        model_new = Model(xnew)
        Broadening = np.array(Broadening)[..., 0]

        #Find the peak
        #maxindex = np.argmax(Broadening)
        maxindex = m // 2

        #Find the first zero-crossings on either size of maxindex
        left = maxindex - self.FindZeroCrossing(Broadening[:maxindex + 1][::-1]) + 1
        right = maxindex + self.FindZeroCrossing(Broadening[maxindex:])
        Broadening[:left] = 0.0
        Broadening[right:] = 0.0
        Broadening /= Broadening.sum()


        #Ensure that the broadening function is appropriate:
        maxindex = Broadening.argmax()
        badfit = False
        if maxindex > m / 2 + m / 10 or maxindex < m / 2 - m / 10:
            #The maximum should be in the middle because we already did wavelength calibration!
            outfilename = "SVD_Error2.log"
            np.savetxt(outfilename, np.transpose((Broadening, )))
            warnings.warn("Warning! SVD Broadening function peaked at the wrong location! "
                          "See SVD_Error2.log for the broadening function")
            badfit = True


        elif Broadening[maxindex] < 3 * np.mean(Broadening[maxindex + int(m / 5):]):
            outfilename = "SVD_Error2.log"
            np.savetxt(outfilename, np.transpose((Broadening, )))
            warnings.warn("Warning! SVD Broadening function is not strongly peaked! "
                          "See SVD_Error2.log for the broadening function")
            badfit = True

        if badfit:
            resolution = self.GetValue('resolution')
            model, resolution = self.FitResolution(data, model, resolution)

            if full_output:
                # Make broadening function from the gaussian
                centralwavelength = (data.x[0] + data.x[-1]) / 2.0
                FWHM = centralwavelength / resolution;
                sigma = FWHM / (2.0 * np.sqrt(2.0 * np.log(2.0)))
                left = 0
                right = np.searchsorted(xnew, 10 * sigma)
                x = np.arange(0, 10 * sigma, xnew[1] - xnew[0])
                gaussian = np.exp(-(x - 5 * sigma) ** 2 / (2 * sigma ** 2))
                return model, [gaussian / gaussian.sum(), xnew]
            else:
                return model


        #If we get here, the broadening function looks okay.
        #Convolve the model with the broadening function
        model = DataStructures.xypoint(x=xnew)
        Broadened = UnivariateSpline(xnew, np.convolve(model_new, Broadening, mode="same"), s=0)
        model.y = Broadened(model.x)

        #Fit the broadening function to a gaussian
        params = [0.0, -Broadening[maxindex], maxindex, 10.0]
        params, success = leastsq(self.GaussianErrorFunction, params, args=(np.arange(Broadening.size), Broadening))
        sigma = params[3] * (xnew[1] - xnew[0])
        FWHM = sigma * 2.0 * np.sqrt(2.0 * np.log(2.0))
        resolution = np.median(data.x) / FWHM
        #idx = self.parnames.index("resolution")
        #self.const_pars[idx] = resolution

        logging.info("Approximate resolution = %g" % resolution)

        #x2 = np.arange(Broadening.size)

        if full_output:
            return FittingUtilities.RebinData(model, data.x), [Broadening, xnew]
        else:
            return FittingUtilities.RebinData(model, data.x)
Fits the broadening profile using singular value decomposition. This function is
dataset/ETHPy150Open kgullikson88/Telluric-Fitter/src/TelluricFitter.py/TelluricFitter.Broaden2