opticspy
LG doughnut mode
Could you provide the mathematical expression for the LG doughnut mode in LightPipes? Thanks in advance!
Dear linmanshuo,
I have never worked with the higher order Gauss beams in LightPipes (or outside), maybe Fred has some more insight. For the moment I can suggest to you to study the source code of the functions to be sure.
The GaussBeam(...,doughnut=True) method will mix (i.e. add the complex fields) two Laguerre Gauss modes of order (n=0,m=1) where one is rotated and phase shifted:
Fout = GaussLaguerre(Fin, w0, p=n, l=m, A = 1.0 ) #first copy
Fout = Interpol( Fout, Fout.siz, Fout.N, 0, 0, 360 / (4 * m), 1) #rotate by 90°/m, e.g. 90° for m=1
Fout = MultPhase(Fout, PI/2 ) #shift phase by pi/2=90°
Fout = BeamMix(GaussLaguerre(Fin, w0, p=n, l=m, A=1 ), Fout) #add complex fields to create spiral phase and doughnut intensity
The mathematical basis used for the GaussLaguerre function is listed in the command reference in addition to the code.
https://opticspy.github.io/lightpipes/command-reference.html#LightPipes.GaussLaguerre
I must admit it looks like there is a discrepancy in the prefactors ($\rho$ instead of $\rho/2$) between the command reference and the actual code. This will lead to a wrong amplitude/intensity for a given $A$, but correct shape.
I worked with laguer modes in my master's research. I tried using your laguerre command but it has some problems with it.
First, which only has the real part of the" iltheta" phase term, it should have the others as well. According to the phase plots it is only correct for a certain region of the simulation. I looked at all the source code and this is partly due to the scipy functions you call to generate the laguerre-gauss functions
See https://en.wikipedia.org/wiki/Gaussian_beam in Laguerre-Gaussian modes.
If you compare the phase of laguerre-gauss doughnut from the page https://opticspy.github.io/lightpipes/DoughnutModes.html with fig 1 of paper Physical meaning of the radial index of Laguerre-Gauss beams and fig 1 page 3 of Coupling of atmospheric perturbed optical beams with guided modes of propagation you will notice some difference .
I have some code that I used in my master's degree that I can try to implement in lights pipes. I used light pipes to do all my research simulations, but that part had problems. I have other suggestions for implementation for the generation of partially coherent beams (subject of my research).
thanks
I propose the following modification of the GaussLaguerre command (in core.py):
def GaussLaguerre(Fin, w0, p = 0, l = 0, A = 1.0, ecs = 1 ):
"""
*Substitutes a Laguerre-Gauss mode (beam waist) in the field.*
:math:`F_{p,l}(x,y,z=0) = A \\left(\\rho\\right)^{\\frac{|l|}{2} }L^p_l(\\rho)e^{-\\frac{\\rho}{2}}\\cos(l\\theta)`,
with: :math:`\\rho=\\frac{2(x^2+y^2)}{w_0^2}`
:math:`\\theta=atan(y/x)`
if :math:`sincos = 0` replace :math:`cos(l\\theta)` by :math:`exp(-il\\theta)`
if :math:`sincos = 2` replace :math:`cos(l\\theta)` by :math:`sin(l\\theta)`
:param Fin: input field
:type Fin: Field
:param w0: Gaussian spot size parameter in the beam waist (1/e amplitude point)
:type w0: int, float
:param p: mode index (default = 0.0)
:param l: mode index (default = 0.0)
:type p: int, float
:type l: int, float
:param A: amplitude (default = 1.0)
:type A: int, float
:param ecs: 0 = exp, 1 = cos, 2 = sin (default = 1)
:type ecs: int, float
:return: output field (N x N square array of complex numbers).
:rtype: `LightPipes.field.Field`
:Example:
>>> F = GaussLaguerre(F, 3*mm) # Fundamental Gauss mode, LG0,0 with a beam radius of 3 mm
>>> F = GaussLaguerre(F, 3*mm, p=3) # Idem, LG3,0
>>> F = GaussLaguerre(F, 3*mm, p=3, l=1, A=2.0) # Idem, LG3,1, amplitude 2.0
>>> F = GaussLaguerre(F, 3*mm, 3, 1, 2.0) # Idem
>>> F = GaussLaguerre(F, 3*mm, p=3, l=2, ecs=0) LG3,1 with exponential phase factor
>>> F = GaussLaguerre(F, 3*mm, p=3, l=2, ecs=2) LG3,1 with sine phase factor
.. seealso::
* :ref:`Examples: Laguerre Gauss modes.<Laguerre Gauss modes.>`
Reference::
A. Siegman, "Lasers", p. 642
"""
# ************* Backward compatibility section ****************
#The general backward_compatible decorator does not work for this command,
#because of the positional argument w0.
#Old style: GaussLaguerre(p, l, A, w0,Fin)
#New style: GaussLaguerre(Fin, w0, p=0, l=0, A=1.0)
_using_oldstyle = False
if not isinstance(Fin, Field):
#first arg is not a field, either backward compat syntax or
# complete usage error -> find out if Field is last, else error
if isinstance(A, Field):
#found field in last arg
_using_oldstyle = True #just in case code wants to know this later
# in function
Fin, w0, p, l, A = A, l, Fin, w0, p
#caution: python can swap the values only if written on single
# line, if split up a temporary assignment is necessary
# (since a=b, b=a would not work, only temp=a, a=b, b=temp)
#-> now all the variables contain what is expected in new style
else:
raise ValueError('GaussLaguerre: Field is neither first nor '
+ 'last parameter (backward compatibility check)'
+ ', please check syntax/usage.')
# ************* end of Backward compatibility section *********
if ecs not in [0, 1, 2]:
raise ValueError('GaussLaguerre: invalid value for ecs. ecs must be: 0, 1 or 2')
Fout = Field.copy(Fin)
R, Phi = Fout.mgrid_polar
w02=w0*w0
la=abs(l)
rho = 2*R*R/w02
if ecs == 1:
Fout.field = A * rho**(la/2) * genlaguerre(p,la)(rho) * _np.exp(-rho/2) * _np.cos(l*Phi)
if ecs == 2 and l == 0:
Fout.field = A * rho**(la/2) * genlaguerre(p,la)(rho) * _np.exp(-rho/2)
if ecs == 2 and l != 0:
Fout.field = A * rho**(la/2) * genlaguerre(p,la)(rho) * _np.exp(-rho/2) * _np.sin(l*Phi)
#################################
# """
# Corrections by: GubioGL 23-8-2023:
# I changed only that part, not calculating the R and PHI with the function,
# because in the part of the phase of the .mgrid_polar function there is a "+ _np.pi",
# which shouldn't be here.
# In order not to change that function, it is simpler to just change it here.
# optional exponential(..), cos( ..) or sin( .. ) (esc = 0,1,2)
# Thank you GubioGL!
# Fred van Goor, 26-8-2023
# """
if ecs == 0:
Y, X = Fout.mgrid_cartesian
R = _np.sqrt(X**2+Y**2)
rho = 2*R**2/w0**2
Fout.field = A * rho**(la/2) * genlaguerre(p,la)(rho) * _np.exp(-rho/2) * _np.exp(-1j*l*_np.arctan2(Y, X))
Fout._IsGauss=False
return Fout
Now you can choose your phase factor: exp, cosine or sine with ecs (default = 1: cosine). All these are eigenmodes. Please comment. Fred
