@@ -0,0 +1,2633 @@
+#!/usr/bin/env python%0A# -*- coding: utf-8 -*-%0A%0Aimport scipy.integrate as integrate%0Aimport numpy as np%0Aimport partials%0A%0Aclass Kernel(object):%0A    %22%22%22 creates kernels %22%22%22%0A    def __init__(self, random_noise=False, noise_level=0.001):%0A        self.random_noise = random_noise%0A        if (random_noise):%0A            self.noise_level = noise_level%0A%0A    def Kg2(self, X, Y, sigma, muX=0.0, muY=0.0):%0A        %22%22%22 return a 2 dimensional Gaussian kernel %22%22%22%0A        kernel = np.zeros(%5BX, Y%5D)%0A        theta = %5B1.0, sigma, muX, muY%5D%0A        for i in xrange(X):%0A            for j in xrange(Y):%0A                kernel%5Bi%5D%5Bj%5D = integrate.dblquad(lambda x, y: %5C%0A                                partials.g2(%5Bx + float(i) - (X-1.0)/2.0, %5C%0A                                y + float(j) - (Y-1.0)/2.0%5D, theta), %5C%0A                                -0.5, 0.5, lambda y: -0.5, lambda y: 0.5)%5B0%5D%0A%0A        if (self.random_noise):%0A            kernel = kernel + np.abs(np.random.randn(X, Y))*self.noise_level%0A            maximum = np.max(kernel)%0A            minimum = np.min(kernel)%0A            return map(lambda arr: (arr-minimum)/(maximum-minimum), kernel)%0A        else:%0A            maximum = np.max(kernel)%0A            minimum = np.min(kernel)%0A            print maximum, minimum%0A            return map(lambda arr: (arr-minimum)/(maximum-minimum), kernel)%0A%0A    @staticmethod%0A    def FWHMg2(theta):%0A        %22%22%22 get the FWHM of the Gaussian distribution %22%22%22%0A        return 2*theta%5B3%5D*sqrt(log(1.0 / %5C%0A            ((pi**2.0)*(theta%5B3%5D)**4.0)))%0A    %0A    def Km2(self, X, Y, alpha, beta, muX=0.0, muY=0.0):%0A        %22%22%22 return a 2 dimensional Moffat kernel %22%22%22%0A        kernel = np.zeros(%5BX, Y%5D)%0A        theta = %5B1.0, alpha, beta, muX, muY%5D%0A        for i in xrange(X):%0A            for j in xrange(Y):%0A                kernel%5Bi%5D%5Bj%5D = integrate.dblquad(lambda x, y: %5C%0A                                partials.m2(%5Bx+ float(i) - (X-1.0)/2.0, %5C%0A                                y + float(j) - (Y-1.0)/2.0%5D, theta), %5C%0A                                -0.5, 0.5, lambda y: -0.5, lambda y: 0.5)%5B0%5D%0A%0A        if (self.random_noise):%0A            kernel = kernel + np.abs(np.random.randn(X, Y))*self.noise_level%0A            maximum = np.max(kernel)%0A            minimum = np.min(kernel)%0A            return map(lambda arr: (arr-minimum)/(maximum-minimum), kernel)%0A        else:%0A            maximum = np.max(kernel)%0A            minimum = np.min(kernel)%0A            return map(lambda arr: (arr-minimum)/(maximum-minimum), kernel)%0A%0A    @staticmethod%0A    def FWHMm2(theta):%0A        %22%22%22 get the FWHM of the Moffat distribution %22%22%22%0A        return 2*theta%5B0%5D*sqrt(pow(2, 1 / theta%5B1%5D) - 1.0)%0A

Create kernel.py