Diffraction Primer - 0      H. Q. Le  copyrighted

ECE 5368/6358                              han q le - copyrighted
Use solely for students registered for UH ECE 6358/5368 during courses - DO NOT distributed (copyrighted materials)

Segment 1

0. Introduction

How does light moves, interacts with objects, structures?

http://www.youtube.com/watch?v=IZgYswtwlT8
http://www.youtube.com/watch?v=4EDr2YY9lyA

http://www.microtechweb.com/2d/2d.htm

What is it? essentially?
An approach to solve optical problems involving free space radiating lightwave interacting with objects and structures.
Why is it useful?
Many optical applications and simple optical effects can be accounted for with diffraction theory, which is an approximation in lieu of exact solutions, that are often just too difficult to solve, nor really necessary when approximations are “good enough”.
Beyond highly simple objects and structures, exact lightwave solutions are too difficult to obtain. If an application requires numerical accuracy, there are numerical algorithms that can be used. However, diffraction theory is more than just a “convenient” approximation, it can give insight and perceptive understanding of light behavior and optical phenenomena.

1. Green theorem approach - summary

1.1. Review

Green's theorem is used as an approach for the Helmholtz equation:
                                                       (∇^2)U + (k^2) U=0                            (1.1.1)
The Green function satisfies:
                                      (∇^2)G[r-r'] + (k^2) G[r-r'] =δ[r-r']         (1.1.2)

Theorem: If we know the field and its gradient on an enclosed surface, then the field at anypoint in the volume within the enclosed surface S is:
                                   U[r]=-(∮_S) (U[r'] ∇G[r-r'] - G[r-r'] ∇U[r'])·d(Subscript[S, in])           (1.1.8b)
where the notation d(Subscript[S, in]) is used to indicate inward normal vector as opposed to outward.

1.2. Green's function for the Helmholtz equation

The Green function for the Helmholtz Eq. is:

G[r-r']= -(1/(4 π ))((e^(i k |r-r'|))/(|r-r'|))                           (1.2.1)

1.3. Kirchhoff-Helmholtz theorem

Applying Eq. (1.2.1) in (1.1.8b):
                         U[r]=(1/(4 π))(∮_S) (U[r'] ∇((e^(i k |r-r'|))/(|r-r'|)) - ((e^(i k |r-r'|))/(|r-r'|)) ∇U[r'])·d(Subscript[S, in])            (1.3.1)
which is known as Kirchhoff-Helmholtz diffraction theorem and the right hand side is called the Kirchhoff integral. It is worth noting that for static field, let k-> 0, then: 
                          U[r]=(1/(4 π))(∮_S) (U[r'] ∇(1/(|r-r'|)) - (1/(|r-r'|)) ∇U[r'])·d(Subscript[S, in])                (1.3.2)
which is the electrostatic potential solution given boundary condition.
What this means is that if we know the field over a complete surface, then we can infer the field everywhere within the volume, but the challenge, as we see is to find the field on the surface.

Exercise: Part 1: Find   (Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|))

Needs["VectorAnalysis`"]
Grad[((E^(I k Sqrt[(x-ξ)^2+(y-η)^2+(z-ζ)^2]))/(Sqrt[(x-ξ)^2+(y-η)^2+(z-ζ)^2])),Cartesian[ξ,η,ζ]]

{(((x-ξ) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)))^(3/2)))-((I k (x-ξ) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2))),(((y-η) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)))^(3/2)))-((I k (y-η) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2))),(((z-ζ) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)))^(3/2)))-((I k (z-ζ) (E^(I k (Sqrt[(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)]))))/((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)))} /.(((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)->  ((Abs[r-r'])^2)   /. ((((((x-ξ))^2)+(((y-η))^2)+(((z-ζ))^2)))^(3/2))-> Abs[r-r'] ((Abs[r-r'])^2)

Hence:        (Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|))=-i k (((r -r'))/(|r-r'|))((e^(i k |r-r'|))/(|r-r'|))   (1 +(i/(k |r-r'|)) ) =-i k (OverHat[u])((e^(i k|r-r'|))/(|r-r'|))   (1 +(i/(k |r-r'|)) ) 
where                               (OverHat[u]) ≡(((r -r'))/(Abs[r-r']))

Exercise: Part 2: what is (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|))? and compare with (Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|)).
Since  ((e^(i k |r-r'|))/(|r-r'|))=((e^(i k |r'-r|))/(|r'-r|)) , why would   (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|)) is opposite to (Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|))?

Since ((e^(i k |r-r'|))/(|r-r'|))=((e^(i k |r'-r|))/(|r'-r|)), 
    (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|))=(Subscript[∇, r])((e^(i k |r'-r|))/(|r'-r|))=-i k (((r'-r ))/(|r'-r|))((e^(i k |r'-r|))/(|r'-r|))   (1 +(i/(k |r'-r|)) ) 
Thus:                            (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|))=-(Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|)) 
Gradient has a direction, hence, a function of 2 variables may have the same value but the gradient for each variable is different: the increase of one variable should have an opposite direction to the increase of the other variable.
Consider the function:    ((e^(i k s))/s)    where  s is a point in 3-D space, hence s is always positive. We have:   ∇((e^(i k s))/s)=(OverHat[s])((i k (e^(i k s)))/s) (1+(i/(k s)))

which is the vector pointing away from the origin.
Hence, (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|)) is the same as ∇((e^(i k s))/s) with r' being the origin in the relative coordinate:
                 (Subscript[∇, r])((e^(i k |r-r'|))/(|r-r'|))=i k (((r -r'))/(Abs[r-r']))((e^(i k (Abs[r-r'])))/(Abs[r-r']))   (1 +(i/(k (Abs[r-r']))) ) 
Thus, consider the function |r-r'| and for simplicity, let y and z be zero. Then:
                                     |r-r'|=|x-x'| 
we see that it matters if we increase x or x'. Let say x> x'. If we increase x, the value increases, hence positive gradient. But if we increase x', the value decreases, hence the gradient is negative. This show that even the function has the same value, its derivative depends differently from each variable.

(continued - after exercise)

Using the result:         (Subscript[∇, r'])((e^(i k |r-r'|))/(|r-r'|))=- i k (OverHat[u])((e^(i k |r-r'|))/(|r-r'|))    (1 +(i/(k |r-r'|)) )                        (1.3.3)
where                                        (OverHat[u])=(((r -r'))/(Abs[r-r']))                                                           (1.3.3a)
We can write Eq. (12) as:
    U[r]=-(1/(4 π))(∮_S) (U[r'] i k (OverHat[u])((e^(i k |r-r'|))/(|r-r'|))  (1 +(i/(k |r-r'|)) ) + ((e^(i k |r-r'|))/(|r-r'|)) ∇U[r'])·d(Subscript[S, in])          (1.3.4a)
    U[r]=-(1/(4 π))(∮_S) ((e^(i k |r-r'|))/(|r-r'|))(U[r'] i k (OverHat[u])  (1 +(i/(k |r-r'|)) ) +  ∇U[r'])·d(Subscript[S, in])                    (1.3.4b)
We see that the term:
                                   (U[r'] i k (OverHat[u])  (1 +(i/(k |r-r'|)) ) +  ∇U[r'])                                  (1.3.4c)
 acts as a surface "source"; each point r' on the surface S contributes this value for a spherical wave ((e^(i k |r-r'|))/(|r-r'|)) to the observation point r.
A very common approximation is to drop the term (i/(k |r-r'|))=((i λ)/(2 π  |r-r'|))  since the observation point is usually several wavelengths away from the surface for farfield application.

1.4. Time-dependent solution (skip)

2. Kirchhoff diffraction theory applied to aperture

2.1. Diffraction by an aperture of a point source

A key application of Kirchhoff-Helmholtz diffraction theorem in early days is to solve the diffraction problem from a finite aperture. But we will see that it doesn’t have to be an aperture in the strict sense. It can be an object or structure.

An exact solution is generally difficult. But Kirchhoff's theorem allows reasonable approximation.

2.1.1

Consider the geometry as shown:

Consider the source at (Subscript[x, s]), the light is blocked by a screen with aperture A, what is the field at observation point (Subscript[x, o])? We apply Eq. (1.12) to the region bounded by the dashed line.
                   U[(Subscript[x, o])]=(1/(4 π))((∮_B) + (∮_A)+(∮_C) ) (U[r'] ∇((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) - ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) ∇U[r'])·d(Subscript[S, in])      (2.1.1)

Do we really know U[r']  for r' on those surface A, B, C? Strictly: not completely. Let B be the perfect absorber so that any radiation from the source falling on it is zero, how do we know that the field on the other side of the screen is zero? what is the field at the aperture? Is it the same as the field without the screen and aperture? Not necessarily. The screen affect the EM waves from the source and hence the field on A with the screen is NOT necessarily the same as the field from the source as if the screen did not exist.

Hence, we see that Eq. (2.1.1) is not as rigorous as it appears. However, we can use it as an approximation.
And the key word here is approximation. We can't use Kirchhoff-Hemholtz theorem to really solvev the exact EM solution for diffraction, but by choosing the right approach, we can obtain approximated solution.
More importantly however, the mathematics of this approach can give us a mental image of the diffraction process and significant insight about the phenomenon.

The approximation is summarized as follow:
1- The field on A is the same as if the screen does not exist (strictly speaking, not true, but reasonable with large aperture)
2- The field on B is zero, and so is ∇U[r’] (mostly true away from a short distance at the aperture boundary - aperture boundary near field)
3- For the field at C, we can reasonably make C infinitely large and the field may have a starting time such that C is at a distance further than the propagation of the light, and hence the field and its gradient vanish everywhere. (very reasonable assumption and true at the limit).

Hence, we have very simple result:
                     U[(Subscript[x, o])]=(1/(4 π)) (∮_A)(U[r'] ∇((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) - ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) ∇U[r'])· (Subscript[(OverHat[n]), in]) dS           (2.1.2)
where (Subscript[(OverHat[n]), in]) is the inward-pointing normal vector of aperture A as shown in the figure.
Let the source field be: A((e^(i k |r-(Subscript[x, s])|))/(|r-(Subscript[x, s])|)) where A is amplitude, not to be confused with the symbol for aperture, we can show that:
            U[(Subscript[x, o])]=-A(i/(2 λ )) (∮_A)((e^(i k (|(Subscript[x, o])-r'|+|r'-(Subscript[x, s])|)))/(|(Subscript[x, o])-r'| |r'-(Subscript[x, s])|)) (Cos[(Subscript[(OverHat[n]), in]),r'-(Subscript[x, s]) ]-Cos[(Subscript[(OverHat[n]), in]),r'-(Subscript[x, o])])  dS           (2.1.3)  
This is Fresnel-Kirchhoff diffraction formula.

Exercise: prove 2.1.3

Consider source field: ((e^(i k |r-(Subscript[x, s])|))/(|r-(Subscript[x, s])|)), then:  U[r']=((e^(i k |r'-(Subscript[x, s])|))/(|r'-(Subscript[x, s])|))  (we drop amplitude term A)
      (Subscript[∇, r'])U[r']=(Subscript[∇, r'])((e^(i k |r'-(Subscript[x, s])|))/(|r'-(Subscript[x, s])|))=i k (((r'-(Subscript[x, s]) ))/(|r'-(Subscript[x, s])|))((e^(i k |r'-(Subscript[x, s])|))/(|r'-(Subscript[x, s])|))   (1 +(i/(k |r'-(Subscript[x, s])|))) 
     (Subscript[∇, r'])((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))=-i k ((((Subscript[x, o])-r' ))/(|(Subscript[x, o])-r'|))((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) (1 +(i/(k |(Subscript[x, o])-r'|)))  
Thus: U[r'] (Subscript[∇, r'])((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) - ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) (Subscript[∇, r'])U[r']=
      -i k ((e^(i k |r'-(Subscript[x, s])|))/(|r'-(Subscript[x, s])|)) ((((Subscript[x, o])-r' ))/(|(Subscript[x, o])-r'|))((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) (1 +(i/(k |(Subscript[x, o])-r'|))) -((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))i k (((r'-(Subscript[x, s]) ))/(|r'-(Subscript[x, s])|))((e^(i k |r'-(Subscript[x, s])|))/(|r'-(Subscript[x, s])|))   (1 +(i/(k |r'-(Subscript[x, s])|)))  
      -i k ((e^(i k (|r'-(Subscript[x, s])|+|(Subscript[x, o])-r'|)))/(|r'-(Subscript[x, s])| |(Subscript[x, o])-r'|))( ((((Subscript[x, o])-r' ))/(|(Subscript[x, o])-r'|)) (1 +(i/(k |(Subscript[x, o])-r'|))) +(((r'-(Subscript[x, s]) ))/(|r'-(Subscript[x, s])|))   (1 +(i/(k |r'-(Subscript[x, s])|))))  
Since the terms (i/(k |(Subscript[x, o])-r'|)) and (i/(k |r'-(Subscript[x, s])|)) are negligible as both |(Subscript[x, o])-r'| and |r'-(Subscript[x, s])| are much larger than λ:
U[r'] (Subscript[∇, r'])((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) - ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) (Subscript[∇, r'])U[r']=
              -i k ((e^(i k (|r'-(Subscript[x, s])|+|(Subscript[x, o])-r'|)))/(|r'-(Subscript[x, s])| |(Subscript[x, o])-r'|))( ((((Subscript[x, o])-r' ))/(|(Subscript[x, o])-r'|))  +(((r'-(Subscript[x, s]) ))/(|r'-(Subscript[x, s])|)) )  
 and hence:
    U[(Subscript[x, o])]=-((i k)/(4 π)) (∮_A)((e^(i k (|r'-(Subscript[x, s])|+|(Subscript[x, o])-r'|)))/(|r'-(Subscript[x, s])| |(Subscript[x, o])-r'|))(((((Subscript[x, o])-r' ))/(|(Subscript[x, o])-r'|))  +(((r'-(Subscript[x, s]) ))/(|r'-(Subscript[x, s])|)))· (Subscript[(OverHat[n]), in]) dS 
    U[(Subscript[x, o])]=-((i k)/(4 π)) (∮_A)((e^(i k (|r'-(Subscript[x, s])|+|(Subscript[x, o])-r'|)))/(|r'-(Subscript[x, s])| |(Subscript[x, o])-r'|))(Cos[(Subscript[(OverHat[n]), in]),(Subscript[x, o])-r']  +Cos[(Subscript[(OverHat[n]), in]),r'-(Subscript[x, s]) ]) dS  
    U[(Subscript[x, o])]=-((i )/(2 λ)) (∮_A)((e^(i k (|r'-(Subscript[x, s])|+|(Subscript[x, o])-r'|)))/(|r'-(Subscript[x, s])| |(Subscript[x, o])-r'|))(Cos[(Subscript[(OverHat[n]), in]),r'-(Subscript[x, s])]-Cos[(Subscript[(OverHat[n]), in]),r'-(Subscript[x, o])]) dS

Exercise: Alternative geometry

An alternative expression is to consider geometry as shown:

Here, we can choose A to be a part of the source wavefront, then the integral becomes:
                        U[(Subscript[x, o])]=-A(i/(2 λ ))((e^(i k (Subscript[r, 0])))/(Subscript[r, 0])) (∮_A)((e^(i k |r'-(Subscript[x, s])|)))/(|r'-(Subscript[x, s])|)) (1+Cos[χ])  dS           (2.1.4)

2.2. Diffraction by an aperture of a plane wave

What if we have a plane wave instead of a point source?

Point source is somewhat an “anachronism” as scientists in 19th century never knew about practical plane-wave like laser. Light source is always thought of as a point, fundamentally, which is a δ-function, and a spatially extended light source is a integral convolution of it. But light field can also be constructed from plane waves as we see and hence, it is worthwhile to consider the K-H diffraction formula for plane wave.

Consider the following simple problem: a plane wave incident upon a perfectly conducting half-plane. What is the diffraction pattern? It may sound simple, but the exact solution, first obtained by Sommerfeld, is far from simple mathematically. (One of those few beautiful exact EM solutions).

We will consider the solution of this problem in a homework that compares the exact solution and the approximated solution. The key point here is that even for some simplest-looking problem, it is not easy to obtain exact solution. But one can use Kirchhoff approximation to obtain a reasonable result. While approximation can be useful, it is just as important to remember the condition for which an approximation is NOT valid. This depends on the specific problem and specific aspects of interest in a problem.

2.2.1 Plane wave diffraction

Let's consider a plane wave as shown above. Then, the result of applying the K-H theory of Eq. 
                     U[(Subscript[x, o])]=(1/(4 π)) (∮_A)(U[r'] ∇((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) - ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) ∇U[r'])· (Subscript[(OverHat[n]), in]) dS           (2.1.2)
with                      U[r']=A (e^(i (Subscript[k, 0]) . r')) 
gives us:                   
           U[(Subscript[x, o])]=-A((i k)/(4 π)) (∮_A) (e^(i (Subscript[k, 0]) · r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))(((i Cos[χ])/(k |(Subscript[x, o])-r'|))+Cos[χ] + Cos[α]) dS           (2.2.4)
where A is the plane wave amplitude, and Cos[χ]=(Subscript[(OverHat[n]), in]) . ((((Subscript[x, o])-r') )/(|(Subscript[x, o])-r'|))  (note the direction of vector  ((((Subscript[x, o])-r') )/(|(Subscript[x, o])-r'|)))

At this point, we may ask, what if the aperture is not just an opening, but some object with a transmission function T[r’].

Strictly speaking, if there is such an object, then the full solution must be obtain for the EM wave that include the part that scatters or reflects backward, the part inside the objects, and the part that radiates away in the FF. It is a very complex problem.

However, in many cases such as imaging application, as we can assume that the EM field after the object is simply modified by a transmission function T[r’] that includes both amplitude and phase modification and the output field is:
                                         T[r'] (e^(i (Subscript[k, 0]) · r'))      (2.2.5)

(We disregard whatever happens inside the object, which, for example, can be a thin film with internal structure, and care only for the transmitted field. Thus T[r’] is already the solution of the EM field at the object output - which can be solved accurately as a separate problem with any method).

Then, (2.2.4) is simplify modified to be:
           U[(Subscript[x, o])]=-A((i k)/(4 π)) (∮_A)T[r'] (e^(i (Subscript[k, 0]) · r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))(((i Cos[χ])/(k |(Subscript[x, o])-r'|))+Cos[χ] + Cos[α]) dS           (2.2.6)
An approximation for farfield is that ((i Cos[χ])/(k |(Subscript[x, o])-r'|)) can be neglected compared to Cos[χ], as k|(Subscript[x, o])-r'|=2π ((|(Subscript[x, o])-r'|)/λ)>>1, and thus: 
           U[(Subscript[x, o])]=-A((i k)/(4 π)) (∮_A) T[r'](e^(i (Subscript[k, 0]) · r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))(Cos[χ] + Cos[α]) dS           (2.2.7a)
Sometimes, (Cos[χ] + Cos[α]) is taken outside, and for paraxial approximation, both of these angles are small and  Cos[χ] + Cos[α]~ 2. (Note that this approximation is NOT crucial, and we can include a small variation of Cos[χ] in the integral, while Cos[α] is a constant - it is just an overall amplitude term).
                               U[(Subscript[x, o])]~A(( k  (Cos[χ]+Cos[α]))/(i 4 π)) (∮_A) T[r'](e^(i (Subscript[k, 0]) ·r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) dS          (2.2.7b)
Or:                          U[(Subscript[x, o])]~A(( k)/(i 2 π)) (∮_A) T[r'] (e^(i (Subscript[k, 0]) ·r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)) dS                            (2.2.7c)
This is a commonly used expression, and the most obvious is that it is a Fourier transform of function T[r'], modified with an extra phase term from ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|)), which is about the phase and amplitude from the object to the observation point.

Review previous HW about Fourier transform of Far Field.

2.3 Fresnel and Fraunhofer approximation

Even expression (2.1.3) for point source or  (2.2.7b or c) for plane wave is still a bit complicated to use directly as a formula. Here, we study 2 other approximations to obtain a more practical formula.

2.3.1 Point source

First, definition for the geometry

In aperture plane

Depending on the various limiting range in the figure above, we divide into three main regimes of approximations:

1- When all the terms must be kept, we call this Kirchhoff diffraction, which is just beyond what we consider as nearfield. (KH theory is not valid for nearfield. This is where numerical method must be used for accurate solutions).

2- Fresnel diffraction: occurs when the quadratic terms (((r')^2)/2)((1/(Subscript[r, s]))+(1/(Subscript[r, o]))) is significant, but all higher order terms can be ignored. Thus:
            U[(Subscript[x, o])]=-A(i/(2 λ )) (∮_A)((e^(i k (|(Subscript[x, o])-r'|+|r'-(Subscript[x, s])|)))/(|(Subscript[x, o])-r'| |r'-(Subscript[x, s])|)) (Cos[(OverHat[n]),r'-(Subscript[x, s]) ]-Cos[(OverHat[n]),r'-(Subscript[x, o])])  dS           (2.1.3)  
             (e^(i k (|(Subscript[x, o])-r'|+|r'-(Subscript[x, s])|)))=(e^(i k  ((Subscript[r, o])+(Subscript[r, s]))))(e^(i k f[r'])) ~ (e^(i k  ((Subscript[r, o])+(Subscript[r, s]))))  (e^(i k (-r'·((( (Subscript[ρ, o]))/(Subscript[r, o]))+((Subscript[ρ, s])/(Subscript[r, s])))+(((r')^2)/2)((1/(Subscript[r, s]))+(1/(Subscript[r, o]))) )))             (2.3.11)
 Let the source be at infinity (approx incident plane wave):
             (e^(i k (|(Subscript[x, o])-r'|+|r'-(Subscript[x, s])|))) ~(e^(i k  (Subscript[r, o])))  (e^(i k (-r'·((( (Subscript[ρ, o]))/(Subscript[r, o]))+ (Subscript[(OverHat[ρ]), s]))+(((r')^2)/(2 (Subscript[r, o]))) )))                                                  (2.3.12a)
 where               (Subscript[(OverHat[ρ]), s])=((Subscript[ρ, s])/(Subscript[r, s]))  and we drop the meaningless phase term (e^(i k  (Subscript[r, s])))                    (2.3.12b) 
 
 Note that we can get this same result for plane wave incident by using (2.2.7a) (without the transmission term T[r’]):             
             U[(Subscript[x, o])]=-A((i k)/(4 π)) (∮_A) (e^(i (Subscript[k, 0]) · r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))(Cos[χ] + Cos[α]) dS           (2.2.7a)
 The term (e^(-i k r'· (Subscript[(OverHat[ρ]), s]))) is the same as (e^(i (Subscript[k, 0]) · r')), and the rest is an approximation for (e^(i k |(Subscript[x, o])-r'|)).
 
3- Fresnel paraxial approximation:
  
There is one subcategory of Fresnel diffraction is paraxial approximation. We notice an alternate expression:
             |(Subscript[x, o])-r'|=(Sqrt[((((Subscript[ρ, o])-r'))^2) +((Subscript[z, o])^2)])                                              (2.3.13a)
             |(Subscript[x, o])-r'|~ (Subscript[z, o]) (1+(1/2)(((((Subscript[ρ, o])-r'))^2)/((Subscript[z, o])^2)) ) =(Subscript[z, o])+(1/(2 (Subscript[z, o]) ))((((Subscript[ρ, o])-r'))^2)            (2.3.13b)
However, this is valid only for:  (Subscript[z, o])>>|(Subscript[ρ, o])-r'|                                (2.3.13c) 
which means it has to be truly paraxial approximation as the observation angle is small and parallel to the axis. With large NA observation, this is NOT necessarily valid.

On the other hand, in the other more general approximation,  ((r' ·(Subscript[ρ, o]))/((Subscript[r, o])^2))  in Eq. (2.13 a) is valid as long as the aperture dimension is << (Subscript[r, o]). 

4- Fraunhofer diffraction when the quadratic and all higher order terms can be ignored and only -r'·((( (Subscript[ρ, o]))/(Subscript[r, o]))+((Subscript[ρ, s])/(Subscript[r, s]))) remains. This is the case of truly far field: beyond the Rayleigh range of the aperture to infinity.
Again, for source at infinity:
               (e^(i k (|(Subscript[x, o])-r'|+|r'-(Subscript[x, s])|))) ~(e^(i k  ((Subscript[r, o])+(Subscript[r, s]))))  (e^(i k (-r'·((( (Subscript[ρ, o]))/(Subscript[r, o]))+ (Subscript[(OverHat[ρ]), s])) )))                             (2.3.14)

2.3.3 Approximation for incident plane wave - Fresnel with paraxial

Note that if the incident angle α for incident plane wave is not small, the axis for Fresnel paraxial approx should be along the plane wave axis to account for the angular tilt. There are 2 configurations shown above. The result is the same, except that we don’t carry an explicit plane wave phase term for the second.

From the above:
                  U[(Subscript[x, o])]=-A((i k)/(4 π)) (∮_A) (e^(i (Subscript[k, 0]) · r')) ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))(Cos[χ] + Cos[α]) dS           (2.2.7a)
Then, as usual (see Eq. 2.3.11) for paraxial approx:
                   | (Subscript[x, o])-r'|=(Sqrt[((Subscript[x, o])^2)-2 (Subscript[x, o]) . r'+((r')^2)])=|(Subscript[x, o])| (Sqrt[1-2 ((Subscript[x, o])/((Subscript[x, o])^2)) . r'+(((r')^2)/((Subscript[x, o])^2))])   (2.3.15a)
                                ≈|(Subscript[x, o])| -((Subscript[x, o])/(|(Subscript[x, o])|)) . r' +(((r')^2)/(2 |(Subscript[x, o])|)) 
Substitute                      (Subscript[x, o])=(Subscript[ζ, 0]) +(Subscript[ρ', o]) 
                   | (Subscript[x, o])-r'|≈|(Subscript[x, o])| -(1/(|(Subscript[x, o])|)) ((Subscript[ζ, 0]) +(Subscript[ρ', o])) . r' +(((r')^2)/(2 |(Subscript[x, o])|))                                
                                 = |(Subscript[x, o])| -((Subscript[ζ, 0])/(|(Subscript[x, o])|)) (Subscript[(OverHat[k]), 0 ]). r'-(1/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r' +(((r')^2)/(2 |(Subscript[x, o])|))                  (2.3.15b)

Hence:           ((e^(i k |(Subscript[x, o])-r'|))/(|(Subscript[x, o])-r'|))=((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|)) (e^(-i ((Subscript[ζ, 0])/(|(Subscript[x, o])|)) (Subscript[k, 0]) . r'))  (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r')) (e^(i (k/(2 |(Subscript[x, o])|))((r')^2)))                    (2.3.15c)
Substitute (2.3.15c) into (2.2.7a):
        U[(Subscript[x, o])]=-A((i k)/(4 π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))(Cos[χ]+Cos[α])(∮_A) (e^(i (Subscript[k, 0]) · r))' (e^(-i ((Subscript[ζ, 0])/(|(Subscript[x, o])|)) (Subscript[k, 0]) . r')) (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r'))(e^(i (k/(2 |(Subscript[x, o])|))((r')^2))) dS  
        U[(Subscript[x, o])]=-A((i k)/(4 π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))(Cos[χ]+Cos[α])(∮_A) (e^(i (1-((Subscript[ζ, 0])/(|(Subscript[x, o])|)))(Subscript[k, 0]) · r')) (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r')) (e^(i (k/(2 |(Subscript[x, o])|))((r')^2))) dS       (2.3.16)
With numerical method, (2.3.16) can actually be computed directly, but further approximation is also useful and accurate with appropriate paraxial condition.
We notice that we can make an approximation for small angle of observation:
                             (1-((Subscript[ζ, 0])/(|(Subscript[x, o])|))) =1-Cos[θ] ~((θ^2)/2)~(((Subscript[ρ', o])^2)/(2 ((Subscript[x, o])^2)))=0          (2.3.17)
then the phase  term (e^(i (Subscript[k, 0]) · r))' at the aperture is canceled as it should, since it doesn’t really add any thing of the aperture. Hence:        
        U[(Subscript[x, o])]=-A((i k)/(4 π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))(Cos[χ]+Cos[α])(∮_A) (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r')) (e^(i (k/(2 |(Subscript[x, o])|))((r')^2))) dS    (2.3.18a)
        U[(Subscript[x, o])]=-A((i k)/(2 π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r')) (e^(i (k/(2 |(Subscript[x, o])|))((r')^2))) dS    (2.3.18b)
by letting Cos[χ] ~ Cos[α].

At this point we can also stop making any distinction between notation (Subscript[ρ, 0]) and (Subscript[ρ', 0]). It is understood heretofore that for an incident plane wave with a substantial angle with respect to the plane aperture normal axis, the reference axis of observation should be along the plane wave direction.
However, note that if we use the tilted plane coordinate as shown below, dot product (Subscript[ρ', o]). r' corresponds to ((Subscript[(OverTilde[ρ]), o])+(Subscript[ρ', o//]) ) . r' which is (Subscript[(OverTilde[ρ]), o]) . r'. They differs by Cos[α] factor in one dimension, but are equal in the orthogonal dimension, which accounts for the fact that the aperture appears to have different effective cross section in different dimensions: a simple tilted view effect. Object appears “squashed” in the tilted view dimension.

If the aperture has a geometry that sipulates its own axes, then the beam angle cannot be arbitrarlily assigned by any choice of coordinate. One can write  (Subscript[ρ', o]). r'=x'(Subscript[ξ, o]) Cos[(Subscript[α, x])]+y'(Subscript[η, o]) Cos[(Subscript[α, y])]  where Cos[(Subscript[α, x])] and Cos[(Subscript[α, y])] are the cosine of angles between unit vectors (OverHat[x]),  (OverHat[y]) and (OverHat[ξ]),  (OverHat[η]) of the two systems of coordinate. Clearly, this choice depends on the specific geometry of a problem.

Now we can put back the aperture transmission function like that of an object in Eq. (2.3.18b):
        U[(Subscript[x, o])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) (T[r'] (e^(i (k/(2 |(Subscript[x, o])|))((r')^2))) )(e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ', o]). r'))  dS    (2.3.19a)       
 we notice that this is a Fourier transform:

U[(Subscript[x, o])]~-A i k ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|)) Cos[α] (ℱ^(-1))[T[r']   (e^(i (k/(2 |(Subscript[x, o])|))((r')^2))),r',(k/(|(Subscript[x, o])|))(Subscript[ρ, 0])]            (2.3.19b)
We have encountered this before in HW (HW 7) of waveguide laser FF output, this is no different, except that T[r’] is the E field.
            E[x,z;t]=(e^(i (k z -ω t)))(e^(i k((x^2)/(2z)) ))(Sqrt[(k/(2 π i z))])(∫_(-∞)^∞)dξ (Subscript[E, 0])[ξ] (e^(-i k ξ (x/z) ))(e^(i k  ((ξ^2)/(2z)) ))   (2.3.20a)
Expressed as a Fourier transform:
              (Sqrt[(1/(2 π))])(∫_(-∞)^∞)dξ ((Subscript[E, 0])[ξ]  (e^(i k  ((ξ^2)/(2z)) )) )(e^(-i k ξ (x/z) ))=(ℱ^(-1))[(Subscript[E, 0])[ξ](e^(i k  ((ξ^2)/(2z)) )),ξ,k(x/z)]   (2.3.20b)
(2.3.19) and (2.3.20) are equivalent with the approximation of |(Subscript[x, o])| for z, the extension to 2 D instead of just x-dimension.

A further approximation for very small paraxial angle is to replace |(Subscript[x, o])| for (Subscript[ζ, o]). In such a case, it can be shown that:
       U[(Subscript[x, o])]~-A((i k)/(2 π)) ((e^(i k (Subscript[ζ, o])))/(Subscript[ζ, o]))(e^(i (k/(2 (Subscript[ζ, o])))((Subscript[ρ, 0])^2)))(∮_A)   (T[r']   (e^(i (k/(2 (Subscript[ζ, o])))((r')^2)))) (e^(-i (k/(Subscript[ζ, o]))(Subscript[ρ, 0]).r'))dS        (2.3.21a) 
       U[(Subscript[x, o])]~-A i k ((e^(i k (Subscript[ζ, o])))/(Subscript[ζ, o]))(e^(i (k/(2 (Subscript[ζ, o])))((Subscript[ρ, 0])^2))) (ℱ^(-1))[T[r']   (e^(i (k/(2 (Subscript[ζ, o])))((r')^2))),r',(k/(Subscript[ζ, o]))(Subscript[ρ, 0])]             (2.3.21b)
This expression (2.3.21b) is often used in discussion of lens imaging system, as the quadratic phase term (e^(i (k/(2 (Subscript[ζ, o])))((Subscript[ρ, 0])^2))) or (e^(i (k/(2 (Subscript[ζ, o])))((r')^2))) are very relevant to any quadratic phase transformation such as the effect of a lens.

2.3.4 Alternative expression for Fresnel diffraction with approximately normal incident

U[(Subscript[x, o])]=k ((e^(i k (Subscript[z, o])))/(i 2 π (Subscript[z, o])))(∮_A)  (e^(i (k/(2 (Subscript[z, o]) ))((((Subscript[ρ, o])-r'))^2)  )) dS                               (2.3.23)

2.3.5 Fraunhofer with plane wave

At very large distance, even the quadratic term can be dropped, we have the farfield condition and obtain the Fraunhofer approximation formula:
               U[(Subscript[x, o])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) T[r']  (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ, o]). r'))  dS        (2.3.22a) 
Converting coordinate ratios ((Subscript[ρ, 0])/(|(Subscript[x, o])|)) into angle for Cartesian coordinate:

U[(Subscript[x, o])]~-A i k ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](ℱ^(-1))[T[r'],r',k{((Tan[(Subscript[θ, x])] Cos[(Subscript[β, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)])),((Tan[(Subscript[θ, y])] Cos[(Subscript[β, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))}]         (2.3.22b) 
We see this as just a Fourier transform of the aperture transmission function.
If we use the aperture coordinate, then:
               U[(Subscript[x, o])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) T[r']  (e^(-i(k/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s) r'))  dS        (2.3.22a) 
where s is the translation vector corresponding to the incident beam axis intersecting with the observation plane. 
                   (1/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s)={((Tan[(Subscript[θ, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-(((Subscript[z, o]) Tan[(Subscript[α, x])])/(|(Subscript[x, o])|)),((Tan[(Subscript[θ, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-(((Subscript[z, o]) Tan[(Subscript[α, y])])/(|(Subscript[x, o])|))} 
                   (1/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s)={((Tan[(Subscript[θ, x])]-Tan[(Subscript[α, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, x])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)])),
                         ((Tan[(Subscript[θ, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, y])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)]))}         (2.3.23)
Note: it is not imperative to use Cartesian. For another problem, one can use spherical or cylindrical as appropriate. See more in HW about coordinate: if we choose the projection plane perpendicular to the aperture axis (parallel with the aperture plane), we have integration variable terms as below. Note that (Subscript[k, 0,∟]) generally does not vanish.

However, a more general point of view is that the diffractive object needs not be an aperture with a specific plane. The natural axis in the farfield is the main direction of the diffracted plane wave.

3. Additional discussion and illustration of Fresnel diffraction

3.1 Concept of free space transfer function and impulse response

3.1.1 Review of free space propagation

We recall in the discussion of basic wave progation that given a scalar field profile at z=0:
                                     U[x,y,z=0]=∫a[(Subscript[k, x]),(Subscript[k, y]) ] (e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) d(Subscript[k, x]) d(Subscript[k, y])                      (3.1.1)
we can construct the wave at arbitrary z:
                                     U[x,y,z]=∫a[(Subscript[k, x]),(Subscript[k, y]) ] (e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) (e^(i (Subscript[k, z]) z))d(Subscript[k, x]) d(Subscript[k, y])                      (3.1.2a)
where:                                         (Subscript[k, z])=(Sqrt[((Subscript[k, 0])^2)-(((Subscript[k, x])^2)+((Subscript[k, y])^2))])                                         (3.1.2b)
For propagating wave, we restrict the Fourier integration to:
                                                      ((Subscript[k, T])^2)≡((Subscript[k, x])^2)+((Subscript[k, y])^2)≤ ((Subscript[k, 0])^2)
Hence:                          U[x,y,z]=(∫_((Subscript[k, T])≤ (Subscript[k, 0])))a[(Subscript[k, x]),(Subscript[k, y]) ] (e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) (e^(i (Subscript[k, z]) z))d(Subscript[k, x]) d(Subscript[k, y])        (3.1.2c)
Alternatively:         U[x,y,z]=∫a[(Subscript[k, x]),(Subscript[k, y]) ] Θ[(Subscript[k, 0])-(Subscript[k, T])](e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) (e^(i (Subscript[k, z]) z))d(Subscript[k, x]) d(Subscript[k, y])        (3.1.2d)
where the unitstep function (Heaviside) Θ[(Subscript[k, 0])-(Subscript[k, T])] indicates the restriction to the circle (Subscript[k, T])≤ (Subscript[k, 0]).

If we look at Eq. (3.1.1) and (3.1.2c), and treat the effect of propagation in z as a linear system output, we can write:
                         U[x,y,z=0]=∫A[(Subscript[k, x]),(Subscript[k, y]) ;0] (e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) d(Subscript[k, x]) d(Subscript[k, y])         (3.1.3a)
                         U[x,y,z]=∫A[(Subscript[k, x]),(Subscript[k, y]) ;z] (e^(i ((Subscript[k, x]) x+(Subscript[k, y]) y ))) d(Subscript[k, x]) d(Subscript[k, y])                (3.1.3b)
then, the transfer function (or spatial frequency response function) connecting between A[(Subscript[k, x]),(Subscript[k, y]) ;z] and A[(Subscript[k, x]),(Subscript[k, y]) ;0] is:
                                  A[(Subscript[k, x]),(Subscript[k, y]) ;z]=H[(Subscript[k, x]),(Subscript[k, y]);z]A[(Subscript[k, x]),(Subscript[k, y]) ;0]                   (3.1.4)
where:                  H[(Subscript[k, x]),(Subscript[k, y]);z]=  (e^(i z(Sqrt[((Subscript[k, 0])^2)-(((Subscript[k, x])^2)+((Subscript[k, y])^2))]))) Θ[(Subscript[k, 0])-(Sqrt[((Subscript[k, x])^2)+((Subscript[k, y])^2)])]   (3.1.4a)           
Then, we can ask, how is this related to Fresnel diffraction, as it is also capable to connect the field at z=0 to arbitrary z?

3.1.2 Comparison of free space propagation to Fresnel diffraction

To see the connection, let's look at the Fresnel integral above for a normal plane wave incident on an aperture: 
                    U[(Subscript[x, o])]=k ((e^(i k (Subscript[z, o])))/(i 2 π (Subscript[z, o])))(∮_A)  (e^(i (k/(2 (Subscript[z, o]) ))((((Subscript[ρ, o])-r'))^2)  )) dS                               (2.28.c)    
Instead of aperture, suppose there is a mask at the aperture with certain transmission function such that the field at the aperture is U[x, y], then, instead of (2.28.c):
             U[(Subscript[x, o])]=k ((e^(i k (Subscript[z, o])))/(i 2 π (Subscript[z, o])))(∮_A) U[x',y'] (e^(i (k/(2 (Subscript[z, o]) ))((((Subscript[ρ, o])-r'))^2)  )) dS                             (3.1.5)                    
We drop the subscript "o" for observation point, and write explicitly:           
            U[x,y,z]=k ((e^(i k z))/(i 2 π z))∫ U[x',y'] (e^(i (k/(2 z)) ((((x-x'))^2)+(((y-y'))^2))  )) dx' dy'            (3.1.6)          
and let the integration extends to infinity with the assumption that U[x', y'] vanishes sufficiently fast at infinity so that the integral exists.
Then the above expression is clearly a convolution:      
                   U[x,y,z]=∫ U[x',y'] h[x-x',y-y';z] dx' dy'            (3.1.6)    
where:                h[x - x', y - y';z]=k ((e^(i k z))/(i 2 π z)) (e^(i (k/(2 z)) ((((x-x'))^2)+(((y-y'))^2))  ))         (3.1.6a)
is the impulse response function
And therefore, if we express in terms of spatial frequency:
                 A[(Subscript[k, x]),(Subscript[k, y]);z]= A[(Subscript[k, x]),(Subscript[k, y]);0]  (Subscript[H, F])[(Subscript[k, x]),(Subscript[k, y]);z]              (3.1.7)  
then:             (Subscript[H, F])[(Subscript[k, x]),(Subscript[k, y]);z]=∫h[x , y ;z] (e^(-i ((Subscript[k, x])x +(Subscript[k, y])y ) ))dx dy           (3.1.7a)
                (Subscript[H, F])[(Subscript[k, x]),(Subscript[k, y]);z]=k ((e^(i k z))/(i 2 π z))∫ (e^(i (k/(2 z)) ((x^2)+(y^2))  )) (e^(-i ((Subscript[k, x])x +(Subscript[k, y])y ) ))dx dy            (3.1.7a)
And:         (Subscript[H, F])[(Subscript[k, x]),(Subscript[k, y]);z]=(e^(i k z)) (e^(-i (z/(2 k)) (((Subscript[k, x])^2)+((Subscript[k, y])^2))))           (3.1.7b)
If we compare (3.1.7b) to (3.1.4a) we expect that:
      (e^(i z(Sqrt[(k^2)-(((Subscript[k, x])^2)+((Subscript[k, y])^2))]))) Θ[k-(Sqrt[((Subscript[k, x])^2)+((Subscript[k, y])^2)])](Underscript[->, Fresnel approximation]) (e^(i k z)) (e^(-i (z/(2 k)) (((Subscript[k, x])^2)+((Subscript[k, y])^2))))   
Indeed, we see that:
                 (e^(i z(Sqrt[(k^2)-(((Subscript[k, x])^2)+((Subscript[k, y])^2))])))->(e^(i k z(Sqrt[1-(((((Subscript[k, x])^2)+((Subscript[k, y])^2)))/(k^2))]))) ~(e^(i k z))(e^(-i ((k z)/(2 (k^2) ))(((Subscript[k, x])^2)+((Subscript[k, y])^2))))  
                                                    =(e^(i k z))(e^(-i (z/(2 k ))(((Subscript[k, x])^2)+((Subscript[k, y])^2)))) 
Hence, Fresnel diffraction transfer function is the paraxial approximation to the propagation transfer function. Notice that the function Θ[k-(Sqrt[((Subscript[k, x])^2)+((Subscript[k, y])^2)])] is needed with (e^(i z(Sqrt[(k^2)-(((Subscript[k, x])^2)+((Subscript[k, y])^2))]))) to avoid evanescent wave, but for the paraxial approximation implicitly assumes Θ[k-(Sqrt[((Subscript[k, x])^2)+((Subscript[k, y])^2)])], as it is valid only for ((Subscript[k, x])^2)+((Subscript[k, y])^2)<< (k^2).

In fact, we will see that (Subscript[H, F])[(Subscript[k, x]),(Subscript[k, y]);z]=(e^(i k z)) (e^(-i (z/(2 k)) (((Subscript[k, x])^2)+((Subscript[k, y])^2)))) Θ[k-(Sqrt[((Subscript[k, x])^2)+((Subscript[k, y])^2)])]  can be useful.
Note the integration:

(∫_0^(2 π ))(E^(-I k r Cos[ϕ]))\[DifferentialD]ϕ

(∫_0^∞)(BesselJ[0, k r])(E^(I s (r^2)  ))r \[DifferentialD]r

Hence:               ∫ (e^(i (k/(2 z)) ((x^2)+(y^2))  )) (e^(-i ((Subscript[k, x])x +(Subscript[k, y])y ) ))dx dy =((i 2 π z)/k)(e^(-i((  z  )/(2 k))(((Subscript[k, x])^2)+((Subscript[k, y])^2))))

3.2 Paraxial phase front revisited

We see that the phase front of a spherical wave is... well spherical!

y=0 ;
ContourPlot[Arg[(E^(I 2 π (Sqrt[x^2+y^2+z^2]))) ]==0,{z,-10,10}, {x,-10,10}]

Now we have very very large z, and we look in the z x plane, we compare 2 wavefronts: spherical and paraxial:

y=0. ;
ContourPlot[{Arg[(E^(I 2 π (Sqrt[x^2+y^2+z^2]))) ]==0, Arg[(E^(I 2 π  (z+((x^2)/(2 z))))) ]==0}
,{z,40,60}, {x,-10,10},ContourStyle-> {{Thick, Red},{Dashed, Thick, Blue}}]

We see that the paraxial approximation works very well. Of course, we don't need contourplot to even plot it out. We know:
                     (z+((x^2)/(2 z)))== m

z= .; x=. ;
Solve[(z+((x^2)/(2 z)))== m ,{x}]

Plot[Flatten[Table[{- (Sqrt[(m^2) -(z^2)]),  (Sqrt[(m^2) -(z^2)]), -(Sqrt[2]) (Sqrt[m z-(z^2)]),(Sqrt[2]) (Sqrt[m z-(z^2)])} ,{m,40,60,2}]], {z,40,60}, PlotRange-> {{40,60},{-40,40}},AspectRatio-> 4]

We will see that even for larger angle where paraxial approx does not appear to be correct, it is actually still fine if equivalent terms are used.

4. Some Fraunhofer illustration

See Fresnel diffraction in HW

4.1 Illustration: Fraunhofer - planewave with rectangle aperture

For this problem, we will use the above result:
               U[(Subscript[x, o])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) T[r']  (e^(-i(k/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s) . r'))  dS        (2.3.22a) 
where s is the translation vector corresponding to the incident beam axis intersecting with the observation plane. 
                   (1/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s)={((Tan[(Subscript[θ, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-(((Subscript[z, o]) Tan[(Subscript[α, x])])/(|(Subscript[x, o])|)),((Tan[(Subscript[θ, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-(((Subscript[z, o]) Tan[(Subscript[α, y])])/(|(Subscript[x, o])|))} 
                   (1/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s)={((Tan[(Subscript[θ, x])]-Tan[(Subscript[α, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, x])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)]))
                      ,((Tan[(Subscript[θ, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, y])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)]))}         (2.3.23)

For a plane wave making an angle to the screen, calculate diffraction pattern for a rectangle 2 a and 2 b

Here, we will use the projection plane to be parallel with the aperture.
      U[(OverHat[Ω])]= -A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A)  (e^(-i(k/(|(Subscript[x, o])|)) ((Subscript[ρ, o])-s) . r'))  dS 
             =-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A) (e^(-i k (x (ξ-p)+y (η-q))))  dS    (4.1.1)
where we define variables: 
       ξ-p=((Tan[(Subscript[θ, x])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, x])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)]))  ;     (4.1.2)
           η-q=((Tan[(Subscript[θ, y])])/(Sqrt[1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)]))-((Tan[(Subscript[α, y])])/(Sqrt[1+((Tan[(Subscript[α, x])])^2)+((Tan[(Subscript[α, y])])^2)]))

(∫_(-a)^a)(E^(-I k x (ξ-p)))\[DifferentialD]x

Hence:         U[(OverHat[Ω])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))((sin(2 π a (p-ξ)))/(π (p-ξ)))((sin(2 π b (q-η)))/(π (q-η)))         (4.1.3)
We drop ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|)) in the limit |(Subscript[x, o])|-> ∞

Source code: Fraunhofer of rectangular aperture

Block[{dat3D,datdensity,D3PLOT, DENSPLOT,plot3d,aperture,plotdens,paraxwindow},
Manipulate[
If[INIT== 1,
{off,on}={Button["OFF",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Red]
,Button["ON",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Hue[0.36,1,0.4]]} ;

FHRect[{θ_,ϕ_},Rect_,λ_]:=Block[{Tθx,Tθy,u,v, k,a,b},
Tθx=Tan[θ] Cos[ϕ] ;Tθy=Tan[θ] Sin[ϕ] ;
{a,b}=Rect ; k=2.π/λ ;
u=Tθx/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
v=Tθy/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
((Sin[a k u])/(k  u)) ((Sin[b k v])/(k v)) ]; 
FHRect2[{θ_,ϕ_},Rect_,incid_,λ_]:=Block[{Tθx,Tθy,u,v,αx,αy, p,q, k,a,b},
Tθx=Tan[θ] Cos[ϕ] ;Tθy=Tan[θ] Sin[ϕ] ;
{a,b}=Rect ; k=2.π/λ ; {αx,αy}=incid  ;
u=Tθx/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
v=Tθy/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
p=Tan[αx]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
q=Tan[αy]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
((Sin[a k (u-p)])/(k  (u-p))) ((Sin[b k (v-q)])/(k (v-q))) ]; 
ContrastF[x_,a_,b_]:=(1./(1.+Exp[-(x-a)/b])) ;
 
INIT=0 
,TextCell["Not ready"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]];
paraxwindow=Graphics[{Hue[0.,1,1.,0.],EdgeForm[{Red,Dashed}]
,Rectangle[{-0.5,-0.5},{0.5,0.5}]}] ;
Manipulate[
If[RUN== 1,
λ=1. ;k=2. π/λ   ;
Manipulate[
aperture=Graphics[{Black,Rectangle[{-10.,-10.},{10.,10.}]
,White,Rectangle[-rect,rect]}] ;
If[TrueQ[D3PLOT==on] , 
dat3D=Table[Chop[FHRect2[{(Sqrt[θRx^2+((θRy)^2)]),Arg[θRx+I θRy]},rect,incd,λ]], {θRy,-1.001,1.,2./60},{θRx,-1.001,1.,2./60}] ;
plot3d=ListPlot3D[dat3D, PlotRange-> All,Mesh-> False,BoxRatios-> {1.,1.,0.5}
,DataRange->{{-1.,1.},{-1.,1.}}],];

{(* Plot3D[ Chop[FHRect2[{(Sqrt[θRx^2+((θRy)^2)]),Arg[θRx+I θRy]},rect,incd,λ]],{θRx,-1,1}, {θRy,-1,1}, PlotRange-> All, Mesh-> False
,PlotPoints-> {40,40}] *)
aperture
,Show[plot3d] 
,
If[TrueQ[DENSPLOT==on] ,
Manipulate[
datdensity=Table[ContrastF[(FHRect2[{(Sqrt[θRx^2+((θRy)^2)]),Arg[θRx+I θRy]},rect,incd,λ])^2,a,b], {θRy,-1.001,1.,2./80.},{θRx,-1.001,1,2./80.}];
plotdens=ListDensityPlot[datdensity, PlotRange-> All,ColorFunction->GrayLevel
,DataRange->{{-1.,1.},{-1.,1.}}] ;
Show[plotdens,paraxwindow]
,{a,0.01,10},{b,1,10}
(*,SynchronousUpdating->False *)],
Show[plotdens,paraxwindow] ] 
}

,Grid[{{Control[{{rect,{2.5,4},"rectangle dimension in λ"},{0.5,0.5},{10,10}}]
,Control[{{incd,{0.,0.},"incident angle"},{-1.,-1},{1.,1.}}]}
,{Control[{{D3PLOT,on,Button["3D Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.12,1,1,1]]},{on,off}
,ControlType-> Toggler}]
,Control[{{DENSPLOT,on,Button["Density Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.52,1,1,1]]},{on,off}
,ControlType-> Toggler}]
}
}]
(*,SynchronousUpdating->False *)] 

,TextCell["Not running"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]]
,Row[{Button["Toggle box to start, and turn off when not running",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Black],Control[{{RUN,0,"->"},{0,1}}]
}]
]
,Button["FRAUNHOFER DIFFRACTION OF RECTANGULAR APERTURE",BaseStyle->{24,Black, Bold,FontFamily-> "Calibri"}
,Background->Hue[0.4,1,1]]
,Row[{Button["Click once to initialize",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Red],Control[{{INIT,0,"->"},{0,1}}]
}]
]
]

Demo: Fraunhofer of rectangular aperture - view from aperture plane axis

Note some distortion from just geometrical linear projection. (see discussion about paraxial approximation and coordinates)

HW: Do the same as above but with the viewpoint parallel with the incident plane wave axis.

Source code: Fraunhofer of rectangular aperture for the view along beam axis

Manipulate[
If[INIT== 1,
FHRect[θx_, θy_,Rect_,λ_]:=Block[{u,v, k,a,b},
{a,b}=Rect ; k=2.π/λ ;
u=Tan[θx]/(Sqrt[1.+((Tan[θx])^2)+((Tan[θy])^2)]) ;
v=Tan[θy]/(Sqrt[1.+((Tan[θx])^2)+((Tan[θy])^2)]) ;
((Sin[a k u])/(k  u)) ((Sin[b k v])/(k v)) ]; 
FHRect2[θx_, θy_,Rect_,incid_,λ_]:=Block[{u,v,αx,αy, p,q, k,a,b},
{a,b}=Rect ; k=2.π/λ ; {αx,αy}=incid  ;
u=Tan[θx]/(Sqrt[1.+((Tan[θx])^2)+((Tan[θy])^2)]) ;
v=Tan[θy]/(Sqrt[1.+((Tan[θx])^2)+((Tan[θy])^2)]) ;
p=Tan[αx]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
q=Tan[αy]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
((Sin[a k (u-p)])/(k  (u-p))) ((Sin[b k (v-q)])/(k (v-q))) ]; 
FHRect3[{θRx_, θRy_},Rect_,incid_,λ_]:=Block[{u,v,αx,αy, p,q, k,a,b,tx,ty},
tx=((Tan[(Sqrt[θRx^2+ θRy^2])]  θRx)/(Sqrt[θRx^2+ θRy^2]))  ; ty= ((Tan[(Sqrt[θRx^2+ θRy^2])] θRy)/(Sqrt[θRx^2+ θRy^2]))   ;
{a,b}=Rect ; k=2.π/λ ; {αx,αy}=incid  ;


u=tx*Cos[αx]/(Sqrt[1.+(tx^2)+(ty^2)]) ;
v=ty*Cos[αy]/(Sqrt[1.+(tx^2)+(ty^2)]) ;
((Sin[a k u])/(k  u)) ((Sin[b k v])/(k v)) ]; 
ContrastF[x_,a_,b_]:=(1./(1.+Exp[-(x-a)/b])) ;
 
INIT=0 
,TextCell["Not ready"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]];
Manipulate[
If[RUN== 1,

λ=1. ;k=2. π/λ   ;
Manipulate[
{Plot3D[FHRect3[{θRx, θRy},rect,incd,λ],{θRx,-1,1}, {θRy,-1,1}, PlotRange-> All, Mesh-> False]
,Manipulate[Row[{Graphics[{Black,Rectangle[{-10.,-10.},{10.,10.}]
,White,Rectangle[-rect,rect]}],
Show[DensityPlot[ ContrastF[(FHRect3[{θRx, θRy},rect,incd,λ])^2,a,b],{θRx,-1,1}, {θRy,-1.,1.}, PlotRange-> All, PlotPoints-> {80,80}, ColorFunction->GrayLevel]
,Graphics[{Hue[0.,1,1.,0.],EdgeForm[{Red,Dashed}]
,Rectangle[{-0.5,-0.5},{0.5,0.5}]}]]
}]
,{a,0.01,10},{b,1,10}]
}

,{{rect,{2.5,4},"rectangle dimension in λ"},{0.5,0.5},{10,10}},{{incd,{0.,0.},"incident angle"},{-1.,-1},{1.,1.}}] 

,TextCell["Not running"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]]
,Row[{Button["Toggle box to start, and turn off when not running",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Black],Control[{{RUN,0,"->"},{0,1}}]
}]
]
,Button["FRAUNHOFER DIFFRACTION OF RECTANGULAR APERTURE",BaseStyle->{24,Black, Bold,FontFamily-> "Calibri"}
,Background->Hue[0.4,1,1]]
,Row[{Button["Click once to initialize",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Red],Control[{{INIT,0,"->"},{0,1}}]
}]
]

Block[{dat3D,datdensity,D3PLOT, DENSPLOT,plot3d,aperture,plotdens,paraxwindow},
Manipulate[
If[INIT== 1,
{off,on}={Button["OFF",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Red]
,Button["ON",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Hue[0.36,1,0.4]]} ;

FHRect[{θ_,ϕ_},Rect_,λ_]:=Block[{Tθx,Tθy,u,v, k,a,b},
Tθx=Tan[θ] Cos[ϕ] ;Tθy=Tan[θ] Sin[ϕ] ;
{a,b}=Rect ; k=2.π/λ ;
u=Tθx/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
v=Tθy/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
((Sin[a k u])/(k  u)) ((Sin[b k v])/(k v)) ]; 
FHRect2[{θ_,ϕ_},Rect_,incid_,λ_]:=Block[{Tθx,Tθy,u,v,αx,αy, p,q, k,a,b},
Tθx=Tan[θ] Cos[ϕ] ;Tθy=Tan[θ] Sin[ϕ] ;
{a,b}=Rect ; k=2.π/λ ; {αx,αy}=incid  ;
u=Tθx/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
v=Tθy/(Sqrt[1.+((Tθx)^2)+((Tθy)^2)]) ;
p=Tan[αx]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
q=Tan[αy]/(Sqrt[1.+((Tan[αx])^2)+((Tan[αy])^2)]) ;
((Sin[a k (u-p)])/(k  (u-p))) ((Sin[b k (v-q)])/(k (v-q))) ]; 
FHRect3[{θRx_, θRy_},Rect_,incid_,λ_]:=Block[{u,v,αx,αy, p,q, k,a,b,tx,ty},
tx=((Tan[(Sqrt[θRx^2+ θRy^2])]  θRx)/(Sqrt[θRx^2+ θRy^2]))  ; ty= ((Tan[(Sqrt[θRx^2+ θRy^2])] θRy)/(Sqrt[θRx^2+ θRy^2]))   ;
{a,b}=Rect ; k=2.π/λ ; {αx,αy}=incid  ;


u=tx*Cos[αx]/(Sqrt[1.+(tx^2)+(ty^2)]) ;
v=ty*Cos[αy]/(Sqrt[1.+(tx^2)+(ty^2)]) ;
((Sin[a k u])/(k  u)) ((Sin[b k v])/(k v)) ];
ContrastF[x_,a_,b_]:=(1./(1.+Exp[-(x-a)/b])) ;
 
INIT=0 
,TextCell["Not ready"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]];
paraxwindow=Graphics[{Hue[0.,1,1.,0.],EdgeForm[{Red,Dashed}]
,Rectangle[{-0.5,-0.5},{0.5,0.5}]}] ;
Manipulate[
If[RUN== 1,
λ=1. ;k=2. π/λ   ;
Manipulate[
aperture=Graphics[{Black,Rectangle[{-10.,-10.},{10.,10.}]
,White,Rectangle[-rect,rect]}] ;
If[TrueQ[D3PLOT==on] , 
dat3D=Table[Chop[FHRect3[{θRx, θRy},rect,incd,λ]], {θRy,-1.001,1.,2./60},{θRx,-1.001,1.,2./60}] ;
plot3d=ListPlot3D[dat3D, PlotRange-> All,Mesh-> False,BoxRatios-> {1.,1.,0.5}
,DataRange->{{-1.,1.},{-1.,1.}}],];

{(* Plot3D[ Chop[FHRect2[{(Sqrt[θRx^2+((θRy)^2)]),Arg[θRx+I θRy]},rect,incd,λ]],{θRx,-1,1}, {θRy,-1,1}, PlotRange-> All, Mesh-> False
,PlotPoints-> {40,40}] *)
aperture
,Show[plot3d] 
,
If[TrueQ[DENSPLOT==on] ,
Manipulate[
datdensity=Table[ContrastF[(FHRect3[{θRx, θRy},rect,incd,λ])^2,a,b], {θRy,-1.001,1.,2./80.},{θRx,-1.001,1,2./80.}];
plotdens=ListDensityPlot[datdensity, PlotRange-> All,ColorFunction->GrayLevel
,DataRange->{{-1.,1.},{-1.,1.}}] ;
Show[plotdens,paraxwindow]
,{a,0.01,10},{b,1,10}
(*,SynchronousUpdating->False *)],
Show[plotdens,paraxwindow] ] 
}

,Grid[{{Control[{{rect,{2.5,4},"rectangle dimension in λ"},{0.5,0.5},{10,10}}]
,Control[{{incd,{0.,0.},"incident angle"},{-1.,-1},{1.,1.}}]}
,{Control[{{D3PLOT,on,Button["3D Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.12,1,1,1]]},{on,off}
,ControlType-> Toggler}]
,Control[{{DENSPLOT,on,Button["Density Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.52,1,1,1]]},{on,off}
,ControlType-> Toggler}]
}
}]
(*,SynchronousUpdating->False *)] 

,TextCell["Not running"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]]
,Row[{Button["Toggle box to start, and turn off when not running",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Black],Control[{{RUN,0,"->"},{0,1}}]
}]
]
,Button["FRAUNHOFER DIFFRACTION OF RECTANGULAR APERTURE",BaseStyle->{24,Black, Bold,FontFamily-> "Calibri"}
,Background->Hue[0.4,1,1]]
,Row[{Button["Click once to initialize",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Red],Control[{{INIT,0,"->"},{0,1}}]
}]
]
]

Demo Fraunhofer of rectangular aperture for the view along beam axis

4.2 Illustration: Fraunhofer - plane wave with circular aperture

We use the diffraction formula:
               U[(Subscript[x, o])]=-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A)  (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ, o]). r'))  dS                 (2.3.22a)

Here, we will allow a general incident plane wave with a tilt angle instead of only normal incident. Hence we use a more general expression of (Subscript[ρ, o]) for a FF view that is along the plane wave axis:

Thus:       U[(Subscript[x, o])]==-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∮_A)  (e^(-i(k/(|(Subscript[x, o])|)) (Subscript[ρ, oT])r'))  dS     
                          =-A((i k)/(2π)) ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∫_0^a)r'dr'(∫_0^(2 π))(e^(-i(k/(|(Subscript[x, o])|)) |(Subscript[ρ, oT])|.r' Cos[ϕ])) dϕ         (4.2.1)

(∫_0^(2 π))(E^(-I q r Cos[ϕ])) \[DifferentialD]ϕ

U[(Subscript[x, o])]=-A i k ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α](∫_0^a)r'dr'(BesselJ[0, (k/(|(Subscript[x, o])|)) |(Subscript[ρ, oT])| r'])

(∫_0^a)(BesselJ[0, q u]) u \[DifferentialD]u

This is also known as the point spread function, which is very important for lens application.

U[(Subscript[x, o])]=-A i k ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))Cos[α] (a/q) (BesselJ[1, a (k/(|(Subscript[x, o])|)) |(Subscript[ρ, oT])|])                (4.2.2)

To convert this in the far field, we see that for normal incident:
                                    (k/(|(Subscript[x, o])|)) |(Subscript[ρ, oT])|=k Sin[θ]                                       (4.2.3)
which gives us the general expression:
                         U[(Subscript[x, o])]=-A i k ((e^(i k |(Subscript[x, o])| ))/(|(Subscript[x, o])|))(a^2) (((BesselJ[1, k a Sin[θ]]) )/(k a Sin[θ]))                      (4.2.4)

For a more general expression with non-normal incident, we note:
            ((|(Subscript[ρ, oT])|)/(|(Subscript[x, o])|))=(Sqrt[(((Subscript[ρ, o;x])^2)/((|(Subscript[x, o])|)^2))((Cos[α])^2)+(((Subscript[ρ, o;y])^2)/((|(Subscript[x, o])|)^2))])
          =(Sqrt[(((Tan[(Subscript[θ, x])])^2)/(1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)))((Cos[α])^2)+(((Tan[(Subscript[θ, y])])^2)/(1+((Tan[(Subscript[θ, x])])^2)+((Tan[(Subscript[θ, y])])^2)))])     (4.2.5)
 where we have chosen the incident direction to be in the x-plane. With circular geometry, we can make this choice.

Source code: Fraunhofer of circular aperture

Block[{λ,k,dat3D,datdensity,D3PLOT, DENSPLOT,aperture,plotdens,paraxwindow},
Manipulate[
If[INIT== 1,
{off,on}={Button["OFF",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Red]
,Button["ON",BaseStyle->{16, Bold,White,FontFamily-> "Arial"}
,Background->Hue[0.36,1,0.4]]} ;

FHCirc[{θ_,ϕ_},a_,λ_,α_]:=Block[{ ki,q,tx,ty},
 ki=2π/λ ; (*  q=k Sin[θ] ;   *)
tx=Tan[θ] Cos[ϕ] ;  ty=Tan[θ] Sin[ϕ] ;
q=k  (Sqrt[(((tx^2) ((Cos[α])^2)+(ty^2))/(1+(tx^2)+(ty^2)))])  ;
2.(( (BesselJ[1, a q]))/(a q)) ];

ContrastF[x_,a_,b_]:=(1./(1.+Exp[-(x-a)/b])) ;
λ=1. ;k=2. π/λ   ;
 
INIT=0 
,TextCell["Not ready"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]];
paraxwindow=Graphics[{Hue[0.,1,1.,0.],EdgeForm[{Red,Dashed}]
,Rectangle[{-0.5,-0.5},{0.5,0.5}]}] ;
Manipulate[
If[RUN== 1,

Manipulate[
aperture=Graphics[{Black,Rectangle[{-10.,-10.},{10.,10.}]
,White,Disk[{0.,0.},a]}] ;
{aperture
,
Manipulate[
If[TrueQ[D3PLOT==on] , 
ParametricPlot3D[ {Sin[θ] Cos[ϕ],Sin[θ] Sin[ϕ],FHCirc[{θ,ϕ},a,λ,α]},{ϕ,0.,2.π},{θ,0.,1.}, PlotRange->{{-1,1},{-1,1},{-0.2,1}}, BoxRatios-> {1,1,1}, Mesh-> False
,PlotPoints-> {60,30}],]
,Control[{{D3PLOT,off,Button["3D Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.12,1,1,1]]},{on,off}
,ControlType-> Toggler}]
]

,Manipulate[
If[TrueQ[DENSPLOT==on] ,
Quiet[Manipulate[
datdensity=Table[ContrastF[FHCirc[{ArcTan[Sqrt[θRx^2+θRy^2]],Arg[θRx+I θRy]},a,λ,α]^2,c1,10.^(0.1c2)], {θRy,-1.001,1.,2./80.},{θRx,-1.001,1,2./80.}];
plotdens=ListDensityPlot[datdensity, PlotRange-> All,ColorFunction->GrayLevel
,DataRange->{{-1.,1.},{-1.,1.}}] ;
Show[plotdens,paraxwindow]
,{c1,0.01,1.},{{c2,-10},-20.,10.}
(*,SynchronousUpdating->False *)]
],
Show[plotdens,paraxwindow] ] 
,Control[{{DENSPLOT,on,Button["Density Plot",BaseStyle->{16, Bold,Black,FontFamily-> "Calibri"},Background->Hue[0.52,1,1,1]]},{on,off}
,ControlType-> Toggler}]
]
}

,Grid[{{Control[{{a,1.,"aperture radius in λ"},0.2,4}]
,Control[{{α,0,"incident angle"},-1.,1.}]}
}]
(*,SynchronousUpdating->False *)] 

,TextCell["Not running"
,Black,20,Italic,FontFamily->"Calibri",Background->GrayLevel[0.5]]]
,Row[{Button["Toggle box to start, and turn off when not running",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Black],Control[{{RUN,0,"->"},{0,1}}]
}]

]
,Button["FRAUNHOFER DIFFRACTION OF CIRCULAR APERTURE",BaseStyle->{24,Black, Bold,FontFamily-> "Calibri"}
,Background->Hue[0.4,1,1]]
,Row[{Button["Click once to initialize",BaseStyle->{16,White, Bold,FontFamily-> "Calibri"}
,Background->Red],Control[{{INIT,0,"->"},{0,1}}]
}]
]
]

Demonstration

Discussion: the shape of the diffraction pattern vs the tilt angle. Near field pattern vs far field diffraction.

4.3 Illustration HW review: Fraunhofer - far field from rectangular waveguides/edge-emitting diode lasers

See more in HW.