/******************************************************************************* * McXtrace, X-ray tracing package * Copyright, All rights reserved * Risoe National Laboratory, Roskilde, Denmark * Institut Laue Langevin, Grenoble, France * * Component: SANSQMonitor * * %I * Written by: Søren Kynde (kynde@nbi.dk) * Rewritten by: Martin Cramer Pedersen (mcpe@nbi.dk) * Date: October 17, 2012 * Origin: KU-Science * * A circular detector measuring the radial average of intensity as a function * of the momentum transform in the sample. * * %D * A simple component simulating the scattering from a box-shaped, thin solution * of monodisperse, spherical particles. * * %P * RadiusDetector: [m] Radius of the detector (in the xy-plane). * DistanceFromSample: [m] Distance from the sample to this component. * LambdaMin: [AA] Max sensitivity in lambda - used to compute the highest possible value of momentum transfer, q. * NumberOfBins: [] Number of bins in the r (and q). * RFilename: [] File used for storing I(r). * qFilename: [] File used for storing I(q). * Lambda0: [] If lambda is given, the momentum transfers of all rays are computed from this value. Otherwise, instrumental effects are neglected. * restore_neutron: [] If set to 1, the component restores the original neutron. * * %E *******************************************************************************/ DEFINE COMPONENT SANSQMonitor DEFINITION PARAMETERS () SETTING PARAMETERS (NumberOfBins = 100, string RFilename = "RDetector", string qFilename = "QDetector", RadiusDetector, DistanceFromSample, LambdaMin = 1.0, Lambda0 = 0.0, restore_neutron = 0) OUTPUT PARAMETERS () DECLARE %{ // Declarations double TwoThetaMax; double qMax; DArray1d Nofq; DArray1d Iofq; DArray1d IofqSquared; DArray1d NofR; DArray1d IofR; DArray1d IofRSquared; %} INITIALIZE %{ // Declarations int i; Nofq = create_darr1d(NumberOfBins); Iofq = create_darr1d(NumberOfBins); IofqSquared = create_darr1d(NumberOfBins); NofR = create_darr1d(NumberOfBins); IofR = create_darr1d(NumberOfBins); IofRSquared = create_darr1d(NumberOfBins); // Initializations for (i = 0; i < NumberOfBins; ++i) { Nofq[i] = 0.0; Iofq[i] = 0.0; IofqSquared[i] = 0.0; } TwoThetaMax = atan(RadiusDetector / DistanceFromSample); qMax = 4 * PI * sin(TwoThetaMax / 2.0) / LambdaMin; %} TRACE %{ // Declarations int i; double TwoTheta; double Lambda; double R; double RLow; double RHigh; double q; double qLow; double qHigh; double TwoThetaLow; double TwoThetaHigh; double AreaOfSlice; PROP_Z0; // Computation of R R = sqrt(pow(x, 2) + pow(y, 2)); // Computation of q if (Lambda0 <= 0.0) { Lambda = 2.0 * PI / (V2K * sqrt(pow(vx, 2) + pow(vy, 2) + pow(vz, 2))); } else { Lambda = Lambda0; } TwoTheta = atan(R / DistanceFromSample); q = 4.0 * PI * sin(TwoTheta / 2.0) / Lambda; // Put neutron in the correct r-bin if (R < RadiusDetector) { i = floor(NumberOfBins * R / RadiusDetector); RLow = RadiusDetector / NumberOfBins * i; RHigh = RadiusDetector / NumberOfBins * (i + 1); TwoThetaLow = atan(RLow / DistanceFromSample); TwoThetaHigh = atan(RHigh / DistanceFromSample); AreaOfSlice = fabs((cos(2.0 * TwoThetaLow) - cos(2.0 * TwoThetaHigh)) * 2.0 * PI); #pragma acc atomic NofR[i] = NofR[i] + 1; double p_A=p/AreaOfSlice; #pragma acc atomic IofR[i] = IofR[i] + p_A; double p2_A2= p*p/(AreaOfSlice*AreaOfSlice); #pragma acc atomic IofRSquared[i] = IofRSquared[i] + p2_A2; } // Put neutron in the correct q-bin if (q < qMax) { i = floor(NumberOfBins * q / qMax); qLow = qMax / NumberOfBins * i; qHigh = qMax / NumberOfBins * (i + 1); TwoThetaLow = asin(qLow * Lambda / (4.0 * PI)); TwoThetaHigh = asin(qHigh * Lambda / (4.0 * PI)); AreaOfSlice = fabs((cos(2.0 * TwoThetaLow) - cos(2.0 * TwoThetaHigh)) * 2.0 * PI); #pragma acc atomic Nofq[i] = Nofq[i] + 1; double p_A=p/AreaOfSlice; #pragma acc atomic Iofq[i] = Iofq[i] + p_A; double p2_A2= p*p/(AreaOfSlice*AreaOfSlice); #pragma acc atomic IofqSquared[i] = IofqSquared[i] + p2_A2; SCATTER; } // Restore neutron if requested if (restore_neutron) { // RESTORE_NEUTRON(INDEX_CURRENT_COMP, x, y, z, vx, vy, vz, t, sx, sy, sz, p); } %} SAVE %{ // Output I(r) DETECTOR_OUT_1D( "QMonitor - Radially averaged distribution", "Radius [m]", "I(r)", "r", 0.0, RadiusDetector, NumberOfBins, &NofR[0], &IofR[0], &IofRSquared[0], RFilename ); // Output I(q) DETECTOR_OUT_1D( "QMonitor - Distribution in q (Radially averaged)", "q [1 / AA]", "I(q)", "q", 0.0, qMax, NumberOfBins, &Nofq[0], &Iofq[0], &IofqSquared[0], qFilename ); %} FINALLY %{ destroy_darr1d(Nofq); destroy_darr1d(Iofq); destroy_darr1d(IofqSquared); destroy_darr1d(NofR); destroy_darr1d(IofR); destroy_darr1d(IofRSquared); %} MCDISPLAY %{ circle("xy", 0, 0, 0, RadiusDetector); %} END