Commit 051096c0 authored by Andreas Freise's avatar Andreas Freise
Browse files

changing HG_beam to HG_mode, maths to math, adding LG2HG function

parent 714cc923
...@@ -10,7 +10,8 @@ import copy ...@@ -10,7 +10,8 @@ import copy
import warnings import warnings
import cmath import cmath
from math import factorial from math import factorial
from scipy.special import hermite, jacobi from scipy.special import hermite, eval_jacobi
from pykat.math.jacobi import jacobi
from pykat.SIfloat import SIfloat from pykat.SIfloat import SIfloat
class gauss_param(object): class gauss_param(object):
...@@ -271,11 +272,11 @@ class beam_param(gauss_param): ...@@ -271,11 +272,11 @@ class beam_param(gauss_param):
pass pass
# Should be renamed to HG_mode? # Should be renamed to HG_mode?
class HG_beam(object): class HG_mode(object):
""" Hermite-Gauss beam profile. Example usage: """ Hermite-Gauss beam profile. Example usage:
import pykat.optics.gaussian_beams as gb import pykat.optics.gaussian_beams as gb
qx=gb.beam_param(w0=1e-3,z=0) qx=gb.beam_param(w0=1e-3,z=0)
beam=gb.HG_beam(qx,n=2,m=0) beam=gb.HG_mode(qx,n=2,m=0)
beam.plot() beam.plot()
""" """
def __init__(self, qx, qy=None, n=0, m=0): def __init__(self, qx, qy=None, n=0, m=0):
...@@ -386,7 +387,7 @@ class HG_beam(object): ...@@ -386,7 +387,7 @@ class HG_beam(object):
return np.outer(_un, _um) return np.outer(_un, _um)
def plot(self, ndx=100, ndy=100, xscale=4, yscale=4): def plot(self, ndx=100, ndy=100, xscale=4, yscale=4):
""" Make a simple plot the HG_beam """ """ Make a simple plot the HG_mode """
import pykat.plotting import pykat.plotting
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
...@@ -405,15 +406,11 @@ class HG_beam(object): ...@@ -405,15 +406,11 @@ class HG_beam(object):
cbar = fig.colorbar(axes) cbar = fig.colorbar(axes)
plt.show() plt.show()
def HG2LG(n,m):
def hg2lg(n,m):
"""A function for Matlab which returns the coefficients and mode indices of """A function for Matlab which returns the coefficients and mode indices of
the LG modes required to create a particular HG mode. the LG modes required to create a particular HG mode.
Usage: coefficients,ps,ls = hg2lg(n,m) Usage: coefficients,ps,ls = HG2LG(n,m)
n,m: Indces of the HG mode to re-create. n,m: Indces of the HG mode to re-create.
coeffcients: Complex coefficients for each order=n+m LG mode required to coeffcients: Complex coefficients for each order=n+m LG mode required to
...@@ -444,10 +441,61 @@ def hg2lg(n,m): ...@@ -444,10 +441,61 @@ def hg2lg(n,m):
# Coefficient # Coefficient
c = (signl*1j)**m * np.sqrt(factorial(N-m)*factorial(m)/(2**N * factorial(np.abs(l)+p)*factorial(p))) c = (signl*1j)**m * np.sqrt(factorial(N-m)*factorial(m)/(2**N * factorial(np.abs(l)+p)*factorial(p)))
c = c * (-1.0)**p * (-2.0)**m * scipy.special.eval_jacobi(m,np.abs(l)+p-m,p-m,0.0) c = c * (-1.0)**p * (-2.0)**m * eval_jacobi(m,np.abs(l)+p-m,p-m,0.0)
coefficients[j] = c coefficients[j] = c
return coefficients, ps, ls return coefficients, ps, ls
def LG2HG(p,l):
""" Function to compute the amplitude coefficients
of Hermite-Gauss modes whose sum yields a Laguerre Gauss mode
of order n,m.
Usage: coefficients, ns, ms = LG2HG(p,l)
p: Radial LG index
l: Azimuthal LG index
The LG mode is written as u_pl with 0<=|l|<=p.
The output is a series of coefficients for u_nm modes,
given as complex numbers and respective n,m indices
coefficients (complex array): field amplitude for mode u_nm
ns (int array): n-index of u_nm
ms (int array): m-index of u_nm
The formula is adpated from M.W. Beijersbergen et al 'Astigmatic
laser mode converters and transfer of orbital angular momentum',
Opt. Comm. 96 123-132 (1993)
We adapted coefficients to be compatible with our
definition of an LG mode, it differs from
Beijersbergen by a (-1)^p factor and has exp(il\phi) rather
than exp(-il\phi). Also adapted for allowing -l.
Andreas Freise, Charlotte Bond 25.03.2007"""
# Mode order
N=2*p+np.abs(l)
# Indices for coefficients
n = np.abs(l)+p
m = p
# create empty vectors
coefficients = np.linspace(0,0,N+1,dtype=np.complex_)
ns = np.linspace(0,0,N+1)
ms = np.linspace(0,0,N+1)
# l positive or negative
signl = np.sign(l)
if (l==0):
signl = 1.0
# Beijersbergen coefficients
for j in np.arange(0,N+1):
ns[j]=N-j
ms[j]=j
c=(-signl*1j)**j * math.sqrt(factorial(N-j)*factorial(j)/(2**N * factorial(n)*factorial(m)))
coefficients[j] = c * (-1.0)**p * (-2)**j * jacobi(j,n-j,m-j,0.0)
return coefficients, ns, ms
from itertools import combinations_with_replacement as combinations from itertools import combinations_with_replacement as combinations
from pykat.optics.gaussian_beams import beam_param, HG_beam from pykat.optics.gaussian_beams import beam_param, HG_mode
from pykat.exceptions import BasePyKatException from pykat.exceptions import BasePyKatException
from pykat.optics.romhom import u_star_u from pykat.optics.romhom import u_star_u
from pykat.external.progressbar import ProgressBar, ETA, Percentage, Bar from pykat.external.progressbar import ProgressBar, ETA, Percentage, Bar
...@@ -7,10 +7,10 @@ from scipy.interpolate import interp2d ...@@ -7,10 +7,10 @@ from scipy.interpolate import interp2d
from scipy.integrate import dblquad from scipy.integrate import dblquad
from pykat.optics.romhom import ROMWeights from pykat.optics.romhom import ROMWeights
from math import factorial from math import factorial
from pykat.maths.hermite import hermite from pykat.math.hermite import hermite
from scipy.misc import comb from scipy.misc import comb
from scipy.integrate import newton_cotes from scipy.integrate import newton_cotes
from pykat.maths import newton_weights from pykat.math import newton_weights
import time import time
import pykat.optics.maps import pykat.optics.maps
...@@ -65,8 +65,8 @@ def adaptive_knm(mode_in, mode_out, q1, q2, q1y=None, q2y=None, smap=None, delta ...@@ -65,8 +65,8 @@ def adaptive_knm(mode_in, mode_out, q1, q2, q1y=None, q2y=None, smap=None, delta
if len(mode_in) != 2 or len(mode_out) != 2: if len(mode_in) != 2 or len(mode_out) != 2:
raise BasePyKatException("Both mode in and out should be a container with modes [n m]") raise BasePyKatException("Both mode in and out should be a container with modes [n m]")
Hg_in = HG_beam(qx=q1, qy=q1y, n=mode_in[0], m=mode_in[1]) Hg_in = HG_mode(qx=q1, qy=q1y, n=mode_in[0], m=mode_in[1])
Hg_out = HG_beam(qx=q2, qy=q2y, n=mode_out[0], m=mode_out[1]) Hg_out = HG_mode(qx=q2, qy=q2y, n=mode_out[0], m=mode_out[1])
Nfuncs = [] Nfuncs = []
Nfuncs.append(0) Nfuncs.append(0)
...@@ -168,8 +168,8 @@ def riemann_HG_knm(x, y, mode_in, mode_out, q1, q2, q1y=None, q2y=None, ...@@ -168,8 +168,8 @@ def riemann_HG_knm(x, y, mode_in, mode_out, q1, q2, q1y=None, q2y=None,
dy = abs(y[1] - y[0]) dy = abs(y[1] - y[0])
if cache is None: if cache is None:
Hg_in = HG_beam(qx=q1, qy=q1y, n=mode_in[0], m=mode_in[1]) Hg_in = HG_mode(qx=q1, qy=q1y, n=mode_in[0], m=mode_in[1])
Hg_out = HG_beam(qx=q2, qy=q2y, n=mode_out[0], m=mode_out[1]) Hg_out = HG_mode(qx=q2, qy=q2y, n=mode_out[0], m=mode_out[1])
U1 = Hg_in.Unm(x+delta[0], y+delta[1]) U1 = Hg_in.Unm(x+delta[0], y+delta[1])
U2 = Hg_out.Unm(x,y).conjugate() U2 = Hg_out.Unm(x,y).conjugate()
...@@ -505,8 +505,8 @@ def square_aperture_HG_knm(mode_in, mode_out, q, R): ...@@ -505,8 +505,8 @@ def square_aperture_HG_knm(mode_in, mode_out, q, R):
m = mode_in[1] m = mode_in[1]
_m = mode_out[1] _m = mode_out[1]
hg1 = HG_beam(q, n=n, m=m) hg1 = HG_mode(q, n=n, m=m)
hg2 = HG_beam(q, n=_n, m=_m) hg2 = HG_mode(q, n=_n, m=_m)
kx = hg1.constant_x * hg2.constant_x.conjugate() kx = hg1.constant_x * hg2.constant_x.conjugate()
ky = hg1.constant_y * hg2.constant_y.conjugate() ky = hg1.constant_y * hg2.constant_y.conjugate()
......
...@@ -13,11 +13,11 @@ import itertools ...@@ -13,11 +13,11 @@ import itertools
from copy import copy from copy import copy
from pykat.external.progressbar import ProgressBar, ETA, Percentage, Bar from pykat.external.progressbar import ProgressBar, ETA, Percentage, Bar
from itertools import combinations_with_replacement as combinations from itertools import combinations_with_replacement as combinations
from pykat.optics.gaussian_beams import beam_param, HG_beam from pykat.optics.gaussian_beams import beam_param
from scipy.linalg import inv from scipy.linalg import inv
from math import factorial from math import factorial
from pykat.maths.hermite import * from pykat.math.hermite import *
from pykat.maths import newton_weights from pykat.math import newton_weights
from scipy.integrate import newton_cotes from scipy.integrate import newton_cotes
from multiprocessing import Process, Queue, Array, Value, Event from multiprocessing import Process, Queue, Array, Value, Event
from pykat.exceptions import BasePyKatException from pykat.exceptions import BasePyKatException
...@@ -658,4 +658,4 @@ def makeWeightsNew(smap, EIxFilename, EIyFilename=None, verbose=True, newtonCote ...@@ -658,4 +658,4 @@ def makeWeightsNew(smap, EIxFilename, EIyFilename=None, verbose=True, newtonCote
if verbose: if verbose:
p.finish() p.finish()
return ROMWeights(w_ij_Q1=w_ij_Q1, w_ij_Q2=w_ij_Q2, w_ij_Q3=w_ij_Q3, w_ij_Q4=w_ij_Q4, EIx=EIx, EIy=EIy, direction=direction) return ROMWeights(w_ij_Q1=w_ij_Q1, w_ij_Q2=w_ij_Q2, w_ij_Q3=w_ij_Q3, w_ij_Q4=w_ij_Q4, EIx=EIx, EIy=EIy, direction=direction)
\ No newline at end of file
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment