-
Notifications
You must be signed in to change notification settings - Fork 0
/
example_mask_design_GSA.py
78 lines (52 loc) · 3.32 KB
/
example_mask_design_GSA.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
from beamshapy.BeamShaper import BeamShaper
from beamshapy.helpers import load_yaml_config, save_input_beam
from LightPipes import mm, um,Power, Intensity
import matplotlib.pyplot as plt
import numpy as np
simulation_config = load_yaml_config("./config/simulation.yml")
input_beam_config = load_yaml_config("./config/input_beam.yml")
optical_system_config = load_yaml_config("./config/optical_system.yml")
results_directory = "experiment_results"
beam_shaper = BeamShaper(simulation_config,input_beam_config, optical_system_config)
# Generate System Input Beam Field
input_field = beam_shaper.generate_input_beam(input_beam_config)
# Save Input Beam Field
# save_input_beam(results_directory, beam_shaper, input_field)
# Fresnel Lens parameters
radius = 2*mm
parabola_coef = 10**6.62
hyper_gauss_order = 12
# Generate Target Intensity Profile
target_intensity = beam_shaper.intensity_generator.generate_target_intensity_profile(profile_type="Fresnel Lens",
radius= radius,
parabola_coef = parabola_coef,
hyper_gauss_order = hyper_gauss_order)
conv_target_intensity = beam_shaper.intensity_generator.convolve_with_gaussian(target_intensity, sigma=1*um,kernel_size=10)
plt.plot(target_intensity[target_intensity.shape[0]//2,:],label="Target Intensity")
plt.plot(conv_target_intensity[conv_target_intensity.shape[0]//2,:],label="Convolved Target Intensity")
plt.legend()
plt.show()
# Generate Target Field from the given Target Intensity Profile
target_field = beam_shaper.generate_target_field_from_intensity(conv_target_intensity)
#Apply Gerchberg-Saxton Algorithm to generate the mask
results_dict = beam_shaper.mask_generator.generate_GSA_mask(input_field, target_field,init_gsa_parabola_coef=10**6)
# get the index of the minimum RMSE
min_rmse_index = np.argmin(results_dict["list rmse"])
best_mask = results_dict["list phase masks"][min_rmse_index]
best_intensity = results_dict["list image plane intensity"][min_rmse_index]
fig, ax = plt.subplots(1,2)
ax[0].imshow(best_mask,extent=[(beam_shaper.x_array_in/mm).min(),(beam_shaper.x_array_in/mm).max(),(beam_shaper.x_array_in/mm).min(),(beam_shaper.x_array_in/mm).max()])
ax[0].set_title("Best Mask")
ax[1].imshow(best_intensity,extent=[(beam_shaper.x_array_out/mm).min(),(beam_shaper.x_array_out/mm).max(),(beam_shaper.x_array_out/mm).min(),(beam_shaper.x_array_out/mm).max()])
ax[1].set_title(f"Best Intensity - radius={radius/mm}mm - coeff={int(parabola_coef)} - order={hyper_gauss_order}")
plt.show()
fig, ax = plt.subplots(1,2)
ax[0].plot(beam_shaper.x_array_out/mm,best_intensity[best_intensity.shape[0]//2,:],label="Achieved Intensity")
ax[0].plot(beam_shaper.x_array_out/mm,Intensity(target_field)[Intensity(target_field).shape[0]//2,:],label="Target Intensity")
ax[0].set_title("Best Intensity - Cut along X (min RMSE)")
ax[0].legend()
ax[1].plot(beam_shaper.x_array_out/mm,best_intensity[:,best_intensity.shape[1]//2],label="Achieved Intensity")
ax[1].plot(beam_shaper.x_array_out/mm,Intensity(target_field)[:,Intensity(target_field).shape[1]//2],label="Target Intensity")
ax[1].set_title("Best Intensity - Cut along Y (min RMSE)")
ax[1].legend()
plt.show()