/*******************************************************************************
*
* McStas, version 1.0, released October 26, 1998
*         Maintained by Kristian Nielsen and Kim Lefmann,
*         Risoe National Laboratory, Roskilde, Denmark
*
* Component: Al-Window
*
* %I
* Written by: S. Roth
* Date: ???
* Version: 0.1
* Origin: McStas 1.5
* Modified by: EF, 2004. Renamed into Al_window
*
// * An Al window plane with exponential absorption. No scattering
*
* %D
* OBSOLETE: rather use optics/Al_window.comp
*
* INPUT PARAMETERS:
* win_thick: Window thickness (m)
*
*
* Variables calculated in the component
*
* Al_my_s     : Attenuation factor due to scattering (m^-1)
* Al_my_a     : Attenuation factor due to absorbtion (m^-1)
*
*
*
*
*******************************************************************************/

DEFINE COMPONENT window
DEFINITION PARAMETERS (win_thick)
SETTING PARAMETERS ()
STATE PARAMETERS (x,y,z,vx,vy,vz,t,s1,s2,p)

SHARE
%{
/* ToDo: Should be component local names. */
#ifndef AL_WINDOW
#define avogadro 6.022 /* 10E23 Atoms per mole (mol-1) */
#define Al_sigma_a .231 /* Absorption cross section per atom (barns) at 2200m/s */
#define Al_sigma_i .0082 /* Incoherent scattering cross section per atom (barns) */
#define Al_rho 2.7 /* density (gcm-3) */
#define Al_mmol 27 /* molar mass Al (gmol-1) */
#define Al_my_s (Al_rho / Al_mmol * Al_sigma_i * avogadro * 10) /* inc. XS (barn) */
#define Al_my_a_v (Al_rho / Al_mmol * Al_sigma_a * avogadro * 10 * 2200 )
/* Define Constants for Polynomial Fit of
  sigma_tot(lambda)=A+B1*X+B2*X^2+B3*X^3+B4*X^4+... */
#define Al_pf_A 1.34722
#define Al_pf_B1 .12409
#define Al_pf_B2 .01078
#define Al_pf_B3 -3.25895e-5
#define Al_pf_B4 3.74731e-6
#define AL_WINDOW
#endif
%}

TRACE
%{
  double t0, t3;                /* Entry/exit time for outer cylinder */
  double t1, t2;                /* Entry/exit time for inner cylinder */
  double v;                     /* Neutron velocity */
  double dt0, dt1, dt2, dt;     /* Flight times through sample */
  double l_full;                /* Flight path length for non-scattered neutron */
  double l_i, l_o;              /* Flight path lenght in/out for scattered neutron */
  double my_a;                  /* Velocity-dependent attenuation factor */
  double solid_angle;           /* Solid angle of target as seen from scattering point */
  double dist;
  double nx, ny, nz, xt, yt, zt,xi,yi,zi;
  double Al_s_tot_lambda,Al_my_tot,Al_my_a ; /* total XS (barn), total scattering length (m-1), absorption scat. length */
  double lambda; /* neutrons wavelength */

  //printf("! win1: x=%lf y%lf z%lf ",x,y,z);

  PROP_Z0;

  //printf("! win2: x=%lf y%lf z%lf ",x,y,z);

  dt0=win_thick/vz;
  v=sqrt(vx*vx+vy*vy+vz*vz);
  PROP_DT(dt0);
  dist=v*dt0;

  //printf("! win3: x=%lf y%lf z%lf ",x,y,z);


  lambda=sqrt(81.81/(VS2E*v*v));
  Al_s_tot_lambda=Al_pf_A+Al_pf_B1*lambda+Al_pf_B2*lambda*lambda+Al_pf_B3*lambda*lambda*lambda;
  Al_s_tot_lambda+=Al_pf_B4*lambda*lambda*lambda*lambda;
  Al_my_tot=Al_rho / Al_mmol * Al_s_tot_lambda * avogadro * 10;
  Al_my_a = Al_my_a_v/v;


  p *=exp(-Al_my_a*dist);/* neutron passes window without any interaction */


  /* TODO: scatter in Debye-Scherrer cone */

%}

MCDISPLAY
%{
  /* A bit ugly; hard-coded dimensions. */
  magnify("");
  line(0,0,0,0.2,0,0);
  line(0,0,0,0,0.2,0);
  line(0,0,0,0,0,0.2);
%}

END