13. Difracción en campo cercano

Autor:Luis Miguel Sánchez Brea
Revisor:José Luis Vilas Prieto

Herramientas utilizadas en este tutorial:

Archivos necesarios para esta sección:


ejemplos de visualización 3D

13.1. Introducción

En esta sección analizaremos la difracción en campo cercano bajo la aproximación de Fresnel.

Según esta ecuación, el campo propagado a una cierta distancia z de un elemento difractor, se calcula a través de la siguiente integral

\mathbf{E}(x,y,z)=\frac{1}{i\lambda z}e^{ikz}\iint\mathbf{E}_{0}(\xi,\eta)e^{i\frac{k}{2z}\left[\left(x-\xi\right)^{2}+\left(y-\eta\right)^{2}\right]}d\xi d\eta

donde \mathbf{E}_{0}(\xi,\eta)=t(\xi,\eta)\mathbf{E}_{inc}(\xi,\eta) es el campo a la salida del elemento difractor que, bajo la aproximación de elemento delgado se obtiene como la multiplicación del campo incidente \mathbf{E}_{inc}(\xi,\eta) por el coeficiente de transmisión de la máscara difractora t(\xi,\eta)

13.2. Programación

Para el desarrollo de los ejemplos de esta sección hemos programado numéricamente el cálculo de esta integral a partir de 'fast-Fourier-transform based direct integration (FFT-DI) method' según el artículo de Shen y Wang [1] . Este algoritmo lo hemos implementado en la clase camposXY (que es heredada por mascarasXY y fuentesXY) con el nombre de RS (a partir de Rayleigh-Sommerfed, que es una versión ligeramente mejorada de la integral de Fresnel).

Por consiguiente, si tenemos un campo, máscara o fuente, tendrá un método RS que podemos aplicar. Si tenemos un campo de tamaño NxM, el resultado de la propagación será también NxM, aunque existe un parámetro de amplificación que nos permite calcular tamaños del plano de observación mayores (jN)x(jM).

La ventaja que tiene es que nos devuelve un parametro de calidad que nos informa si la aproximación se ha realizado debidamente. Si el parámetro de calidad (almacenado en self.calidad) tiene un valor superior a 1, la propagación se ha realizado convenientemente.

	def RS(self, z = 10 * mm, newField = True, tipo='z', xout = None, yout = None):
		"""Metodo de Fast-Fourier-Transform para la integracion numerica
		de la Formula de la difraccion de Rayleigh-Sommerfeld.
		Extraido de Applied Optics vol 45 num 6 (1102-1110) - 2006
		- z es la distancia (si es <0 se realiza la propagacion inversa
		- newfield, si es False el calculo se queda en la instancia, si es True se genera nueva
		- xout, yout se utilizan para decir las coordenadas de salida (para la amplificacion)
  
 		el tamaño del pixel debe ser igual, aunque la máscara puede tener distinto tamaño
		"""
		
  

		xin = self.x; 		yin = self.y
		x1 = self.x[0];		y1 = self.y[0]
		xin1 = self.x[0];	yin1 = self.y[0]

		if xout == None:
			xout = self.x
		if yout == None:
			yout = self.y

		nx = len(xout);				ny = len(yout)
		dx = xout[1] - xout[0];		dy = yout[1] - yout[0]

		#parametro de calidad
		dr_real = sp.sqrt(dx ** 2 + dy ** 2)
		rmax = sp.sqrt((xout ** 2).max() + (yout ** 2).max())
		dr_ideal = sp.sqrt(self.wavelength ** 2 + rmax ** 2 + 2 * self.wavelength * sp.sqrt(rmax ** 2 + z ** 2)) - rmax
		self.calidad = dr_ideal / dr_real
		#siempre que se hace un cálculo se informa por consola de la calidad de este.
		if(self.calidad > 1):
			print 'Buena aproximacion: factor ', self.calidad
		else:
			print 'Necesita un muestreo mayor: factor', self.calidad

		#matriz computada
		W = 1
		U = sp.zeros((2 * ny - 1, 2 * nx - 1), dtype = complex)
		U[0:ny, 0:nx] = W*self.u   #el transpose lo he cambiado porque daba problemas para matrices no cuadradas
		xext = x1 - xin[::-1]  #da la vuelta
		xext = xext[0:-1]
		xext = sp.concatenate((xext, self.x - xin1))
		yext = y1 - yin[::-1]
		yext = yext[0:-1]
		yext = sp.concatenate((yext, self.y - yin1))
		Xext, Yext = sp.meshgrid(xext, yext)

		#permite calcula la propagacion y la propagacion inversa, cuando z<0.
		if z > 0:
			H = kernelRS(Xext, Yext, self.wavelength, z, tipo=tipo)
		else:
			H = kernelRSinversa(Xext, Yext, self.wavelength, z, tipo=tipo)

		#calculo de la transformada de Fourier
		S = ifft2(fft2(U) * fft2(H)) * dx * dy  #el transpose lo he cambiado porque daba problemas para matrices no cuadradas
		Usalida = S[ny - 1:, nx - 1:] #hasta el final

		#los calculos se pueden dejar en la instancia o crear un nuevo campo (más usual).
		if newField == True:
			campoSalida = campoXY(self.x, self.y, self.wavelength)
			campoSalida.u = Usalida / z
			campoSalida.calidad=self.calidad
			return campoSalida
		else:
			self.u = Usalida / z

	

Todos los cálculos realizados en esta sección serán númericos, lo cual nos da una gran ventaja computacional: solamente tenemos que describir la fuente y la máscara, y tendremos el campo difractado a una distancia z.

[1]Fabin Shen and Anbo Wang &quot;Fast-Fourier-transform based numerical integration method for the Rayleigh–Sommerfeld diffraction formula&quot; Applied Optics vol. 45 num 6. pp. 1102-1110 (2006).

13.3. Difracción por una abertura cuadrada

Como primer ejemplo realizaremos la propagación en campo cercano de una abertura bidimensional cuadrada. Utilizaremos este ejemplo para explicar el proceso.

  • Se definen los parámetros dimensionales de las máscaras y los parámetros ópticos (longitud de onda).
  • Se crea una fuente de iluminación a través de la clase fuenteXY().
  • Se crea una máscara con las dimensiones apropiadas.
  • Se realiza la Aproximación de elemento delgado: u2 = u1 * t1
  • se realiza la propagación: u3 = u2.RS(z = z_difraccion, newField = True)
  • Se dibuja el resultado.

El código empleado y las gráficas resultantes son las siguientes.

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
#-------------------------------------

import sys
sys.path.append('../')
sys.path.append('../general/')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY


def difraccion_cuadrado():
	#tamano de area de visualizacion
	tamano = 100 * um
	ndatos = 256
	x0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	y0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	#longitud de onda
	lambda0 = 0.6238 * um

	#fuente de iluminacion
	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	#mascara
	radio = 20 * um
	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0)
	t1.cuadrado(r0 = (0 * um, 0 * um), size = (60 * um, 60 * um), angulo = 0 * grados)
	u2 = u1 * t1

	#plano de observacion
	z_difraccion = 500*um

	#propagacion y acondicionamiento grafica
	u3 = u2.RS(z = z_difraccion, newField = True)
	texto = "z=%d $\mu m$" % (z_difraccion)

	dibujar_varios_campos(campos=(u2,u3),titulos=('cuadrado', texto))


difraccion_cuadrado()
plt.show()

(Source code)

13.4. Difracción por una abertura cuadrada: varias distancias

En el ejemplo anterior se vislumbra el comportamiento general en campo cercano: la distribución de intensidad se asemeja algo al objeto, pero los efectos difractivos producen fluctuaciones de intensidad que aumentan a medida que nos separamos del objeto difractor. Para ver esto de forma más rigurosa, en el siguiente ejemplo se muestra la distribución de intensidad en varios planos. Para ello generamos dos funciones, mascara_cuadrado() que calcula el campo inicial

def mascara_cuadrado():
	tamano = 25 * um
	ndatos = 256
	x0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	y0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	lambda0 = 0.6238 * um

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0)
	t1.cuadrado(r0 = (0 * um, 0 * um), size = (20 * um, 20 * um), angulo = 0 * grados)

	u2 = u1 * t1
	return u2
	

y otra que recibe el campo u2 y genera un bucle con difrerentes distancias

def difraccion_varias_posiciones(u2, distancias):

	plt.figure(figsize=(12,6))

	num_figuras = len(distancias)

	filas = 	[1, 1, 1, 2, 2, 2, 3, 3, 3, 2, 4, 4, 4, 4, 5, 4]
	columnas = 	[1, 2, 3, 2, 3, 3, 3, 3, 3, 5, 3, 3, 4, 4, 3, 4]


	extension=[u2.x[0],u2.x[-1],u2.y[0],u2.y[-1]]
	plt.subplot(filas[num_figuras], columnas[num_figuras], 1)
	h1 = plt.imshow(abs(u2.u) ** 2, interpolation = 'bilinear', 
					origin = 'lower', extent = extension)	
	h1.set_cmap("gist_heat")
	plt.title(u'mascara')

	for z, i in zip(distancias, range(num_figuras)):
		#difraccion a la distancia estipulada
		u3 = u2.RS(z = z, newField = True)
		texto = "z=%d $\mu m$" % (z)
		#carga de los dibujos
		plt.subplot(filas[num_figuras], columnas[num_figuras], i + 2)
		plt.axis('off')
		plt.title(texto)
		h1 = plt.imshow(abs(u3.u) ** 2, interpolation = 'bilinear', 
					origin = 'lower', extent = extension)
		h1.set_cmap("gist_heat")
	plt.tight_layout()

El campo total se genera llamando a las dos funciones secuencialmente, transmitiendo el valor u2 de una a la otra función:

u2=mascara_cuadrado()
difraccion_varias_posiciones(u2, distancias = sp.linspace(25 * um, 2* mm, 5))
plt.show()

(Source code)

Eventualmente, si la distancia de propagación fuera mucho mayor que la distancia de Fresnel, z >> a^{2}/\lambda, el campo difractado sería la transformada de Fourier del objeto (pues estamos iluminando con onda plana), en este caso una función sinc.

13.5. Difracción por una abertura circular

Una máscara de gran importancia en la óptica es la abertura circular. Según hemos definido la clase mascarasXY, para determinar el campo generado por la máscara circular, los ejemplos anteriores son completamente válidos, y solamente es necesaio intercambiar la función de creación de la máscara

t1.circulo(r0 = (0 * um, 0 * um), radius = (radio, radio))

Se muestran los resultados, pero no el código que es casi idéntico

(Source code)

Para este ejemplo se ha elegido que la distancia sea exactamente la mitad de la distancia de Fresnel

t1.circulo(r0 = (0 * um, 0 * um), radius = (radio, radio))
#distancia de fresnel
z_fresnel = radio ** 2 / lambda0

#plano de observacion
z_difraccion = .5 * z_fresnel

puesto que a esta distancia aparece un anillo bien definido.

En el archivo difraccion_circulo_varias_distancias.py se muestra la difracción en campo cercano a varias distancias. Nótese que en la función difraccion_varias_posiciones las variables filas y columnas permiten reordenar las figuras (hasta 16) según su número.

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
#-------------------------------------


import sys
sys.path.append('../')
sys.path.append('../general/')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY



def difraccion_varias_posiciones(u2, distancias):

	plt.figure(figsize=(12,8))

	num_figuras = len(distancias)

	filas = 	[1, 1, 1, 2, 2, 2, 3, 3, 3, 2, 4, 4, 4, 4, 5, 4]
	columnas = 	[1, 2, 3, 2, 3, 3, 3, 3, 3, 5, 3, 3, 4, 4, 3, 4]

	extension=[u2.x[0],u2.x[-1],u2.y[0],u2.y[-1]]
	plt.subplot(filas[num_figuras], columnas[num_figuras], 1)
	h1 = plt.imshow(abs(u2.u) ** 2, interpolation = 'bilinear',
					origin = 'lower', extent = extension)
	h1.set_cmap("gist_heat")

	plt.title(u'circle')

	for z, i in zip(distancias, range(num_figuras)):
		#difraccion a la distancia estipulada
		u3 = u2.RS(z = z, newField = True)
		texto = "z=%d $\mu m$" % (z)
		#carga de los dibujos
		plt.subplot(filas[num_figuras], columnas[num_figuras], i + 2)
		plt.axis('off')
		plt.title(texto)
		h1 = plt.imshow(abs(u3.u) ** 2, interpolation = 'bilinear',
					origin = 'lower', extent = extension)
		h1.set_cmap("gist_heat")
	plt.tight_layout()



def mascara_circulo():
	tamano = 100 * um
	ndatos = 256
	x0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	y0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	lambda0 = 0.6238 * um

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0)
	t1.circulo(r0 = (0 * um, 0 * um), radius = (40*um, 40*um))

	num_figuras = 3
	z_fresnel = (tamano) ** 2 / lambda0
	distancias_difraccion = sp.linspace(25 * um, .25 * z_fresnel, num_figuras)

	u2 = u1 * t1
	return u2


u2=mascara_circulo()
difraccion_varias_posiciones(u2, distancias = sp.linspace(25 * um, 1* mm, 3))
plt.show()

(Source code)

13.6. Difracción por un borde

En este epígrafe se considera como elemento difractor un borde. Así se considera como máscara una máscara binaria, y se calcula la propagacióna una distancia de 25 um de la máscara. El proceder es análogo a una máscara cuadrada.

#!/usr/local/bin/python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
# Objetivo:	Ejemplos de difraccion en campo lejano
#-------------------------------------


import sys
sys.path.append('../')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY


def difraccion_borde():
	#tamano de area de visualizacion
	tamano = 100 * um
	ndatos = 256
	x0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	y0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	#longitud de onda
	lambda0 = 0.6238 * um

	#fuente de iluminacion
	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	#mascara
	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0)
	t1.dosNiveles(nivel1 = 0, nivel2 = 1, xcorte = 0 * um)
	u2 = u1 * t1


	#plano de observacion
	z_difraccion = 25*um

	u3 = u2.RS(z = z_difraccion, newField = True)
	texto = "z=%d $\mu m$" % (z_difraccion)
	dibujar_varios_campos(campos=(u2,u3),titulos=('circulo', texto))
	u3.dibujar_perfil(punto1=(-15, 0),punto2=(15,0), tipo='intensidad')
	h,perfil,p1,p2=t1.perfil(punto1=(-15, 0),punto2=(15,0), tipo='intensidad',order=0)
	plt.hold(True)
	plt.plot(h,perfil,'r',lw=2)

difraccion_borde()
plt.show()

(Source code)

Debido a la difracción, la luz al propagarse se introduce en la sombra, como se observa en el perfil.

Veamos que ocurre a varias distancias. A medida que nos separamos del borde, más luz entra en la zona de sombra y mas perceptibles son los efectos difractivos

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
#-------------------------------------


import sys
sys.path.append('../')
sys.path.append('../general/')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY



def difraccion_varias_posiciones(u2, distancias):

	plt.figure(figsize=(12,6))

	num_figuras = len(distancias)

	filas = 	[1, 1, 1, 2, 2, 2, 3, 3, 3, 2, 4, 4, 4, 4, 5, 4]
	columnas = 	[1, 2, 3, 2, 3, 3, 3, 3, 3, 5, 3, 3, 4, 4, 3, 4]

	plt.subplot(filas[num_figuras], columnas[num_figuras], 1)
	h1 = plt.imshow(abs(u2.u) ** 2)
	h1.set_cmap("gist_heat")

	plt.title(u'mascara')

	for z, i in zip(distancias, range(num_figuras)):
		#difraccion a la distancia estipulada
		u3 = u2.RS(z = z, newField = True)
		texto = "z=%d $\mu m$" % (z)
		#carga de los dibujos
		plt.subplot(filas[num_figuras], columnas[num_figuras], i + 2)
		plt.axis('off')
		plt.title(texto)
		h1 = plt.imshow(abs(u3.u) ** 2)
		h1.set_cmap("gist_heat")



def mascara_borde():
	tamano = 100 * um
	ndatos = 256
	x0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	y0 = sp.linspace(-tamano / 2, tamano / 2, ndatos)
	lambda0 = 0.6238 * um

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0)
	t1.dosNiveles(nivel1 = 0, nivel2 = 1, xcorte = 0 * um)

	num_figuras = 8
	z_fresnel = (tamano) ** 2 / lambda0
	distancias_difraccion = sp.linspace(25 * um, .25 * z_fresnel, num_figuras)

	u2 = u1 * t1
	return u2


u2=mascara_borde()
difraccion_varias_posiciones(u2, distancias = sp.linspace(10 * um, 100* um, 2))
plt.tight_layout()
plt.show()

(Source code)

13.7. Difracción por una doble rendija

En este caso el elemento difractor es una doble rendija. Se calcula la propagación a varias disntancias talesque la distribución de intensidad no está en la posición de la zona clara, sino en la zona opaca, en contradicción con la teoría geométrica (cosas que tienen los efectos difractivos), es decir entre las dos rendijas.

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
# Objetivo:	difraccion por doble rendija
#-------------------------------------

import sys
sys.path.append('../')
sys.path.append('../general/')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY


def difraccionDobleRendija():
	tamano = 250 * um
	numdatos = 128
	x0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	y0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	lambda0 = 0.6238 * um
	z = 4 * mm

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana(A = 1, theta = 0 * grados, phi = 0 * grados)

	t = mascaraXY(x = x0, y = y0, wavelength = lambda0);
	t.dobleRendija(x0 = 0 * um, size = 25 * um, separacion = 50 * um, angulo = 0 * grados)

	u2 = u1 * t
	u3 = u2.RS(z=z, newField = True)

	texto = "z=%d $\mu m$" % (z)
	dibujar_varios_campos(campos=(u2,u3),titulos=('doble rendija', texto))

difraccionDobleRendija()
plt.show()

(Source code)

Nótese como aunque la óptica geométrica predice la no existencia de luz en dicha región, experimentalmente (en este caso simulación) se observa un patrón de intensidad en dicha región. Este muestra como la difracción es la capacidad de la luz para sortear bordes y agujeros.

13.8. Difracción por una lente

Otro ejemplo es la difracción por una lente, de la cual mostramos la intensidad y la fase. Debido a que la lente se ha suspuesto rotacionalmente simétrica, el patrón de difracción es perfectamente simétrico. Así se muestra un perfil de intensidad.

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
#-------------------------------------

import sys
sys.path.append('../')
sys.path.append('../general/')

from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY

def difraccion_lente():
	tamano = 200 * um
	numdatos = 256
	x0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	y0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	lambda0 = 0.6238 * um

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana();

	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0);
	t1.circulo(r0 = (0 * um, 0 * um), radius = (90 * um, 90 * um), angulo = 45 * grados)

	t2 = mascaraXY(x = x0, y = y0, wavelength = lambda0);
	t2.lente(r0 = (0 * um, 0 * um), radius = (90 * um, 90 * um), focal = (2 * mm, 2 * mm), angulo = 45 * grados)

	t = t1 * t2
	t.dibujar(tipo = 'campo')

	u3 = u1 * t
	u3 = u3.RS(z = 2 * mm, newField = True)
	u3.dibujar(tipo = 'intensidad')

	u3.dibujar_perfil(punto1=(-30*um, 0),punto2=(30*um,0), tipo='intensidad')


difraccion_lente()
plt.show()

(Source code)

13.9. Difracción por una lente binaria

En este último ejemplo se supone una lente de Fresnel. Este tipo de lentes es puramente difractiva. Así se obtiene

#!/usr/local/bin/python
#!/usr/bin/env python
# -*- coding: utf-8 -*-
#------------------------------------
# Autor:	Luis Miguel Sanchez Brea
# Fecha		2013/10/24 (version 1.0)
# Licencia:	GPL
#-------------------------------------

import sys
sys.path.append('../')
from camposXY import * #@UnusedWildImport
from fuentesXY import fuenteXY
from mascarasXY import mascaraXY


def difraccion_lente_binaria():
	tamano = 200 * um
	numdatos = 256
	x0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	y0 = sp.linspace(-tamano / 2 , tamano / 2, numdatos)
	lambda0 = 0.6238 * um

	u1 = fuenteXY(x = x0, y = y0, wavelength = lambda0);
	u1.onda_plana();

	t1 = mascaraXY(x = x0, y = y0, wavelength = lambda0);
	t1.circulo(r0 = (0 * um, 0 * um), radius = (90 * um, 90 * um), angulo = 45 * grados)

	t2 = mascaraXY(x = x0, y = y0, wavelength = lambda0);
	t2.lente(r0 = (0 * um, 0 * um), radius = (90 * um, 90 * um), focal = (2 * mm, 2 * mm), angulo = 45 * grados)
	t2.discretizar(tipo = 'fase', numNiveles = 2, factor = 1, campoNuevo = False, matriz = False)

	t = t1 * t2
	t.dibujar(tipo = 'campo')
	u3 = u1 * t
	u3 = u3.RS(z = 2 * mm, newField = True)
	u3.dibujar(tipo = 'intensidad')
	u3.dibujar_perfil(punto1=(-30*um, 0),punto2=(30*um,0), tipo='intensidad')


difraccion_lente_binaria()
plt.show()

(Source code)