Commit c6252746 authored by Sean Leavey's avatar Sean Leavey
Browse files

Scale FFT propagation functionality added to oscar.py (with example)

parent b0861468
from __future__ import absolute_import
from __future__ import division
from __future__ import print_function
from __future__ import unicode_literals
import pylab as pl
import pykat.oscar as oscar
import numpy as np
def main():
grid1 = oscar.grid(256, 256, 1e-2, 1e-2, xoffset=0, yoffset=0)
grid2 = oscar.grid(256, 256, 0.5e-2, 0.5e-2)
field1 = oscar.field(grid1, w=1e-3, power=1, mode='HG 0 0')
field2 = oscar.field(grid1, w=1e-3, power=1, mode='HG 0 0')
field1.propagate(10)
field2.scalePropagate(10, 1, grid2)
Z1 = np.abs(field1.amplitude) ** 2
Z2 = np.abs(field2.amplitude) ** 2
fig = pl.figure()
pl.subplot(2, 1, 1)
ax = pl.gca()
xLim = (field1.grid.xaxis.min(), field1.grid.xaxis.max())
yLim = (field1.grid.yaxis.min(), field1.grid.yaxis.max())
ax.set_xlim(*xLim)
ax.set_ylim(*yLim)
im = ax.imshow(Z1, extent=[xLim[0], xLim[1], yLim[0], yLim[1]])
cb = fig.colorbar(im, ax=ax)
pl.subplot(2, 1, 2)
ax = pl.gca()
xLim = (field2.grid.xaxis.min(), field2.grid.xaxis.max())
yLim = (field2.grid.yaxis.min(), field2.grid.yaxis.max())
ax.set_xlim(*xLim)
ax.set_ylim(*yLim)
im = ax.imshow(Z2, extent=[xLim[0], xLim[1], yLim[0], yLim[1]])
cb = fig.colorbar(im, ax=ax)
pl.show()
if __name__ == '__main__':
main()
......@@ -233,7 +233,7 @@ class field(object):
self.grid = grid
self.ref_index = 1 # reffractive index
self.ref_index = 1 # refractive index
self.wavelength = 1064E-9 # [m]
# TODO does this need to include ref_index??
......@@ -424,7 +424,37 @@ class field(object):
tmpfield = tmpfield * propMatrix
self.amplitude = np.fft.ifft2(tmpfield)
return self
def scalePropagate(self, distance, scale, newGrid):
# FFT propagation code with an adaptive grid size.
# Propagates to a scaled coordinate system, see Virgo Book of
# Physics pages 179-184, the scaling factor is given
# as w1/w0 with w0 the beam size at the start of propagation
# and w1 the expected beam size at the end of propatation.
# NOT YET TESTED
# scale = w0 / w1
if scale <= 0:
raise Exception('Specified scale factor must be > 0')
if newGrid.xpoints != self.grid.xpoints or newGrid.ypoints != self.grid.ypoints:
raise Exception('New grid must have same number of points as existing grid')
plD = np.pi * self.wavelength * distance * scale / self.ref_index
invz0 = (1.0 / scale - 1.0) / distance
# initial scaling
field = self.amplitude * np.exp(-1j * self.k_prop * self.grid.r_squared * invz0 / 2.0)
field = np.fft.fft2(field)
# scaled propagator
field = field * np.exp(-1j * self.k_prop * distance) * np.exp(1j * plD * self.grid.fft_ir_squared)
field = np.fft.ifft2(field)
# final scaling
self.amplitude = field * scale * np.exp(1j * newGrid.r_squared * (invz0 + distance * invz0 * invz0) / 2.0)
# update grid
self.grid = newGrid
def passAperture(self , diam = 0, shape = 'round'):
# Pass through a round aperture of a fixed diameter
......
Markdown is supported
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