IN14(Ki,Q,E)
/* A template IN14 simulation for McStas. E. Farhi <farhi@ill.fr> Dec 99
* input parameters :
* Ki : incident n-wavevector [Angs-1]
* Q : transfert wavevector [Angs-1] (modulus)
* E : transfert energy [meV]
* building : mcstas in14_tut.instr ; gcc -O -o in14_tut in14_tut.c -lm
* execute : in14_tut -n 1e6 Ki=2.662 Q=2 E=2
*
* Emmanuel FARHI, <farhi@ill.fr>.CS-group (c) ILL. Oct 2000
*/
%{
/* all variables defined here can be used, and modified, in the whole instrument */
/* variables used in instrument definition and for user computations */
double L1 = 16.68; /* source-mono distance (m) */
double L2 = 2.12; /* mono-sample */
double L3 = 1.35; /* sample-ana */
double L4 = 0.70; /* ana-det */
double DM = 3.355; /* mono d-spacing (Angs) */
double DA = 3.355; /* ana d-spacing */
double A1,A2,A3,A4,A5,A6; /* geometrical angles (deg) */
double mono_q; /* Q vector for bragg scattering with monochromator */
double ana_q; /* Q vector for bragg scattering with analyzer */
double Ei = 0; /* incident n-energy*/
double Ef = 0; /* out-going n-energy*/
double Kf = 0; /*out-going n-wavevector */
%}
/* end of DECLARE */
%{
/* variables defined here are local (only in INIT section) */
/* can be computed / displayed (see below) */
double Vi = 0; /* incident n-velocity (m/s) */
double Vf = 0; /* out-going n-velocity */
int SM,SS,SA; /* rotation sign*/
double arg = 0; /* tempo var */
int ORDER = 1;
/* dealing with ORDERth Bragg one-neutron diffraction on mono/ana */
/* SM : scattering at mono to the right (-1)/left(+1) */
/* SS : scattering at sample to the right (-1)/left(+1) */
/* SA : scattering at analyser to the right (-1)/left(+1) */
SM = 1; SS = -1; SA = 1;
/* Monocrhomator characteristics (A1,A2 from Ki,DM) */
if (DM <= 0.1)
{
printf("TAS: Mono DM = %f Angs, too small\n", DM);
exit(-1);
}
mono_q = 2*PI*ORDER/DM; /* Q mono in Angs-1 */
if ((Ki == 0) || (fabs(mono_q/2/Ki) > 1))
{
printf("TAS: Can not reach this incident wavevector\n");
printf("with such a monochromator (A2)\n");
printf("Mono DM = %f Angs, requested Ki = %f Angs-1\n", DM, Ki);
exit(-1);
}
Vi = K2V*fabs(Ki);
Ei = VS2E*Vi*Vi;
A2 = asin(mono_q/2/Ki)*RAD2DEG*2;
A2 *= SM; /*A1 : mono theta (crystal) */
A1 = A2/2; /*A2 : mono 2 theta (arm to sample) */
Ef = Ei - E;
if (Ef < 0)
{
printf("TAS: Out-going neutrons have negative energy !!\n");
printf("Incident Ei = %f meV, Exchange E = %f meV\n",Ei,E);
exit(-1);
}
/* Analyzer characteristics (A5,A6 from Kf, DA) */
Vf = SE2V*sqrt(Ef); /*final energy from energy transfert */
Kf = V2K*Vf;
if (DA <= 0.1)
{
printf("TAS: Ana DA = %f Angs, too small\n", DM);
exit(-1);
}
ana_q = 2*PI*ORDER/DA; /* Q ana in Angs-1 */
if ((Kf == 0) || (fabs(mono_q/2/Kf) > 1))
{
printf("TAS: Can not reach this out-going wavevector\n");
printf("with such an analyzer (A6)\n");
printf("Mono DA = %f Angs, requested Kf = %f Angs-1\n",DA, Kf);
exit(-1);
}
A6 = asin(ana_q/2/Kf)*RAD2DEG*2;
A6 *= SA; /* A5 : analyser theta (crystal) */
A5 = A6/2; /* A6 : analyser 2 theta (arm to detector) */
/* Sample characteristics */
arg = (Ki*Ki + Kf*Kf - Q*Q)/2/fabs(Ki*Kf);
if (fabs(arg) > 1)
{
printf("TAS: Can not close wavevector triangle (A4)\n");
printf("Incident Ki = %f Angs-1, Exchange Q=%f Angs-1, ",Ki,Q);
printf("Final Kf = %f Angs-1\n",Kf);
exit(-1);
}
A4 = acos(arg)*RAD2DEG;
A4 *= SS; /*A4 : sample 2 theta (arm to analyser) */
/* computing of A3 requires sample reciprocal axes characteristics */
A3 = A4/2; /* A3 : sample theta */
printf("IN14@ILL 3-axis spectrometer\n");
printf("Angles (deg) A2 = %f ; A4 = %f ; A6 = %f\n",A2,A4,A6);
%}
/* end of INITIALIZE */
/* Source description */
COMPONENT source = Source_flux( /* a flat constant source */
radius = 0.20,
dist = 2.16,
xw = 0.06, yh = 0.12,
E0 = Ei,
dE = 0.5,
flux=1e13)
AT (0,0,0) ABSOLUTE
COMPONENT cold_guide = Guide(
w1 = 0.06, h1 = 0.12,
w2 = 0.06, h2 = 0.12,
l = 14.5, /*guide between source and mono */
R0 = 1,
Qc = 0.021,
alpha = 2.33,
m = 1.2, /*1.2 Ni 58, 3 Super mirror */
W = 2e-3)
AT (0,0,2.15) RELATIVE source
/* Monochromator description */
COMPONENT mono_craddle = Arm() /*this is the monochromator craddle */
AT (0, 0, L1) RELATIVE source ROTATED (0, A1, 0) RELATIVE source
COMPONENT mono = Monochromator(
zmin = -0.075,
zmax = 0.075,
ymin = -0.06,
ymax = 0.06,
mosaich=30,
mosaicv=30,
r0=.9, Q=mono_q)
AT (0,0,0) RELATIVE mono_craddle
COMPONENT out_mono = Arm() /*this is the mono-sample axis */
AT (0,0,0) RELATIVE mono_craddle ROTATED (0, A2, 0) RELATIVE source
COMPONENT mono_soller = Soller(
xmin = -0.02, xmax = 0.02, ymin = -0.04, ymax = 0.04,
len = 0.4, divergence = 40)
AT (0, 0, .9) RELATIVE out_mono
/* Sample description */
COMPONENT focus_sample = Arm()
AT (0, 0,L2) RELATIVE out_mono ROTATED (0, A3, 0) RELATIVE out_mono
/* to be added in the section 'A-9: Looking at results : mcplot' of tutorial */
COMPONENT PSD_sample = PSD_monitor(
xmin = -.1, xmax = .1,
ymin = -.1, ymax = .1,
nx = 10, ny = 10,
filename="sample.psd")
AT (0, 0, 0) RELATIVE focus_sample
COMPONENT out_sample = Arm() /* this is the sample-ana axis */
AT (0,0,0) RELATIVE focus_sample ROTATED (0, A4, 0) RELATIVE out_mono
COMPONENT sample = V_sample(
radius_i = 0.008, radius_o = 0.012, h = 0.015,
focus_r = L3, pack = 1,
target_x = 01, target_y = 0, target_z = 1000)
AT (0,0,0) RELATIVE out_sample
/* Analyzer description */
COMPONENT focus_ana = Arm()
AT (0,0, L3) RELATIVE out_sample ROTATED (0, A5, 0) RELATIVE out_sample
COMPONENT ana = Monochromator(
zmin = -0.05,
zmax = 0.05,
ymin = -0.07,
ymax = 0.07,
mosaich=30,
mosaicv=30,
r0=.9, Q=ana_q)
AT (0,0,0)
RELATIVE focus_ana
COMPONENT out_ana = Arm() /* this is the ana-det axis */
AT (0,0,0) RELATIVE focus_ana ROTATED (0, A6, 0) RELATIVE out_sample
/* Detector description */
COMPONENT focus_det = Arm()
AT (0, 0, L4) RELATIVE out_ana
COMPONENT Detector = Monitor(
xmin = -0.02, xmax = 0.02,
ymin = -0.05, ymax = 0.05)
AT(0,0,0) RELATIVE focus_det