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
|
/**************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright 1997-2006, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* Component: Pol_monitor
*
* %I
* Written by: Peter Christiansen
* Modified by Erik B Knudsen
* Date: July 2006
* Origin: Risoe
*
* Polarisation sensitive monitor.
*
* %D A square single monitor that measures the projection of the
* polarisation along a given normalized m-vector (mx, my, mz).
* The measured quantity is: sx*mx+sy*my+mz*sz
*
* Example: Pol_monitor(xwidth=0.1, yheight=0.1, nchan=11, mx=0, my=1, mz=0)
*
* %P
* INPUT PARAMETERS:
*
* xwidth: [m] Width of detector
* yheight: [m] Height of detector
* mx: [1] X-component of monitor vector (can be negative)
* my: [1] Y-component of monitor vector (can be negative)
* mz: [1] Z-component of monitor vector (can be negative)
* restore_neutron: [1] If set, the monitor does not influence the neutron state
* nowritefile: [1] If set, monitor will skip writing to disk
*
* CALCULATED PARAMETERS:
*
* Pol_N: [] Array of neutron counts
* Pol_p: [] Array of neutron weight counts
* Pol_p2: [] Array of second moments
*
* %E
*************************************************************************/
DEFINE COMPONENT Pol_monitor
SETTING PARAMETERS (xwidth=0.1, yheight=0.1, int restore_neutron=0, mx=0, my=0, mz=0, int nowritefile=0)
/* Neutron parameters: (x,y,z,vx,vy,vz,t,sx,sy,sz,p) */
DECLARE
%{
double Pol_N;
double Pol_p;
double Pol_p2;
double Psum;
char titlestring[128];
%}
INITIALIZE
%{
// Check that input parameteters makes sense
if (mx==0 && my==0 && mz==0) {
fprintf(stderr, "Pol_monitor: %s: NULL vector defined!\n"
"ERROR (mx, my, mz). Exiting",
NAME_CURRENT_COMP);
exit(1);
}
if ((xwidth<=0) || (yheight <= 0)) {
fprintf(stderr, "Pol_monitor: %s: Null detection area !\n"
"ERROR (xwidth,yheight). Exiting",
NAME_CURRENT_COMP);
exit(1);
}
sprintf(titlestring, "Polarisation monitor m=(%g %g %g) %s", mx, my, mz, NAME_CURRENT_COMP);
// Initialize variables
NORM(mx, my, mz);
Psum=0;
Pol_N = 0;
Pol_p = 0;
Pol_p2 = 0;
%}
TRACE
%{
double pol_proj;
PROP_Z0;
if (inside_rectangle(x, y, xwidth, yheight)){
pol_proj = scalar_prod(mx, my, mz, sx, sy, sz);
if(fabs(pol_proj)>1) {
if (fabs(pol_proj)<1+FLT_EPSILON){
pol_proj=1;
}else{
ABSORB;
}
}
double p2 = p*p;
double pp = pol_proj*p;
double p2p = pol_proj*pol_proj*p;
#pragma acc atomic
Pol_N = Pol_N+1;
#pragma acc atomic
Psum = Psum+p;
#pragma acc atomic
Pol_p = Pol_p+pp;
#pragma acc atomic
Pol_p2 = Pol_p2+p2p;
SCATTER;
}
if (restore_neutron) {
RESTORE_NEUTRON(INDEX_CURRENT_COMP, x, y, z, vx, vy, vz, t, sx, sy, sz, p);
}
%}
SAVE
%{
if (!nowritefile) {
if (Psum && Pol_N){
Pol_p /= Psum;
#ifdef USE_MPI
Pol_p /= mpi_node_count;
#endif
Pol_p2 /= Psum;
Pol_p2 -= Pol_p*Pol_p;
Pol_p2 /= Pol_N;
}
DETECTOR_OUT_0D(titlestring, Pol_N, Pol_p, Pol_p2);
}
%}
MCDISPLAY
%{
rectangle("xy", 0, 0, 0, xwidth, yheight);
%}
END
|
