Calculating the profile


Functions

double calcProfile (double z, double *a, Parratt *par)
 Calculate the profile.


Function Documentation

double calcProfile ( double  z,
double *  a,
Parratt par 
)

Parameters:
z The value of $z$ to evaluate $\rho(z)$
a The set of all adjustable parameters in one array.
p The Parratt struct containing all the other information necessary for calculation of $\rho(z)$.
Returns:
$\rho(z)$
For a set of $N$ layers, that run from bulk ($i=0$) to substrate ($i=N+1$), the scattering length density profile can be given as

\begin{eqnarray*} \rho(z) &=& \rho_0 + \sum_{i=1}^{N+1} \frac{(\rho_i-\rho_{i-1})}{2} \left( 1+\mathrm{erf}\left(\frac{z-z_i}{\sqrt{2}\sigma_i}\right) \right) \\ z_i &=& \sum_{j=1}^{i} t_j \end{eqnarray*}

where $t_i$ is the thickness, $\rho_i$ is the real part of the scattering length density, $\sigma_i$ the roughness, each of layer $i$. The location of the error function, $z_i$, is calculated by a running total of the thicknesses.

In actuality, the total number of layers, which is given by num, includes the substrate. For practical purposes of C arrays, the layer numbering runs as $i=0\rightarrow$ num-1. Therefore the parameters (6 for each layer), can be found in the following manner:

for each layer i = 0 to num-1
	index = idx[i]
	par[index+0] = Thickness
	par[index+1] = Roughness
	par[index+2] = Real part of nuclear SLD
	par[index+3] = Imaginary part of nuclear SLD
	par[index+4] = Real part of magnetic SLD
	par[index+5] = Imaginary part of magnetic SLD

When the roughness is zero, the error function becomes a Heavyside step function.

00023 {
00024     int     j, indx, indx1;
00025     double  thk, tmp, tmp1, tmp2;
00026     
00027     /********************
00028      *  The first layer is the layer 
00029      *  the neutrons enter from.
00030      */
00031     indx = par->l.idx[0];
00032     if (par->g.plz == POLARIZED_UP) tmp = a[indx+2] + a[indx+4];
00033     else                            tmp = a[indx+2] - a[indx+4];
00034 
00035     /********************
00036      *  Loop over the layers
00037      */
00038     thk = 0.0;
00039     for (j=1; j<par->l.num; j++)
00040     {
00041         indx  = par->l.idx[j];
00042         indx1 = par->l.idx[j-1];
00043 
00044         if (par->g.plz == POLARIZED_UP) 
00045         {
00046             tmp2 = 0.5*((a[indx+2] + a[indx+4]) - (a[indx1+2] + a[indx1+4]));
00047         }
00048         else
00049         {
00050             tmp2 = 0.5*((a[indx+2] - a[indx+4]) - (a[indx1+2] - a[indx1+4]));
00051         }
00052         
00053         /* zero roughness : erf becomes step function */
00054         if (a[indx+1] == 0)
00055         {
00056             if (z>=thk) tmp += 2.0 * tmp2;
00057         }
00058         else
00059         {
00060             tmp1 = (z-thk)/(M_SQRT2*fabs(a[indx+1]));
00061             tmp += tmp2 * (1.0 + gsl_sf_erf(tmp1));
00062         }
00063         thk += a[indx+0];
00064     }
00065         
00066     return tmp;
00067 }


Generated on Wed Mar 14 13:24:56 2007 for Parratt by  doxygen 1.4.7