help me figure out refraction equation, please! its easy!

[unxzst]

This is my function, it calculates the normal reflectivity (%) of a non-absorbing multi-layer film on a substrate.
Function Rpc(Lambda As Single, SubstIndex As Single, n As Single, FilmProperty As Variant)
'FilmProperty(2N) made of R(1), T(1), R(2), T(2), ..., R(N), T(N).
'R(I) = FilmIndex
'T(I) = FilmThickness
'The layers are numbered starting from the one most away from the substrate.
'MEDIUM INDEX = R(0)
ReDim R(1 To n), T(1 To n), X(1 To n), C(1 To n), S(1 To n) As Variant
AA = 1 / SubstIndex
BB = 1
CC = 1 / SubstIndex
For i = 1 To n
R(i) = FilmProperty(2 * i - 1)
T(i) = FilmProperty(2 * i)
X(i) = (2 * 3.141592654 * R(i) * T(i)) / Lambda
'CALCULATE SINE AND COSINE TERMS FOR MATRIX ELEMENTS
C(i) = Cos(X(i))
S(i) = Sin(X(i))
Next i
'SETUP MATRIX ELEMENTS
B11 = C(1)
B12 = S(1) / R(1)
B21 = S(1) * R(1)
B22 = B11
'CHECK FOR ONE FILM
If n = 1 Then GoTo 500
For i = 1 To n
C11 = C(i)
C12 = S(i) / R(i)
C21 = S(i) * R(i)
C22 = C11
A11 = B11 * C11 - B12 * C21 'C11
A12 = B11 * C12 + B12 * C22 'C12
A21 = B21 * C11 + B22 * C21 'C21
A22 = -B21 * C12 + B22 * C22 'C22
B11 = A11
B12 = A12
B21 = A21
B22 = A22
Next i
500 CA = (AA * B11 - B22) ^ 2 + (BB * B12 - CC * B21) ^ 2
CB = (AA * B11 + B22) ^ 2 + (BB * B12 + CC * B21) ^ 2
Rpc = (CA / CB)
End Function


i will also have data being gathered from an instrument, which is monitoring the deposition of each thin film, by mearusing the transmission of light through quartz. trasferring between transmission and reflection is easy:
T=(1 - R) * 100

what I'm having trouble with is the crux of my program. I need to figure out at what thickness I am, and adjust that as deposition progresses.
I found some code online that will do what Solver does in Excel, plug in variables, get an answer, and iterate until the answer matches prediction, outputting the neccesary thickness values.
here, the example, uses a simple equation p1+p2*cx+p3*cx^2+p4*cx^3+p5*sin(p6*cx)
Private Function Fitter_OnGetLMFunction(ByVal X As Variant, ByVal Parameters As Variant) As Double
Fitter_OnGetLMFunction = Parameters(0) + Parameters(1) * X + _
Parameters(2) * X ^ 2 + Parameters(3) * X ^ 3 + _
Parameters(4) * Sin(Parameters(5) * X)
End Function
Sub Fitter_OnProgress()
Dim I, X, Parameters
Workbooks(1).Worksheets(1).Cells(39, 7) = Fitter.Iterations
For I = 1 To NumPoints
X = Workbooks(1).Worksheets(2).Cells(I, 1) ' X in 1st column
Parameters = Fitter.Parameters
Workbooks(1).Worksheets(2).Cells(I, 3) = Fitter_OnGetLMFunction(X, Parameters)
Next
DoEvents ' refresh chart
End Sub

Private Sub FitCommandButton_Click()
Dim I, datasheet, X, Y, Parameters
Set datasheet = Workbooks(1).Worksheets(2)
ReDim X(1000)
ReDim Y(1000)
I = 1
While datasheet.Cells(I, 1) <> ""
X(I - 1) = datasheet.Cells(I, 1)
Y(I - 1) = datasheet.Cells(I, 2)
I = I + 1
Wend
ReDim Preserve X(I - 2)
ReDim Preserve Y(I - 2)
NumPoints = UBound(X) + 1 ' save to use in Fitter_OnProgress

ReDim Parameters(5)
For I = 0 To UBound(Parameters)
If Workbooks(1).Worksheets(1).Cells(34 + I, 4) <> "" Then
Parameters(I) = Workbooks(1).Worksheets(1).Cells(34 + I, 4)
Else
Parameters(I) = 1
End If
Next

Set Fitter = CreateObject("Fitter.DMFitter")
Fitter.X = X
Fitter.Y = Y
Fitter.Expression = -2 ' perform differentiation numerically
Fitter.ParamCount = 6
' Fitter.Parameters = Array(1, -1, -0.1, 0, 50, 0.192)
Fitter.Parameters = Parameters
Fitter.Sigmas = Array(0, 0, 0, 0, 0, 0)
Fitter.WeightType = 0
Fitter.Options = Array(0, 0.0000001, 0)
Fitter.Iterations = 60
If Fitter.LMFit Then
Workbooks(1).Worksheets(1).Cells(33, 4) = Fitter.Deviation
Parameters = Fitter.Parameters
For I = 0 To UBound(Parameters)
Workbooks(1).Worksheets(1).Cells(34 + I, 4) = Parameters(I)
Next
Else
MsgBox "Fitter error " & Fitter.ResultCode
End If
End Sub

but i need to adapt it to use my matrix equation. I am doing this in VB6, and I have no idea how I can begin to adopt this Levenberg-Marquardt least-squares fitter (ActiveX from http://fitting.datamaster2003.com/index.htm )to use my matrix function...but if anybody could help point me in the right direction for curve-fitting, it'll be greatly appreciated.


I found this Levenberg-Marquardt algorithm for multivariate optimization, sadly this is not something covered in school (or maybe i was sleeping?)
I cannot make head or tail of this example http://www.alglib.net/optimization/levenbergmarquardt.php
maybe someone can point me to whatever is being optimized? I THINK i have to extrace a Jacobian determinant from the function in form of a matrix, and then plug it in...but where are the dependant variables??

[zaza]

Is this a fair summary?

You're depositing multilayers of some film on a substrate and measuring the reflectivity of the film. From this, you want to know how many multilayers you've deposited so far, so that you can adjust the deposition rate.

My question is, what are you trying to fit to with the latter routine? Do you have a function that relates the thickness of the film to the reflectivity? If so, why not just calculate an array for each integer number of possible multilayers and check your reflectivity against that?


zaza

[unxzst]

Is this a fair summary?

You're depositing multilayers of some film on a substrate and measuring the reflectivity of the film. From this, you want to know how many multilayers you've deposited so far, so that you can adjust the deposition rate.

My question is, what are you trying to fit to with the latter routine? Do you have a function that relates the thickness of the film to the reflectivity? If so, why not just calculate an array for each integer number of possible multilayers and check your reflectivity against that?


zaza
youre almost right, layers get deposited one by one, and i'm trying to determine the thickness of each layer as it is being deposited by measuring the transmission.

from this transmission data, i need to backpedal through the matrix algorithm to figure out exactly what thickness i have. im trying to accomplish this by fitting the data to the original matrix, and with help of optimization find out the thickness of the film.

the matrix solves for a stack of films, i can set thicknesses of all the not-yet-deposited films to 0, and start building the first layer by incrementing its thickness to the desired limit, then start the next layer.

this is the code i use to output transmission as a function of thickness, simulating the deposition process, and what the graph looks like of the data (more or less, i took out all the stuff that wouldnt be clear, just showing the process)
Private Sub Command1_Click()
ReDim FilmP(1 To (Nl * 2)) 'film property array
Dim z, Stp, mB, sI, sB, ssB, Lambda As Single 'InLayer Thickness, step size, substrate index, Backside reflection, 1-BS refl.
sI = 1.52141 'index of substrate
sB = ((1 - sI) ^ 2) / ((1 + sI) ^ 2) 'substrate backside reflection
ssB = 1 - sB

For i = 1 To Nl
FilmP(i * 2) = 0 'each thickness set to 0 for a clean start
Next i

nn = 0
LayerNumber = Nl
tcur = 0
TStart = 0
resolution = 10
CurvEnd = resolution

Do While LayerNumber > 0
Step = Ta(LayerNumber) / resolution 'number of sample points per layer
ReDim Preserve curve(0 To CurvEnd) 'make room for new data
LayerStart = tcur
Lambda = La(LayerNumber) 'each layer has its own monitor wavelength
FilmP((LayerNumber * 2)) = 0 'thickness start at 0
curve(tcur) = RFun((La(LayerNumber)), (sI), (Nl), (ssB), (sB), (FilmP)) 'get value
tcur = tcur + 1 'move up one

For Counter = 1 To resolution
z = Round(Step * Counter, 4)
FilmP((LayerNumber * 2)) = z ' Thickness_a
curve(tcur) = RFun((La(LayerNumber)), (sI), (Nl), (ssB), (sB), (FilmP)) 'resolution-samples

tcur = tcur + 1

Next Counter

CurvEnd = CurvEnd + resolution + 1
cuts(nn) = curve(tcur - 1)
nn = nn + 1
LayerNumber = LayerNumber - 1
Loop

does this make sense? thank you so much for helping me with this

the problem is that this is supposed to be dynamic, and calculations are done as the layer is being deposited, in effect, i need to know where i am, and how much has been deposited, so I can find the rate of deposition, and model how long it will take, and when I need to stop. As soon as deposition starts, every 2 seconds I get a data point, from this i need to figure out rate (thickness as a function of time) so i can figure out when to stop with the precision of 0.1 seconds.:afrog:

[zaza]

I am still not clear on exactly what it is that you are doing to work out your current thickness from the data that you have. Surely you must have either a calibration curve or some function that allows you to work out the thickness from the reflectivity.

You can't use experimental data to work out an exact value, because the data points have errors.
If you don't have either a function or a curve, you will need to generate a calibration curve. For that, you'll need a different way of measuring the thickness of your layer (eg conductivity?), then deposit for a variety of different times and measure the thicknesses and reflectivities of each.


I am not yet convinced by your explanation of pushing your data back through some matrix calc and fitting something to it to try to get a thickness.


Are you saying that you are planning to bung a variety of thicknesses into the first bit of code until you get a reflectivity out that matches what you measure?


zaza

[unxzst]

the matrix calc is THE source of transmission/reflection values, by plugging in thickness, index of individual layers, index of substrate, number of layers and monitoring wavelength. getting the values into it is sort of complex (i need to use Cauchy's equation, reflections from substrate that are not accounted) but I can reliably say that its accurate with comparison to few other sources of info.

if i really wanted to, i could generate a curve of thickness vs. transmission, and then see how far my data is from that. but then i run into problems with the dynamic part of the program, because i need to know at each point, what the thicknes is, to determine rate. if it wasnt for that, i can just generate a table T(thickness) vs %T(transmission), and at each %T data point from the table get T value.

i sort of need to guesstimate the thickness, pushing it back in to see how close it is to my data. by getting a few data points, even with noise (oh the noise...) i can accurately predict where the process is.

this can be simulated in Excel using solver. it picks values for transmission to minimize the squared difference. but i can't use solver since its proprietary

[unxzst]

Are you saying that you are planning to bung a variety of thicknesses into the first bit of code until you get a reflectivity out that matches what you measure?
zaza
basically yes.
except thickness will be determined by difference in predicted thickness and measured, then calculated thickness. yeah...over and over again, getting the best fit curve to the data, from the matrix, at each data point, into the future.

[zaza]

if i really wanted to, i could generate a curve of thickness vs. transmission, and then see how far my data is from that. but then i run into problems with the dynamic part of the program, because i need to know at each point, what the thicknes is, to determine rate. if it wasnt for that, i can just generate a table T(thickness) vs %T(transmission), and at each %T data point from the table get T value.



You will have a lot less computational trouble with looking a value up from an array than by completing an optimization calc before the next datapoint arrives. I am surprised that the Excel solver would be able to solve anything remotely complex in 2 secs.
And if it is as difficult as you say to even get the data into your function, I would have thought you'd want to avoid using it during data collection if at all possible.

Why not compute the thickness in advance for each monolayer and then look it up. Your table of thicknesses won't change because the only thing you will vary is the deposition rate and the data should be independent of that. It sounds to me like the most important thing is getting some idea of how many layers you have got deposited before the next datapoint comes in, so that you can stop the process where necessary.

Possibly a compromise between the two would be the best option if you wanted good accuracy and some speed - use a look up table to get pretty close to where you want to go and then slow the deposition rate down as you approach your required value.

[unxzst]

You will have a lot less computational trouble with looking a value up from an array than by completing an optimization calc before the next datapoint arrives. I am surprised that the Excel solver would be able to solve anything remotely complex in 2 secs.
And if it is as difficult as you say to even get the data into your function, I would have thought you'd want to avoid using it during data collection if at all possible.

Why not compute the thickness in advance for each monolayer and then look it up. Your table of thicknesses won't change because the only thing you will vary is the deposition rate and the data should be independent of that. It sounds to me like the most important thing is getting some idea of how many layers you have got deposited before the next datapoint comes in, so that you can stop the process where necessary.

Possibly a compromise between the two would be the best option if you wanted good accuracy and some speed - use a look up table to get pretty close to where you want to go and then slow the deposition rate down as you approach your required value.
right. i was trying to automate it to obsoletion. but a lookup table will still need to use something like solver to determine the rate, that is, current thickness. getting the data into the function is not hard, but it takes some calculation, thats all (this part is already taken care of, just need to work in a little fidge coefficient)

deposition rate stays almost constant, but it can change slightly during the same layer (these take about 50-150 seconds each to deposit, this is where rate comes from)
even the index can change, or be different than prediction.
so a lookup table will need to be updated with each new data point any way.

also there needs to be compensation for compensation of the measuring instrument...:eek2:

so here's how the process looks like in theory

measure transmission of bare quartz, verify its index(as a function of wavelength), calcualte index of substrate to use later
start deposition, get first data point
attempt best fit
get next data point, make best fit...
determine rate, calculate time to finish
using rate determine current thickness, when to stop

[zaza]

So are you saying that the reflectivity measurement will give you different values for the thickness depending on the history of the sample? i.e. if the deposition rate is fast the reflectivity measurement will give you a certain value for the thickness but if it is slow then you will get a different value for the thickness?

I don't see this as a feature in your initial function. And to be honest, I'm not really sure why this would be the case. I don't profess to be an expert in thin films, but I would have thought that issues with transmission / reflectivity would largely come about due to quality of the sample, which I wouldn't have thought to be a significant issue from this perspective for thin films.

[unxzst]

So are you saying that the reflectivity measurement will give you different values for the thickness depending on the history of the sample? i.e. if the deposition rate is fast the reflectivity measurement will give you a certain value for the thickness but if it is slow then you will get a different value for the thickness?

I don't see this as a feature in your initial function. And to be honest, I'm not really sure why this would be the case. I don't profess to be an expert in thin films, but I would have thought that issues with transmission / reflectivity would largely come about due to quality of the sample, which I wouldn't have thought to be a significant issue from this perspective for thin films.
deposition rate may affect the index of the layer, which also needs to be calculated via solver-like approximation. rate does not change drastically. typically it is 1.5-2.5 Angstroms, and change only in the second or third decimal, however it does chage.

Rate gets calculated during the first few points of data collected, after that it remains pretty much the same. So optimization has to occur in the first few seconds, thats the thing, afterwards its just kind of like keeping an eye on it.

furthermore, index can change with change in temperature, quality of coating, etc...again, in the second or third decimal, but still significant enough and needs to be tracked.

it is a significant issue, even if the quartz witness chip is off by a fraction of a degree in the fixture, or if the floor is vibrating, or if the heater is on and emitting infrared radiation...the measuring instrument will compensate, this is where fudging(coefficients ;) ) come in.

[zaza]

What actually is the index of a layer?

[unxzst]

What actually is the index of a layer?
sort of like index of refraction of the material, but its different because layer can be porous, contain traces of other materials etc...
usually two contrasting indexes are used in thin films, high and low, ~1.7 and ~2.3...something like that.
this is approximately known, and is further refined during the process.

[zaza]

It sounds like a dangerous path down which you are embarking. If you need to do a fit then you'll need at least 3 data points, and probably more before your fit has a reasonable chance of working. If you are relying on the results of these early fits to do the later ones, then you will rapidly veer off track, I fear, because the noise will just build up and be aggregated into subsequent points.

Your alternative is to re-fit to all of the data as each new data point comes in, meaning that you will be doing an increasing number of calculations as the sample gets thicker.

Also, once you get into the realm of multi-variant fitting, you may find that the fit doesn't come up with a physically reasonable answer although it does successfully minimise. Often you need to fit multiple times and use acceptable boundary conditions to arrive at a justifiable answer. Doing this bit computationally is not really an option, because the computer can't make the judgement that, say fitting an index of 200 and thickness of 1cm is worse than an index of 5 and thickness of 100A.

But anyway. Fitting. How's your C? I reckon a good place to go would be the Numerical Recipes site: nr.com. Chapter 10 covers minimization of functions, including the downhill simplex method, and it has a load of code in it as well. And you can download it free.

Take a look and see what you think.

[unxzst]

noise is not an issue, check out this pic.
http://i51.photobucket.com/albums/f356/Unxzst/rpcntm.jpg
thats just plain number-crunching in excel, doing the same thing i'm trying to do in VB.

i'm afraid you guys are too smart, and are seeing a snake instead of an elephant.

all i'm trying to do, is to replicate solver in excel. If i can fit Rpc into excel, and have solver optimize for thickness and index, like in the picture, why can't i do the same thing in VB? all it is, is trying to avoid having to use excel.

here's the plan: run Rpc with approximate index and thickness, store values into array, get data point, re-optimise, get another data point, re-optimize again. getting least squares fit each time would leave me pretty close, noise or no noise.

[zaza]

Fair do's. I'd hate to be accused of being too smart. I'd prefer the expression "trying to avoid wasting time on something that may not work". Or "trying to understand the problem before going off in the wrong direction".

Anyhow, you seem to understand the idea; you need to have the following:

1) One or more functions, which take a set of input parameters and return a value

2) A comparison set of data to which you are fitting

3) An iteration routine to optimise your parameters


Your functions you have got. One for reflectivity and one for index.

Your comparison dataset you will be generating as you go.

Your optimisation routine; that's a matter of choice. I've used downhill simplex with much success in the past; it's in the link I gave you earlier.


You choose a set of starting parameters and pass them to the optimisation routine. The optimisation routine calls the function with said parameters and calculates the predicted datapoint for each input X value. The optimisation routine then calls the least-squares fit routine (of which many examples on these forums, if you don't have one already) which returns some goodness-of-fit parameter. The optimisation routine then tweaks the parameters and recalculates the goodness-of-fit. Loops until goodness-of-fit is less than a specified tolerance.

If you look at the optimiser in the module in the link you sent, you'll see that the user is intended to specify FuncVecJac, which is then plumbed into the routine and optimised. Have you tried just importing the module, putting your function in and seeing what happens?

The only way to get somewhere with these sorts of problems is to try them out. As I said before, for a more wordy guide to how these things work, look up Numerical Recipes.


zaza

[unxzst]

1) One or more functions, which take a set of input parameters and return a value

2) A comparison set of data to which you are fitting

3) An iteration routine to optimise your parameters

zaza

Beer and thanks to you for the good word mate, you've got it figured out!

This is not much easier to understand than it is to explain, for an intern...
I've been fiddling my way thought the LM non-linear optimizator, and one thing I've learned is that I could use more VB experience. Although the language is not hard, formulating the logic and then expressing it is kinda difficult, and thats kind of what I've been doing for the past two weeks. I am going to rework my code, using Excel solver as a proof-of-concept try and work my way back from it.

I dont see a link for the downhill simplex and I'm curious to find out what it is, could you give it again, please?

[zaza]

Numerical Recipes, online versions (http://www.nrbook.com/a/bookcpdf.php)

Chapters 10 and 15 will probably be of interest to you.