1 /* This file is part of the hkl library.
3 * The hkl library is free software: you can redistribute it and/or modify
4 * it under the terms of the GNU General Public License as published by
5 * the Free Software Foundation, either version 3 of the License, or
6 * (at your option) any later version.
8 * The hkl library is distributed in the hope that it will be useful,
9 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 * GNU General Public License for more details.
13 * You should have received a copy of the GNU General Public License
14 * along with the hkl library. If not, see <http://www.gnu.org/licenses/>.
16 * Copyright (C) 2003-2015 Synchrotron SOLEIL
17 * L'Orme des Merisiers Saint-Aubin
18 * BP 48 91192 GIF-sur-YVETTE CEDEX
20 * Authors: Picca Frédéric-Emmanuel <picca@synchrotron-soleil.fr>
22 #include <gsl/gsl_multiroots.h>
23 #include "hkl-factory-private.h" // for autodata_factories_, etc
24 #include "hkl-axis-private.h" // for HklAxis
25 #include "hkl-pseudoaxis-common-hkl-private.h" // for hkl_engine_hkl_new, etc
26 #include "hkl-pseudoaxis-common-q-private.h" // for hkl_engine_q2_new, etc
27 #include "hkl-pseudoaxis-common-tth-private.h" // for hkl_engine_tth2_new, etc
31 /*********************/
32 /* MED 2+2 HklEngine */
33 /*********************/
35 static int _reflectivity_func(const gsl_vector
*x
, void *params
, gsl_vector
*f
)
37 const double mu
= x
->data
[0];
38 const double gamma
= x
->data
[2];
40 CHECK_NAN(x
->data
, x
->size
);
42 RUBh_minus_Q(x
->data
, params
, f
->data
);
43 f
->data
[3] = gamma
- 2 * mu
;
48 static const HklFunction reflectivity_func
= {
49 .function
= _reflectivity_func
,
53 static HklMode
* mu_fixed_2_2()
55 static const char* axes_r
[] = {"beta", "mu", "omega", "gamma", "delta"};
56 static const char* axes_w
[] = {"omega", "gamma", "delta"};
57 static const HklFunction
*functions
[] = {&RUBh_minus_Q_func
};
58 static const HklModeAutoInfo info
= {
59 HKL_MODE_AUTO_INFO("mu_fixed", axes_r
, axes_w
, functions
),
62 return hkl_mode_auto_new(&info
,
63 &hkl_full_mode_operations
,
67 static HklMode
* reflectivity_2_2()
69 static const char* axes_r
[] = {"beta", "mu", "omega", "gamma", "delta"};
70 static const char* axes_w
[] = {"mu", "omega", "gamma", "delta"};
71 static const HklFunction
*functions
[] = {&reflectivity_func
};
72 static const HklModeAutoInfo info
= {
73 HKL_MODE_AUTO_INFO("reflectivity", axes_r
, axes_w
, functions
),
76 return hkl_mode_auto_new(&info
,
77 &hkl_full_mode_operations
,
81 static HklEngine
*hkl_engine_soleil_sixs_med_2_2_hkl_new(HklEngineList
*engines
)
84 HklMode
*default_mode
;
86 self
= hkl_engine_hkl_new(engines
);
88 default_mode
= mu_fixed_2_2();
89 hkl_engine_add_mode(self
, default_mode
);
90 hkl_engine_mode_set(self
, default_mode
);
92 hkl_engine_add_mode(self
, reflectivity_2_2());
97 /*********************/
98 /* MED 1+2 HklEngine */
99 /*********************/
101 static HklMode
* pitch_fixed()
103 static const char *axes_r
[] = {"pitch", "mu", "gamma", "delta"};
104 static const char* axes_w
[] = {"mu", "gamma", "delta"};
105 static const HklFunction
*functions
[] = {&RUBh_minus_Q_func
};
106 static const HklModeAutoInfo info
= {
107 HKL_MODE_AUTO_INFO(__func__
, axes_r
, axes_w
, functions
),
110 return hkl_mode_auto_new(&info
,
111 &hkl_full_mode_operations
,
115 static HklEngine
*hkl_engine_soleil_sixs_med_1_2_hkl_new(HklEngineList
*engines
)
118 HklMode
*default_mode
;
120 self
= hkl_engine_hkl_new(engines
);
122 default_mode
= pitch_fixed();
123 hkl_engine_add_mode(self
, default_mode
);
124 hkl_engine_mode_set(self
, default_mode
);
129 /*********************/
130 /* MED 2+3 HklEngine */
131 /*********************/
133 typedef struct _HklSlitsFit HklSlitsFit
;
136 HklGeometry
*geometry
;
138 unsigned int slits_id
;
143 static int slits_func(const gsl_vector
*x
, void *params
, gsl_vector
*f
)
145 double const *x_data
= gsl_vector_const_ptr(x
, 0);
146 double *f_data
= gsl_vector_ptr(f
, 0);
147 HklVector n_slits
= {{0, 0, 1}};
148 HklSlitsFit
*parameters
= params
;
150 hkl_parameter_value_set(parameters
->axis
, x_data
[0], HKL_UNIT_DEFAULT
, NULL
);
151 hkl_geometry_update(parameters
->geometry
);
153 /* compute the orientation of the slits */
154 hkl_vector_rotated_quaternion(&n_slits
,
155 &darray_item(parameters
->geometry
->holders
, 1)->q
);
157 /* both directions must be perpendicular */
158 f_data
[0] = hkl_vector_scalar_product(¶meters
->surface
, &n_slits
);
163 static int fit_slits_orientation(HklSlitsFit
*params
)
166 gsl_multiroot_fsolver_type
const *T
;
167 gsl_multiroot_fsolver
*s
;
168 gsl_multiroot_function f
;
175 /* now solve the system */
176 /* Initialize method */
177 T
= gsl_multiroot_fsolver_hybrid
;
178 s
= gsl_multiroot_fsolver_alloc (T
, params
->len
);
179 x
= gsl_vector_alloc(params
->len
);
180 x_data
= gsl_vector_ptr(x
, 0);
182 /* initialize x with the right values */
183 x_data
[0] = params
->axis
->_value
;
188 gsl_multiroot_fsolver_set (s
, &f
, x
);
190 /* iterate to find the solution */
194 status
= gsl_multiroot_fsolver_iterate(s
);
195 if (status
|| iter
% 100 == 0) {
196 /* Restart from another point. */
197 for(i
=0; i
<params
->len
; ++i
)
198 x_data
[i
] = (double)rand() / RAND_MAX
* 180. / M_PI
;
199 gsl_multiroot_fsolver_set(s
, &f
, x
);
200 gsl_multiroot_fsolver_iterate(s
);
202 status
= gsl_multiroot_test_residual (s
->f
, HKL_EPSILON
);
203 } while (status
== GSL_CONTINUE
&& iter
< 1000);
206 fprintf(stdout
, "\n fitting the detector position using thoses axes :");
207 for(i
=0; i
<params
->len
; ++i
)
208 fprintf(stdout
, " \"%s\"", params
->axis
->name
);
209 fprintf(stdout
, " status : %d iter : %d", status
, iter
);
210 fprintf(stdout
, " x: [");
211 for(i
=0; i
<params
->len
; ++i
)
212 fprintf(stdout
, " %.7f", s
->x
->data
[i
]);
213 fprintf(stdout
, "] f: [");
214 for(i
=0; i
<params
->len
; ++i
)
215 fprintf(stdout
, " %.7f", s
->f
->data
[i
]);
216 fprintf(stdout
, "]\n");
217 hkl_geometry_fprintf(stdout
, params
->geometry
);
219 if(status
!= GSL_CONTINUE
){
221 /* put the axes in the -pi, pi range. */
222 gsl_sf_angle_restrict_pos_e(¶ms
->axis
->_value
);
226 gsl_multiroot_fsolver_free(s
);
231 static void hkl_geometry_list_multiply_soleil_sixs_med_2_3(HklGeometryList
*self
,
232 HklGeometryListItem
*item
)
237 HklGeometry
*geometry
;
238 double slits_position
;
239 HklHolder
*sample_holder
;
240 HklHolder
*detector_holder
;
242 /* For each solution already found we will generate another one */
243 /* we will set the right slit orientation for a given detector arm position */
244 geometry
= item
->geometry
;
245 sample_holder
= darray_item(geometry
->holders
, 0);
246 detector_holder
= darray_item(geometry
->holders
, 1);
248 /* get the index of the axis corresponding to the slits */
249 /* for now the last holder is the detector one */
250 params
.slits_id
= detector_holder
->config
->idx
[detector_holder
->config
->len
-1];
251 params
.len
= 1; /* only one axis to fit */
252 params
.geometry
= geometry
;
253 params
.axis
= darray_item(params
.geometry
->axes
, params
.slits_id
);
255 /* compute the surface orientation fixed during the fit */
256 /* use the last sample axis as sample surface normal */
257 params
.surface
= container_of(darray_item(geometry
->axes
,
258 sample_holder
->config
->idx
[sample_holder
->config
->len
- 1]),
259 HklAxis
, parameter
)->axis_v
;
260 hkl_vector_rotated_quaternion(¶ms
.surface
,
264 /* we just need to fit the slits orientation */
265 /* save it's value before */
266 slits_position
= hkl_parameter_value_get(params
.axis
, HKL_UNIT_DEFAULT
);
267 if (fit_slits_orientation(¶ms
) != TRUE
)
268 hkl_parameter_value_set(params
.axis
, slits_position
, HKL_UNIT_DEFAULT
, NULL
);
271 static HklMode
* mu_fixed_2_3()
273 static const char *axes_r
[] = {"beta", "mu", "omega", "gamma", "delta", "eta_a"};
274 static const char* axes_w
[] = {"omega", "gamma", "delta"};
275 static const HklFunction
*functions
[] = {&RUBh_minus_Q_func
};
276 static const HklModeAutoInfo info
= {
277 HKL_MODE_AUTO_INFO("mu_fixed", axes_r
, axes_w
, functions
),
280 return hkl_mode_auto_new(&info
,
281 &hkl_full_mode_operations
,
285 static HklEngine
*hkl_engine_soleil_sixs_med_2_3_hkl_new(HklEngineList
*engines
)
288 HklMode
*default_mode
;
290 self
= hkl_engine_hkl_new(engines
);
292 default_mode
= mu_fixed_2_3();
293 hkl_engine_add_mode(self
, default_mode
);
294 hkl_engine_mode_set(self
, default_mode
);
299 /***********************/
300 /* SOLEIL SIXS MED 2+2 */
301 /***********************/
303 #define HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_2_2_DESCRIPTION \
304 "+ xrays source fix allong the :math:`\\vec{x}` direction (1, 0, 0)\n" \
305 "+ 3 axes for the sample\n" \
307 " + **beta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
308 " + **mu** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
309 " + **omega** : rotating around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
311 "+ 3 axis for the detector\n" \
313 " + **beta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
314 " + **gamma** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
315 " + **delta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n"
317 static const char* hkl_geometry_soleil_sixs_med_2_2_axes
[] = {"beta", "mu", "omega", "gamma", "delta"};
319 static HklGeometry
*hkl_geometry_new_soleil_sixs_med_2_2(const HklFactory
*factory
)
321 HklGeometry
*self
= hkl_geometry_new(factory
);
324 h
= hkl_geometry_add_holder(self
);
325 hkl_holder_add_rotation_axis(h
, "beta", 0, -1, 0);
326 hkl_holder_add_rotation_axis(h
, "mu", 0, 0, 1);
327 hkl_holder_add_rotation_axis(h
, "omega", 0, -1, 0);
329 h
= hkl_geometry_add_holder(self
);
330 hkl_holder_add_rotation_axis(h
, "beta", 0, -1, 0);
331 hkl_holder_add_rotation_axis(h
, "gamma", 0, 0, 1);
332 hkl_holder_add_rotation_axis(h
, "delta", 0, -1, 0);
337 static HklEngineList
*hkl_engine_list_new_soleil_sixs_med_2_2(const HklFactory
*factory
)
339 HklEngineList
*self
= hkl_engine_list_new();
341 hkl_engine_soleil_sixs_med_2_2_hkl_new(self
);
342 hkl_engine_q2_new(self
);
343 hkl_engine_qper_qpar_new(self
);
344 hkl_engine_tth2_new(self
);
349 REGISTER_DIFFRACTOMETER(soleil_sixs_med_2_2
,"SOLEIL SIXS MED2+2", HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_2_2_DESCRIPTION
);
351 /***********************/
352 /* SOLEIL SIXS MED 1+2 */
353 /***********************/
355 #define HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_1_2_DESCRIPTION \
356 "+ xrays source fix allong the :math:`\\vec{x}` direction (1, 0, 0)\n" \
357 "+ 2 axes for the sample\n" \
359 " + **pitch** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
360 " + **mu** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
362 "+ 3 axis for the detector\n" \
364 " + **pitch** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
365 " + **gamma** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
366 " + **delta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n"
368 static const char* hkl_geometry_soleil_sixs_med_1_2_axes
[] = {"pitch", "mu", "gamma", "delta"};
370 static HklGeometry
*hkl_geometry_new_soleil_sixs_med_1_2(const HklFactory
*factory
)
372 HklGeometry
*self
= hkl_geometry_new(factory
);
375 h
= hkl_geometry_add_holder(self
);
376 hkl_holder_add_rotation_axis(h
, "pitch", 0, -1, 0);
377 hkl_holder_add_rotation_axis(h
, "mu", 0, 0, 1);
379 h
= hkl_geometry_add_holder(self
);
380 hkl_holder_add_rotation_axis(h
, "pitch", 0, -1, 0);
381 hkl_holder_add_rotation_axis(h
, "gamma", 0, 0, 1);
382 hkl_holder_add_rotation_axis(h
, "delta", 0, -1, 0);
387 static HklEngineList
*hkl_engine_list_new_soleil_sixs_med_1_2(const HklFactory
*factory
)
389 HklEngineList
*self
= hkl_engine_list_new();
391 hkl_engine_soleil_sixs_med_1_2_hkl_new(self
);
392 hkl_engine_q2_new(self
);
393 hkl_engine_qper_qpar_new(self
);
394 hkl_engine_tth2_new(self
);
399 REGISTER_DIFFRACTOMETER(soleil_sixs_med_1_2
, "SOLEIL SIXS MED1+2", HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_1_2_DESCRIPTION
);
402 /***********************/
403 /* SOLEIL SIXS MED 2+3 */
404 /***********************/
406 #define HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_2_3_DESCRIPTION \
407 "+ xrays source fix allong the :math:`\\vec{x}` direction (1, 0, 0)\n" \
408 "+ 3 axes for the sample\n" \
410 " + **beta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
411 " + **mu** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
412 " + **omega** : rotating around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
414 "+ 4 axis for the detector\n" \
416 " + **beta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
417 " + **gamma** : rotation around the :math:`\\vec{z}` direction (0, 0, 1)\n" \
418 " + **delta** : rotation around the :math:`-\\vec{y}` direction (0, -1, 0)\n" \
419 " + **eta_a** : rotation around the :math:`-\\vec{x}` direction (-1, 0, 0)\n"
421 static const char* hkl_geometry_soleil_sixs_med_2_3_axes
[] = {"beta", "mu", "omega", "gamma", "delta", "eta_a"};
423 static HklGeometry
*hkl_geometry_new_soleil_sixs_med_2_3(const HklFactory
*factory
)
425 HklGeometry
*self
= hkl_geometry_new(factory
);
428 h
= hkl_geometry_add_holder(self
);
429 hkl_holder_add_rotation_axis(h
, "beta", 0, -1, 0);
430 hkl_holder_add_rotation_axis(h
, "mu", 0, 0, 1);
431 hkl_holder_add_rotation_axis(h
, "omega", 0, -1, 0);
433 h
= hkl_geometry_add_holder(self
);
434 hkl_holder_add_rotation_axis(h
, "beta", 0, -1, 0);
435 hkl_holder_add_rotation_axis(h
, "gamma", 0, 0, 1);
436 hkl_holder_add_rotation_axis(h
, "delta", 0, -1, 0);
437 hkl_holder_add_rotation_axis(h
, "eta_a", -1, 0, 0);
442 static HklEngineList
*hkl_engine_list_new_soleil_sixs_med_2_3(const HklFactory
*factory
)
444 HklEngineList
*self
= hkl_engine_list_new();
446 self
->geometries
->multiply
= hkl_geometry_list_multiply_soleil_sixs_med_2_3
;
447 hkl_engine_soleil_sixs_med_2_3_hkl_new(self
);
448 hkl_engine_q2_new(self
);
449 hkl_engine_qper_qpar_new(self
);
450 hkl_engine_tth2_new(self
);
455 REGISTER_DIFFRACTOMETER(soleil_sixs_med_2_3
, "SOLEIL SIXS MED2+3", HKL_GEOMETRY_TYPE_SOLEIL_SIXS_MED_2_3_DESCRIPTION
);