/*******************************************************************************
* Instrument: Test_Mirror_conic
*
* %I
* Written by: Peter Willendrup
* Date: October 2026
* Origin: DTU
* %INSTRUMENT_SITE: Tests_optics
*
* Compares the conic.h reflection routines with Mirror.comp
*
* %D
* A collimated, polychromatic beam hits a single flat mirror at grazing angle theta.
* The mirror is either Mirror.comp or Mirror_conic.comp, which uses the reflection
* routines from conic.h (as used by Conics_*, FlatEllipse_finite_mirror and NMO).
* Reflected and transmitted beams are recorded separately as function of wavelength.
* For a given theta, pairs of COMP settings should give identical results:
*
* COMP=0: Mirror, transmit=1           (reference for COMP=1)
* COMP=1: Mirror_conic, conic=0        (FlatSurf, reflectNeutronFlat)
* COMP=2: Mirror, transmit=0           (reference for COMP=3)
* COMP=3: Mirror_conic, conic=1        (ConicSurf, reflectNeutronConic)
*
* COMP=0 and COMP=1 should agree exactly, and so should COMP=2 and COMP=3 when
* no figure error is applied (sigma_m=sigma_s=0).
* With sigma_m comparable to the grazing angle (e.g. theta=0.1 sigma_m=1.745e-3),
* COMP=3 uses the waviness model of Hansen et al., arXiv:1612.03763 (see conic.h).
*
* %Example: COMP=0 Detector: L_refl_I=8.98208e-07
* %Example: COMP=1 Detector: L_refl_I=8.98208e-07
* %Example: COMP=0 Detector: L_trans_I=7.01789e-07
* %Example: COMP=1 Detector: L_trans_I=7.01789e-07
* %Example: COMP=2 Detector: L_refl_I=8.9763e-07
* %Example: COMP=3 Detector: L_refl_I=8.9763e-07
* %Example: COMP=3 theta=0.1 sigma_m=1.745e-3 Detector: L_refl_I=1.35886e-06
*
* %P
* COMP:   [1]     Mirror implementation, see description
* theta:  [deg]   Grazing angle of the beam on the mirror
* length: [m]     Mirror length
* m:      [1]     Mirror m-value
* R0:     [1]     Low-angle reflectivity
* Qc:     [AA^-1] Critical scattering vector
* alpha:  [AA]    Slope of reflectivity
* W:      [AA^-1] Width of supermirror cut-off
* Lmin:   [AA]    Minimum source wavelength
* Lmax:   [AA]    Maximum source wavelength
* sigma_m: [rad]  RMS meridional slope (figure) error, COMP=3 only
* sigma_s: [rad]  RMS sagittal slope (figure) error, COMP=3 only
*
* %E
*******************************************************************************/
DEFINE INSTRUMENT Test_Mirror_conic(int COMP=0, theta=1, length=1, m=2, R0=0.99, Qc=0.0217, alpha=6.07, W=0.003, Lmin=1, Lmax=11, sigma_m=0, sigma_s=0)

DECLARE
%{
  double v_refl;
  #pragma acc declare create(v_refl)
%}

USERVARS
%{
  int reflected;
%}

INITIALIZE
%{
  /* Transverse velocity threshold separating reflected from transmitted neutrons */
  v_refl = sin (DEG2RAD * theta);
  #pragma acc update device(v_refl)
%}

TRACE

COMPONENT origin = Progress_bar()
AT (0, 0, 0) RELATIVE ABSOLUTE

COMPONENT source = Source_simple(
    xwidth=0.002, yheight=0.01, dist=5, focus_xw=0.002, focus_yh=0.01,
    lambda0=(Lmin+Lmax)/2, dlambda=(Lmax-Lmin)/2, flux=1)
AT (0, 0, 0) RELATIVE origin

COMPONENT L_in = L_monitor(nL=100, filename="L_in", xwidth=0.01, yheight=0.02,
    Lmin=Lmin, Lmax=Lmax, restore_neutron=1)
AT (0, 0, 4.9) RELATIVE origin

/* Mirror frame: mirror surface is the local y-z plane, from z=0 to z=length,
   reflecting side x<0. The beam axis hits the mirror centre. */
COMPONENT mirror_arm = Arm()
AT (length/2*sin(DEG2RAD*theta), 0, 5-length/2*cos(DEG2RAD*theta)) RELATIVE origin
ROTATED (0, -theta, 0) RELATIVE origin

COMPONENT mirror_transmit = Mirror(xwidth=length, yheight=0.1, center=1, transmit=1,
    R0=R0, Qc=Qc, alpha=alpha, m=m, W=W)
WHEN (COMP==0)
AT (0, 0, length/2) RELATIVE mirror_arm
ROTATED (0, 90, 0) RELATIVE mirror_arm

COMPONENT mirror_conic_flat = Mirror_conic(length=length, yheight=0.1, conic=0,
    R0=R0, Qc=Qc, alpha=alpha, m=m, W=W)
WHEN (COMP==1)
AT (0, 0, 0) RELATIVE mirror_arm

COMPONENT mirror_weight = Mirror(xwidth=length, yheight=0.1, center=1, transmit=0,
    R0=R0, Qc=Qc, alpha=alpha, m=m, W=W)
WHEN (COMP==2)
AT (0, 0, length/2) RELATIVE mirror_arm
ROTATED (0, 90, 0) RELATIVE mirror_arm

COMPONENT mirror_conic_conic = Mirror_conic(length=length, yheight=0.1, conic=1,
    R0=R0, Qc=Qc, alpha=alpha, m=m, W=W, sigma_m=sigma_m, sigma_s=sigma_s)
WHEN (COMP==3)
AT (0, 0, 0) RELATIVE mirror_arm

COMPONENT sorter = Arm()
AT (0, 0, 5+length) RELATIVE origin
EXTEND
%{
  reflected = (fabs (vx) > v_refl * sqrt (vx * vx + vy * vy + vz * vz));
%}

COMPONENT L_refl = L_monitor(nL=100, filename="L_refl", xwidth=1, yheight=0.2,
    Lmin=Lmin, Lmax=Lmax, restore_neutron=1)
WHEN (reflected)
AT (0, 0, 0.01) RELATIVE sorter

COMPONENT L_trans = L_monitor(nL=100, filename="L_trans", xwidth=1, yheight=0.2,
    Lmin=Lmin, Lmax=Lmax, restore_neutron=1)
WHEN (!reflected)
AT (0, 0, 0.01) RELATIVE sorter

END