def Broaden(self, data, model, oversampling=5, m=101, dimension=20, full_output=False):
        """
        Fits the broadening profile using singular value decomposition. This function is
        called by GenerateModel, and is not meant to be called directly!

        -oversampling is the oversampling factor to use before doing the SVD
        -m is the size of the broadening function, in oversampled units
        -dimension is the number of eigenvalues to keep in the broadening function. (Keeping too many starts fitting noise)

        -NOTE: This function works well when there are strong telluric lines and a flat continuum.
               If there are weak telluric lines, it's hard to not fit noise.
               If the continuum is not very flat (i.e. from the spectrum of the actual
                 object you are trying to telluric correct), the broadening function
                 can become multiply-peaked and oscillatory. Use with care!
        """
        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)

        #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]

        #Ensure that the broadening function is appropriate:
        maxindex = Broadening.argmax()
        badfit = False
        if maxindex > m / 2.0 + m / 10.0 or maxindex < m / 2.0 - m / 10.0:
            #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 np.mean(Broadening[maxindex - int(m / 10.0):maxindex + int(m / 10.0)]) < 3 * np.mean(
                Broadening[int(m / 5.0):]):
            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)
dataset/ETHPy150Open kgullikson88/Telluric-Fitter/src/TelluricFitter.py/TelluricFitter.Broaden/original