/*******************************************************************************
*
* McStas, neutron ray-tracing package
*         Copyright (C) 1997-2011, All rights reserved
*         Risoe National Laboratory, Roskilde, Denmark
*         Institut Laue Langevin, Grenoble, France
*
* Instrument: Crisp3 test instrument for the Multilayer_NR_Sample component.
*
* %Identification
* Written by: <a href="robert.dalgliesh@stfc.ac.uk">Robert Dalgliesh</a>
* Date: 2010
* Origin: <a href="http://www.isis.stfc.ac.uk/">ISIS (UK)</a>
* Release: McStas 1.12cb
* Version: 1.0
* %INSTRUMENT_SITE: ISIS
*
* Model of the ISIS CRISP reflectometer, including the Multilayer_Sample reflectivity sample.
*
* %Description
* This model of the ISIS CRISP reflectometer demonstrates the use of the Multilayer_Sample component.
*
* The algorithm of the component requires complex numbers, facilitated by the 
* <a href="http://www.gnu.org/software/gsl/">GNU Scientific Library (GSL)</a>. To link
* your the instrument with an installed GSL, you should use MCSTAS_CFLAGS like
*
* MCSTAS_CFLAGS = -g -O2 -lm -lgsl -lgslcblas
*
* Example: mcrun ISIS_CRISP.instr glen=1.4 flen=0.4 w1=0.05 vm=3.0 FRAC=0.1
*
* %Parameters
* INPUT PARAMETERS:
* glen: Length of the elliptical guides [m]
* flen: Focal-length of elliptical guides [m]
* w1: Entry-width of elliptical guides [m]
* vm: m-value of material for left and right vertical guide mirrors [m]
* FRAC: Fraction of statistics used to model incoherent scattering from sample [1]
*
* %Link
* <a href="http://www.isis.stfc.ac.uk/instruments/crisp/">Website of the CRISP instrument</a>
*
* %End
*******************************************************************************/
DEFINE INSTRUMENT ISIS_CRISP(glen=1.4,flen=0.4,w1=0.05,vm=3.0, FRAC=0)


DECLARE
%{
  double glen,flen,spos,tend,t15,fp1; 
  double linw,loutw,l,u1,u2,div1,b_ell_q,w1,w12,a_ell_q,lbw,w2;
%}

INITIALIZE
%{
spos=10.2055;
l=glen;
loutw=flen;
linw=10.0;
w12=w1/2.0;

// calculate the width of the guide exit
    lbw = l + linw + loutw;
    u1 = sqrt((linw*linw)+(w12*w12));
    u2 = sqrt((w12*w12) + ((l+loutw)*(l+loutw)));
    a_ell_q = ((u1 + u2)/2)*((u1 + u2)/2);
    b_ell_q = a_ell_q - ((lbw/2)*(lbw/2));
    /* calculate width of guide exit (w2) */
    div1 = ((lbw/2-loutw)*(lbw/2-loutw))/a_ell_q;
    w2 = sqrt(b_ell_q*(1-div1));
    w2 = w2*2;

tend=w2;
t15=tan(-1.5*DEG2RAD);
%}
TRACE

COMPONENT a1 = Arm()
  AT (0,0,0) ABSOLUTE

COMPONENT isis_source = ISIS_moderator(
    Face = "crisp", E0 = -6.5, E1 = -0.55, dist = 7.2695,
    xw = 0.045, yh = 0.0045, modXsize = -1, modYsize = -1,
    CAngle = 0.0, SAC = 1)
  AT (0.0, 0.0, 0.00001) RELATIVE a1
  ROTATED (1.5, 0.0, 0.0) ABSOLUTE

COMPONENT defslit1 = Slit(
    xmin = -0.0345, xmax = 0.0345, ymin = -0.01731,
    ymax = 0.01731)
  AT (0.0, -0.09829, 3.7537) ABSOLUTE 

COMPONENT defslit2 = Slit(
    xmin = -0.0289, xmax = 0.0289, ymin = -0.0123,
    ymax = 0.0123)
  AT (0.0, -0.15976, 6.101) ABSOLUTE 

COMPONENT coarseslit1 = Slit(
    xmin = -0.03, xmax = 0.03, ymin = -0.1,
    ymax = 0.1)
  AT (0.0, -0.1757, 6.7096) ABSOLUTE 

COMPONENT coarseslit2 = Slit(
    xmin = -0.1, xmax = 0.1, ymin = -0.01,
    ymax = 0.01)
  AT (0.0, -0.17832, 6.8096) ABSOLUTE

COMPONENT lmon1 = L_monitor(
    nchan = 500, filename = "lmon1.dat", xmin = -0.06,
    xmax = 0.06, ymin = -0.04, ymax = 0.04, Lmin = 0.0,
    Lmax = 10.0)
  AT (0.0, -0.1833, 7.0) ABSOLUTE

COMPONENT PSDmon1 = PSD_monitor(
    nx = 100, ny = 100, filename = "PSD1.dat",
    xmin = -0.05, xmax = 0.05, ymin = -0.01, ymax = 0.01)
  AT (0.0, -0.1833, 7.01) ABSOLUTE 

COMPONENT slit1 = Slit(
    width=0.04, height=0.004)
  AT (0.0, -0.1904, 7.2695) ABSOLUTE

COMPONENT PSDmons1 = PSD_monitor(
    nx = 100, ny = 100, filename = "PSDs1.dat",
    xmin = -0.05, xmax = 0.05, ymin = -0.01, ymax = 0.01)
  AT (0.0, 7.27*t15, 7.27) ABSOLUTE 

COMPONENT eguide1 = Guide_tapering(
    segno = 10, w1 = w1, h1 = 0.1, l = glen, option="elliptical",
    linw=10.0, loutw = flen, mx = vm, my = 0)
  AT (0.0, (spos-flen-glen)*t15, spos-flen-glen) ABSOLUTE
  ROTATED (1.5, 0.0, 0.0) ABSOLUTE

COMPONENT slit2 = Slit(
    width=tend, height=0.0025)
  AT (0.0, -0.2581, 9.8555) ABSOLUTE

COMPONENT lmon2 = L_monitor(
    nchan = 500, filename = "lmon2.dat", xmin = -0.06,
    xmax = 0.06, ymin = -0.04, ymax = 0.04, Lmin = 0.0,
    Lmax = 10.0)
  AT (0.0, -0.2608, 9.96) ABSOLUTE

COMPONENT samp1 = Multilayer_Sample(
    xwidth = 0.05, zlength = 0.15, mu_inc=0.138, ythick=0.01,
    nlayer = 0, sldPar = {0.0,6.35e-6}, dPar = {0.0},
    sigmaPar = {5.0}, target_index=1, focus_xw=2*tend, focus_yh=0.01,
    frac_inc=FRAC)
  AT (0.0, spos*t15, spos) ABSOLUTE

COMPONENT slit3 = Slit(
    width= tend, height=0.003)
  AT (0.0, (spos+flen-0.01)*t15-2.0*(flen-0.01)*t15, spos+flen-0.01) ABSOLUTE

COMPONENT eguide2 = Guide_tapering(
    segno = 10, w1 = tend, h1 = 0.1, l = glen, option="elliptical",
    linw=flen, loutw = 10.0, mx = vm, my = 0)
  AT (0.0, (spos+flen)*t15-2.0*flen*t15, spos+flen ) ABSOLUTE
  ROTATED(-1.5,0.0,0.0) ABSOLUTE

COMPONENT lmon3 = L_monitor(
    nchan = 500, filename = "lmon3.dat", xmin = -0.05,
    xmax = 0.05, ymin = -0.05, ymax = 0.05, Lmin = 0.0,
    Lmax = 10.0)
  AT (0.0, 12.11*t15-2.0*(12.11-spos)*t15, 12.11) ABSOLUTE

COMPONENT PSDdet1 = PSD_monitor(
    nx = 100, ny = 100, filename = "PSD3.dat",
    xmin = -0.05, xmax = 0.05, ymin = -0.05, ymax = 0.05)
  AT (0.0, 12.111*t15-2.0*(12.111-spos)*t15, 12.111) ABSOLUTE 

FINALLY
%{
%}
END