1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
|
/*******************************************************************************
*
* McXtrace, x-ray tracing package
* Copyright, All rights reserved
* DTU Physics, Kgs. Lyngby, Denmark
* Synchrotron SOLEIL, Saint-Aubin, France
*
* Component: Grating_reflect
*
* %Identification
*
* Written by: Erik B Knudsen (erkn@fysik.dtu.dk), Kristian Sorensen and Philip Smith
* Date: June 2021
* Version: 1.6
* Origin: DTU
*
* A reflective grating.
*
* %Description
* A reflective grating that diffracts incident photons.
* The grating is in the XZ-plane. It then reflects the incoming photon using a MC picked angle,
* where the angle is picked from a uniform distribution of width d_phi, i.e. U[-d_phi/2,d_phi/2]
* The Monte Carlo wight of the ray is then adjusted wrt. to the grating interference pattern, and
* the diffraction pattern associated with each grating line. All lines are considered equal.
* For more efficient sampling of a particular direction the centre of the d_phi may be shifted
* using the parameters order or phi0. In the latter case a set angle is chosen as the centre of the
* sampled interval, in the former the centre angle is computed from the specified grating order.
*
* In an upcoming release this grating model will also include a blazed grating.
*
* Example: Grating_reflect(
* d_phi=1,order=0,rho_l=100,zdepth=102e-3,xwidth=102e-3)
*
* %Parameters
* Input parameters:
* xwidth: [m] The width of the grating.
* zdepth: [m] The length of the grating.
* R0: [0-1] Constant reflecticity of the grating [0;1].
* rho_l: [l/mm] Number of lines pr mm of the grating.
* b: [AA] Width of the spacing in Angstrom. If zero, default is found using rho_l/3.
* d: [AA] Width of the slits in Angstrom. If zero, default is found using rho_l.
* N: [1] Number of slits. If Zero, default is found using rho_l.
* d_phi: [deg] Range of diffraction angle that is to be simulated -d_phi/2 ; d_phi/2.
* order: [1] The target order of the grating. If non-zero d_phi will be centered around this scattering line.
* phi0: [deg] Target angle to center d_phi. If this is set to 0 the 0th (or any other chosen by the parameter order) order line will be used.
* verbose:[0/1] If non-zero, more information will be displayed. Nb. generates much output.
*
*
* %End
*******************************************************************************/
DEFINE COMPONENT Grating_reflect
SETTING PARAMETERS (d_phi=1, R0=1, rho_l=800,
order=0, phi0=0,
zdepth=0.015, xwidth=0.136, int verbose=0)
SHARE
%{
#include <complex.h>
%include "read_table-lib"
%}
/********************************************************************************************
In DECLARE section, small functions or variables can be defined and used in the entire grating.
********************************************************************************************/
DECLARE
%{
double pdir;
double d;
double b;
int nslits;
%}
/********************************************************************************************
The INITIALIZE section run each statement once.
********************************************************************************************/
INITIALIZE
%{
int status;
if (R0<0 || R0>1){
fprintf(stderr,"ERROR: (%s) reflectivity (%f) is specified but is not in [0:1]. \n",NAME_CURRENT_COMP,R0);
exit(-1);
}
d=1e7/rho_l;/*grating pitch in Angstrom*/
nslits=(zdepth*1e3)*rho_l; /*the relative width of openings*/
if(!b){
/* Approximated slit-width in Angstrom*/
b=d/3;
}
if(verbose){
printf("INFO: (%s): Line width, d=%f [AA]. Number of slits, N=%f. Slit width, b=%f [AA].\n \n",NAME_CURRENT_COMP,d,nslits,b);
}
/*We are not sampling the full 2pi range*/
pdir=(((d_phi*DEG2RAD))/(2*M_PI));
%}
TRACE
%{
/********************************************************************************************
Initializing by simulating a perfectly reflecting mirror:
********************************************************************************************/
/*Placing the grating along the Y-axis. */
PROP_Y0;
/*If the photon is passing the grating, it shouldn't be reflected. Instead, for book-keeping, it is restored to previoulsy state before the component using RESTORE_XRAY */
if(x<-xwidth/2.0|| x>xwidth/2.0 || z<-zdepth/2.0 || z>zdepth/2.0){
RESTORE_XRAY(INDEX_CURRENT_COMP, x,y,z, kx,ky,kz, phi,t, Ex,Ey,Ez, p);
}else{
double nx,ny,nz;
double Gamma,gamma,psi_i,psi_o,spsi_i,spsi_o;
double kx_notouch,ky_notouch,kz_notouch;
double delta_theta,pdiff;
/* Normal vector for the grating*/
nx=0;
ny=1;
nz=0;
/* scalar-product, s and length of k.*/
double s=scalar_prod(kx,ky,kz,nx,ny,nz);
double k=sqrt(scalar_prod(kx,ky,kz,kx,ky,kz));
/* Outgoing grazing angle for a perfectly reflecting mirror using angle between two vectors.*/
psi_i = acos(s/k)-M_PI_2;
/* outgoing vector for a perfectly reflecting mirror. Order 0.*/
kx=kx-2*s*nx;
ky=ky-2*s*ny;
kz=kz-2*s*nz;
// MC picked angle going from +- d_phi
delta_theta=DEG2RAD*(rand01()*d_phi - (d_phi*0.5));
double alpha,beta;
alpha=M_PI_2-psi_i;
/********************************************************************************************
The grating is lamellar:
********************************************************************************************/
/* Lamellar grating is used.
Thus finding outgoing angle only from ingoing angle.
using old k and n and rotation matrix.*/
beta=asin(-order*(2*M_PI/k)/d+sin(alpha));
//Could have also used an angle convention where the angle between the grating's norm and the ray of order m is positive if it is on the samer side as the incident angle, negative otherwise.
if(phi0){
/*a target angle is specified*/
psi_o=phi0*DEG2RAD;
}else{
psi_o=M_PI_2-fabs(beta);/*the fabs is necessary due to the sign convention of beta*/
}
kx_notouch=kx;
ky_notouch=ky;
kz_notouch=kz;
kx=kx_notouch;
if(!isnan(psi_o) && !isnan(psi_i)){
//Still need to take into account the case where the incident ray does not hit perpendicular to the grating
if(beta>=0){
if(kz>=0){
ky=(ky_notouch*cos(-(psi_o-psi_i) + delta_theta) - kz_notouch*sin(-(psi_o-psi_i) + delta_theta));
kz=(ky_notouch*sin(-(psi_o-psi_i) + delta_theta) + kz_notouch*cos(-(psi_o-psi_i) + delta_theta));
}
else{
ky=(ky_notouch*cos((psi_o-psi_i) + delta_theta) - kz_notouch*sin((psi_o-psi_i) + delta_theta));
kz=(ky_notouch*sin((psi_o-psi_i) + delta_theta) + kz_notouch*cos((psi_o-psi_i) + delta_theta));
}
}
else{
if(kz>=0){
ky=(ky_notouch*cos((-2*alpha+(psi_o-psi_i)) + delta_theta) - kz_notouch*sin((-2*alpha+(psi_o-psi_i)) + delta_theta));
kz=(ky_notouch*sin((-2*alpha+(psi_o-psi_i)) + delta_theta) + kz_notouch*cos((-2*alpha+(psi_o-psi_i)) + delta_theta));
}
else{
ky=(ky_notouch*cos((2*alpha-(psi_o-psi_i)) + delta_theta) - kz_notouch*sin((2*alpha-(psi_o-psi_i)) + delta_theta));
kz=(ky_notouch*sin((2*alpha-(psi_o-psi_i)) + delta_theta) + kz_notouch*cos((2*alpha-(psi_o-psi_i)) + delta_theta));
}
}
/*Finding the weight using diffraction theory.*/
s=scalar_prod(kx,ky,kz,nx,ny,nz);
psi_o = acos(s/k)-M_PI_2;
/*pay attention to sign conventions here. We might cross to the "incoming" side of the normal.*/
spsi_i = sin(M_PI_2-psi_i);
spsi_o = sin(M_PI_2-psi_o);//+delta_theta)); // asin angle out
/* Phase for interference pattern: k = wave vector, d = lines/Angstrom, b=width in Angstrom */
gamma = k*d*(spsi_o-spsi_i);
/* Phase for diffraction pattern */
Gamma = k*b*(spsi_o-spsi_i);
pdiff=((sin(Gamma/2)/(Gamma/2))*(sin(Gamma/2)/(Gamma/2)))*((sin(nslits*gamma/2)/sin(gamma/2))*(sin(nslits*gamma/2)/sin(gamma/2)));
p=p*pdir*(pdiff/(pow(nslits,2)));
}else{
p=0;
}
}
%}
/************************************************************************************************************
Making lines to illustrace with the TRACE option.
*************************************************************************************************************/
MCDISPLAY
%{
magnify("");
line(-xwidth/2.0,0,-zdepth/2.0, xwidth/2.0,0,-zdepth/2.0);
line(-xwidth/2.0,0, zdepth/2.0, xwidth/2.0,0, zdepth/2.0);
line(-xwidth/2.0,0,-zdepth/2.0,-xwidth/2.0,0, zdepth/2.0);
line( xwidth/2.0,0,-zdepth/2.0, xwidth/2.0,0, zdepth/2.0);
%}
END
|
