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
|
/*******************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright (C) 1997-2010, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* Component: Lens_simple
*
* %I
* Written by: Henrich Frielinghaus
* Date: June 2010
* Origin: FZ Juelich
*
* Rectangular/circular slit with parabolic/spherical LENS.
*
* %D
* A simple rectangular or circular slit. You may either
* specify the radius (circular shape), or the rectangular bounds.
* No transmission around the slit is allowed.
*
* At the slit position there is/are also a LENS or multiple LENSES present.
* Spherical/parabolic abberation is taken into account in a simplified
* manner. The focal point becomes a function of the distance ss from the
* origin of the lens where the ray hits the lens (plane).
* With this focal distance the refraction is calculated based on simple
* geometrical optics (as tought already in school).
*
* The approximation gets wrong when the distance ss is too big, and
* total reflection becomes possible. Then the corrections of the focal
* distance are strong (see foc*= ... lines. When this correction factor
* is too close to zero, the corrections become highly wrong).
*
* The advantage of this routine is the simplicity which makes the
* simulation very fast - especially for multiple lenses. The argument
* for this way of treating the lenses is that the collimation and
* detector distance are large (say 20m) compared to the lens thickness
* (say 2-5cm) and the lens array thickness (say max. 1m).
*
* For lens arrays one still can simulate several of these simplified
* lenses instead of an exact treatment (because 5cm<<1m). The author
* believes that this medium detailed treatment meets usually the
* desired accuracy.
*
* Example: Lens_simple(rho=5.16e10,Rc=0.0254,Nl=7,xmin=-0.01,xmax=0.01,ymin=-0.01,ymax=0.01)
* Lens_simple(rho=5.16e10,Rc=0.0254,Nl=7,radius=0.01)
* Lens_simple(rho=5.16e10,Rc=0.0254,Nl=7,radius=0.035,SigmaAL=0.141, d0=0.002)
* ^^^^^^^^^^^ last example includes absorption of MgF2-lens.
*
* %P
* INPUT PARAMETERS
*
* rho: [m-2] Scattering length density
* Rc: [m] Radius (concave: Rc>0) of biconcave lens
* Nl: [1] Number of single lenses
* parab: [] Switch (not 0 -> parabolic, for 0 spherical)
* radius: [m] Radius of slit in the z=0 plane, centered at Origin
* xmin: [m] Lower x bound
* xmax: [m] Upper x bound
* ymin: [m] Lower y bound
* ymax: [m] Upper y bound
* SigmaAL: Absorption cross section per lambda (wavelength) [1/(m*AA)]
* d0: [m] Minimum thickness in the centre of the lens
*
* %E
*******************************************************************************/
DEFINE COMPONENT Lens_simple
SETTING PARAMETERS (rho=5.16e14, Rc=0.02, Nl=7.0, parab=1.1,
xmin=0, xmax=0, ymin=0, ymax=0, radius=0,
SigmaAL=0.141, d0=0.002)
// (x,y,z,vx,vy,vz,t,s1,s2,p)
INITIALIZE
%{
if ( (xmin >= xmax || ymin >= ymax) && radius == 0)
{ fprintf(stderr,"Lens_simple: %s: Error: give geometry\n", NAME_CURRENT_COMP); exit(-1); }
%}
TRACE
%{
double vvv = vx*vx + vy*vy + vz*vz;
double Xi = rho/PI *(4e-20*PI*PI)/(V2Q*V2Q*vvv);
double foc = Rc/Xi/Nl;
if (z>0e0 || vvv==0e0) ABSORB;
PROP_Z0;
double ss = x*x + y*y;
if (((radius == 0) && (x<xmin || x>xmax || y<ymin || y>ymax))
|| ((radius != 0) && (ss > radius*radius)))
ABSORB;
else
SCATTER;
if (parab==0e0 && ss >= Rc*Rc) ABSORB;
if (parab!=0e0)
foc*= 1e0 - 0.5*Nl*Xi*(1e0-0.5*ss/(Rc*Rc));
else
foc*= (1e0 - 0.5*Nl*Xi)*sqrt(1e0-ss/(Rc*Rc));
double tt2 = (-0.7*fabs(foc))/vz;
double xx2 = x + vx*tt2;
double yy2 = y + vy*tt2;
double zz2 = -0.7*fabs(foc);
double ll1 = -zz2;
double ll2 = 1.0/(1.0/foc-1.0/ll1);
double llr = -ll2/ll1;
xx2*= llr;
yy2*= llr;
zz2*= llr;
xx2 = xx2 - x;
yy2 = yy2 - y;
double zdir = zz2/fabs(zz2);
double xyzlen = xx2*xx2 + yy2*yy2 + zz2*zz2;
vx = xx2*zdir*sqrt(vvv/xyzlen);
vy = yy2*zdir*sqrt(vvv/xyzlen);
vz = zz2*zdir*sqrt(vvv/xyzlen);
double thck;
if (parab!=0e0)
thck = ss/Rc + d0;
else
thck = 2.0*(Rc-sqrt(Rc*Rc-ss)) + d0;
p*=exp(-Nl*SigmaAL*(2.0*PI/(V2Q*sqrt(vvv)))*thck); // Transmission //
%}
MCDISPLAY
%{
double xw, yh;
xw = (xmax - xmin)/2.0;
yh = (ymax - ymin)/2.0;
multiline(3, xmin-xw, (double)ymax, 0.0,
(double)xmin, (double)ymax, 0.0,
(double)xmin, ymax+yh, 0.0);
multiline(3, xmax+xw, (double)ymax, 0.0,
(double)xmax, (double)ymax, 0.0,
(double)xmax, ymax+yh, 0.0);
multiline(3, xmin-xw, (double)ymin, 0.0,
(double)xmin, (double)ymin, 0.0,
(double)xmin, ymin-yh, 0.0);
multiline(3, xmax+xw, (double)ymin, 0.0,
(double)xmax, (double)ymin, 0.0,
(double)xmax, ymin-yh, 0.0);
%}
END
|
