4225.sagews

var('x y')
density_plot(sin(x^2 + y^2), (x, -5, 5), (y, -5, 5), plot_points = 100)

x, y = var('x y')
lm = 555 * 10**-9
len = 1
d = 0.5 * 10**-2

def norm(t):
    return sqrt(t * t)

def path_difference(t, z):
    return norm(vector([t, z - d / 2, len])) - norm(vector([t, z + d / 2, len]))

density_plot((len**2 / (x**2 + y**2 + len**2)) * sin(path_difference(x, y) / lm)**2, (x, -1, 1), (y, -1, 1), plot_points = 400)