Skip to content

Commit 62c94b2

Browse files
ebknudsenebknudsen
committed
union system updated for McXtrace 3
Including new Compton and Powder scattering process
1 parent ba683bc commit 62c94b2

18 files changed

+14764
-30
lines changed

‎mcxtrace-comps/share/Geometry_functions.c‎

Lines changed: 3844 additions & 0 deletions
Large diffs are not rendered by default.

‎mcxtrace-comps/share/Union_functions.c‎

Lines changed: 4574 additions & 0 deletions
Large diffs are not rendered by default.

‎mcxtrace-comps/share/Union_initialization.c‎

Lines changed: 70 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,70 @@
1+
/*******************************************************************************
2+
*
3+
* McXtrace, photon ray-tracing package
4+
* Copyright(C) 2007 Risoe National Laboratory.
5+
*
6+
* %I
7+
* Written by: Mads Bertelsen
8+
* Date: 20.08.15
9+
* Version: $Revision: 0.1 $
10+
* Origin: Svanevej 19
11+
*
12+
* A sample component to separate geometry and phsysics
13+
*
14+
* %D
15+
* Alpha version, no input system yet
16+
* Hardcode input to geometry engine
17+
* Allows complicated geometry by combination of simple shapes
18+
*
19+
* Algorithm:
20+
* Described elsewhere
21+
*
22+
* %P
23+
* INPUT PARAMETERS:
24+
* radius: [m] Outer radius of sample in (x,z) plane
25+
*
26+
* OUTPUT PARAMETERS:
27+
* V_rho: [AA^-3] Atomic density
28+
*
29+
* %L
30+
* The test/example instrument <a href="../examples/Test_Phonon.instr">Test_Phonon.instr</a>.
31+
*
32+
* %E
33+
******************************************************************************/
34+
35+
36+
// This file sets up global lists needed for Union components to communicate with each other
37+
// These all have dynamically allocated memory somewhere in the structure, which is deallocated
38+
// by the last Union_master.
39+
40+
41+
// Initialize global positions / rotations to transform lists
42+
// These are lists of pointers to positons / rotations, that will be updated from global frame
43+
// to the frame of the master component that uses them in that masters initialize section.
44+
struct global_positions_to_transform_list_struct global_positions_to_transform_list = {0,NULL};
45+
struct global_rotations_to_transform_list_struct global_rotations_to_transform_list = {0,NULL};
46+
47+
// Initialize global_process_list
48+
// Used to facilitate communication between processes and the other types of Union components
49+
struct pointer_to_global_process_list global_process_list = {0,NULL};
50+
51+
// Initialize global_material_list
52+
// Used to facilitate communication between materials and the other types of Union components
53+
struct pointer_to_global_material_list global_material_list = {0,NULL};
54+
55+
// Initialize global_geometry_list
56+
// Used to facilitate communication between geometries and the other types of Union components
57+
struct pointer_to_global_geometry_list global_geometry_list = {0,NULL};
58+
59+
// Initialize global_logger_lists
60+
// Used to facilitate communication between loggers and the other types of Union components
61+
struct pointer_to_global_logger_list global_all_volume_logger_list = {0,NULL};
62+
struct pointer_to_global_logger_list global_specific_volumes_logger_list = {0,NULL};
63+
64+
// Initialize global_tagging_conditional_list
65+
// Used to facilitate communication between conditionals and the other types of Union components
66+
struct global_tagging_conditional_list_struct global_tagging_conditional_list = {0,0,NULL};
67+
68+
// Initialize global_master_list
69+
// Used to facilitate communication between Master components (mainly for deallocation)
70+
struct pointer_to_global_master_list global_master_list = {0,NULL};

‎mcxtrace-comps/share/Union_last_functions.c‎

Lines changed: 113 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,113 @@
1+
/*
2+
TODO
3+
4+
update Union master to use the flexible functions
5+
6+
deal with loggers and conditionals
7+
*/
8+
9+
int physics_my(enum process choice, double *my,double *k_initial, union data_transfer_union data_transfer, struct focus_data_struct *focus_data, _class_particle *_particle) {
10+
11+
int output = 0; // Error return value
12+
#ifdef PROCESS_DETECTOR
13+
switch(choice) {
14+
#ifdef PROCESS_INCOHERENT_DETECTOR
15+
case Incoherent:
16+
output = Incoherent_physics_my(my, k_initial, data_transfer, focus_data, _particle);
17+
break;
18+
#endif
19+
#ifdef PROCESS_COMPTON_XRL_DETECTOR
20+
case Compton_xrl:
21+
output = Compton_xrl_physics_my(my, k_initial, data_transfer, focus_data, _particle);
22+
break;
23+
#endif
24+
#ifdef PROCESS_KN_XRL_DETECTOR
25+
case KN_xrl:
26+
output = KN_xrl_physics_my(my, k_initial, data_transfer, focus_data, _particle);
27+
break;
28+
#endif
29+
#ifdef PROCESS_RAYLEIGH_XRL_DETECTOR
30+
case Rayleigh_xrl:
31+
output = Rayleigh_xrl_physics_my(my, k_initial, data_transfer, focus_data, _particle);
32+
break;
33+
#endif
34+
#ifdef PROCESS_POWDER_DETECTOR
35+
case Powder:
36+
output = Powder_physics_my(my, k_initial, data_transfer, focus_data, _particle);
37+
break;
38+
#endif
39+
#ifdef PROCESS_SINGLE_CRYSTAL_DETECTOR
40+
case Single_crystal:
41+
output = Single_crystal_physics_my(my, k_initial, data_transfer, focus_data, _particle);
42+
break;
43+
#endif
44+
#ifdef PROCESS_NCRYSTAL_DETECTOR
45+
case NCrystal:
46+
output = NCrystal_physics_my(my, k_initial, data_transfer, focus_data, _particle);
47+
break;
48+
#endif
49+
#ifdef PROCESS_TEMPLATE_DETECTOR
50+
case Template:
51+
output = Template_physics_my(my, k_initial, data_transfer, focus_data, _particle);
52+
break;
53+
#endif
54+
default:
55+
printf("physics_my: No scattering process matches input!\n");
56+
break;
57+
}
58+
#endif
59+
return output;
60+
}
61+
62+
63+
int physics_scattering(enum process choice, double *k_final, double *k_initial, double *weight, union data_transfer_union data_transfer, struct focus_data_struct *focus_data, _class_particle *_particle) {
64+
65+
int output = 0; // Error return value
66+
#ifdef PROCESS_DETECTOR
67+
switch(choice) {
68+
#ifdef PROCESS_INCOHERENT_DETECTOR
69+
case Incoherent:
70+
output = Incoherent_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
71+
break;
72+
#endif
73+
#ifdef PROCESS_COMPTON_XRL_DETECTOR
74+
case Compton_xrl:
75+
output = Compton_xrl_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
76+
break;
77+
#endif
78+
#ifdef PROCESS_KN_XRL_DETECTOR
79+
case KN_xrl:
80+
output = KN_xrl_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
81+
break;
82+
#endif
83+
#ifdef PROCESS_RAYLEIGH_XRL_DETECTOR
84+
case Rayleigh_xrl:
85+
output = Rayleigh_xrl_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
86+
break;
87+
#endif
88+
#ifdef PROCESS_POWDER_DETECTOR
89+
case Powder:
90+
output = Powder_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
91+
break;
92+
#endif
93+
#ifdef PROCESS_SINGLE_CRYSTAL_DETECTOR
94+
case Single_crystal:
95+
output = Single_crystal_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
96+
break;
97+
#endif
98+
#ifdef PROCESS_NCRYSTAL_DETECTOR
99+
case NCrystal:
100+
output = NCrystal_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
101+
break;
102+
#endif
103+
#ifdef PROCESS_TEMPLATE_DETECTOR
104+
case Template:
105+
output = Incoherent_physics_scattering(k_final, k_initial, weight, data_transfer, focus_data, _particle);
106+
break;
107+
#endif
108+
default: printf("physics_scattering: No scattering process matches input!\n");
109+
break;
110+
}
111+
#endif
112+
return output;
113+
}

‎mcxtrace-comps/union/Compton_xrl_process.comp‎

Lines changed: 205 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,205 @@
1+
/*******************************************************************************
2+
*
3+
* McXtrace, x-ray ray-tracing package
4+
* Copyright(C) 2007 Risoe National Laboratory.
5+
*
6+
* %I
7+
* Written by: Mads Bertelsen and Erik B Knudsen
8+
* Date: 20.08.15
9+
* Version: $Revision: 0.1 $
10+
* Origin: ESS DMSC & DTU Physics
11+
*
12+
* A sample component to separate geometry and phsysics
13+
*
14+
* %D
15+
*
16+
* This Union_process is based on the Incoherent.comp component originally written
17+
* by Kim Lefmann and Kristian Nielsen
18+
*
19+
* Part of the Union components, a set of components that work together and thus
20+
* sperates geometry and physics within McXtrace.
21+
* The use of this component requires other components to be used.
22+
*
23+
* 1) One specifies a number of processes using process components like this one
24+
* 2) These are gathered into material definitions using Union_make_material
25+
* 3) Geometries are placed using Union_box / Union_cylinder, assigned a material
26+
* 4) A Union_master component placed after all of the above
27+
*
28+
* Only in step 4 will any simulation happen, and per default all geometries
29+
* defined before the master, but after the previous will be simulated here.
30+
*
31+
* There is a dedicated manual available for the Union_components
32+
*
33+
* Algorithm:
34+
* Described elsewhere
35+
*
36+
* %P
37+
* INPUT PARAMETERS:
38+
* density: [g/cm^3]   Nominal density of material.
39+
* element: [str] The element (symbol) of the material. Overrides Z.
40+
* Z: [1] Atomic number.
41+
* OUTPUT PARAMETERS:
42+
*
43+
* %L
44+
* The test/example instrument <a href="../examples/Test_Phonon.instr">Test_Phonon.instr</a>.
45+
*
46+
* %E
47+
******************************************************************************/
48+
49+
DEFINE COMPONENT Compton_xrl_process
50+
DEFINITION PARAMETERS ()
51+
SETTING PARAMETERS (density=0,atomno=14, string element="", string init="init")
52+
OUTPUT PARAMETERS (This_process,Compton_xrl_storage)
53+
DEPENDENCY "-lxrl"
54+
/* X-ray parameters: (x,y,z,kx,ky,kz,phi,t,Ex,Ey,Ez,p) */
55+
56+
SHARE
57+
%{
58+
#ifndef Union
59+
#define Union $Revision: 0.8 $
60+
61+
%include "Union_functions.c"
62+
%include "Union_initialization.c"
63+
64+
#endif
65+
66+
/* If on windows this should likely be set like this
67+
* since xraylib.h is not (as a standard) installed
68+
* into an xraylib subdir*/
69+
#ifdef _WIN32
70+
#include <xraylib.h>
71+
#else
72+
#include <xraylib/xraylib.h>
73+
#endif
74+
75+
struct Compton_xrl_physics_storage_struct{
76+
char name[13];
77+
int Z;
78+
double rho;
79+
double sigma;
80+
};
81+
82+
int Compton_xrl_physics_my(double *my, double *k_i, union data_transfer_union data_transfer, struct focus_data_struct *focus_data, _class_particle *_particle){
83+
/*get the partial cross section using xraylib*/
84+
double k_length = sqrt(k_i[0]*k_i[0]+k_i[1]*k_i[1]+k_i[2]*k_i[2]);
85+
double E=K2E*k_length;
86+
87+
/*how to handle the fact that we may be dealing with polarized radiation - goes in data_transfer union?*/
88+
/*we have to somehow transfer input data to this function - is this also in data_transfer union? yes. It contains pointers to all the possible physics_storage structs.*/
89+
struct Compton_xrl_physics_storage_struct *p = data_transfer.pointer_to_a_Compton_xrl_physics_storage_struct;
90+
91+
/*call xraylib to get cross_section in cm^2/g - scale by density to get mu in cm^-1 - scale by 10^2 to get mu in m^-1*/
92+
p->sigma=CS_Compt(p->Z,E,NULL);
93+
*my=p->rho*p->sigma*100;
94+
return 1;
95+
}
96+
97+
int Compton_xrl_physics_scattering(double *k_f, double *k_i, double *weight, union data_transfer_union data_transfer, struct focus_data_struct *focus_data, _class_particle *_particle) {
98+
double k_length = sqrt(k_i[0]*k_i[0]+k_i[1]*k_i[1]+k_i[2]*k_i[2]);
99+
double E=K2E*k_length;
100+
double E0;
101+
102+
struct Compton_xrl_physics_storage_struct *p = data_transfer.pointer_to_a_Compton_xrl_physics_storage_struct;
103+
104+
/*pick a direction randomly*/
105+
Coords k_out;
106+
// Here is the focusing system in action, get a vector
107+
double solid_angle;
108+
focus_data->focusing_function(&k_out,&solid_angle,focus_data);
109+
NORM(k_out.x,k_out.y,k_out.z);
110+
*weight *= solid_angle*0.25/PI;
111+
112+
double theta,dsigma;
113+
theta = acos(scalar_prod(k_out.x,k_out.y, k_out.z,k_i[0], k_i[1],k_i[2])/k_length);
114+
115+
/*make the call to xraylib to get the partial cross section for the particular scattering triangle defined by k_i, k_f.
116+
We rely on a previous call to CS_Compt by _mu-function*/
117+
dsigma = DCS_Compt(p->Z,E,theta,NULL);
118+
*weight *= dsigma/p->sigma;
119+
E0=E;
120+
E=ComptonEnergy(E0,theta,NULL);
121+
k_length=E*E2K;
122+
k_f[0]=k_out.x*k_length;
123+
k_f[1]=k_out.y*k_length;
124+
k_f[2]=k_out.z*k_length;
125+
/*get */
126+
return 1;
127+
}
128+
129+
#ifndef PROCESS_DETECTOR
130+
#define PROCESS_DETECTOR dummy
131+
#endif
132+
133+
#ifndef PROCESS_COMPTON_XRL_DETECTOR
134+
#define PROCESS_COMPTON_XRL_DETECTOR dummy
135+
#endif
136+
%}
137+
138+
DECLARE
139+
%{
140+
// Needed for transport to the main component
141+
struct global_process_element_struct global_process_element;
142+
struct scattering_process_struct This_process;
143+
144+
// Declare for this component, to do calculations on the input / store in the transported data
145+
struct Compton_xrl_physics_storage_struct Compton_xrl_storage;
146+
147+
%}
148+
149+
INITIALIZE
150+
%{
151+
// Initialize done in the component
152+
xrl_error *xe;
153+
if(atomno==0 && strlen(element)!=0){
154+
Compton_xrl_storage.Z=SymbolToAtomicNumber(element,&xe);
155+
snprintf(Compton_xrl_storage.name,3,element);
156+
} else if (atomno>0 || atomno<116){
157+
Compton_xrl_storage.Z=atomno;
158+
snprintf(Compton_xrl_storage.name,3,AtomicNumberToSymbol(atomno,&xe));
159+
}else{
160+
fprintf(stderr,"ERROR: (%s): Must set either an element by symbol or by atomic number Got (%s,%d).\n",NAME_CURRENT_COMP,element,atomno);
161+
exit(-1);
162+
}
163+
Compton_xrl_storage.rho=density;
164+
165+
// The type of the process must be saved in the global enum process
166+
This_process.eProcess = Compton_xrl;
167+
168+
// Need to specify if this process is isotropic
169+
This_process.non_isotropic_rot_index = -1; // Yes (powder)
170+
//This_process.non_isotropic_rot_index = 1; // No (single crystal)
171+
172+
// Packing the data into a structure that is transported to the main component
173+
sprintf(This_process.name,NAME_CURRENT_COMP);
174+
This_process.process_p_interact = -1;
175+
This_process.data_transfer.pointer_to_a_Compton_xrl_physics_storage_struct = &Compton_xrl_storage;
176+
//This_process.data_transfer.pointer_to_a_Incoherent_physics_storage_struct->my_scattering = effective_my_scattering;
177+
This_process.probability_for_scattering_function = &Compton_xrl_physics_my;
178+
This_process.scattering_function = &Compton_xrl_physics_scattering;
179+
180+
// This will be the same for all process's, and can thus be moved to an include.
181+
sprintf(global_process_element.name,NAME_CURRENT_COMP);
182+
global_process_element.component_index = INDEX_CURRENT_COMP;
183+
global_process_element.p_scattering_process = &This_process;
184+
185+
if (_getcomp_index(init) < 0) {
186+
fprintf(stderr,"Compton_xrl_process:%s: Error identifying Union_init component, %s is not a known component name.\n",
187+
NAME_CURRENT_COMP, init);
188+
exit(-1);
189+
}
190+
191+
struct pointer_to_global_process_list *global_process_list = COMP_GETPAR3(Union_init, init, global_process_list);
192+
add_element_to_process_list(global_process_list, global_process_element);
193+
%}
194+
195+
TRACE
196+
%{
197+
%}
198+
199+
FINALLY
200+
%{
201+
// Since the process and it's storage is a static allocation, there is nothing to deallocate
202+
203+
%}
204+
205+
END

0 commit comments

Comments
 (0)