/*******************************************************************************
*
* McStas, neutron ray-tracing package
*         Copyright 1997-2011, All rights reserved
*         Risoe National Laboratory, Roskilde, Denmark
*         Institut Laue Langevin, Grenoble, France
*
* Instrument: RITA-II @ PSI with hkl calculator
*
* %Identification
* Written by: <a href="mailto:udby@fys.ku.dk">Linda Udby</a>  and <a href="mailto:pkwi@fysik.dtu.dk">Peter Willendrup</a>
* Date: 2012
* Origin: <a href="http://www.risoe.dk">Ris&oslash; (Denmark)</a>
* %INSTRUMENT_SITE: PSI
*
* RITA type triple-axis spectrometer (TAS)
*
* %Description
* This instrument is model of RITA-II at the RNR13 guide at SINQ, PSI as of 2009. The energy and q-resolution as been verified against measured data in 2- and 3-axis modes. See Udby et al., Nuclear Instruments and Methods 634 (2011) s138-.<br>
* The model is based on the following components in downstream order
* <ul>
* <li> The SINQ source in 2009 - tested up to Ei=15 meV versus measurements. The spectrum of genereated rays is controlled by the BPL and BPH parameters.</li>
* <li> The RNR13 guide with average m-value based on color-codes on guide pieces and gaps as measured/from technical drawings.</li>
* <li> The Druckal monochromator: 5 slabs of PG(002) with vertically focusing option. The horizontal rotation angle (A1) and  scattering angle (A2) can either be set by user input or calculated automatically from EI/KI. </li>
* <li> Optional filter after the monochromator, before the sample. </li>
* <li> Variable linear collimator after the monochromator, before the sample. </li>
* <li> Optional perspex attenuation. </li>
* <li> Inbound Diaphagm slit. </li>
* <li> Sample: Incoherent scatterer, powder or single crystal from reflection list. </li>
* <li> Outbound Diaphragm slit. </li>
* <li> Optional filter after the sample, with radial collimator.</li>
* <li> 9 - blade analyser: Base-rotation (AA5) and individual blade angles (C1-C9) can either be set by user or calculated automatically from EF/KF. </li>
* <li> Coarse collimator: Removes cross-talk by allowing only passage of scattered neutrons from specified blade in each channel.</li>
* <li> Detector: PSD detector with 100% efficiency and specified point-spread function from measurements.</li>
* <li> Windows: Monitors the neutrons in each electronically specified window [XWinMin,XWinMax;YWinMin,YWinMax].</li>
* </ul>
* EXAMPLES
* <ul>
* <li>Generating a virtual source </li>
* mcrun RITA-II.instr -n3e6  BPL=0.97 BPH=1.03 EI=5 EF=5 COLL_MS=40 SAMPLE=1 OUTFILTER=0 REP=10 VIRTUALOUT=1 -d  RITA-II_Vsource40 <br>
* <li>Vanadium sample, no filters, seeing 1st-3rd order neutrons </li>
* mcrun RITA-II.instr -n3e6 BPL=0.30 BPH=1.03 EI=5 EF=5 SAMPLE=1 INFILTER=0 OUTFILTER=0  REP=10 -d RITA-II_Vanadium <br>
* <li>Powder sample (default sapphire), simulation 1) through entire instrument or  2)from virtual source 3) in 2-axis mode, with particular slit settings scaling the source yield to 2008 values </li>
* 1) mcrun RITA-II.instr -n3e6 EF=5 EN=0 SAMPLE=2 COLL_MS=40 SAMPLEFILE="Al2O3_sapphire.lau" BARNS=0 SOURCEFILE="RITA-II_Vsource40/Vin_default.dat" REP=10 -d RITA-II_40_powder<br>
* 2) mcrun RITA-II.instr EI=5 EF=5 SAMPLE=2 SAMPLEFILE="Al2O3_sapphire.lau" BARNS=0 SOURCEFILE="RITA-II_Vsource40/Vin_default.dat" VIRTUALIN=1 REP=3 -d RITA-II_VINsource40_powder <br>
* 3) mcrun RITA-II.instr -n 1e7 ITAR=0.71 COLL_MS=19.6 BPL=0.97 BPH=1.03 EI=5 EF=5 QH=0 QK=0 QL=0 QM=1 A3=0.001 A4=-71 AA5=-90 A6=1e-3 MST=25 MSB=25 MSL=10 MSR=10 SAMPLE=2 SAMPLEFILE="Al2O3_sapphire.lau" SST=25 SSB=25 SSL=10 SSR=10 OUTFILTER=0 OUTFILTERFILE=Be.trm COARSE=1 REP=5 LC=6 RC=4 REP=1 BARNS=0 ANAFORCE=1 -d RITA-II_Al2O3_2axis<br>
* <li>Single crystal sample. Define first lattice vector (a) along -x (to the right) , second lattice vector (b) along z (forward), third lattice vector (c) along y (up). </li>
*  mcrun RITA-II.instr -n1e6 BPL=0.97 BPH=1.03 EI=5 EF=5 SAMPLE=3 REP=10 QH=2 QM=0 AS=4.95 BS=4.95 CS=4.95 AAX=-4.95 BBZ=4.95 CCY=4.95 AH=0 AK=1 AL=0  BH=1 BK=0 BL=0 SAMPLEFILE="Pb.laz" -d RITA-II_SXPb200 <br>
* <li> Phonon sample, inelastic scattering, no Bragg peaks in sample </li>
* mcrun RITA-II.instr -n1e6 L0=3.418 BPL=0.97 BPH=1.03 EI=7 EF=5 SAMPLE=4 REP=10 QH=2 QK=0.16 QM=0 AS=4.95 BS=4.95 CS=4.95 AAX=-4.95 BBZ=4.95 CCY=4.95 AH=0 AK=1 AL=0  BH=1 BK=0 BL=0 SAMPLEFILE="Pb.laz" -d RITA-II_SXPb200_2meV <br>
* </ul>
* %Parameters
* ITAR: [mA]               Relative neutron yield from the spallation target. Value relative to 2009
* EI: [meV]                Incoming neutron energy. Also used to set range of energy monitors
* EF: [meV]                Outgoing neutron energy. Also used to set range of energy monitors
* QH: [rlu]                Measurement QH position in crystal
* QK: [rlu]                Measurement QK position in crystal
* QL: [rlu]                Measurement QL position in crystal
* EN: [meV]                Energy transferred in crystal
* QM: [Angs-1]             Wavevector transferred in sample, use QM=0 if (QH,QK,QL) is specified
* A1: [deg]                Monohromator rotation angle
* A2: [deg]                Monohromator take-off angle
* A3: [deg]                Sample rotation angle
* A4: [deg]                Sample take-off angle
* AA5: [deg]               Analyzer rack rotation angle
* A6: [deg]                Analyzer take-off angle
* C1: [deg]                Analyzer blade 1 position
* C2: [deg]                Analyzer blade 2 position
* C3: [deg]                Analyzer blade 3 position
* C4: [deg]                Analyzer blade 4 position
* C5: [deg]                Analyzer blade 5 position
* C6: [deg]                Analyzer blade 6 position
* C7: [deg]                Analyzer blade 7 position
* C8: [deg]                Analyzer blade 8 position
* C9: [deg]                Analyzer blade 9 position
* COLL_MS: [deg]           Primary collimator max divergence
* MONO_N: [1]              Order of reflection used on mono
* MONOFORCE: [1]           If set, monochromator is flat
* L0: [Angs]               Centre of generated wavelength distribution from source
* BPL:  Band Pass Low factor, multiplied on source wavelength L0 to allow neutrons with wavelengths with lambda {BPL*L0,BPH*L0} to be traced from the source.
* BPH:  Band Pass High factor, multiplied on source wavelength L0 to allow neutrons with wavelengths with lambda {BPL*L0,BPH*L0} to be traced from the source.
* MST: [m]                 Monochromator TOP slit
* MSB: [m]                 Monochromator BOTTOM slit
* MSR: [m]                 Monochromator RIGHT slit
* MSL: [m]                 Monochromator LEFT slit
* SST: [m]                 Sample TOP slit
* SSB: [m]                 Sample BOTTOM slit
* SSR: [m]                 Sample RIGHT slit
* SSL: [m]                 Sample LEFT slit
* SM: [1]
* SS: [1]                  Scattering configuration signs. 'W' is SM=1,SS=-1,SA=1
* SA: [1]
* OUTFILTER: [1]           Flag to indicate if the outbound filter is in or out
* OUTFILTERFILE: [string]  An input file of wavevector vs transmission to use with outgoing filter
* INFILTER: [1]            Flag to indicate if mono-block filter is in or out
* INFILTERFILE: [string]   An input file of wavevector vs transmission to use with incoming filter
* COARSE: [1]              Flag to indicate if Detector collimator is in or out
* REP: [1]                 Repetition factor of virtual_input or Monochromator/Sample/Analyzer SPLITs
* ANAFORCE: [1]            If set the analyzer is forced flat
* MONO_MOS_H: [min]        Monochromator, horizontal mosaicity
* MONO_MOS_V: [min]        Monochromator, vertical mosaicity
* LC: [deg]                Detector-collimator angle of the leftmost blade
* RC: [deg]                Detector-collimator angle of the rightmost blade
* PERSPEX: [1]             Flag to indicate if perspex attenuator is in or out
* PTHICK: [m]              Thickness of perspex attenuator
* VIRTUALOUT: [1]          If set, flag indicates that a Virtual source is created after the monochromator collimator
* VIRTUALIN: [1]           If set, flag to indicates that Virtual source is used after the monochromator collimator
* SOURCEFILE: [string]     Name of file to be used/genereated as virtual input/output
* verbose: [1]             print TAS configuration. 0 to be quiet
* SAMPLE: [1]              1 is incoherent scatterer, 2 is powder, 3 is single crystal.
* SAMPLEFILE: [string]     Name of samplefile (with reflectionlist etc)
* SAMPLESIZE: [m]          Length, height and width of single crystal sample, or radius and height of phonon sample
* MOS: []                 Isotropic 'mosaicity' of single crystal
* DD_D: []                spead of lattice parameter
* AS: [Angs]               Sample lattice parameter A
* BS: [Angs]               Sample lattice parameter B
* CS: [Angs]               Sample lattice parameter C
* AA: [deg]                Angle between lattice vectors B,C
* BB: [deg]                Angle between lattice vectors C,A
* CC: [deg]                Angle between lattice vectors A,B
* AH: [rlu]                First reciprocal lattice vector in scattering plane, X
* AK: [rlu]                First reciprocal lattice vector in scattering plane, Y
* AL: [rlu]                First reciprocal lattice vector in scattering plane, Z
* BH: [rlu]                Second reciprocal lattice vector in scattering plane, X
* BK: [rlu]                Second reciprocal lattice vector in scattering plane, Y
* BL: [rlu]                Second reciprocal lattice vector in scattering plane, Z
* AAX: []
* AAY: []                 Orientation vector of unit cell, single_crystal
* AAZ: []
* BBX: []
* BBY: []                 Orientation vector of unit cell, single_crystal
* BBZ: []
* CCX: []
* CCY: []                 Orientation vector of unit cell, single_crystal
* CCZ: []
* TILT: [deg]              Tilt angle out of plane of the sample scattering vector out of plane (rotation around z of the A3 arm)
* BARNS: [1]               If set the flag indicates that reflection list structure factors are in units of barns, otherwise fm^2
* BPL: []
* BPH: []
* MOS: []
* DD_D: []
* AAX: []
* AAY: []
* AAZ: []
* BBX: []
* BBY: []
* BBZ: []
* CCX: []
* CCY: []
* CCZ: []
*
* %Link
* Rescal for Matlab at http://www.ill.fr/tas/matlab
* %Link
* Restrax at http://omega.ujf.cas.cz/restrax/
* %End
*******************************************************************************/
DEFINE INSTRUMENT RITA_II( ITAR=1.0,L0=4.045,BPL=0.97,BPH=1.03, MONO_N=1,MONOFORCE=0,MONO_MOS_H=37,
  MONO_MOS_V=37, EI=0,EF=0,EN=0, SM=1,SS=-1,SA=1, QH=0, QK=0, QL=0, QM=1.8051,
  AS=4.95, BS=4.95, CS=4.95, AA=90, BB=90, CC=90, AH=0, AK=1, AL=0, BH=1, BK=0, BL=0,
  INFILTER=0, string INFILTERFILE="Be.trm", COLL_MS=80, MST=40, MSB=40, MSL=40, MSR=40,
  PTHICK=1e-3,PERSPEX=0, SAMPLE=1, string SAMPLEFILE="default",MOS=100,DD_D=1e-3,
  SAMPLESIZE=0.01,BARNS=1,
  AAX=-4.95, AAY=0, AAZ=0, BBX=0, BBY=0, BBZ=4.95, CCX=0, CCY=4.95, CCZ=0,
  A1=0,A2=0,A3=0,A4=0,A6=0,TILT=0,  SST=100, SSB=100, SSL=100, SSR=100,
  OUTFILTER=1, string OUTFILTERFILE="Be.trm", ANAFORCE=0,
  AA5=0,C1=0,C2=0,C3=0,C4=0,C5=0,C6=0,C7=0,C8=0,C9=0, COARSE=1,LC=4, RC=3, REP=1,
  VIRTUALOUT=0, VIRTUALIN=0, string SOURCEFILE="Vin_default.dat",verbose=1 )

DECLARE
%{
/* The following is from RITA2front, Kim Lefmann / Linda Udby */


/* Static values... */
double l0,lmin,lmax;
double emini,emaxi;
double eminf,emaxf;
double nefchan;
double neichan;
double lI,KI,KF;

/* Internal variable for SPLIT repetitions */
int SPLITREP;
int SPLITMREP;
int SPLITAREP;
#pragma acc declare create (SPLITREP,SPLITMREP,SPLITAREP)
/* Guide element parameters*/
double angleGuideCurved;
double R = 0.88;
double R0 = 0.995;
double Qc = 0.0217;
double W = 0;
double  M= 2.15;
double ALPHA;

/* Monochromator material parameters */
double mono_q = 1.87325;
double mono_r0 = 0.8;
double DM;      /*d-spacing monochromator*/
double mono_mosaic_h;
double mono_mosaic_v;
/* Monochromator curvature parmeters */
//double u,v;
double dms;  /* Target vector for focusing */
double tx,tz;    /* Target vector for focusing */
double sintm,sinta;
double rv; /* Mono focusing parameter */
/* Monochromator geometrical parameters */
/* Size of monochromator blades are hard-coded in the component*/
double mono_d = 0.026;  /* Distance between mono blades. From drawings */
double dmc; /* Distance monochromator to front of collimator. Was 0.32 from Stine */
double rmh = 0.58; /* Radius of monohousing. Measured 2008/11/05 */
double lc = 0.198; /*Length of monocollimator. Measured 2008/11/05*/

/* RITA Monitor efficiency 1.7e-5 */
double EFF = 1.7e-5;


/* Sample parameters */
double d_sample_slit = 0.35; /* Measured 2008/11/05, was 25 cm from Stine */
double d_sample_filter = 0.51; /*To centre of filter. Measured 2008/11/05, was to filter front 25 cm from Stine */
double dsa = 1.195;  /* distance sample-analyzer (m). Was 1.256 from Stine */
double dsb ; /*distance sample to blade, depends on analyser rack postion*/
double dmv = 1.00; /* distance monochromator to virtual out */
//double dvs = 0.67; /* distance virtual in to sample. dmv+dvs=1.67=rmh+1.09. Ma15 setting */
double dvs = 0.54; /* distance virtual in to sample. dmv+dvs=1.54 */
double is_incoh;/* Random number pr. neutron event for incoherent V scattering */
/* Filenames for the sample comps: */
char *PowderFile;
char *SingleXFile;

/* Analyser material parameters*/
double ana_mosaic_h;
double ana_mosaic_v;
double ana_q = 1.87325;
double ana_r0 = 0.8;
double DA;    /* d-spacing analyser*/
/*Analyser geometrical parameters */
double ana_d = 0.025;  /* Width of analyser blades. From drawings */
double ana_h = 0.15;   /* Height of analyser blades. From drawings */
double dad = 0.338;  /* distance analyzer to detector front (m) */
double dadw = 0.3595; //= 0.3595;  /* distance analyzer to detector wire (m) */
double wan = 0.024;  /* width of analyzer blades , by ruler (m) */
double RR; //  = 0.3008; /* ratio of dadw and dsa*/
//double DA4; /*angle between blades*/
double A5; /* scattering angle for flat analyser*/
/* Declarations for 'Coarse Collimator' at the PSD detector surface */
int EntrySlit;
int ExitSlit;
double BladeThickness = 0.007;// detector coll after 2006, from drawings
double WindowSize = 0.025;
double BladeLength = 0.179;// detector coll after 2006, from drawings
double BladeHeight = 0.272;// detector coll after 2006, from drawings
//double BladeThickness = 0.002;// detector coll before 2006
//double BladeLength = 0.180;// detector coll before 2006
//double BladeHeight = 0.250;// detector coll before 2006
double FirstWindowSizeL;
double FirstWindowSizeR;
double deltaL;
int coarse;

/* Detector parameters */
//double det_width = 0.2735;/* was 0.3 from Stine*/
double det_width = 0.275;
double det_height =0.5;
double PSF = 0.0074/2.35;// FWHM=0.0074m, measured by C. Bahl NIMB 246, 452.
/* Electronic Window positions in the PSD */
int XwinMin[9];
int YwinMin[9];
int XwinMax[9];
int YwinMax[9];
#pragma acc declare create (XwinMin,YwinMin,XwinMax,YwinMax)




/* BEGIN declaring stuff for the HKL calculator adapted from templateTAS */
  struct sample_struct {
    double as, bs, cs; // Lattice parameters
    double aa, bb, cc; // Lattice angles
    double ax, ay, az; // First scattering plane vector
    double bx, by, bz; // Second scattering plane vector
  } sample;

  struct machine_hkl_struct {
    double dm, da;     // Mono and ana d-spacings
    double sm, ss, sa; // Mono, sample, ana angle signs
    double ki, kf, ei, ef; // Initial and Final wavevectors and energies
    double qh, qk, ql, en; // Momentum transfer and energy transfer in sample
  } machine_hkl;

  struct machine_real_struct {
    double a1,a2,a3,a4,a5,a6;
    double da4,aa5;
    double c1,c2,c3,c4,c5,c6,c7,c8,c9;
    double rmh, rmv, rah, rav;
    double qm, qs, qt[3];
    char   message[256];
  } machine_real;

  struct machine_real_struct qhkl2angles(
    struct sample_struct      sample,
    struct machine_hkl_struct machine_hkl,
    struct machine_real_struct machine_real) {

      /* code from TASMAD/t_rlp.F:SETRLP */
      double qhkl[3];
      double alpha[3];
      double a[3];
      double aspv[3][2];
      double cosa[3], sina[3];
      double cosb[3], sinb[3];
      double b[3], c[3], s[4][4];
      double vv[3][3], bb[3][3];
      double arg, cc;
      int i,j,k,l,m,n;
      char liquid_case=1;
      /* transfered parameters to local arrays */
      qhkl[0]   = machine_hkl.qh; /* HKL target */
      qhkl[1]   = machine_hkl.qk;
      qhkl[2]   = machine_hkl.ql;
      alpha[0]  = sample.aa; /* cell angles */
      alpha[1]  = sample.bb;
      alpha[2]  = sample.cc;
      a[0]      = sample.as; /* cell parameters */
      a[1]      = sample.bs;
      a[2]      = sample.cs;
      aspv[0][0]= sample.ax; /* cell axis A */
      aspv[1][0]= sample.ay;
      aspv[2][0]= sample.az;
      aspv[0][1]= sample.bx; /* cell axis B */
      aspv[1][1]= sample.by;
      aspv[2][1]= sample.bz;

      /* default return values */
      strcpy(machine_real.message, "");
      machine_real.a3 = machine_real.a4 = 0;
      machine_real.a1 = machine_real.a5 = 0;

      /* if using HKL positioning in crystal (QM = 0) */
      if (machine_real.qm <= 0) {
        liquid_case = 0;
        /* compute reciprocal cell */
        for (i=0; i< 3; i++)
          if (a[i] <=0) sprintf(machine_real.message, "Lattice parameters a[%i]=%g", i, a[i]);
          else {
            a[i]    /= 2*PI;
            alpha[i]*= DEG2RAD;
            cosa[i]  = cos(alpha[i]);
            sina[i]  = sin(alpha[i]);
          }
        cc = cosa[0]*cosa[0]+cosa[1]*cosa[1]+cosa[2]*cosa[2]; /* nprm */
        cc = 1 + 2*cosa[0]*cosa[1]*cosa[2] - cc;
        if (cc <= 0) sprintf(machine_real.message, "Lattice angles (AA,BB,CC) cc=%g", cc);
        else cc = sqrt(cc);

        if (strlen(machine_real.message)) return machine_real;

        /* compute bb */
        j=1; k=2;
        for (i=0; i<3; i++) {
          b[i] = sina[i]/(a[i]*cc);
          cosb[i] = (cosa[j]*cosa[k] - cosa[i])/(sina[j]*sina[k]);
          sinb[i] = sqrt(1 - cosb[i]*cosb[i]);
          j=k; k=i;
        }

        bb[0][0] = b[0];
        bb[1][0] = 0;
        bb[2][0] = 0;
        bb[0][1] = b[1]*cosb[2];
        bb[1][1] = b[1]*sinb[2];
        bb[2][1] = 0;
        bb[0][2] = b[2]*cosb[1];
        bb[1][2] =-b[2]*sinb[1]*cosa[0];
        bb[2][2] = 1/a[2];

        /* compute vv */
        for (k=0; k< 3; k++)
          for (i=0; i< 3; i++) vv[k][i] = 0;

        for (k=0; k< 2; k++)
          for (i=0; i< 3; i++)
            for (j=0; j< 3; j++)
              vv[k][i] += bb[i][j]*aspv[j][k];

        for (m=2; m>=1; m--)
          for (n=0; n<3; n++) {
            i = (int)fmod(m+1,3); j= (int)fmod(m+2,3);
            k = (int)fmod(n+1,3); l= (int)fmod(n+2,3);
            vv[m][n]=vv[i][k]*vv[j][l]-vv[i][l]*vv[j][k];
          }

        for (i=0; i< 3; i++) { /* compute norm(vv) */
          c[i]=0;
          for (j=0; j< 3; j++)
            c[i] += vv[i][j]*vv[i][j];
          if (c[i]>0) c[i] =  sqrt(c[i]);
          else {
            sprintf(machine_real.message, "Vectors A and B, c[%i]=%g", i, c[i]);
            return machine_real;
          }
        }

        for (i=0; i< 3; i++) /* normalize vv */
          for (j=0; j< 3; j++)
            vv[j][i] /= c[j];

        for (i=0; i< 3; i++) /* compute S */
          for (j=0; j< 3; j++) {
            s[i][j] = 0;
            for (k=0; k< 3; k++)
              s[i][j] += vv[i][k]*bb[k][j];
          }
        s[3][3]=1;
        for (i=0;  i< 3;  i++)  s[3][i]=s[i][3]=0;

        /* compute q modulus and transverse component */
        machine_real.qs = 0;
        for (i=0; i< 3; i++) {
          machine_real.qt[i] = 0;
          for (j=0; j< 3; j++) machine_real.qt[i] += qhkl[j]*s[i][j];
          machine_real.qs += machine_real.qt[i]*machine_real.qt[i];
        }
        if (machine_real.qs > 0) machine_real.qm = sqrt(machine_real.qs);
        else sprintf(machine_real.message, "Q modulus too small QM^2=%g", machine_real.qs);
      } else {
        machine_real.qs = machine_real.qm*machine_real.qm;
      }
      /* end if  qm <= 0 ********************************************* */

      /* positioning of monochromator and analyser */
      arg = PI/machine_hkl.dm/machine_hkl.ki;
      if (fabs(arg > 1))
        sprintf(machine_real.message, "Monochromator can not reach this KI. arg=%g", arg);
      else {
        if (machine_hkl.dm <= 0 || machine_hkl.ki <= 0)
          strcpy(machine_real.message, "Monochromator DM=0 or KI=0.");
        else
          machine_real.a1 = asin(arg)*RAD2DEG;
        machine_real.a1 *= machine_hkl.sm;
      }
      machine_real.a2=2*machine_real.a1;

      arg = PI/machine_hkl.da/machine_hkl.kf;
      if (fabs(arg > 1))
        sprintf(machine_real.message, "Analyzer can not reach this KF. arg=%g",arg);
      else {
        if (machine_hkl.da <= 0 || machine_hkl.kf <= 0)
          strcpy(machine_real.message, "Analyzer DA=0 or KF=0.");
        else
          machine_real.a5 = asin(arg)*RAD2DEG;
        machine_real.a5 *= machine_hkl.sa;
      }
      machine_real.a6=2*machine_real.a5;
      if (strlen(machine_real.message)) return machine_real;


      /* code from TASMAD/t_conv.F:SAM_CASE */
      arg = (machine_hkl.ki*machine_hkl.ki + machine_hkl.kf*machine_hkl.kf - machine_real.qs)
          / (2*machine_hkl.ki*machine_hkl.kf);
      if (fabs(arg) < 1)
        machine_real.a4 = RAD2DEG*acos(arg);
      else
        sprintf(machine_real.message, "Q modulus too big. Can not close triangle. arg=%g", arg);
      machine_real.a4 *= machine_hkl.ss;

      if (!liquid_case) { /* compute a3 in crystals */
        machine_real.a3 =
            -atan2(machine_real.qt[1],machine_real.qt[0])
            -acos( (machine_hkl.kf*machine_hkl.kf-machine_real.qs-machine_hkl.ki*machine_hkl.ki)
                  /(-2*machine_real.qm*machine_hkl.ki) );
        machine_real.a3 *= RAD2DEG*(machine_real.a4 > 0 ? 1 : -1 );
  //machine_real.a3 = machine_real.a3 -90;// Add by PW & LU


      }

/* Analyser angles */
/* Angle of the analyser rack aa5 */
RR = dadw/dsa;
machine_real.aa5 = -RAD2DEG*(acos((sin(DEG2RAD*machine_real.a6)-(cos(DEG2RAD*machine_real.a6)+RR)*sqrt(RR*RR+2*RR*cos(DEG2RAD*machine_real.a6)))/(1+RR*RR+2*RR*cos(DEG2RAD*machine_real.a6)))); // From Bahl et al. , NIMB 226 (2004)
/*printf("(RITA Analyzer Rack Angle: AA5=%.4g[deg])\n", machine_real.aa5);*/
/*A4 angle between blades */
dsb = sqrt(dsa*dsa + ana_d*ana_d - 2*dsa*ana_d*cos(DEG2RAD*machine_real.aa5));
 machine_real.da4 = RAD2DEG*asin(ana_d*sin(DEG2RAD*abs(machine_real.aa5))/dsb);
//machine_real.da4 = RAD2DEG*atan(ana_d*cos(DEG2RAD*machine_real.aa5)/dsa);
printf("(RITA Analyzer Angle between blades: DA4=%.4g[deg])\n", machine_real.da4);
/*Blade angle rotation */
 machine_real.c1= -machine_real.aa5+machine_real.a5-(1-5)*machine_real.da4;
 machine_real.c2= -machine_real.aa5+machine_real.a5-(2-5)*machine_real.da4;
 machine_real.c3= -machine_real.aa5+machine_real.a5-(3-5)*machine_real.da4;
 machine_real.c4= -machine_real.aa5+machine_real.a5-(4-5)*machine_real.da4;
 machine_real.c5= -machine_real.aa5+machine_real.a5-(5-5)*machine_real.da4;
 machine_real.c6= -machine_real.aa5+machine_real.a5-(6-5)*machine_real.da4;
 machine_real.c7= -machine_real.aa5+machine_real.a5-(7-5)*machine_real.da4;
 machine_real.c8= -machine_real.aa5+machine_real.a5-(8-5)*machine_real.da4;
 machine_real.c9= -machine_real.aa5+machine_real.a5-(9-5)*machine_real.da4;

    return machine_real;
  }
/* END declaring stuff for the HKL calculator adapted from templateTAS */



%}
/* end of DECLARE */

USERVARS
%{
  /* Order of scattering from the monochromator */
  double Mono_order;
  int AnaBlade;
  int BinX;
  int BinY;
%}


INITIALIZE %{

/* Source parameters */
//l0=9.045/sqrt(EI);
l0=L0;// fix 20100222 by LU since we want source independent of mono
lmin=BPL*l0;/* MONO_N is the order of the reflection */
lmax=BPH*l0;
printf("Source wavelength interval=%.4g - %.4g[Ang] \n",lmin, lmax);
/* Calculate min and max energy for source to use with energy monitors before sample */
emini=(9.045*9.045)/(lmax*lmax);
emaxi=(9.045*9.045)/(lmin*lmin);
//neichan=floor((emaxi-emini)/0.01);
printf("Source energy interval=%.4g - %.4g[meV] \n",emini, emaxi);
/* Calculate min and max energy for the detectors after the sample */
eminf=EF*0.95;
emaxf=EF*1.05;
//nefchan=floor((emaxf-eminf)/0.005);
printf("Detector energy interval=%.4g - %.4g[meV] \n",eminf, emaxf);

/* Only SPLIT if we are not running with VIRTUAL parms */
if (VIRTUALIN || VIRTUALOUT) {
  SPLITREP = 1;
  if (!VIRTUALOUT){
    SPLITMREP = 1;
  }else{
    SPLITMREP = REP;
  }
  if (!VIRTUALIN){
    SPLITAREP = 1;
  }else{
    SPLITAREP = REP;
  }
 } else {
  SPLITREP = REP;
  SPLITMREP = REP;
  SPLITAREP = REP;
}


/* Guide parameters */
angleGuideCurved=20.0/2408.0;
/* calculate mirror reflectivity slope */
 ALPHA=0;//(R0-R)/(Qc*(M-1));
printf("* GUIDE MIRROR ALPHA=%g [AA]\n", ALPHA);

/* Monochromator parameters */
double Vi, Vf;
double tmp=0;
char Qmode = 0;
DM = 2*PI/mono_q;
DA = 2*PI/ana_q;
mono_mosaic_h = MONO_MOS_H; /* MON_MOSAIC; */
mono_mosaic_v = MONO_MOS_V; /* MON_MOSAIC; */
dms = dmv+dvs;/* Distance between mono and sample*/
dmc = rmh-lc; /*distance monochromator to front of monocollimator*/
/* fix 08/02/2009 by LU */
/* Forcing monochromator to set value */
 if (MONOFORCE) {/* Flat mono */
   rv=0;
   printf("* MONOCHROMATOR IS VERTICALLY FOCUSING: RV=%g \n",rv);
 }
 else {rv=dms; /* curved mono */
 printf("* MONOCHROMATOR IS FLAT: RV=%g \n", rv);
 }


/* BEGIN HKL calculator adapted from templateTAS */

machine_real.a1 = A1;
machine_real.a2 = A2;
machine_real.a3 = A3;
machine_real.a4 = A4;
machine_real.a5 = A5;
machine_real.a6 = A6;

machine_real.c1 = -C1;
machine_real.c2 = -C2;
machine_real.c3 = -C3;
machine_real.c4 = -C4;
machine_real.c5 = -C5;
machine_real.c6 = -C6;
machine_real.c7 = -C7;
machine_real.c8 = -C8;
machine_real.c9 = -C9;

/* energy conservation */
 if (EI && EF) {
   EN = EI - EF;
   fprintf(stderr,"%s WARNING: EN is now set to %g since you provided both EI=%g (KI) and EF=%g (KF)\n", NAME_INSTRUMENT, EN, EI, EF);
 } else if (EI && !EF){
   EF = EI - EN;
   fprintf(stderr,"%s WARNING: EF is now set to %g since you provided both EI=%g (KI) and EN=%g\n", NAME_INSTRUMENT, EF, EI, EN);
   fprintf(stderr,"%s POSSIBLE ERROR: On RITA, EF is normally SET, EI calculated - are you sure about this?\n", NAME_INSTRUMENT);
 } else if (EF && !EI) {
   EI = EF + EN;
   fprintf(stderr,"%s WARNING: EI is now set to %g since you provided both EF=%g (KF) and EN=%g\n", NAME_INSTRUMENT, EI, EF, EN);
 } else {
   fprintf(stderr,"%s WARNING: Neither EI, EF nor EN defined: Energies are set from user angle input:\n", NAME_INSTRUMENT);
   l0 = 2 * DM * sin(DEG2RAD*A2/2)/MONO_N;
   EI = 9.045/l0;
   EI = EI*EI;
   fprintf(stderr,"%s: WARNING: EI has been adjusted to %g[meV]  (A2 = %g[deg])\n", NAME_INSTRUMENT, EI, A2);
   l0 = 2 * DA * sin(DEG2RAD*A6/2)/MONO_N;
   EF= 9.045/l0;
   EF = EF*EF;
   fprintf(stderr,"%s: WARNING: EF has been adjusted to %g[meV]  (A6 = %g[deg])\n", NAME_INSTRUMENT, EF, A2);
 }
 /* determine remaining neutron energies */
 Vi = SE2V*sqrt(EI);
 KI = V2K*Vi;
 Vf = SE2V*sqrt(EF);
 KF = V2K*Vf;

/* transfered sample parameters */
sample.aa = AA;
sample.bb = BB;
sample.cc = CC;
sample.as = AS;
sample.bs = BS;
sample.cs = CS;
sample.ax = AH;
sample.ay = AK;
sample.az = AL;
sample.bx = BH;
sample.by = BK;
sample.bz = BL;

/* transfered target parameters */
machine_hkl.ki = KI;
machine_hkl.kf = KF;
machine_hkl.ei = EI;
machine_hkl.ef = EF;
machine_hkl.qh = QH;
machine_hkl.qk = QK;
machine_hkl.ql = QL;
machine_hkl.en = EN;
machine_real.qm = QM;

if (QM || QH || QK || QL) {
  Qmode=1;
  fprintf(stderr,"%s: Running in HKL mode\n", NAME_INSTRUMENT);
} else {
  fprintf(stderr,"%s: Running in angle mode\n", NAME_INSTRUMENT);
}

if (verbose && Qmode) {
  printf("%s: Detailed TAS configuration\n", NAME_INSTRUMENT);
  printf("* Incoming beam: EI=%.4g [meV] KI=%.4g [Angs-1] Vi=%g [m/s]\n", EI, KI, Vi);
  printf("* Outgoing beam: EF=%.4g [meV] KF=%.4g [Angs-1] Vf=%g [m/s]\n", EF, KF, Vf);
}

/* transfered machine parameters */
/* For W configuartion of TAS: */
machine_hkl.sm = SM;
machine_hkl.ss = SS;
machine_hkl.sa = SA;
/* These two are actually constants, see top of INITIALIZE */
machine_hkl.dm = DM;
machine_hkl.da = DA;

if (Qmode) {
  machine_real = qhkl2angles(sample, machine_hkl, machine_real);
  if (strlen(machine_real.message)) {
    exit(fprintf(stderr, "%s: ERROR: %s [qhkl2angles]\n", NAME_INSTRUMENT, machine_real.message));
  }
}



if (A1 && (machine_real.a1 != A1)) {
   printf("Warning, resetting A1 from calculated %g to user provided %g\n",machine_real.a1,A1);
   machine_real.a1 = A1;
}
if (A2 && (machine_real.a2 != A2)) {
   printf("Warning, resetting A2 from calculated %g to user provided %g\n",machine_real.a2,A2);
   machine_real.a2 = A2;
}
if (A3 && (machine_real.a3 != A3)) {
   printf("Warning, resetting A3 from calculated %g to user provided %g\n",machine_real.a3,A3);
   machine_real.a3 = A3;
}
if (A4 && (machine_real.a4 != A4)) {
   printf("Warning, resetting A4 from calculated %g to user provided %g\n",machine_real.a4,A4);
   machine_real.a4 = A4;
}
if (A5 && (machine_real.a5 != A5)) {
   printf("Warning, resetting A5 from calculated %g to user provided %g\n",machine_real.a5,A5);
   machine_real.a5 = A5;
}
if (A6 && (machine_real.a6 != A6)) {
   printf("Warning, resetting A6 from calculated %g to user provided %g\n",machine_real.a6,A6);
   machine_real.a6 = A6;
}


if (AA5 && (machine_real.aa5 != AA5) && !ANAFORCE) {
   printf("Warning, resetting AA5 from calculated %g to user provided %g\n",machine_real.aa5,AA5);
   machine_real.aa5 = AA5;
   /*A4 angle between blades */
machine_real.da4 = RAD2DEG*atan(ana_d*cos(DEG2RAD*machine_real.aa5)/dsa);
/*Blade angle rotation */
 machine_real.c1= -machine_real.aa5+machine_real.a5-(1-5)*machine_real.da4;
 machine_real.c2= -machine_real.aa5+machine_real.a5-(2-5)*machine_real.da4;
 machine_real.c3= -machine_real.aa5+machine_real.a5-(3-5)*machine_real.da4;
 machine_real.c4= -machine_real.aa5+machine_real.a5-(4-5)*machine_real.da4;
 machine_real.c5= -machine_real.aa5+machine_real.a5-(5-5)*machine_real.da4;
 machine_real.c6= -machine_real.aa5+machine_real.a5-(6-5)*machine_real.da4;
 machine_real.c7= -machine_real.aa5+machine_real.a5-(7-5)*machine_real.da4;
 machine_real.c8= -machine_real.aa5+machine_real.a5-(8-5)*machine_real.da4;
 machine_real.c9= -machine_real.aa5+machine_real.a5-(9-5)*machine_real.da4;
 printf("Warning, recalculating C's to %g %g %g %g %g %g %g %g %g \n",machine_real.c1,machine_real.c2,machine_real.c3,machine_real.c4,machine_real.c5,machine_real.c6,machine_real.c7,machine_real.c8,machine_real.c9);
}


if (C1 && (machine_real.c1 != -C1)) {
   printf("Warning, resetting C1 from calculated %g to user provided %g\n",machine_real.c1,-C1);
   machine_real.c1 = -C1;
}
if (C2 && (machine_real.c2 != -C2)) {
   printf("Warning, resetting C2 from calculated %g to user provided %g\n",machine_real.c2,-C2);
   machine_real.c2 = -C2;
}
if (C3 && (machine_real.c3 != -C3)) {
   printf("Warning, resetting C3 from calculated %g to user provided %g\n",machine_real.c3,-C3);
   machine_real.c3 = -C3;
}
if (C4 && (machine_real.c4 != -C4)) {
   printf("Warning, resetting C4 from calculated %g to user provided %g\n",machine_real.c4,-C4);
   machine_real.c4 = -C4;
}
if (C5 && (machine_real.c5 != -C5)) {
   printf("Warning, resetting C5 from calculated %g to user provided %g\n",machine_real.c5,-C5);
   machine_real.c5 = -C5;
}
if (C6 && (machine_real.c6 != -C6)) {
   printf("Warning, resetting C6 from calculated %g to user provided %g\n",machine_real.c6,-C6);
   machine_real.c6 = -C6;
}
if (C7 && (machine_real.c7 != -C7)) {
   printf("Warning, resetting C7 from calculated %g to user provided %g\n",machine_real.c7,-C7);
   machine_real.c7 = -C7;
}
if (C8 && (machine_real.c8 != -C8)) {
   printf("Warning, resetting C8 from calculated %g to user provided %g\n",machine_real.c8,-C8);
   machine_real.c8 = -C8;
}
if (C9 && (machine_real.c9 != -C9)) {
   printf("Warning, resetting C9 from calculated %g to user provided %g\n",machine_real.c9,-C9);
   machine_real.c9 = -C9;
}

/* Forcing ana to set value */
 if (ANAFORCE) {/* Flat ana mode */
   if (AA5 && machine_real.a5 != AA5){
     machine_real.aa5=AA5;
     printf("Warning, resetting AA5 from calculated %g to user provided %g\n",machine_real.a5,AA5);
   } else{
     machine_real.aa5=machine_real.a5;
   }
   machine_real.c1=0;
   machine_real.c2=0;
   machine_real.c3=0;
   machine_real.c4=0;
   machine_real.c5=0;
   machine_real.c6=0;
   machine_real.c7=0;
   machine_real.c8=0;
   machine_real.c9=0;
   C1 = machine_real.c1;
   C2 = machine_real.c2;
   C3 = machine_real.c3;
   C4 = machine_real.c4;
   C5 = machine_real.c5;
   C6 = machine_real.c6;
   C7 = machine_real.c7;
   C8 = machine_real.c8;
   C9 = machine_real.c9;
   printf("* ANALYSER IN FLAT MODE \n");
 }
 else {/* else imaging mode ana */
   printf("* ANALYSER IN IMAGING MODE \n");
 }


if (verbose) {
  printf("* Transfered:     EN=%g [meV] QM=%g [Angs-1]\n", EN, machine_real.qm);
  printf("Angles: A1=%.4g A2=%.4g A3=%.4g A4=%.4g A5=%.4g A6=%.4g [deg]\n",
    machine_real.a1, machine_real.a2,
    machine_real.a3, machine_real.a4,
    machine_real.a5, machine_real.a6);
    printf("(RITA Analyzer Angles: AA5=%.4g C1=%.4g C2=%.4g C3=%.4g C4=%.4g C5=%.4g C6=%.4g C7=%.4g C8=%.4g C9=%.4g [deg])\n",
   machine_real.aa5,
   machine_real.c1,machine_real.c2,machine_real.c3,
   machine_real.c4,machine_real.c5,machine_real.c6,
   machine_real.c7,machine_real.c8,machine_real.c9);
}

/* END HKL calculator adapted from templateTAS */



/* Sample paramaters */
if (SAMPLE == 2){
  /* Powder sample, relevant WHEN at sample position */
  if (!strcmp(SAMPLEFILE,"default")) {
    PowderFile="Al2O3_sapphire.lau";
  } else {
    PowderFile=SAMPLEFILE;
  }
  SingleXFile="";
} else if (SAMPLE == 3) {
  /* Single crystal sample, relevant WHEN at sample position */
  if (!strcmp(SAMPLEFILE,"default")) {
    SingleXFile="Pb.laz";
  } else {
    SingleXFile=SAMPLEFILE;
  }
  PowderFile="";
} else if (SAMPLE==4){
  /* Phonon sample */
 printf("Phonon sample used - no elastic contribution. \n");
}else {
  /* Purely incoherent scatterer, relevant WHEN at sample position */
  SingleXFile="";
  PowderFile="";
}


/* Coarse collimator */
if (COARSE) {coarse = 1;}
else {coarse = 0;}

/* Opening slits of coarse collimator depends on analyzer settings: */
deltaL = 2 * WindowSize * cos(DEG2RAD*machine_real.a6 + fabs(DEG2RAD*machine_real.a5));
FirstWindowSizeL = WindowSize * (dad-BladeLength-deltaL+BladeLength*sin(DEG2RAD*machine_real.a6 + fabs(DEG2RAD*machine_real.a5)))/(dad-deltaL);
FirstWindowSizeR = WindowSize * (dad-BladeLength+deltaL+BladeLength*sin(DEG2RAD*machine_real.a6 + fabs(DEG2RAD*machine_real.a5)))/(dad+deltaL);

/* Window positions */
XwinMin[0] = 0;
XwinMin[1] = 12;
XwinMin[2] = 24;
XwinMin[3] = 36;
XwinMin[4] = 48;
XwinMin[5] = 60;
XwinMin[6] = 72;
XwinMin[7] = 84;
XwinMin[8] = 95;
XwinMin[9] = 106;

XwinMax[0] = 0;
XwinMax[1] = 21;
XwinMax[2] = 33;
XwinMax[3] = 45;
XwinMax[4] = 57;
XwinMax[5] = 69;
XwinMax[6] = 81;
XwinMax[7] = 93;
XwinMax[8] = 104;
XwinMax[9] = 115;

YwinMin[0] = 0;
YwinMin[1] = 39;
YwinMin[2] = 39;
YwinMin[3] = 39;
YwinMin[4] = 39;
YwinMin[5] = 39;
YwinMin[6] = 39;
YwinMin[7] = 39;
YwinMin[8] = 39;
YwinMin[9] = 39;

YwinMax[0] = 0;
YwinMax[1] = 91;
YwinMax[2] = 91;
YwinMax[3] = 91;
YwinMax[4] = 91;
YwinMax[5] = 91;
YwinMax[6] = 91;
YwinMax[7] = 91;
YwinMax[8] = 91;
YwinMax[9] = 91;


%}

/* end of INITIALIZE */

TRACE

/* Source description */
COMPONENT armSource = Progress_bar()
  AT (0,0,0) ABSOLUTE


COMPONENT source = Source_gen4 (
    h = 0.135, w = 0.08, xw = 0.03, yh = 0.12,
    dist = 1.465, Lmin=lmin, Lmax=lmax, /*Lmin=0.1 LMax=10*/
    T1=301.287, I1=ITAR*(1.27e13/4/PI),
    T2=105.655,I2=ITAR*(3.818e12/4/PI),
    T3=25.379,I3=ITAR*(2.331e12/4/PI),
    HEtailA=ITAR*8.306e11/4/PI, HEtailL0=-0.398)
  WHEN (!VIRTUALIN)  AT (0,0,0) RELATIVE armSource  ROTATED (0,0,0) RELATIVE armSource


COMPONENT slitGuideBegin = Slit(
  xmin = -0.015, xmax = 0.015,
  ymin = -0.06, ymax = 0.06)
  WHEN (!VIRTUALIN)  AT (0,0,1.464999) RELATIVE armSource


COMPONENT lmon_guide_start = L_monitor(
    nL = 100, filename = "lmon_guide_start.dat", xmin = -0.02, xmax = 0.02,
    ymin = -0.075, ymax = 0.075, Lmin = lmin, Lmax = lmax, restore_neutron = 1)
WHEN (!VIRTUALIN)  AT (0, 0, 1.4649992 ) RELATIVE armSource

COMPONENT guideStraight = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 4.628,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 1.465) RELATIVE armSource

 /* 0.035m gap*/

COMPONENT guideCurved1 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 4.663) RELATIVE guideStraight
  ROTATED (0,angleGuideCurved,0) RELATIVE guideStraight

COMPONENT guideCurved2 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved1
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved1

COMPONENT guideCurved3 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved2
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved2

COMPONENT guideCurved4 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved3
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved3

COMPONENT guideCurved5 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved4
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved4

COMPONENT guideCurved6 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved5
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved5

COMPONENT guideCurved7 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved6
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved6

COMPONENT guideCurved8 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved7
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved7

COMPONENT guideCurved9 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN) AT (0, 0, 0.5) RELATIVE guideCurved8
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved8

COMPONENT guideCurved10 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved9
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved9

COMPONENT guideCurved11 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved10
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved10

COMPONENT guideCurved12 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved11
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved11

COMPONENT guideCurved13 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved12
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved12

COMPONENT guideCurved14 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN) AT (0, 0, 0.5) RELATIVE guideCurved13
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved13

COMPONENT guideCurved15 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved14
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved14

COMPONENT guideCurved16 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved15
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved15

COMPONENT guideCurved17 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved16
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved16

COMPONENT guideCurved18 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved17
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved17

COMPONENT guideCurved19 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved18
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved18

COMPONENT guideCurved20 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved19
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved19

COMPONENT guideCurved21 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved20
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved20

COMPONENT guideCurved22 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved21
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved21

COMPONENT guideCurved23 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved22
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved22

COMPONENT guideCurved24 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved23
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved23

COMPONENT guideCurved25 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved24
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved24

COMPONENT guideCurved26 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved25
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved25

COMPONENT guideCurved27 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN) AT (0, 0, 0.5) RELATIVE guideCurved26
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved26

COMPONENT guideCurved28 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved27
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved27

COMPONENT guideCurved29 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved28
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved28

COMPONENT guideCurved30 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved29
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved29

COMPONENT guideCurved31 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved30
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved30

COMPONENT guideCurved32 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved31
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved31

COMPONENT guideCurved33 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved32
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved32

COMPONENT guideCurved34 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved33
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved33

COMPONENT guideCurved35 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved34
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved34

COMPONENT guideCurved36 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved35
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved35

COMPONENT guideCurved37 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved36
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved36

COMPONENT guideCurved38 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved37
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved37

COMPONENT guideCurved39 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved38
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved38

COMPONENT guideCurved40 = Guide(
        w1 = 0.03, h1 = 0.12, w2 = 0.03, h2 = 0.12, l = 0.499995,
         R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)  AT (0, 0, 0.5) RELATIVE guideCurved39
  ROTATED (0,angleGuideCurved,0) RELATIVE guideCurved39


/* bunker wall, m=2, 3.0 m */
 COMPONENT bunker = Guide(w1= 0.03, h1=0.12, w2=0.03, h2=0.12,
        l=3.45, R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)    AT (0,0,0.5) RELATIVE guideCurved40

/*15cm space for main shutter*/

/* guide segment 3, m=2, 1.3 m */
 COMPONENT guide3 = Guide(w1=0.03, h1=0.12, w2=0.03, h2=0.12,
        l=5.2, R0 = R0, Qc = Qc, alpha = ALPHA, m = M, W = W)
WHEN (!VIRTUALIN)    AT (0,0,3.6) RELATIVE bunker


COMPONENT slitGuideEnd = Slit(xmin=-0.016,xmax=0.016, ymin=-0.061,ymax=0.061)
WHEN (!VIRTUALIN)  AT (0,0,5.2001) RELATIVE guide3

COMPONENT psd_guide_end = PSD_monitor(
  xwidth=0.04, yheight=0.15,
  nx=64, ny=64, filename="psd_guide_end.dat", restore_neutron = 1)
WHEN (!VIRTUALIN)  AT (0, 0, 0.04) RELATIVE slitGuideEnd

COMPONENT emon_guide_end = E_monitor(
 nE = 100,filename = "emon_guide_end.dat",xwidth=0.04, yheight=0.15,  Emin = emini, Emax = emaxi, restore_neutron = 1)
  AT (0, 0, 0.05) RELATIVE slitGuideEnd

COMPONENT lmon_guide_end = L_monitor(
  nL = 100, filename = "lmon_guide_end.dat", xwidth=0.04, yheight=0.15,Lmin = lmin, Lmax = lmax, restore_neutron = 1)
WHEN (!VIRTUALIN)  AT (0, 0, 0.06) RELATIVE slitGuideEnd

COMPONENT divmon_guide_end = Divergence_monitor(
    nh = 128, nv = 128, filename = "divmon_guide_end.dat",
    restore_neutron = 1, xwidth = 0.04, yheight = 0.15,
    maxdiv_h = 2, maxdiv_v = 2, restore_neutron = 1)
WHEN (!VIRTUALIN)  AT (0, 0, 0.07) RELATIVE slitGuideEnd


COMPONENT focus_mono = Arm()
WHEN (!VIRTUALIN)  AT (0, 0, 0.15) RELATIVE slitGuideEnd
  ROTATED (0, machine_real.a1, 0) RELATIVE slitGuideEnd
EXTEND %{
  Mono_order=-1;
%}

SPLIT SPLITMREP COMPONENT monochromator_curved = Monochromator_curved(
    reflect = "HOPG.rfl", zwidth = 0.15, yheight = 0.025,
    gap = 0.001, NH = 1, NV = 5, mosaich = mono_mosaic_h,
    mosaicv = mono_mosaic_v, r0 = 1, Q = 1.8734, RV = rv,
    RH = 0)
WHEN (!VIRTUALIN)  AT (0, 0, 0) RELATIVE focus_mono
EXTEND %{
  if(SCATTERED) Mono_order=order;
%}

COMPONENT a2 = Arm()
WHEN (!VIRTUALIN)  AT (0,0,0) RELATIVE focus_mono
  ROTATED (0, machine_real.a2, 0) RELATIVE slitGuideEnd


COMPONENT slitShutter = Slit(
  xmin = -0.02, xmax = 0.02,
  ymin = -0.075, ymax = 0.075)
WHEN (!VIRTUALIN)  AT (0, 0, 0.2) RELATIVE a2
  ROTATED (0,0,0) RELATIVE a2

COMPONENT MSCollimator = Collimator_linear(
     xmin = -0.02, xmax = 0.02, ymin = -0.075, ymax = 0.075, length = 0.20 ,
     divergence = COLL_MS )
WHEN (!VIRTUALIN) AT (0, 0, dmc) RELATIVE a2


COMPONENT infilter = Filter_gen(
xmin = -0.02, xmax = 0.02, ymin = -0.075, ymax = 0.075, options="wavevector multiply", filename=INFILTERFILE)
WHEN (INFILTER>0 && !VIRTUALIN) AT (0, 0,0.201) RELATIVE MSCollimator


COMPONENT psd_virt = PSD_monitor(
    nx=128, ny=128, filename="psd_virt.dat", xwidth = 0.10, yheight = 0.10 )
WHEN (!VIRTUALIN)  AT (0, 0, dmv-0.01) RELATIVE a2


COMPONENT lmon_virt = L_monitor(
    nL = 100, filename = "lmon_virt.dat", xwidth=0.06, yheight=0.10,//the size of the beam
    Lmin = lmin, Lmax = lmax, restore_neutron=1)
WHEN (!VIRTUALIN)  AT (0, 0, dmv) RELATIVE a2 // Needs to be here in order to propagate virtual_out to dmv

COMPONENT aa2 = Arm()
WHEN (!VIRTUALIN)  AT (0,0,dmv) RELATIVE a2


/* If the virtual output is generated the secondary spectrometer IS NOT simulated */
/*  COMPONENT virtualout = Virtual_output(filename=SOURCEFILE, bufsize=0)*/
/*WHEN (!VIRTUALIN && VIRTUALOUT)  AT (0, 0, 0) RELATIVE aa2*/
/**/
/**/
/*COMPONENT virtualsource = Virtual_input(filename=SOURCEFILE, verbose=1, repeat_count=REP)*/
/*WHEN (VIRTUALIN && !VIRTUALOUT) AT (0,0,0) RELATIVE aa2*/


COMPONENT OrderMon = Monitor_nD(
      xwidth = 0.04, yheight = 0.08,
      user1="Mono_order", username1="Intensity of multiples",
      options="user1 bins=5 limits=[0.5 5.5], user1 bins=5 limits=[0.5 5.5]", filename="OrderMonNtest.dat", restore_neutron=1)
 WHEN (Mono_order && !VIRTUALOUT)    AT (0, 0, 0.170-4e-3) RELATIVE aa2

COMPONENT kMoni = PSD_monitor_psf_eff(eff=EFF, k0=1.553, xwidth = 0.04, yheight = 0.08, restore_neutron=1,filename="kMoni.dat")
  WHEN(!VIRTUALOUT)     AT (0, 0, 0.170-3e-3) RELATIVE aa2

COMPONENT kMoni1st = PSD_monitor_psf_eff(eff=EFF, k0=1.553, xwidth = 0.04, yheight = 0.08, restore_neutron=1,filename="kMoni1st.dat")
WHEN (Mono_order==1 && !VIRTUALOUT) AT (0, 0, 0.170-2e-3) RELATIVE aa2

COMPONENT kMoni2nd = PSD_monitor_psf_eff(eff=EFF, k0=1.553, xwidth = 0.04, yheight = 0.08, restore_neutron=1,filename="kMoni2nd.dat")
WHEN (Mono_order==2 && !VIRTUALOUT) AT (0, 0, 0.170-1e-3) RELATIVE aa2


COMPONENT kMoni3rd = PSD_monitor_psf_eff(eff=EFF,  k0=1.553, xwidth = 0.04, yheight = 0.08, restore_neutron=1,filename="kMoni3rd.dat")
WHEN (Mono_order==3 && !VIRTUALOUT) AT (0, 0, 0.170) RELATIVE aa2


COMPONENT slitMonochromator = Slit(
  xmin = -MSL/1000.0, xmax = MSR/1000.0,
  ymin = -MSB/1000.0, ymax = MST/1000.0)
WHEN(!VIRTUALOUT) AT (0, 0, 0.285) RELATIVE aa2

COMPONENT Perspex = Incoherent(Vc=1,sigma_abs=0.019,sigma_inc=4.7, xwidth = 0.1, yheight = 0.1, zdepth = PTHICK, p_interact=1e-2)
 WHEN (PERSPEX>0 && !VIRTUALOUT) AT (0, 0, 0.3) RELATIVE aa2

COMPONENT psd_samplepos_1cm2 = PSD_monitor( xwidth = 0.01, yheight = 0.01, restore_neutron=1,filename="psd_samplepos_1cm2.dat")
WHEN(!VIRTUALOUT)    AT (0, 0, dvs) RELATIVE aa2

COMPONENT emon_samplepos_1cm2 = E_monitor(
  nE=100, filename="emon_samplepos_1cm2.dat",xwidth = 0.01, yheight = 0.01, Emin=emini, Emax=emaxi, restore_neutron=1)
WHEN(!VIRTUALOUT)   AT (0, 0, dvs) RELATIVE aa2

COMPONENT divmon_samplepos_1cm2 = Divergence_monitor(
    nh = 128, nv = 128, filename = "divmon_samplepos_1cm2.dat",
    xwidth = 0.01, yheight = 0.01, maxdiv_h = 3, maxdiv_v = 3, restore_neutron=1)
WHEN(!VIRTUALOUT)  AT (0, 0, dvs) RELATIVE aa2

COMPONENT psd_samplepos_large = PSD_monitor( xwidth = 0.06, yheight = 0.06, restore_neutron=1,filename="psd_samplepos_large.dat")
WHEN(!VIRTUALOUT)    AT (0, 0, dvs) RELATIVE aa2

SPLIT SPLITREP COMPONENT a3 = Arm()
WHEN(!VIRTUALOUT) AT (0,0,dvs) RELATIVE aa2
 ROTATED (0, machine_real.a3, 0) RELATIVE aa2

COMPONENT aa3 = Arm()
WHEN(!VIRTUALOUT) AT (0,0,0) RELATIVE a3
 ROTATED (0, 0, TILT) RELATIVE a3


COMPONENT incohSample = Incoherent( // The size is the standard  Vanadium sample at RITA
/*      V0=1,sig_a=0,sig_i=4.70,radius_i=0, radius_o=0.0025, h=0.068,target_index=10, // perspex sample  */
    Vc=13.827,sigma_abs=5.08,sigma_inc=5.08,thickness=0.0201/2-0.0149/2, radius=0.0201/2, yheight=0.0187,target_index=10, // vanadium sample
    focus_xw=0.16,focus_yh=0.16) // focus on ana_slit1.
  WHEN (SAMPLE==1 && !VIRTUALOUT)
  AT (0, 0, 0) RELATIVE aa3

COMPONENT powderSample = PowderN( // The size is the standard sapphire sample at RITA
    reflections=PowderFile, d_phi = 12 , radius = 0.0068 , yheight = 0.015, pack = 1 , Vc = 254.52, sigma_abs = 0.4625 ,
    sigma_inc = 0.0188, p_inc = 0, barns=BARNS)
  WHEN (SAMPLE==2 && !VIRTUALOUT)
  AT (0, 0, 0) RELATIVE aa3

COMPONENT crystalSample = Single_crystal(//order=1, p_transmit=0,
     reflections = SingleXFile, mosaic=MOS, xwidth = SAMPLESIZE, yheight = SAMPLESIZE,zdepth = SAMPLESIZE, delta_d_d=DD_D,
     ax =AAX , ay = AAY, az =AAZ , bx =BBX , by =BBY, bz =BBZ , cx =CCX , cy =CCY , cz =CCZ, barns=BARNS)
   WHEN (SAMPLE==3 && !VIRTUALOUT)
  AT (0, 0, 0) RELATIVE aa3
EXTEND %{
  if(!SCATTERED) ABSORB;
  %}


COMPONENT phononSample = Phonon_simple(
    radius = SAMPLESIZE/2, yheight = SAMPLESIZE, focus_xw =0.16 ,
    focus_yh = 0.16, sigma_abs = 0.171, sigma_inc = 0.003,
    a = 4.95, b = 4.95, c = 10, M=207.2, DW = 1, T = 300,
    target_index = 7)// Focus on the ana_slit1
WHEN (SAMPLE==4 && !VIRTUALOUT)  AT (0, 0, 0) RELATIVE aa3

COMPONENT psd_4pi = PSD_monitor_4PI(
    nx = 1000, ny = 1000, filename = "psd_4pi",
    restore_neutron = 1, radius = 0.1)
 WHEN (!VIRTUALOUT) AT (0, 0, dvs) RELATIVE aa2

COMPONENT a4 = Arm()
 WHEN (!VIRTUALOUT) AT (0, 0, 0) RELATIVE a3
  ROTATED (0, machine_real.a4, 0) RELATIVE aa2


COMPONENT slitSample = Slit(
  xmin = -SSL/1000.0, xmax = SSR/1000.0,
  ymin = -SSB/1000.0, ymax = SST/1000.0)
  WHEN (!VIRTUALOUT) AT (0, 0, d_sample_slit) RELATIVE a4


COMPONENT filter_coll=Exact_radial_coll(
  radius=0.4525, nslit=9, d=0.000125,verbose=1,h_in=0.2,h_out=0.2,length=0.0988,theta_min=-10.26/2,theta_max=10.26/2
)
WHEN (OUTFILTER>0 && !VIRTUALOUT)
AT (0, 0, d_sample_filter-0.4525-0.0988/2) RELATIVE a4

COMPONENT filter = Filter_gen(
    xmin = -0.1, xmax = 0.1, ymin = -0.1, ymax = 0.1, options="wavevector multiply", filename=OUTFILTERFILE)
  WHEN (OUTFILTER>0 && !VIRTUALOUT)
  AT (0, 0,0.4525+0.0988+0.0001) RELATIVE filter_coll



COMPONENT ana_slit1 = Slit(
  xmin = -0.158/2, xmax = 0.158/2,
  ymin = -0.08/2, ymax = 0.08/2)
WHEN (!VIRTUALOUT)  AT (0, 0, d_sample_filter+0.13) RELATIVE a4
  ROTATED (0,0,0) RELATIVE a4

COMPONENT ana_slit2 = Slit(
  xmin = -0.158/2, xmax = 0.158/2,
  ymin = -0.103/2, ymax = 0.103/2)
WHEN (!VIRTUALOUT)  AT (0, 0, d_sample_filter+0.33) RELATIVE a4
  ROTATED (0,0,0) RELATIVE a4


COMPONENT emon_before_ana = E_monitor(
  xmin = -0.060, xmax = 0.060, ymin = -0.085, ymax = 0.085,
  Emin=eminf, Emax=emaxi, nE=100, filename="emon_before_ana.dat", restore_neutron=1)
WHEN (!VIRTUALOUT)  AT (0, 0, dsa-5*ana_d) RELATIVE a4

COMPONENT psd_before_ana = PSD_monitor(
    nx = 128, ny = 128, filename = "psd_before_ana.dat",
    //xmin = -0.1, xmax = 0.1, ymin = -0.1, ymax = 0.1)
    xmin = -0.060, xmax = 0.060, ymin = -0.085, ymax = 0.085)
WHEN (!VIRTUALOUT)  AT (0, 0, dsa-5*ana_d+0.001) RELATIVE a4

COMPONENT divmon_before_ana = Divergence_monitor(
    nh = 128, nv = 128, filename = "divmon_before_ana",
    xwidth = 0.02, yheight = 0.17, maxdiv_h = 1, maxdiv_v = 1, restore_neutron = 1)
WHEN (!VIRTUALOUT)  AT (0, 0, dsa-5*ana_d+0.002) RELATIVE a4

SPLIT SPLITAREP COMPONENT focus_ana = Arm()
WHEN (!VIRTUALOUT)  AT (0, 0, dsa) RELATIVE a4
ROTATED (0, machine_real.aa5, 0) RELATIVE a4
EXTEND %{
  AnaBlade=0;
%}

COMPONENT an1l= Monochromator_flat(// The mosaic values for this blade have not been checked!
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=40, mosaicv=40,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -4*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c1, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 1;
%}

COMPONENT an1u= Monochromator_flat( // The mosaic values for this blade have not been checked!
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=40, mosaicv=40,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -4*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c1, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 1;
%}


COMPONENT an2l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=38.7, mosaicv=38.7,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -3*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c2, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 2;
%}

COMPONENT an2u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=43.0, mosaicv=43.0,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -3*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c2, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 2;
%}

COMPONENT an3l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=31.1, mosaicv=31.1,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -2*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c3, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 3;
%}

COMPONENT an3u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=35.5, mosaicv=35.5,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -2*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c3, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 3;
%}


COMPONENT an4l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=27.2, mosaicv=27.2,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, -1*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c4, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 4;
%}

COMPONENT an4u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=30.4, mosaicv=30.4,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT) AT (0, 0, -1*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c4, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 4;
%}

COMPONENT an5l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=36.6, mosaicv=36.6,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 0) RELATIVE focus_ana
  ROTATED (0, machine_real.c5, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 5;
%}

COMPONENT an5u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=35.9, mosaicv=35.9,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 0) RELATIVE focus_ana
  ROTATED (0, machine_real.c5, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 5;
%}

COMPONENT an6l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=31.5, mosaicv=31.5,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c6, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 6;
%}

COMPONENT an6u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=36.1, mosaicv=36.1,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c6, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 6;
%}


COMPONENT an7l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=33.1, mosaicv=33.1,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 2*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c7, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 7;
%}

COMPONENT an7u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=37.2, mosaicv=37.2,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 2*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c7, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 7;
%}

COMPONENT an8l= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=46.8, mosaicv=46.8,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 3*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c8, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 8;
%}

COMPONENT an8u= Monochromator_flat(
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=51.3, mosaicv=51.3,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 3*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c8, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 8;
%}

COMPONENT an9l= Monochromator_flat(// The mosaic values for this blade have not been checked!
  zmin=-wan/2.0, zmax=wan/2.0, ymin=-ana_h/2.0, ymax=0,
  mosaich=40, mosaicv=40,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 4*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c9, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 9;
%}

COMPONENT an9u= Monochromator_flat( // The mosaic values for this blade have not been checked!
  zmin=-wan/2.0, zmax=wan/2.0, ymin=0, ymax=ana_h/2.0,
  mosaich=40, mosaicv=40,
  r0=ana_r0, Q=ana_q)
WHEN (!VIRTUALOUT)  AT (0, 0, 4*ana_d) RELATIVE focus_ana
  ROTATED (0, machine_real.c9, 0) RELATIVE focus_ana
//GROUP ANA
EXTEND %{
  if(SCATTERED) AnaBlade = 9;
%}


COMPONENT a6 = Arm()
WHEN (!VIRTUALOUT)  AT (0,0,0) RELATIVE focus_ana
  ROTATED (0, machine_real.a6, 0) RELATIVE a4

COMPONENT emon_before_coarse = E_monitor(
    nE=100, filename = "emon_before_coarse.dat", Emin=eminf, Emax=emaxf,
    xmin = -0.01, xmax = 0.01, ymin = -0.15, ymax = 0.15, restore_neutron=1)
    //xmin = -0.060, xmax = 0.060, ymin = -0.085, ymax = 0.085, restore_neutron=1)
WHEN (!VIRTUALOUT)  AT (0, 0, dad-BladeLength-0.02) RELATIVE a6

COMPONENT psd_before_coarse = PSD_monitor(
    nx = 128, ny = 128, filename = "psd_before_coarse.dat",
    xmin = -0.15, xmax = 0.15, ymin = -0.15, ymax = 0.15, restore_neutron=1)
    //xmin = -0.060, xmax = 0.060, ymin = -0.085, ymax = 0.085, restore_neutron=1)
WHEN (!VIRTUALOUT)  AT (0, 0, dad-BladeLength-0.01) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6


COMPONENT ArmR1 = Arm()
WHEN (!VIRTUALOUT) AT (-WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,-RC/9,0) RELATIVE a6

COMPONENT BladeR1 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0 && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmR1
EXTEND %{
  if (SCATTERED) printf("Absorption in R1\n");
%}

COMPONENT ArmR2 = Arm()
 AT (-3*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,-3*RC/9,0) RELATIVE a6

COMPONENT BladeR2 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0 && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmR2

COMPONENT ArmR3 = Arm()
WHEN (!VIRTUALOUT) AT (-5*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,-5*RC/9,0) RELATIVE a6

COMPONENT BladeR3 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
     WHEN (COARSE>0) AT (0,0,0) RELATIVE ArmR3

COMPONENT ArmR4 = Arm()
WHEN (!VIRTUALOUT) AT (-7*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,-7*RC/9,0) RELATIVE a6

COMPONENT BladeR4 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
     WHEN (COARSE>0) AT (0,0,0) RELATIVE ArmR4

COMPONENT ArmR5 = Arm()
WHEN (!VIRTUALOUT) AT (-9*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,-9*RC/9,0) RELATIVE a6

COMPONENT BladeR5 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
     WHEN (COARSE>0) AT (0,0,0) RELATIVE ArmR5



COMPONENT ArmL1 = Arm()
WHEN (!VIRTUALOUT) AT (WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,LC/9,0) RELATIVE a6

COMPONENT BladeL1 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0 && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmL1


COMPONENT ArmL2 = Arm()
 AT (3*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,3*LC/9,0) RELATIVE a6

COMPONENT BladeL2 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0  && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmL2

COMPONENT ArmL3 = Arm()
 AT (5*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,5*LC/9,0) RELATIVE a6

COMPONENT BladeL3 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0  && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmL3

COMPONENT ArmL4 = Arm()
WHEN (!VIRTUALOUT) AT (7*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,7*LC/9,0) RELATIVE a6

COMPONENT BladeL4 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0 && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmL4

COMPONENT ArmL5 = Arm()
WHEN (!VIRTUALOUT) AT (9*WindowSize/2,0,dad-0.005) RELATIVE a6
ROTATED (0,9*LC/9,0) RELATIVE a6
EXTEND %{
  BinX = 0; BinY=0;
%}

COMPONENT BladeL5 = Absorber(xmin=-BladeThickness/2, xmax=BladeThickness/2,
           ymin=-BladeHeight/2, ymax=BladeHeight/2,
           zmin=-BladeLength,zmax=0)
WHEN (COARSE>0 && !VIRTUALOUT) AT (0,0,0) RELATIVE ArmL5

  COMPONENT psd_detector = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
   ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_detector.dat",restore_neutron=0) // restore_neutron needs to be 0 to propagate the netrons in use for EXTEND section
 WHEN (!VIRTUALOUT) AT (0, 0, dad+0.0215) RELATIVE a6
 ROTATED (0,180,0) RELATIVE a6
EXTEND %{
  BinX = floor((x - xmin)*nx/(xmax - xmin)); BinY = floor((y - ymin)*ny/(ymax - ymin));
%}


 COMPONENT emon_detector = E_monitor(
     nE=100, filename="detector.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (!VIRTUALOUT)    AT (0, 0, dad+0.0215) RELATIVE a6
ROTATED (0,180,0) RELATIVE a6

  COMPONENT psd_window1 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
    xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
    ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window1.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[1] && BinX <= XwinMax[1] && BinY >= YwinMin[1] && BinY <= YwinMax[1] && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6

  COMPONENT emon_window1 = E_monitor(
       nE = 100, filename="window1.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[1] && BinX <= XwinMax[1] && BinY >= YwinMin[1] && BinY <= YwinMax[1]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

  COMPONENT psd_window2 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
    xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
    ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window2.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[2] && BinX <= XwinMax[2] && BinY >= YwinMin[2] && BinY <= YwinMax[2] && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


  COMPONENT emon_window2 = E_monitor(
       nE = 100, filename="window2.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[2] && BinX <= XwinMax[2] && BinY >= YwinMin[2] && BinY <= YwinMax[2]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

  COMPONENT psd_window3 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
    xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
    ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window3.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[3] && BinX <= XwinMax[3] && BinY >= YwinMin[3] && BinY <= YwinMax[3]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


  COMPONENT emon_window3 = E_monitor(
       nE = 100, filename="window3.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[3] && BinX <= XwinMax[3] && BinY >= YwinMin[3] && BinY <= YwinMax[3]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215+0.0008) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

 COMPONENT psd_window4 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
    ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window4.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[4] && BinX <= XwinMax[4] && BinY >= YwinMin[4] && BinY <= YwinMax[4] &&  !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


  COMPONENT emon_window4 = E_monitor(
     nE = 100, filename="window4.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[4] && BinX <= XwinMax[4] && BinY >= YwinMin[4] && BinY <= YwinMax[4]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

  COMPONENT psd_window5 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
    ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window5.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[5] && BinX <= XwinMax[5] && BinY >= YwinMin[5] && BinY <= YwinMax[5] &&  !VIRTUALOUT)    AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


  COMPONENT emon_window5 = E_monitor(
       nE = 100, filename="window5.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[5] && BinX <= XwinMax[5] && BinY >= YwinMin[5] && BinY <= YwinMax[5]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6


   COMPONENT psd_window6 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
     xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
     ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window6.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[6] && BinX <= XwinMax[6] && BinY >= YwinMin[6] && BinY <= YwinMax[6]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


  COMPONENT emon_window6 = E_monitor(
      nE = 100, filename="window6.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[6] && BinX <= XwinMax[6] && BinY >= YwinMin[6] && BinY <= YwinMax[6]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

 COMPONENT psd_window7 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
   ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window7.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[7] && BinX <= XwinMax[7] && BinY >= YwinMin[7] && BinY <= YwinMax[7]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


 COMPONENT emon_window7 = E_monitor(
       nE = 100, filename="window7.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[7] && BinX <= XwinMax[7] && BinY >= YwinMin[7] && BinY <= YwinMax[7]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

 COMPONENT psd_window8 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
   ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window8.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[8] && BinX <= XwinMax[8] && BinY >= YwinMin[8] && BinY <= YwinMax[8]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


 COMPONENT emon_window8 = E_monitor(
       nE = 100, filename="window8.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[8] && BinX <= XwinMax[8] && BinY >= YwinMin[8] && BinY <= YwinMax[8]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6

 COMPONENT psd_window9 = PSD_monitor_psf_eff(eff=0.8, k0=1.553,
   xmin = -det_width/2.0, xmax = det_width/2.0, ymin = -det_height/2.0,
   ymax = det_height/2.0, nx = 128, ny = 128, psf = PSF, filename = "psd_window9.dat",restore_neutron=1)
 WHEN (BinX >= XwinMin[9] && BinX <= XwinMax[9] && BinY >= YwinMin[9] && BinY <= YwinMax[9]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
   ROTATED (0,180,0) RELATIVE a6


 COMPONENT emon_window9 = E_monitor(
       nE = 128, filename="window9.dat", xwidth =det_width , yheight=det_height, Emin=eminf, Emax=emaxf,restore_neutron=1)
 WHEN (BinX >= XwinMin[9] && BinX <= XwinMax[9] && BinY >= YwinMin[9] && BinY <= YwinMax[9]  && !VIRTUALOUT)  AT (0, 0, dad+0.0215) RELATIVE a6
  ROTATED (0,180,0) RELATIVE a6



END