/******************************************************************************* * * McXtrace, X-ray tracing package * Copyright, All rights reserved * Risoe National Laboratory, Roskilde, Denmark * Institut Laue Langevin, Grenoble, France * University of Copenhagen, Copenhagen, Denmark * * Component: PerfectCrystal * * %I * * Written by: Anette Vickery, Andrea Prodi, Erik Knudsen * Date: April 2011 * Version: 1.0 * Release: McXtrace 1.0 * Origin: NBI * * Perfect crystal with diamond or zincblende structure * * %D * Reads atomic formfactors from a data input file. * The PerfectCrystal code reflects ray in an ideal geometry, does not include surface imperfections or mosaicity * * The crystal is positioned such that the long axis of the crystal surface coincides with * z-axis. The angle between the Bragg planes and the crystal surface is alpha * * %D * The algorithm: * Incoming photon's coordinates and direction (k-vector) are transformed into an elliptical reference frame * (elliptical parameters are calculated according to the mirror's position and its focusing distances and the * incident angle), the intersection point is then defined. * A new, reflected photon is then starting at the point of intersection. * Notation follows Tadashi Matsushita and Hiro-O Hashizume, X-RAY MONOCHROMATORS. Handbook on Synchrotron Radiation,North-Holland Publishing Company, 1:263–274, 1983. * * %P * Input parameters: * length [m] length of the crystal (along z-axis) * width [m] width of the crystal (along x-axis) * Material Si, Ge (maybe also GaAs?) * V [AA^3] unit cell volum * h [ ] Miller index of reflection * k [ ] Miller index of reflection * l [ ] Miller index of reflection * alpha [rad] asymmetry angle (alpha=0 for symmetric reflection, ie the Bragg planes are parallel to the crystal surface) * R0 [ ] Reflectivity. Overrides the computed Darwin reflectivity. Probably only useful for debugging. * * * %E *******************************************************************************/ DEFINE COMPONENT Perfect_crystal SETTING PARAMETERS (string form_factors="FormFactors.txt", string material="Si.txt",R0=0, length=0.05, width=0.02, V=160.1826, h=1, k=1, l=1, alpha=0.0) /* X-ray parameters: (x,y,z,kx,ky,kz,phi,t,Ex,Ey,Ez,p) */ SHARE %{ %include "read_table-lib" #include /* something that would be relevant for ALL crystals */ void DarwinReflectivity(double *R, double *Thetah, double *Theta0, double *DeltaTheta0, double f00, double f0h, double fp, double fpp, double V, double alpha, double h, double k, double l, double M, double E, double Thetain, int pol ) { double r0,lambda,theta,theta0,DeltaThetas,a,d,b,C,W,kappa,g,L; double F0r,F0i,Fhr,Fhi,psi0r,psi0i,psihr,psihi; r0 = 2.82e-5; /* Thomson scattering length in AA */ lambda = 12.398/E; /* wavelength in AA, E in keV */ a = cbrt(V); /* side length of unit cubic cell (AA)*/ d = a/sqrt(h*h + k*k + l*l); /* d-spacing (AA)*/ theta = asin(lambda/(2*d)); /* kinematical bragg angle (rad) */ b = sin(theta - alpha)/sin(theta + alpha); /* asymmetry factor */ *Theta0 = Thetain - alpha; /* (rad) angle between Bragg planes and incident ray */ *Thetah = b*(*Theta0 - theta) + theta; /* (rad) Angle betweeb Bragg planes and reflected ray */ /*check if Bragg angle is less than alpha. If so return 0 reflectivity*/ if (theta=0){ /* check whether our intersection lies within the boundaries of the crystal*/ x_int=x+kxu*tin; y_int=y+kyu*tin;V, z_int=z+kzu*tin; if (fabs(x_int)<=width/2 && fabs(z_int)<=length/2){ dist=sqrt(SQR(x-x_int)+SQR(y-y_int)+SQR(z-z_int)); PROP_DL(dist); /* now the photon is on the crystal surface, ready to be reflected... */ SCATTER; /*Check which quadrant the k vector is in to determine sense of alpha. This to allow for hitting the crystal from behind.*/ int quadrant; Thetain=fabs(atan(ky/kz)); double d=cbrt(V)/(sqrt(h*h+k*k+l*l));/*this is valid only for cubic structures*/ f00 = Z; f0h = Table_Value(f0_t,1/(2*d),Z); fp = Table_Value(m_t,E,1)-Z; fpp = Table_Value(m_t,E,2); if (ky<0 && kz>0){ /*4th quadrant - forward hit from above. No change.*/ DarwinReflectivity(&R, &Thetah, &Theta0, &DeltaTheta0, f00, f0h, fp, fpp, V, alpha, h, k, l, M, E, Thetain,pol); }else if(ky<0 && kz<0){ /*3rd quadrant - backward hit from above. Sense of alpha is reversed.*/ DarwinReflectivity(&R, &Thetah, &Theta0, &DeltaTheta0, f00, f0h, fp, fpp, V, -alpha, h, k, l, M, E, Thetain,pol); }else if(ky>0 && kz>0){ /*1st quadrant - forward hit from below. Sense of alpha is reversed.*/ DarwinReflectivity(&R, &Thetah, &Theta0, &DeltaTheta0, f00, f0h, fp, fpp, V, -alpha, h, k, l, M, E, Thetain,pol); }else if(ky>0 && kz<0){ /*2nd quadrant - backward hit from below. No change.*/ DarwinReflectivity(&R, &Thetah, &Theta0, &DeltaTheta0, f00, f0h, fp, fpp, V, alpha, h, k, l, M, E, Thetain,pol); } Thetaout = Thetah - alpha; /* (rad) the angle between the crystal surface and the reflected ray */ /*reflect ray in normal to crystal planes*/ do { double ax,ay,az,vnx,vny,vnz; ax=0;ay=-sin(alpha);az=cos(alpha); vnx=kx-scalar_prod(kx,ky,kz,ax,ay,az)*ax; vny=ky-scalar_prod(kx,ky,kz,ax,ay,az)*ay; vnz=kz-scalar_prod(kx,ky,kz,ax,ay,az)*az; //printf("a=[ %g %g %g ] vn=[ %g %g %g ]\n",ax,ay,az,vnx,vny,vnz); kx=kx-2*vnx; ky=ky-2*vny; kz=kz-2*vnz; }while(0); /* apply Darwin reflectivity if not is supplied from outside*/ if (!R0){ p*=R; }else{ p*=R0; } /*catch dead rays*/ if (p==0) ABSORB; } else { RESTORE_XRAY(INDEX_CURRENT_COMP, x, y, z, kx, ky, kz, phi, t, Ex, Ey, Ez, p); } } %} MCDISPLAY %{ rectangle("xz",0,0,0,width,length); %} END