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
|
/************************************************
* McXtrace, x-ray tracing package
* Copyright, All rights reserved
* DTU Physics, Kgs. Lyngby, Denmark
* Synchrotron SOLEIL, Saint-Aubin, France
*
* %Identification
* Written by: Erik Knudsen
* Date: September 25th, 2009
* Origin: Risoe
*
* A cylindrically curved mirror (in YZ)
*
* %Description
* mirror is in the YZ-plane curved towards positive X if radius is positive
*
* Example: Mirror_curved( radius=2, length=20e-3, width=40e-3, coating="AlMgF2_disco.dat")
*
* %Parameters
* radius: [m] Radius of curvature, along Z.
* yheight:[m] Width of the mirror along Y.
* zdepth: [m] Length of the unbent mirror along Z.
* coating:[str] Name of file containing the material data (i.e. f1 and f2) for the coating
* R0: [1] Constant reflectivity value (mostly relevant for debugging, when coating="").
* width: [m] Width of the mirror along Y = yheight
* length: [m] Length of the unbent mirror along Z = zdepth
*
* %End
***********************************************************************/
DEFINE COMPONENT Mirror_curved
SETTING PARAMETERS (radius=1, R0=1, string coating="Be.txt",
zdepth=0.2, yheight=0.2,
length=0, width=0)
/* X-ray parameters: (x,y,z,kx,ky,kz,phi,t,Ex,Ey,Ez,p) */
SHARE
%{
#include <complex.h>
%include "read_table-lib"
%include "reflectivity-lib"
%}
DECLARE
%{
int Z;
double At;
double rho;
t_Table T;
t_Reflec re;
%}
INITIALIZE
%{
int status;
if (coating && strlen(coating)) {
status=reflec_Init(&re,COATING_UNDEFINED,coating,NULL);
}else{
/*assume a constant reflectivity*/
status=reflec_Init(&re,CONSTANT,NULL, &(R0));
}
if (!radius) radius=1e6;
// compatibility with older use
if (length) zdepth =length;
if (width) yheight=width;
%}
TRACE
%{
double l0,l1,dl,alpha,n,nx,nz,s,k,knx,knz;
/*do we hit the mirror within range?*/
//PROP_Z0;
k=sqrt(scalar_prod(kx,0,kz,kx,0,kz));
knx=kx/k;
knz=kz/k;
if (solve_2nd_order(&l0,&l1, knx*knx+knz*knz,2*(x*knx+z*knz-radius*knx),x*x-2*x*radius+z*z))
{
if (l0<0 && l1>0){
dl=l1;
}else if(l1<0 && l0>0){
dl=l0;
}else if( l0>0 && l1>0){
if (fabs(z+knz*l0) <zdepth/2.0 && fabs(x+knx*l0)<fabs(radius)) {
/*this means that the first positive match is on the mirror z-wise
- this should be picked. This would work even if mirror is hit from behind*/
dl=l0;
}else{
dl=l1;
}
}else{
/*particle misses the mirror completely (l0<0 and l1<0)*/
RESTORE_XRAY(INDEX_CURRENT_COMP, x, y, z, kx, ky, kz, phi, t, Ex, Ey, Ez, p);
}
/*correction for only solving in xz-plane*/
dl*=sqrt(scalar_prod(kx,ky,kz,kx,ky,kz)/scalar_prod(kx,0,kz,kx,0,kz));
PROP_DL(dl);
alpha=asin(z/radius);
if ( (zdepth<radius*alpha) || (fabs(y)>yheight/2.0)) {
RESTORE_XRAY(INDEX_CURRENT_COMP, x, y, z, kx, ky, kz, phi, t, Ex, Ey, Ez, p);
}else{
SCATTER;
/*reflection*/
nx=radius-x;
nz=-z;
n=sqrt(nx*nx+nz*nz);
nx/=n;
nz/=n;
s=scalar_prod(kx,0,kz,nx,0,nz);
kx=kx-2*s*nx;
kz=kz-2*s*nz;
double Q;
/*length of wavevector transfer may be compute from s and n_ above*/
Q=fabs(2*s*sqrt(nx*nx+nz*nz));
double complex R=refleccq(re,Q,0,k,fabs(90-acos(s/k)*RAD2DEG));
p*=sqrt(creal(R*conj(R)));
phi+=atan2(cimag(R),creal(R));
}
}else if(fabs(kx)<DBL_EPSILON && fabs(kz)<DBL_EPSILON && ky!=0){
/*This to catch the extreme case where k//y.*/
if (y==0){
/*We hit the "side" of the mirror.*/
ABSORB;
}
}
%}
MCDISPLAY
%{
double x0,y0,z0,x1,y1,z1,alpha;
int N=20;
y0=yheight/2.0;
y1=-y0;
alpha=-(zdepth/2.0)/fabs(radius);
x0=radius*(1-cos(alpha));
z0=radius*sin(alpha);
line(x0,y0,z0,x0,y1,z0);
while (alpha<(zdepth/2.0/fabs(radius))){
alpha+=zdepth/N/fabs(radius);
x1=radius*(1-cos(alpha));
z1=radius*sin(alpha);
line(x0,y0,z0,x1,y0,z1);
line(x0,y1,z0,x1,y1,z1);
x0=x1;z0=z1;
line(x0,y0,z0,x0,y1,z0);
}
%}
END
|
