Documentation for optical modelling program
Contents
Introduction and overview
This is a program for optical modelling using the Fresnel equations through a series of stratified layers. The permittivities of the materials can be calculated using a variety of models, including the Bruggeman EMT (effective medium theory) for non-homogenous systems. This README document attempts to describe the program from a users perspective and it does not explain the source code in detail. I've made some attempt to comment the source code, and parts of my thesis describe some of the math.
If you are not familliar with optical modelling using the Transfer Matrix Method, you can read part of the Supporting Information for Andvaag, I., Morhart, T., Clarke, O., Burgess, I. ACS Appl. Nano Mater. 2019 2 (3), 1274-1284. It should give a reasonable primer on the various permittivity models and how to use the Fresnel equations.
Here's an rough outline of how I expect the program might be used:
- Create a materials() object, thereby initializing the domain mode, the fixed AOI or frequency, the low and high limits of the domain of interest, and the number of wavenumbers between plotted points
- Assign the material types using the methods of the materials object.
- Create a phaseSys()object
- Define a stratified system on which to perform calculations by calling the setLayers() and setThicknesses() methods of the phaseSys() object
- Perform the calculation using the calcR(), calcT(), calcA(), or calcE() methods of phaseSys().
- Plot up the computed data in the notebook, or export an ASCII file of the computed values.
Instructions for use
There are two main classes that do the heavy lifting: the materials() class, and the phaseSys() class. There is also a createFit() class for fitting to the Drude model. To make use of these classes, you must import them with:
from opticalc_330 import *
Using the materials() class
An object that defines the basic parameters of the system and has various methods to define permittivities of materials.
Using the phaseSys() class
An object used to setup the specific system and calculate reflectivity, transmissivity, and absorptivity.
Using the createFit() class
An object used to visually fit experimental reflectivity data to the Drude model.
This example use case also exists in a jupyter notebook form titled ATR-SEIRAS_example.ipynb
State of the code
Please use this code at your own risk and do not trust its accuracy. I have attempted to implement the mathematics correctly, but there are surely many errors and bugs. In order to verify that my implementation is correct, I have tried to reproduce various plots found in the literature. The notebook <test_cases_for_evaluation.ipynb> contains some of the test cases that I have successfully reproduced, as well as some cases that I have failed to reproduce.
Known Working
- getGeometricFactor(a, b, c, forceGeneral = False) Consistent. Any of Bohren and Huffman, Osborne, Fedotov are all acceptable and in agreement. The general formula is given in Bohren and Huffman, page 145. See DEMO_001 for replication of Fig 5.6 on page 147 in Bohren and Huffman. See DEMO_002 for a check for the consistency of the general formula and the analytical expressions for the ellipsoids of revolution. Depolarization factors for ellipsoids are explained in the following two papers: Stoner, E. C. The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science 1945, 36, 803-821. and Osborn, J. A. Physical Review 1945, 67, 351-357.
- setBasicMaxwellGarnett(self, matName, inclusionMat, matrixMat, f) Working. See DEMO_003 for replication of Fig 3 in Bohren and Battan 1980.
- setBasicBruggeman(self, matName, inclusionMat, matrixMat, f) Working. See DEMO_003 for replication of Fig 3 in Bohren and Battan 1980. Uses the getPhysicalRoot() function. 2-phase, 3-dimensional, no self consistency requirement.
- setCoatedEllipsoid(matName, coreMat, shellMat, a_core, b_core, c_core, rotMat=None, lamda=False) Working. See DEMO_006 which replicates figures 2c, 3d and 4d in: 2017_zhao_mei_J._Phys._D _Appl._Phys._50_505001. TO DO: add rotMat implementation so that the ellipsoid axes need not align with the multilayer stack axes. Could add an edge case to warn the user if they try and call the function using lamda = 1, or other invalid range (lamda < 0; lambda >= 1). xi2 = np.real(rutz[2]) # TODO: Verify that this is the case. -> # Select the root which yields purely real a2, b2, c2. The second root always seems to be correct.
- calcR(), calcT(), calcA(), calcE(), and propMatrix() Broadly, the Transfer Matrix Method which employs the Fresnel equations is working correctly. See DEMO_007 for reproduction of some figures in Ohta, K.; Ishida, H. Appl. Opt. 1990, 29, 1952-1959.. One issue I have not been able to resolve is my inability to reproduce figure 5 from Hansen, W. J. Opt. Soc. Am. 1968, 58, 380-390.. The frequency of the standing waves in the incident semi-infinite phase is double what Hansen reports. I'm not sure why this is the case. The exponent is i*2*pi*nu*eta*cos(theta)*delta(z). If the exponent is halved, the frequency of the standing waves will be halved, but I cannot find any physical justification why any of these terms should be halved.
Unknown Status (as of yet unverified)
- setBruggeman(self, matName, coreMat, shellMat, matrixMat, thickness, molec, ratio1, F, plotResult=False, iterations=1) Bruggeman formula for 3 phase, coated core-shell ellipsoids, NO self consistency condition. Uses equation 6 from Granqvist, C. G.; Hunderi, O.. Phys. Rev. B.1978, 18, 2897
- getPhysicalRoot(wienerZ_ARR, roots_ARR) This function uses the Wiener bounds to select the physically correct solution of the n-phase Bruggeman effective medium approximation. See 1993 Jansson & Arwin. Optics communication 106. p 133-138. They provide a fantastic explanation of the algorithm to select the physically correct root, however, no examples are given of the final permittivities calculated.
- setMultiBruggeman(self, matName, constituentMat_ARR, fill_ARR)The math used by this function is described in my MSc. thesis. I have verified that my implementation produces correct output in the limiting case of a two-phase composite (see DEMO_003, that the black dotted lines produced by the setMultiBruggeman() function match those produced by the setBasicBruggeman() function. I do not have test cases for 3 or more constituents to verify my implementation. This function makes us of getPhysicalRoot().
Known Issues
- setOsawaBruggeman(matName, ellipsoidMat, adsorbateMat, hostMat, thickness, molec, ratio1, F) Unable to quantitatively reproduce the figures in Osawa, M.; Ataka, K.; Yochii, K.; Yotsuyanagi, T. J. Electron. Spectrosc. Relat. Phenom. 1993, 64-65, 371-379.. The general lineshapes are approximately correct, but the magnitudes are incorrect. See DEMO_010.
References for experimental refractive indices
|
Reference |
Materials |
1. |
Babar, S. and Weaver, J. H. 2015. Appl. Opt. |
Ag, Au, Cu |
2. |
Hale, G. M. and Querry, M. R..1973. Appl. Opt.12, 555-563. |
H2O |
3. |
Li, H. H. 1980. J. Phys. Chem. Ref. Data 9, 161-289. |
CaF2 |
4. |
Li, H. H.1993. J. Phys. Chem. Ref. Data 9, 561-658. |
Ge, Si |
5. |
Li, H. H.1976. J. Phys. Chem. Ref. Data 5, 329, 528 |
KBr |
6. |
Ordal, M. A., Bell, R. J., Alexander, R. W., Long, L. L., Querry, R. R.1987. Appl. Opt. 26, 744-752. |
Ni |
7. |
Philip, H.R. and Taft, E. A. 1964. Phys. Rev.136, A1445. |
C (diamond) |
8. |
Querry, M. R. 1987. Contractor ReportCRDEC-CR-88009 |
ZnSe |
Ian Andvaag is the primary author of this program. He can be reached at: ira044@usask.ca
License
This project is licensed under the GNU GPL v2.0 license. If you make something cool, please reach out to me -- I'd love to hear about it!