```python
from sympy.physics.optics import RayTransferMatrix as RTM, FreeSpace as FS, ThinLens as TL
from sympy import symbols, Eq, init_printing, S, sqrt
init_printing()
import numpy as np
import matplotlib.pyplot as plt
from scipy.constants import c
```


```python
# Defining functions
ghz_cm = lambda x: c * 1e2 / x * 1e-9
cm_ghz = lambda x: c*1e-9*1e2/x
```


```python
d_out, d_in, w0_in, w0_out = symbols("d_out, d_in, w0_in, w0_out")
w, lam, w0, z, Z_r = symbols("w, lam, w0, z, Z_r")
M, E, f, f0 = symbols('M, E, f, f0')
f_p, f_m = symbols("f_p, f_m")
f_e, R1, R2 = symbols("f_e, R1, R2")
```


```python
# Beam radius
Wz = Eq(w, w0 * sqrt(1 + (z / Z_r) ** 2))
Wz
```


```python
# Input distance
eq_din_m = Eq(d_in, f - sqrt(f ** 2 - f0 ** 2) / M) # Negative sign
eq_din_p = Eq(d_in, f + sqrt(f ** 2 - f0 ** 2) / M) # Positive sign

# output distance
eq_dout_m = Eq(d_out, f - sqrt(f ** 2 - f0 ** 2) * M) # Negative sign
eq_dout_p = Eq(d_out, f + sqrt(f ** 2 - f0 ** 2) * M) # Positive sign
# NOTE: use same sign for output distance as chosed for input distane

eq_din_m, eq_din_p, eq_dout_m, eq_dout_p
```


```python
# Focal length
E = 1-M**-2
fp = Eq(f_p, Z_r*(d_in/(E*Z_r))*(1+sqrt(1-E*(1+(d_in/Z_r)**-2)))) # positive sign but output is negative focus
fm = Eq(f_m, Z_r*(d_in/(E*Z_r))*(1-sqrt(1-E*(1+(d_in/Z_r)**-2)))) # negative sign but output is positive focus

fp, fm
```


```python
magnification = 0.15/0.34

w0_ = 0.34
confocalLength = np.pi*w0_**2/ghz_cm(325)
d1 = 20

f1 = float(fp.subs({M:magnification, Z_r:confocalLength, d_in:d1}).rhs)
f2 = float(fm.subs({M:magnification, Z_r:confocalLength, d_in:d1}).rhs)

f1, f2
```


```python
fe = Eq(f_e, (R1*R2)/(R1+R2))
fe
```


```python
M_ellps = TL(fe.rhs)
M_ellps
```


```python
dout_direct = (
    (M_ellps.A * d_in + M_ellps.B) * (M_ellps.C * d_in + M_ellps.D)
    + (M_ellps.A * M_ellps.C * Z_r)
) / ((M_ellps.C * d_in + M_ellps.D) ** 2 + (M_ellps.C * Z_r) ** 2)

dout_direct
```


```python
dout_direct.subs({R1:f1, R2:f2, d_in:d1, Z_r:confocalLength})
```


```python
wout_direct = w0_in/((M_ellps.C * d_in + M_ellps.D) ** 2 + (M_ellps.C * Z_r) ** 2)
wout_direct
```


```python
wout_direct.subs({R1:f1, R2:f2, d_in:d1, Z_r:confocalLength, w0_in:w0_})
```


```python

```


```python

```


```python

```


```python

```


```python

```

18b855d2f13dc88fc10d4a88a922958c1145f098
ipynb
Jupyter Notebook
Simulations/THz/gaussian_beam_formula_check.ipynb
aravindhnivas/FELion-Spectrum-Analyser
430f16884482089b2f717ea7dd50625078971e48
MIT
Simulations/THz/gaussian_beam_formula_check.ipynb
aravindhnivas/FELion-Spectrum-Analyser
430f16884482089b2f717ea7dd50625078971e48
MIT
Simulations/THz/gaussian_beam_formula_check.ipynb
aravindhnivas/FELion-Spectrum-Analyser
430f16884482089b2f717ea7dd50625078971e48
MIT
2019-01-25T20:37:57.000Z
2019-01-25T20:37:57.000Z
Qwen/Qwen-72B
1. YES
2. YES