1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240
|
/*******************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright (C) 1997-2011, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* %I
* Written by: Peter Willendrup, derived from TOF_lambda_monitor.comp
* Date: May 23, 2012
* Version: $Revision: 4587 $
* Origin: DTU Physics
* Release: McStas 1.12c
*
* Special "Brilliance" monitor of FIXED size 1x1cm.
*
* %D
* If used in the right setting, will output "instantaneous" and "mean" brilliances in units of Neutrons/cm^2/ster/AA/s. Conditions for proper units:
* <ul>
* <li>Use a with a source of area 1x1cm
* <li>The source must illuminate/focus to an area of 1x1cm a 1m distance
* <li>Parametrise the Brilliance_monitor with the frequency of the source
* <li>To not change the source TOF distribution, place the Brilliance monitor close to the source!
* </ul>
*
* with a source of area 1x1cm illuminating/focusing to an area of 1x1cm a 1m distance, this monitor will output "instantaneous" and "mean" brilliances in units of Neutrons/cm^2/ster/AA/s
*
* Here is an example of the use of the component. Note how the mentioned Unit conditions are implemented in instrument code.
*
*COMPONENT Source = ESS_moderator_long(
* l_low = lambdamin, l_high = lambdamax, dist = 1, xw = 0.01, yh = 0.01,
* freq = 14, T=50, tau=287e-6, tau1=0, tau2=20e-6,
* n=20, n2=5, d=0.00286, chi2=0.9, I0=6.9e11, I2=27.6e10,
* branch1=0, branch2=0.5, twopulses=0, size=0.01)
* AT (0, 0, 0) RELATIVE Origin
*
*COMPONENT BRIL = Brilliance_monitor(nlam=196,nt=401,filename="bril.sim",
* t_0=0,t_1=4000,lambda_0=lambdamin,
* lambda_1=lambdamax, Freq=14)
*AT (0,0,0.000001) RELATIVE Source
*
* %P
* INPUT PARAMETERS:
*
* nlam: Number of bins in wavelength (1)
* nt: Number of bins in TOF (1)
* t_0: Minimum time (us)
* t_1: Maximum time (us)
* lambda_0: Minimum wavelength detected (AA)
* lambda_1: Maximum wavelength detected (AA)
* filename: Defines filenames for the detector images. Stored as:<br>Peak_<filename> and Mean_<filename> (string)
* restore_neutron: If set, the monitor does not influence the neutron state (1)
* Freq: Source frequency. Use freq=1 for reactor source (Hz)
* srcarea: Source area (cm^2)
* tofcuts: Flag to generate TOF-distributions as function of wavelength (1)
* toflambda: Flag to generate TOF-lambda distribution output ยด (1)
*
* CALCULATED PARAMETERS:
*
* Div_N: Array of neutron counts
* Div_p: Array of neutron weight counts
* Div_p2: Array of second moments
*
* %E
*******************************************************************************/
DEFINE COMPONENT Brilliance_monitor
DEFINITION PARAMETERS (nlam=101, nt=1001, string filename, t_0=0, t_1=20000, srcarea=1)
SETTING PARAMETERS (lambda_0=0, lambda_1=20, restore_neutron=0, Freq, int tofcuts=0, int toflambda=0)
OUTPUT PARAMETERS (tt_0, tt_1, BRIL_N, BRIL_p, BRIL_p2, BRIL_mean, BRIL_peak, BRIL_meanN, BRIL_peakN, BRIL_meanE, BRIL_peakE)
// STATE PARAMETERS (x,y,z,vx,vy,vz,t,s1,s2,p)
// POLARISATION PARAMETERS (sx,sy,sz)
DECLARE
%{
double BRIL_N[nt][nlam];
double BRIL_p[nt][nlam];
double BRIL_p2[nt][nlam];
double BRIL_mean[nlam];
double BRIL_meanN[nlam];
double BRIL_meanE[nlam];
double BRIL_peak[nlam];
double BRIL_peakN[nlam];
double BRIL_peakE[nlam];
double BRIL_shape[nt];
double BRIL_shapeN[nt];
double BRIL_shapeE[nt];
double tt_0, tt_1;
double xmin=-0.005, xmax=0.005, ymin=-0.005, ymax=0.005;
double ster=1e-4;
double prsec=1e-6;
double dlam;
double dt;
%}
INITIALIZE
%{
int i,j;
tt_0 = t_0*prsec;
tt_1 = t_1*prsec;
dt=(t_1-t_0)*prsec/nt;
dlam=(lambda_1-lambda_0)/(nlam-1);
for (i=0; i<nlam; i++)
{
BRIL_mean[i] = 0;
BRIL_peak[i] = 0;
BRIL_meanN[i] = 0;
BRIL_peakN[i] = 0;
BRIL_meanE[i] = 0;
BRIL_peakE[i] = 0;
for (j=0; j<nt; j++)
{
BRIL_N[j][i] = 0;
BRIL_p[j][i] = 0;
BRIL_p2[j][i] = 0;
if (i==0) {
BRIL_shape[j] = 0;
BRIL_shapeN[j] = 0;
BRIL_shapeE[j] = 0;
}
}
}
%}
TRACE
%{
int i,j;
double div;
double lambda;
double Pnorm;
PROP_Z0;
lambda = (2*PI/V2K)/sqrt(vx*vx + vy*vy + vz*vz);
if (x>xmin && x<xmax && y>ymin && y<ymax &&
lambda > lambda_0 && lambda < lambda_1)
{
if (t < tt_1 && t > tt_0)
{
i = floor((lambda - lambda_0)*nlam/(lambda_1 - lambda_0));
j = floor((t-tt_0)*nt/(tt_1-tt_0));
Pnorm=p/dlam/ster/srcarea;
BRIL_meanN[i]++;
BRIL_mean[i] += Pnorm;
BRIL_meanE[i] += Pnorm*Pnorm;
Pnorm=Pnorm/Freq/dt;
BRIL_N[j][i]++;
BRIL_p[j][i] += Pnorm;
BRIL_p2[j][i] += Pnorm*Pnorm;
}
}
if (restore_neutron) {
RESTORE_NEUTRON(INDEX_CURRENT_COMP, x, y, z, vx, vy, vz, t, sx, sy, sz, p);
}
%}
SAVE
%{
/* First, dump the 2D monitor */
/* For each Wavelength channel, find peak brilliance */
int i,j,jmax;
double Pnorm;
char ff[256];
char tt[256];
for (i=0; i<nlam; i++)
{
Pnorm = -1;
jmax = -1;
for (j=0; j<nt; j++)
{
if (BRIL_p[j][i]>=Pnorm)
{
Pnorm = BRIL_p[j][i];
jmax=j;
}
BRIL_shape[j] = BRIL_p[j][i];
BRIL_shapeN[j] = BRIL_N[j][i];
BRIL_shapeE[j] = BRIL_p2[j][i];
}
if (tofcuts == 1) {
sprintf(ff, "Shape_%s_%g",filename,lambda_0+i*dlam);
sprintf(tt, "Peak shape at %g AA",lambda_0+i*dlam);
DETECTOR_OUT_1D(
tt,
"TOF [us]",
"Peak Brilliance",
"Shape", t_0, t_1, nt,
&BRIL_shapeN[0],&BRIL_shape[0],&BRIL_shapeE[0],
ff);
}
BRIL_peakN[i] = BRIL_N[jmax][i];
BRIL_peak[i] = BRIL_p[jmax][i];
BRIL_peakE[i] = BRIL_p2[jmax][i];
}
sprintf(ff, "Mean_%s",filename);
DETECTOR_OUT_1D(
"Mean brilliance",
"Wavelength [AA]",
"Mean Brilliance",
"Mean", lambda_0, lambda_1, nlam,
&BRIL_meanN[0],&BRIL_mean[0],&BRIL_meanE[0],
ff);
sprintf(ff, "Peak_%s",filename);
DETECTOR_OUT_1D(
"Peak brilliance",
"Wavelength [AA]",
"Peak Brilliance",
"Peak", lambda_0, lambda_1, nlam,
&BRIL_peakN[0],&BRIL_peak[0],&BRIL_peakE[0],
ff);
/* MPI related NOTE: Order is important here! The 2D-data used to generate wavelength-slices and calculate
the peak brilliance should be done LAST, otherwise we will get a factor of MPI_node_count too much as
scatter/gather has been performed on the arrays... */
if (toflambda == 1) {
sprintf(ff, "TOFL_%s",filename);
DETECTOR_OUT_2D(
"TOF-wavelength brilliance",
"Time-of-flight [\\gms]", "Wavelength [AA]",
t_0, t_1, lambda_0, lambda_1,
nt, nlam,
&BRIL_N[0][0],&BRIL_p[0][0],&BRIL_p2[0][0],
filename);
}
%}
MCDISPLAY
%{
multiline(5, (double)xmin, (double)ymin, 0.0,
(double)xmax, (double)ymin, 0.0,
(double)xmax, (double)ymax, 0.0,
(double)xmin, (double)ymax, 0.0,
(double)xmin, (double)ymin, 0.0);
%}
END
|