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
|
/*******************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright 1997-2002, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* Component: Scatter_log_iterator.comp
*
* %I
*
* Written by: Erik B Knudsen
* Date: November 2012
* Version: $Revision: 1.21 $
* Release: McStas 2.1
* Origin: DTU Physics
*
* Iteration element for a Scatter_log
*
* %D
*
* This component marks the beginning of the region in trace in which pseudo-neutrons
* are to be propagated. Pseudo-neutrons are neutrons which are generated by the function
* <i>compute_func</i> from before and after SCATTER neutron states which have been logged by
* a set of Scatter_logger/Scatter_logger_stop components.
*
* N.B. This component should be immediately preceeded by an Arm. Any components between this one
* and a subsequent Scatter_log_iterator_stop component will be visited by the set of pseudo-neutrons
* as if they were regular neutrons in the classical McStas-manner.
*
* %P
* Input parameters:
*
* compute_func: [] Address of the function that computes a psuedo neutron from before and after scatter states.
*
* %E
*******************************************************************************/
DEFINE COMPONENT Scatter_log_iterator
DEFINITION PARAMETERS (compute_func=NULL)
SETTING PARAMETERS ()
/* Neutron parameters: (x,y,z,vx,vy,vz,t,sx,sy,sz,p) */
SHARE
%{
/*This is the specialized pseudo-neutron function that computes
an escaping neutron from logged before and after SCATTER neutron states*/
int exit_neutron(double *ns_tilde, struct Generalized_State_t *S0, struct Generalized_State_t *S1){
/*!!Note that the transformation into global coordinate system must be done while logging
as we do not have access to neither the component name nor can get the component rotation by index.*/
Coords c1,c2;
Rotation R1,R2;
/*so now compute the pseudo neutron state and possibly user variables*/
/*position comes from "new" state*/
ns_tilde[0]=S1->_x;ns_tilde[1]=S1->_y;ns_tilde[2]=S1->_z;
/*velocity is the "old" state*/
ns_tilde[3]=S0->_vx;ns_tilde[4]=S0->_vy;ns_tilde[5]=S0->_vz;
/*time from new*/
ns_tilde[6]=S1->_t;
/*spin comes from "new" state*/
ns_tilde[7]=S1->_sx;ns_tilde[8]=S1->_sy;ns_tilde[9]=S1->_sz;
/*weight is difference old-new to mean the neutrons "deposited" in the guide wall*/
ns_tilde[10]=S0->_p-S1->_p;
return 0;
}
#define NOABS \
do {/* Nothing*/} while(0)
%}
DECLARE
%{
int (*pseudo_neutron_state_function) (double *, struct Generalized_State_t *, struct Generalized_State_t *);
struct Generalized_State_t *s1,*s0;
double *nstate_initial;
/*need a pointer to the structure set up by the logger*/
%}
INITIALIZE
%{
if (compute_func) {
pseudo_neutron_state_function=compute_func;
}else{
pseudo_neutron_state_function=exit_neutron;
}
nstate_initial=NULL;
%}
TRACE
%{
/*I am the start of the pseudo neutron iterator*/
if (nstate_initial==NULL){
double *ns=nstate_initial=calloc(11,sizeof(double));
ns[0]=x;ns[1]=y; ns[2]=z;
ns[3]=vx;ns[4]=vy;ns[5]=vz;
ns[6]=t;
ns[7]=sx;ns[8]=sy;ns[9]=sz;
ns[10]=p;
s0=Bounce_store;
s1=Bounce_store+1;
/* Remove std. ABSORB to avoid breaking analysis loop */
#undef mcabsorb
#define mcabsorb scatter_iterator_stop
}
if (s1->_p!=-1){
/*a neutron weight of -1 is nonsensical. I.e. if s1->p>=0, it means there are two states in the log from which we can compute a pseudo-neutron*/
double nstate[11];
if ( pseudo_neutron_state_function(nstate,s0,s1) ){
printf("Warning: (%s): error reported when computing pseudo neutron\n",NAME_CURRENT_COMP);
}
/*set neutron state for subsequent components*/
x=nstate[0];y=nstate[1];z=nstate[2];
vx=nstate[3];vy=nstate[4];vz=nstate[5];
t=nstate[6];
sx=nstate[7];sy=nstate[8];sz=nstate[9];
p=nstate[10];
s0++;
s1++;
}else if (Bounce_store[1]._p==-1){
x=s0->_x;y=s0->_y;z=s0->_z;
vx=s0->_vx;vy=s0->_vy;vz=s0->_vz;
t=s0->_t;
sx=s0->_sx;sy=s0->_sy;sz=s0->_sz;
p=s0->_p;
} else {
fprintf(stderr,"This should not happen. Period.\n");
exit(1);
}
%}
MCDISPLAY
%{
/* A bit ugly; hard-coded dimensions. */
line(0,0,0,0.2,0,0);
line(0,0,0,0,0.2,0);
line(0,0,0,0,0,0.2);
%}
END
|
