import numpy as np


def generate_kernel(radius, sigma=None):
    """
    Generate a 2D Gaussian kernel with a given radius.

    Parameters:
    - radius: int, the radius of the kernel (size will be 2*radius + 1)
    - sigma: float (optional), standard deviation of the Gaussian. If None, sigma = radius / 3

    Returns:
    - kernel: 2D numpy array of shape (2*radius+1, 2*radius+1)
    """
    size = 2 * radius + 1
    if sigma is None:
        sigma = radius / 3.0  # Common default choice

    print(f"radius: {radius}, sigma: {sigma}")

    # Create a grid of (x,y) coordinates
    ax = np.arange(-radius, radius + 1)
    xx, yy = np.meshgrid(ax, ax)

    # Apply the 2D Gaussian formula
    kernel = np.exp(-(xx**2 + yy**2) / (2 * sigma**2))
    kernel /= 2 * np.pi * sigma**2  # Normalize based on Gaussian PDF
    kernel /= kernel.sum()          # Normalize to sum to 1

    return kernel


if __name__ == "__main__":
    kernel = generate_kernel(radius=10, sigma=1)
    print(kernel)