/*******************************************************************************
*
* McStas, neutron ray-tracing package
*         Copyright 1997-2002, All rights reserved
*         Risoe National Laboratory, Roskilde, Denmark
*         Institut Laue Langevin, Grenoble, France
*
* Component: SiC
*
* %I
* Written by: <a href="mailto:S.RYCROFT@IRI.TUDELFT.NL">S. Rycroft</a>
* Date: Jan 2000
* Origin: IRI.
*
* SiC layer sample
*
* %D
* SiC layer sample
*
* Example: SiC()
*
* %P
* INPUT PARAMETERS:
*
* xlength: [m]  length of SiC plate
* yheight: [m]  height of SiC plate
*
* %E
*******************************************************************************/

DEFINE COMPONENT SiC

SETTING PARAMETERS (xlength, yheight)

/* Neutron parameters: (x,y,z,vx,vy,vz,t,sx,sy,sz,p) */

DECLARE
%{
  double ScatGrad[85];
  double xlength;
  double yheight;
%}

INITIALIZE
%{
  int i;

  for (i = 0; i < 85; ++i)
    ScatGrad[i] = 0;
  ScatGrad[1] = 7.3106e-6;
  ScatGrad[81] = 2.073e-6 - 7.3106e-6;
%}

TRACE
%{
  double dt, q;

  /* Variables added for Rayleigh Appproximation. */
  double realZ, imagZ;
  double csZ, snZ, Qc, R0;
  int i, Z, q_count;

  Qc = 0.010208;
  R0 = 1.0;

  /* First check if neutron has the right direction. */
  if (vz != 0.0 && (dt = -z / vz) >= 0) {
    double old_x = x, old_y = y;

    x += vx * dt;
    y += vy * dt;
    /* Now check if neutron intersects mirror. */
    if (x >= 0 && x <= xlength && y >= 0 && y <= yheight) {
      z = 0;
      t += dt;
      q = fabs (2 * vz * V2Q);
      vz = -vz;
      /* Reflectivity (see component Guide).
             Changed to calculate from real sample. */

      if (q > Qc) {
        realZ = 0;
        imagZ = 0;
        for (Z = 0; Z < 85; Z++) {
          csZ = cos (-1. * q * Z);
          snZ = sin (-1. * q * Z);
          realZ = realZ + ScatGrad[Z] * csZ;
          imagZ = imagZ + ScatGrad[Z] * snZ;
        }
        p *= 16. * PI * PI * ((realZ * realZ) + (imagZ * imagZ)) / (q * q * q * q);
      }
      p *= R0;
      SCATTER;
    } else {
      x = old_x;
      y = old_y;
    }
  }
%}

MCDISPLAY
%{

%}

END